在matlab使用矩阵的方法计算DFT

在matlab使用矩阵的方法计算DFT

dft矩阵算法原理学习笔记
2021年10月22日

编写一个函数，在matlab使用矩阵的方法计算DFT，输入序列x(n)和N，输出序列x(n)的N点DFT变换结果X(k)

函数的定义

%这是matlab代码，若要使用此函数，请复制此代码保存为MDFT.m文件到你的工作路径
function [Xk]=MDFT(xn,N)
%此函数使用DFT矩阵方法计算序列x(n)的N点DFT
M=length(xn);%记录序列x(n)初始长度
    if M<N   %如果序列长度小于N，补零到N，否则，截取前N项 
        xn=[xn,zeros(1,N-M)];
    else
        xn=xn(:,1:N);  %截取（全部行）前N列 
    end
n=0:N-1; k=0:N-1;      %匹配矩阵维度n,k
W_N=exp(-1i*2*pi/N);   %生成N点DFT对应的旋转因子W(常数)
nk=n'*k;               %生成旋转因子W的幂次系数矩阵(N*N)
W=W_N.^nk;             %生成N点DFT矩阵(N*N),数幂运算
Xk=W*xn';              %[X(k)]=[W]*[x(n)]T,矩阵乘法
Xk=Xk';                %以行向量的方式输出[X(k)]
end
%%%%%October 22, 2021 by Linwenzhe

函数的调用方法：

Xk=MDFT(xn,N);

此函数使用DFT矩阵方法计算序列x(n)的N点DFT。
xn和N是入口参数，其中，xn是一个行向量。
函数的输出Xk也是一个行向量。

调用示例：
绘制矩形序列R_4(n)={1,1,1,1}的16点DFT和32点DFT的频谱图

xn=[1 1 1 1];        %输入时域序列向量xn=R4(n)
Xk16=MDFT(xn,16);    %计算xn的16点DFT
Xk32=MDFT(xn,32);    %计算xn的32点DFT
%以下为绘图部分
k=0:15;wk=2*k/16;       %产生16点DFT对应的采样点频率(关于π归一化值)
subplot(2,2,1);stem(wk,abs(Xk16),'.'); %绘制16点DFT的幅频特性图
title('(a)16点DFT的幅频特性图');xlabel('ω/π');ylabel('幅度');
subplot(2,2,3);stem(wk,angle(Xk16),'.'); %绘制16点DFT的相频特性图
line([0,2],[0,0]);title('(b)16点DFT的相频特性图')
xlabel('ω/π');ylabel('相位');axis([0,2,-3.5,3.5]);
k=0:31;wk=2*k/32;       %产生32点DFT对应的采样点频率(关于π归一化值)
subplot(2,2,2);stem(wk,abs(Xk32),'.'); %绘制32点DFT的幅频特性图
title('(c)32点DFT的幅频特性图');xlabel('ω/π');ylabel('幅度');
subplot(2,2,4);stem(wk,angle(Xk32),'.'); %绘制32点DFT的相频特性图
line([0,2],[0,0]);title('(d)32点DFT的相频特性图');
xlabel('ω/π');ylabel('相位');axis([0,2,-3.5,3.5]);

输出图形：

谢谢阅读，如有错误，请多多指教！