循环谱估计_时间平滑FFT积累方法

function [dalfaNNpf0alfa0SCD] = SCD_FAM（IfsdfM）
% 输入
% I：信号实部
% fs:信号采样率
% df：频率分辨率
% M ：不确定量

% 输出
% dalfa：循环频率分辨率
% N ：观测时间内离散取样总数
% Np：第一个FFT的点数N‘是靠近且大于fs/df的2的幂次方
% f0：频率轴
% alfa0：循环频率轴
% SCD：循环谱
dalfa = df/M;                       % 循环频率分辨率
x = I;                              % 实信号序列
Np = pow2（nextpow2（fs/df））;         % 第一个FFT的点数N‘是靠近且大于fs/df的2的幂次方？                                   
L = Np/4;                           % 滑动因子，为了避免估计中的混叠现象和循环泄漏应该小于等于四分之一的时间窗口Np。                                   
P = pow2（nextpow2（fs/dalfa/L））;     % 用一个P点FFT的傅里叶变换取代计算复解调序列乘积的平均
N = P*L;                            % 观测时间内离散取样总数

if length（x）    x（N） = 0;
elseif length （x）>N
    x = x（1:N）;                    %观测时间内离散取样的总数
end

NN = （P-1）*L+Np;
xx = x;
xx（NN） = 0;
xx = xx（:）;
X = zeros（NpP）;
for k = 0:P-1
    X（:k+1） = xx（k*L+1:k*L+Np）;    %将数据按时域滑动L点保存大小为窗口的大小Np
end

a = hamming （Np）;
XW = diag（a）*X; %数据削度

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 复解调的计算——第一次Np点FFT的计算

XF1 = fft（XW）;                        %滑动Np点FFT
XF1 = fftshift（XF1）;                  %将fft的主值区间变换至[-pipi]
XF1 = [XF1（:P/2+1:P） XF1（:1:P/2）];  %第9列-16列和第1列-8列交换——f

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

E = zeros（NpP）;
for k = -Np/2:Np/2-1                             %分别左右搬移循环频率的一半
    for m = 0:P-1
        E（k+Np/2+1 m+1） = exp（-i*2*pi*k*m*L/Np）;%频移量
    end
end
XD = XF1.*E;                                     %对Np点FFT的输出进行频移
XD = conj（XD‘）;                                  %计算每一小段谱的共轭

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