%%%%%%% -------------------------------------------------------------------
%%%%% -----相關法用于頻率估計
%%%%% -----復相關、頻率校正算法
% 
function fe = correlate2(x,fs,N)
% 
% x:輸入數據
% fs:數據采樣率
% N:FFT的點數
% fe:估計出的頻率

delta_f = fs/N;
n = 0:N-1;   %%% ?

Z = fft(x,N);
P = abs(Z);
[Pm mi] = max(P);
fo = fs*(mi-1)/N

x = hilbert(x);   %%% x的解析信號

y = exp(-j*2*pi*fo*n./fs);   %%% 創建標準信號

z = x.*y;   %%% 相關

m = floor(N/2);   %%%???
g = sum(conj(z(1:N-m)).*z(m+1:N));

phase_g = angle(g);
delta_f = phase_g/(2*pi*m/fs)

fe = fo+delta_f;