前不久團隊貼出一些常用信號的循環(huán)譜仿真圖（見【GreenSim原創(chuàng)】常規(guī)信號的循環(huán)譜），由于篇幅所限這個帖子把剩余的幾幅圖和程序貼出來，供大家交流。

            64QAM信號循環(huán)譜圖

                        16QAM信號循環(huán)譜圖

                         256QAM信號循環(huán)譜圖

計算信號循環(huán)譜的函數(shù)：

這個函數(shù)用來計算數(shù)字調(diào)制信號的循環(huán)譜圖，結(jié)果存放在sxf中，surf和shading interp這兩個函數(shù)配合使用就可以得到信號的循環(huán)譜圖，M一般取64，N取2048。需要注意的是，數(shù)字信號要保留一定的載頻，不能是基帶信號。

function sxf=MyCS_density(s,M,N)

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

%%  譜相關(guān)的計算函數(shù)

%%  參數(shù)說明：M是平滑窗,必須為2的冪

%%            N是信號序列長度,必須為2的冪

%%            s是時域序列

%%            sxf是返回參數(shù)

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

spe_s = fft(s,N);

k = N/2-M;

sxf = zeros(2*k+1,4*k+1);                          %因為菱形區(qū)之外的數(shù)為零

for ff = -k:k;

    for alpha= -k+ff:1:2* ff                        %（菱形區(qū)）

         for v = -M:M

             temp1 = ff+v+alpha/2;

             temp2 = ff+v-alpha/2;

              if(temp1<1)             

               temp1 = temp1+N;

              end 

              if(temp2<1)

               temp2 =temp2+N;

              end

              cospe(v+M+1) = spe_s(temp1)* spe_s(temp2);

          end            

          j=ff+k+1;

          u=alpha+k+1;

          sxf(j,u) =(abs( cospe))/(2*M*N);   

      end

  end

計算信號循環(huán)譜alpha截面的函數(shù)：

這個函數(shù)用來計算數(shù)字調(diào)制信號的alpha截面圖，結(jié)果存放在sxf_amp中，參數(shù)的設置和上一個函數(shù)一樣，需要注意的是ff參數(shù)必須是數(shù)字頻率。這個函數(shù)非常簡潔也非常有用，尤其是在設計基于循環(huán)譜的參數(shù)估計和調(diào)制識別算法中，一般用這個函數(shù)。

function sxf_amp = MyCSalpha_imag_density(s,M,N,ff)

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

%%  譜相關(guān)alpha截面圖的計算函數(shù)

%%  參數(shù)說明：M是平滑窗

%%            N是信號序列長度

%%            s是信號序列

%%            sxf_amp是幅度譜

%%            ff是頻率

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

spe_s = fft(s,N);

spe_s=spe_s./2;

k = N/2-M;

sxf_amp = zeros(1,2*k+1);                         

for alpha=0:1:k-abs(ff))  

    for v =-M:M

        temp1 = ff+v+alpha/2;

        temp2 = ff+v-alpha/2;

        if(temp1<1)             

            temp1 = abs(temp1);

            A(temp1+1)=conj(spe_s(temp1+1));

        else

            A(temp1+1)=spe_s(temp1+1);

        end

        if(temp2<1)

            temp2 = abs(temp2);

            B(temp2+1)=conj(spe_s(temp2+1));

        else

            B(temp2+1)=spe_s(temp2+1);

        end

        cospe(v+M+1) = A(temp1+1)*B(temp2+1);

    end            

    u=alpha+1;

    sxf_amp(u) =(abs(cospe))/(2*M*N);   

end

計算信號循環(huán)譜f截面的函數(shù)：

這個函數(shù)和上一個函數(shù)用法差不多，主要在設計參數(shù)估計和調(diào)制識別算法中使用。

function sxf_amp = MyCSf_imag_density(s,M,N,alpha)

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

%%  譜相關(guān)f截面圖的計算函數(shù)

%%  參數(shù)說明：M是平滑窗

%%            N是信號序列長度

%%            s是復數(shù)序列

%%            sxf_amp是幅度譜

%%            alpha是數(shù)字頻率

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

spe_s = fft(s,N);

spe_s=spe_s./2;

k = N/2-M;

sxf_amp = zeros(k+1,1);                         

for ff = 0:k-abs(alpha/2);

    for v =-M:M

        temp1 = ff+v+ alpha/2;

        temp2 = ff+v- alpha/2;

        if(temp1<1)             

            temp1 = abs(temp1);

            A(temp1+1)=conj(spe_s(temp1+1));

        else

            A(temp1+1)=spe_s(temp1+1);

        end

        if(temp2<1)

            temp2 = abs(temp2);

            B(temp2+1)=conj(spe_s(temp2+1));

        else

            B(temp2+1)=spe_s(temp2+1);

        end

        cospe(v+M+1) = A(temp1+1)* B(temp2+1);

    end            

    j=ff+1;

    sxf_amp(j) =(abs (cospe))/(2*M*N);   

end