PAST+MUSIC算法进行离散时间信号频率估计程序仿真

clear;
A1=sqrt（10^1.5）;
A2=sqrt（10^1.5）;
A3=sqrt（10）;
DIM=6;
N=15;
beta=0.99;
FLAG = 1;
RNG=10;
phi1=2*pi*rand;
phi2=2*pi*rand;
phi3=2*pi*rand;
for t=1:1:1000
    switch FLAG
        case 1
            f1=0.2;
            f2=0.3;
            f3=0.1;
            realf（1t）=0.2;
            realf（2t）=0.3;
            realf（3t）=0.1;
        case 2
            f1=0.2+（5e-5）*t;
            f2=0.3-（5e-5）*t;
            f3=0.1;
            
            a=0.2;
            b=5e-5;
            prd=（sqrt（（2*b*t+a）^2+4*b）-（2*b*t+a））/（2*b）;
            realf（1t）=1/prd;
            
            a=0.3;
            b=-（5e-5）;
            prd=（sqrt（（2*b*t+a）^2+4*b）-（2*b*t+a））/（2*b）;
            realf（2t）=1/prd;
            
            realf（3t）=0.1;
        case 3
            if t<500
                f1=0.2;
                f2=0.3;
                f3=0.1;
                realf（1t）=0.2;
                realf（2t）=0.3;
                realf（3t）=0.1;
            else
                f1=0.2+（5e-5）*（t-500）;
                f2=0.4;
                f3=0.1;
                
                a=0.175;
                b=5e-5;
                prd=（sqrt（（2*b*t+a）^2+4*b）-（2*b*t+a））/（2*b）;
                realf（1t）=1/prd;
                
                realf（2t）=0.4;
                realf（3t）=0.1;
            end
    end
    s（t）=sqrt（2）*（A1*sin（2*pi*f1*t+phi1）+A2*sin（2*pi*f2*t+phi2）+A3*sin（2*pi*f3*t+phi3））+randn;
end


W=rand（NDIM）;
I15=eye（15）;
c1=I15（:1）;
c2=I15（:2）;
c3=I15（:3）;
c4=I15（:4）;
c5=I15（:5）;
c6=I15（:6）;
%W=[c1 c2 c3 c4 c