cyclic_function.m

对数字调制信号计算循环谱相关

function G=cyclic_function(x,N,M,K,alfa,fs,fsc)   %谱相关频率平滑
%N=2048;
%M=32;
X=fft(x,N);
for k=1:5*N
   if k>=1&&k<=N
      XX(k)=X(k);
   else if mod(k,N)~=0
           XX(k)=X(k-floor(k/N)*N);
       else XX(k)=X(N);
       end
   end
end
    
%K=round(1.5*N*fc/fsc);          %频率轴
%alfa=round(1.5*N*fc/fsc);       %循环频率轴  
% K=round(2*N*fc/fsc);   nj
% alfa=round(4*N*fc/fsc);
i=1;
for k=-K:16:K
    j=1;
    for m=-alfa:16:alfa
        sum=0;
        for v=-M:M
            a=k+v+m;b=k+v-m;
            if (k+v+m>=0)&&(k+v-m>=0)
                sum=sum+XX(k+v+m+1)*conj(XX(k+v-m+1));
            else while(a<0||b<0)
                a=a+N;
                b=b+N;
                end
                sum=sum+XX(a+1)*conj(XX(b+1));
            end
        end
       G(i,j)=sum/fs/2/M/N;
        %g_en(i,j)=sum/fs/2/M/N;
        j=j+1;
    end
    i=i+1;
end

p=(-K:16:K)*fsc/N;
q=(-alfa:16:alfa)*2*fsc/N;
figure(3)
%mesh(abs(g_en));
%subplot(2,2,1)
mesh(p,q,abs(G));

%Sm=max(max(abs(G)));
i=K/16+1;j=alfa/16+1;
%subplot(2,2,2)
figure(4)
plot(q,abs(G(i,:)));
xlabel('alpha')
figure(5)
plot(p,abs(G(:,j)));
xlabel('f')