function[f,s]=sy(data,fdata,win)

% Modified periodogram Sy estimator - provides frequency axis
% Returns freq and S axix from 0 - fdata/2
%
% win=0: Rect,  win=1: Hamming,  win=2: Blackman,  
% win=3: BHH, win=4: Hanning, win=5: Hann, win=6: Bartlett.   
%
% syntax: [f,s]=sy(inndata, freq_of_inndata, window) ;
%
% Window loss: BHH-4.2dB, Blackman-2.4dB, Hamming-1.6dB, Rect-0dB


Z=size(data);
if Z(2)~=1  data=data';end;
Z=size(data);
N=Z(1);

%f=0:fdata/N:(fdata-fdata/N);
f=fdata*(0:floor(N/2))/N;

if nargin==4
   f=0:fdata/N:bw;  
end;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

if win==6  
U   = mean(bartlett(N).^2);
s = ((abs( fft(data.*(bartlett(N))') )).^2) / (N*U);
end;

if win==5  
U   = mean(hann(N).^2);
s = ((abs( fft(data.*(hann(N))') )).^2) / (N*U);
end;


if win==4  
U   = mean(hanning(N).^2);
s = ((abs( fft(data.*hanning(N)) )).^2) / (N*U);
end;

if win==3  
U   = mean(blackmanharris(N).^2);
s = ((abs( fft(data.*blackmanharris(N)) )).^2) / (N*U);
end;

if win==2  
U   = mean(blackman(N).^2);
s = ((abs( fft(data.*blackman(N)) )).^2) / (N*U);
end;

if win==1  
U   = mean(hamming(N).^2);
s = ((abs( fft(data.*hamming(N)) )).^2) / (N*U);
end;

if win==0  
s = (abs( fft(data) ).^2) / N;
end;


%s=10.*log10(s);
s=10.*log10(s(1:length(f)));
%s=s(1:length(f));
