|
离散傅里叶变换的简单Python与C++实现
- # -*- coding: utf-8 -*-
- __author__ = u'东方耀 微信:dfy_88888'
- __date__ = '2021/9/1 上午9:57'
- __product__ = 'PyCharm'
- __filename__ = '21_信号处理demo'
- from scipy.fftpack import fft, ifft, fftshift
- from scipy.signal import findfreqs, find_peaks
- import numpy as np
- import scipy
- from scipy.signal import get_window
- import matplotlib.pyplot as plt
- import matplotlib.ticker as mticker
- from matplotlib.font_manager import FontProperties
- font_fname = "/usr/share/fonts/wps-office/simfang.ttf"
- font = FontProperties(fname=font_fname)
- # plt.rcParams['font.sans-serif']=['simfang']#设置作图中文显示
- # 离散傅里叶变换的简单Python与C++实现
- a = np.array([1, 2, 4, 7]) # a是时域的 数字信号 采样后没有时间的概念了
- print("直接调包的结果:", fft(a), ifft(fft(a)))
- a_fft = np.zeros_like(a, dtype=np.complex)
- # k是频域的index n是时域的索引index
- for k in range(len(a_fft)):
- # a_fft[k]
- for n in range(len(a)):
- a_fft[k] += a[n] * np.exp(-1j * 2 * np.pi * k * n / len(a))
- print("a做fft的结果:", a_fft)
- a_fft_ifft = np.zeros_like(a_fft, dtype=np.complex)
- for n in range(len(a_fft_ifft)):
- # a_fft_ifft[n]
- for k in range(len(a_fft)):
- a_fft_ifft[n] += a_fft[k] * np.exp(1j * 2 * np.pi * k * n / len(a_fft))
- a_fft_ifft[n] /= len(a_fft)
- print("a做fft再ifft的结果:", a_fft_ifft)
- # 算法复杂度为 O(n**2)
复制代码
- #include <fftw3.h>
- #include <complex>
- #include <iostream>
- #include <cassert>
- using namespace std;
- typedef std::complex<float> gr_complex; //这个float精度与fftwf一致 fftwf_complex
- void print_complex_array(gr_complex *p, int num){
- for(int i =0; i< num;i++){
- cout << "value of complex: real=" << (*(p + i)).real() << ",imag=" << (*(p + i)).imag() << endl;
- }
-
- }
- //FFTW编程案例:fft与ifft
- // 注意:target_link_libraries(main -lfftw3f)
- //fftw_plan一次创建可以重复使用
- //傅里叶变换后需要回收内存 fftw_free() 不能用 free() delete
- //复数的多维离散傅里叶变换 二维DFT变换 三维DFT变换
- int main(){
- auto N = 4;
- gr_complex* a; // 数组指针
- // 这个malloc还得跟数据类型对应 float
- a = (gr_complex*)fftwf_malloc(sizeof(gr_complex) * N);
- // 输入的数据 赋值 a = np.array([1, 2, 4, 7])
- a[0] = gr_complex(1,0);
- a[1] = gr_complex(2,0);
- a[2] = gr_complex(4,0);
- a[3] = gr_complex(7,0);
- //print_complex_array(a, N);
- // 变量名与Python一致
- gr_complex* a_fft;
- gr_complex* a_fft_ifft;
- a_fft = (gr_complex*)fftwf_malloc(sizeof(gr_complex) * N);
- a_fft_ifft = (gr_complex*)fftwf_malloc(sizeof(gr_complex) * N);
- fftwf_plan d_fft_plan, d_ifft_plan;
- cout << "C++实现离散傅里叶变换的例子[都是float类型的]!" << endl;
- d_fft_plan = fftwf_plan_dft_1d(N,
- reinterpret_cast<fftwf_complex*>(a),
- reinterpret_cast<fftwf_complex*>(a_fft),
- FFTW_FORWARD,
- FFTW_ESTIMATE);
- d_ifft_plan = fftwf_plan_dft_1d(N,
- reinterpret_cast<fftwf_complex*>(a_fft),
- reinterpret_cast<fftwf_complex*>(a_fft_ifft),
- FFTW_BACKWARD, // 逆变换
- FFTW_ESTIMATE);
- // 真正开始执行
- fftwf_execute(d_fft_plan);
- //fftw_execute(p);
- cout << "fftwf变换的结果:" << endl;
- print_complex_array(a_fft, N);
- fftwf_execute(d_ifft_plan);
- cout << "fftwf逆变换的结果:" << endl;
- cout << "需注意FFTW的逆变换没有除以N,即数据正变换再反变换后是原始数据的N倍。" << endl;
- print_complex_array(a_fft_ifft, N);
- //收尾
- fftwf_destroy_plan(d_fft_plan);
- fftwf_destroy_plan(d_ifft_plan);
- assert(a);
- fftwf_free(a);
- assert(a_fft);
- fftwf_free(a_fft);
- assert(a_fft_ifft);
- fftwf_free(a_fft_ifft);
- //assert(a == nullptr);
- //先用fftwf_malloc分配输入输出内存,然后输入数据赋值,
- //然后创建变换方案(fftwf_plan),然后执行变换(fftwf_execute),最后释放资源
-
- return 0;
- }
复制代码
|
-
|