function [fi,ai] = interpolate_PSD(f,a,s,df)
 %   
	fi(1)=f(1);
	ai(1)=a(1);
%
    m=length(f);
%
    MAX = 5000;
%   
    for  i=2:MAX 
 %       
		fi(i)=fi(i-1)+df; 
%
        for j=1:(m-1)
%
			if( ( fi(i) >= f(j) ) & ( fi(i) <= f(j+1) )  )
				ai(i)=a(j)*( ( fi(i) / f(j) )^ s(j) );
				break;
            end
        end
        if( fi(i) >= f(m) )
            break;
        end                
    end
    nn=length(fi);
    ai(nn)=a(m);