fudaopu.m

来自「一个用matlab进行数字信号处理的例程」· M 代码 · 共 24 行

% 最小相位序列复倒谱分析
N=100;
x=[1,1,2,3,4,5,6,7];n=1:8;
k=1:N;w=exp(j*2*pi/N)*k;
Xk=x*(exp(-j*2*pi/N).^(n'*(k-1)));
FXk=log(Xk);
FXrk=log(abs(Xk));
COR=zeros(1,length(n));
for i=2:N
    if angle(Xk(i))-angle(Xk(i-1))>pi
       COR(i)=COR(i-1)-2*pi;
    elseif angle(Xk(i-1))-angle(Xk(i))>pi
       COR(i)=COR(i-1)+2*pi; 
    else 
       COR(i)=COR(i-1);
    end
end   
r=(1/pi)*Xk(N/2+1);sign=Xk(1)/abs(Xk(1));
argX=angle(Xk)+COR-(2*pi/N)*r.*ones(1,N);
fxn=sign*(1/N)*(FXrk+j*argX)*(exp(j*2*pi/N).^(k'*n));
subplot(2,2,1);stem(n,x);title('x(n)');
subplot(2,2,2);stem(n,fxn);title('fx(n)');
subplot(2,2,3);stem(k,FXrk);title('复对数实部');
subplot(2,2,4);stem(k,argX);title('复对数虚部');