信号处理——信号频域变换
来源:互联网 发布:小米床垫怎么样 知乎 编辑:程序博客网 时间:2024/06/05 07:48
作者:桂。
时间:2017-03-07 18:40:05
链接:
声明:欢迎转载,不过记得注明出处哦~
前言
本文主要介绍信号频域变换的常用方式,主要介绍基本的三种方法,实现时域数字信号的傅里叶变换:
1)基本DFT变换;
2)FFT变换
3)逆序级联FFT
三种方法实现方式略有差别,但本质完全相同。为了方便大家理解,本文给出基本MATLAB实现。内容的理论部分,多有借鉴他人,相应链接在最后一并给出。
〇、写在前面
关于信号时域连续、离散以及各自频域的对应关系,可以参考之前的一篇博文。为了便于理解,此处给一个DTFT—>DFT的解释。
MATLAB代码:
%DTFTclc;clear all;close allN=8; %原离散信号有8点n=[0:1:N-1] ; %原信号是1行8列的矩阵xn=0.5.^n; %构建原始信号,为指数信号w0=[-2*N:.1:2*N]*2*pi/N; %方便观察,此处取4个周期-4*pi~4*pi X0=xn*exp(-j*(n'*w0)); %求dtft变换,采用原始定义的方法,对复指数分量求和而得subplot(221)plot(n,xn,'r');hold on;xlabel('时间')stem(n,xn,'k');grid on;title('时域采样后的信号(指数信号)');ylabel('幅值')subplot(222);plot(w0,abs(X0),'k');ylabel('幅值');xlabel('角频率');%逼近连续,其实这么处理不正确,仅仅是方便理解title('DTFT变换');xlim([-4*pi,4*pi]);grid on;%DFTw=[-2*N:1:2*N]*2*pi/8; X=xn*exp(-j*(n'*w)); %对DTFT进行频域采样subplot(223)stem(n,xn,'k');ylabel('幅值');xlabel('时间');grid on;title('频域采样后的信号(指数信号)');subplot(224);plot(w0,abs(X0),'r');hold on;stem(w,abs(X),'k');ylabel('幅值');xlabel('角频率')title('DFT变换');xlim([-4*pi,4*pi]);grid on;
对应结果图:
一、基本DFT变换
DFT变换有:
$F(k) = \sum^{N-1}_{n=0} f(n) e^{\frac{-j2\pi kn}{N}}$
写成矩阵的形式:
对应的转换为VanderMonde矩阵,对应的MATLAB代码:
clc;clear all;close all;fs = 200;f0 = 30;t = 0:1/fs:2;x = sin(2*pi*f0*t');N = length(x);%DFTX_DFT = exp(-j*2*pi/N).^([0:N-1]'*[0:N-1])*x;%绘图xf = linspace(0,fs,N);plot(xf,abs(X_DFT),'k');grid on;xlabel('频率(Hz)');ylabel('幅度');title('DFT结果');
也可以将DFT的代码换为Vander矩阵形式,只需修改一句:
% X_DFT1 = exp(-j*2*pi/N).^([0:N-1]'*[0:N-1])*x;%修改为下句X_DFT = fliplr(vander(exp(-j*2*pi/N).^([0:N-1]))).'*x;
结果图:
二、FFT变换
DFT到FFT常用的有奇偶抽取/前后抽取,两种实现方式,原理此处不再堆砌,可以参考维基百科:
- FFT:https://en.wikipedia.org/wiki/Fast_Fourier_transform
- Cooley-Tukey FFT:https://en.wikipedia.org/wiki/Cooley%E2%80%93Tukey_FFT_algorithm
给出对应的MATLAB代码:
clc;clear all;close all;fs = 200;f0 = 30;t = 0:1/fs:2;x = sin(2*pi*f0*t');N = length(x);%FFTX_FFT = fft(x);%绘图xf = linspace(0,fs,N);plot(xf,abs(X_FFT),'k');grid on;xlabel('频率(Hz)');ylabel('幅度');title('FFT结果');
结果图:
如果严格按照FFT定义,则应该将序列补为$2^N$长度(事实上,内置fft自动完成此操作),代码修改为:
clc;clear all;close all;fs = 200;f0 = 30;t = 0:1/fs:2;x = sin(2*pi*f0*t');N = length(x);%FFTLen = 2^nextpow2(length(x));X_FFT = fft(x,Len);%绘图xf = linspace(0,fs,Len);plot(xf,abs(X_FFT),'k');grid on;xlabel('频率(Hz)');ylabel('幅度');title('FFT结果');
对应结果图:
三、逆序级联FFT
为了方便,假设信号$f(n)$长度为2的整次幂。对于数据量过大的情况,可以通过串/并的不断转化,对数据进行分流,从而利用多个芯片实现一个大芯片的任务(此处为个人理解,仅供参考)。给出对应的流图:
以及具体的推导公式:
步骤可以简化为:分段FFT,即公式[.]部分—>因子补偿,即公式{.}部分—>分段FFT,即公式最外层部分。
给出具体的实现代码:
clc;clear all;close all;fs = 200;f0 = 30;t = 0:1/fs:2;x = sin(2*pi*f0*t');Na = 2^nextpow2(length(x));sig = [x;zeros(Na-length(x),1)];N=16;M=Na/N;sr=zeros(M,N);for m=1:M sr(m,:)=sig((m-1)*N+1:m*N,1);endfor n=1:N sr(:,n)=fft(sr(:,n));endP=zeros(M,N);for p=0:M-1 for n=0:N-1 P(p+1,n+1)=exp(-j*2*pi*n*p/M/N); endendsr=sr.*P;for m=1:M sr(m,:)=fft(sr(m,:));end%save resultX_CaFFT = zeros(length(sig),1);for n=1:N X_CaFFT((n-1)*M+1:n*M,1)=sr(:,n);end%绘图xf = linspace(0,fs,Na);plot(xf,abs(X_CaFFT),'k');grid on;hold on;xlabel('频率(Hz)');ylabel('幅度');title('CaFFT结果');
对应结果图:
四、三种算法对比分析
三种结果效果完全等价,分析其各自性能。
对于快速运算中,一个蝶形运算,对应一个乘法,两次加法。
DFTFFTCaFFT加法器$N^2$$N log_2N$$N log_2N+N$乘法器$N(N-1)$$2Nlog_2N$$2Nlog_2N$可见FFT所需资源最少,但如果数据量过大导致直接FFT,硬件无法满足条件,则CaFFT是最优的选择。
参考:
郑君里:《信号与系统》第二版;
张大炜:一种新的级联FFT算法.
- 信号处理——信号频域变换
- 信号处理——信号频域变换
- 信号处理——信号频域变换
- 信号处理——Hilbert变换及谱分析
- Linux — 信号 信号处理和信号处理函数详解
- 频域信号处理
- 信号处理—卷积
- 信号处理——sigaction
- 信号处理——傅里叶变换
- 信号处理——滤波器
- Houdini学习 —— 使用音频驱动几何体变换之信号处理
- 小波变换和motion信号处理
- 小波变换和motion信号处理
- 小波变换和motion信号处理
- 小波变换和motion信号处理
- Linux&C ——信号以及信号处理
- (三十二)信号——信号处理函数
- (三十五)信号——SIGCHLD信号处理
- 蓝桥杯——算法训练 求指数(Vip试题)
- VS开发Application生硬古怪的问题及原因汇总
- JS隐藏div占用与不占用页面位置方法
- go reflect
- Apache Commons Email邮件发送
- 信号处理——信号频域变换
- HTML5中 audio标签的样式修改
- Postman用法简介-Http请求模拟工具
- Python3.6安装以及numpy库、matplotlib库的安装方法(Win7)
- 快速傅立叶变换(FFT)的一种推导
- 学生信息管理系统(jdbc操作MySQL)
- 文章标题
- linux下svn使用小结 创建 添加仓库 版本管理
- iBET Slot Games Rebate 1% Unlimited Bonus(iBET, iBET Promotion, Online Casino Malaysia, Online Casin