2

嵌入式算法16:傅里叶变换算法

 2 years ago
source link: https://www.eefocus.com/embedded/514797
Go to the source link to view the article. You can view the picture content, updated content and better typesetting reading experience. If the link is broken, please click the button below to view the snapshot at that time.

摘要:傅里叶变换的核心在于,“任何连续周期信号可以由一组适当的正弦曲线组合而成”,在这个基础上对信号的中特定频率的正弦波进行分解或者重组,基于频率方面分析波形。

1、傅里叶变换的意义

近似周期性的方波(橙色),可采用6组正弦波(蓝色)合成,这是傅里叶的基础。

forward?url=https%3A%2F%2Fmmbiz.qpic.cn%2Fmmbiz_png%2F4W1T4tmuLNxf2hYEwUWG05UicV8rzHU5iaFIaCHgHLRcwFIsa7eryiaeHPFbd3Kbudy9DMZfhkxAWc117bZgMAjPA%2F640%3Fwx_fmt%3Dpng&s=504eba

对数字信号处理或者工程数学有一定基础,就明白傅里叶变换的价值。一般情况下的信号或者波形随时间变化,称为时域信号,时域(Time domain)是描述数学函数或物理信号对时间的关系。而以频率和幅度表示信号称为频域,频域(frequency domain)是描述信号在频率方面特性时用到的一种坐标系。其数学理论上暂且不理,针对嵌入式系统开发,只探讨其物理意义或者应用场景。

时域和频域是信号的基本性质,时域的表示较为形象与直观,比较符合一般认知,而频域分析则更为简练,剖析问题更为深刻和方便。

例如下图的左侧时域信号,其可以分解为2路正弦波的叠加效果,右侧为频域信号,表示2路正弦波的频率和幅度。傅里叶变化可简单理解为求解一段信号或波形,由哪些正弦波组成,也可以反向推导多路正弦波合并后的效果。

forward?url=https%3A%2F%2Fmmbiz.qpic.cn%2Fmmbiz_jpg%2F4W1T4tmuLNxf2hYEwUWG05UicV8rzHU5iaJ5gkib9kwXVL5ASVoFXicUolCEMusn0iazlT4Er0uPZR195r5ib6v1LRBg%2F640%3Fwx_fmt%3Djpeg&s=91dad6

基于动画形式表现如下:

forward?url=https%3A%2F%2Fmmbiz.qpic.cn%2Fmmbiz_gif%2F4W1T4tmuLNxf2hYEwUWG05UicV8rzHU5iampyiaGAHsLPXybvsRukAaXUgBYric0DxCqicADLRicquZBumsVicMRzbT0A%2F640%3Fwx_fmt%3Dgif&s=617259

傅里叶变换是一种信号分析方法,让我们对信号的构成和特点进行深入的、定量的研究。把信号通过频谱的方式进行准确的、定量的描述。将原来难以处理的时域信号转换成了易于分析的频域信号,即傅里叶变换的核心是从时域到频域的变换。

2、变换方式

数字信号属于离散值,对应的称为离散傅里叶变换(DFT),是傅里叶变换在时域和频域上都呈现离散的形式,将时域信号的采样变换为在离散时间傅里叶变换(DTFT)频域的采样。在形式上,变换两端(时域和频域上)的序列是有限长的,而实际上这两组序列都应当被认为是离散周期信号的主值序列。即使对有限长的离散信号作DFT,也应当将其看作经过周期延拓成为周期信号再作变换。

在实际应用中通常采用快速傅里叶变换以高效计算效率,快速傅里叶变换 FFT(Fast Fourier Transformation)是离散傅里叶变换 DFT(Discrete Fourier Transform)的快速算法。采用这种算法能使计算机计算离散傅里叶变换所需要的乘法次数大为减少,特别是被变换的抽样点数N越多,FFT算法计算量的节省就越显著。同理,从频域到时域的变换,称为逆变换,快速傅里叶逆变换 IFFT和离散傅里叶逆变换IDFT。

一般嵌入式系统使用快速傅里叶变换是分析某段信号中混合的杂波干扰,或者剔除某个频段后再逆变换。有些高端示波器可以对信号快速傅里叶变换 FFT,直接显示高频杂波的频率或者较大幅度干扰源的频率,以便由针对性的检查电路窜扰。本文主要是讲解使用fft/ifft进行干扰信号的过滤。

1、python生成正弦波,以及混合波形,提取数值作为c语言FFT/IFFT的数据源 

2、先进行FFT,输出幅频数据,导入Excel看看效果 

3、对数据进行简单过滤

4、过滤后的数据进行IFFT运算,再导入Excel看还原的效果 

5、关注微信公众号  嵌入式系统

3.1 基于python生成数据源

# This is a Python script.

import numpy as np
import matplotlib.pyplot as plt
def sin_wave(A, f, fs, phi, t):
    '''
    :params A:    振幅
    :params f:    信号频率
    :params fs:   采样频率
    :params phi:  相位
    :params t:    时间长度
    '''
    # 若时间序列长度为 t=1s,
    # 采样频率 fs=1000 Hz, 则采样时间间隔 Ts=1/fs=0.001s
    # 采样点个数为 n=t/Ts=1/0.001=1000, 即有1000个点,每个点间隔为 Ts
    Ts = 1/fs
    n = t / Ts
    n = np.arange(n)
    y = A*np.sin(2*np.pi*f*n*Ts + phi*(np.pi/180))
    return y

fs = 360*40
my_sin1 = sin_wave(100, 100, fs=fs, phi=0, t=0.08)
my_sin2 = sin_wave(20, 500, fs=fs, phi=0, t=0.08)
my_sin3 = sin_wave(10, 1100, fs=fs, phi=0, t=0.08)

x = np.arange(0, 0.08, 1/fs)
plt.xlabel('x-t (samRate=14400)')
plt.ylabel('y-A')
plt.grid()
plt.plot(x, my_sin1, 'k--',label="f=100")
plt.plot(x, my_sin2, 'b--',label="f=500")
plt.plot(x, my_sin3, 'g--',label="f=1100")
plt.plot(x, my_sin1+my_sin2+my_sin3, 'r',label="mix")
plt.legend()
plt.show()

np.savetxt("sin.txt",my_sin1+my_sin2+my_sin3, fmt='%.06f')

def print_name(name):
    print(f'Hi, {name}')


if __name__ == '__main__':
    print_name('fft test')

生成的波形图如下:

forward?url=https%3A%2F%2Fmmbiz.qpic.cn%2Fmmbiz_png%2F4W1T4tmuLNxf2hYEwUWG05UicV8rzHU5iaYkbaRFSy7dD4pxKjia63qMSD2aJGf7yNJcxYdcJ4icekwuISO8YPNZhQ%2F640%3Fwx_fmt%3Dpng&s=386f2f

三种频率的sin以及合成后的红色曲线,变形的正弦波;同时也将数据存入文件sin.txt备用。

3.2 基于C验证算法

将前面生成的数据使用c语言的FFT/IFFT进行处理,输出结果,并导入Excels生成图表展示效果。

#include <math.h>
#include <stdio.h>
#include "string.h"

typedef struct
{
    float real;
    float imag;
} complex_t;

#ifndef PI
#define PI             (3.14159265)
#endif

#define TYPE_FFT_E     float    /* Type is the same with complex_t member */

typedef complex_t TYPE_FFT;  /* Define complex_t in Config.h */
typedef unsigned int                      uint32_t;

#define     SAMPLE_NODES    (1024)
complex_t   fft_buff[SAMPLE_NODES];

//python生成的3个sin混合的波形数组
float my_sin_wave_table[] =
{
    0.000000,13.308217,25.359460,35.142296,42.082633,46.140080,47.788611,47.887148,47.470293,47.507139,48.682920,51.252751,\
    55.000000,59.307883,63.326941,66.199152,67.286436,66.350198,63.639610,59.866963,56.074110,53.418823,52.928233,55.274323,\
    60.621778,68.582540,78.287693,88.561119,98.156743,106.007043,111.428222,114.237346,114.756729,113.706291,112.009861,\
    110.560634,110.000000,110.560634,112.009861,113.706291,114.756729,114.237346,111.428222,106.007043,98.156743,88.561119,\
    78.287693,68.582540,60.621778,55.274323,52.928233,53.418823,56.074110,59.866963,63.639610,66.350198,67.286436,66.199152,\
    63.326941,59.307883,55.000000,51.252751,48.682920,47.507139,47.470293,47.887148,47.788611,46.140080,42.082633,35.142296,\
    25.359460,13.308217,0.000000,-13.308217,-25.359460,-35.142296,-42.082633,-46.140080,-47.788611,-47.887148,-47.470293,\
    -47.507139,-48.682920,-51.252751,-55.000000,-59.307883,-63.326941,-66.199152,-67.286436,-66.350198,-63.639610,-59.866963,\
    -56.074110,-53.418823,-52.928233,-55.274323,-60.621778,-68.582540,-78.287693,-88.561119,-98.156743,-106.007043,-111.428222,\
    -114.237346,-114.756729,-113.706291,-112.009861,-110.560634,-110.000000,-110.560634,-112.009861,-113.706291,-114.756729,\
    -114.237346,-111.428222,-106.007043,-98.156743,-88.561119,-78.287693,-68.582540,-60.621778,-55.274323,-52.928233,\
    -53.418823,-56.074110,-59.866963,-63.639610,-66.350198,-67.286436,-66.199152,-63.326941,-59.307883,-55.000000,\
    -51.252751,-48.682920,-47.507139,-47.470293,-47.887148,-47.788611,-46.140080,-42.082633,-35.142296,-25.359460,\
    -13.308217,-0.000000,13.308217,25.359460,35.142296,42.082633,46.140080,47.788611,47.887148,47.470293,47.507139,\
    48.682920,51.252751,55.000000,59.307883,63.326941,66.199152,67.286436,66.350198,63.639610,59.866963,56.074110,\
    53.418823,52.928233,55.274323,60.621778,68.582540,78.287693,88.561119,98.156743,106.007043,111.428222,114.237346,\
    114.756729,113.706291,112.009861,110.560634,110.000000,110.560634,112.009861,113.706291,114.756729,114.237346,\
    111.428222,106.007043,98.156743,88.561119,78.287693,68.582540,60.621778,55.274323,52.928233,53.418823,56.074110,\
    59.866963,63.639610,66.350198,67.286436,66.199152,63.326941,59.307883,55.000000,51.252751,48.682920,47.507139,\
    47.470293,47.887148,47.788611,46.140080,42.082633,35.142296,25.359460,13.308217,0.000000,-13.308217,-25.359460,\
    -35.142296,-42.082633,-46.140080,-47.788611,-47.887148,-47.470293,-47.507139,-48.682920,-51.252751,-55.000000,\
    -59.307883,-63.326941,-66.199152,-67.286436,-66.350198,-63.639610,-59.866963,-56.074110,-53.418823,-52.928233,\
    -55.274323,-60.621778,-68.582540,-78.287693,-88.561119,-98.156743,-106.007043,-111.428222,-114.237346,-114.756729,\
    -113.706291,-112.009861,-110.560634,-110.000000,-110.560634,-112.009861,-113.706291,-114.756729,-114.237346,-111.428222,\
    -106.007043,-98.156743,-88.561119,-78.287693,-68.582540,-60.621778,-55.274323,-52.928233,-53.418823,-56.074110,-59.866963,\
    -63.639610,-66.350198,-67.286436,-66.199152,-63.326941,-59.307883,-55.000000,-51.252751,-48.682920,-47.507139,-47.470293,\
    -47.887148,-47.788611,-46.140080,-42.082633,-35.142296,-25.359460,-13.308217,-0.000000,13.308217,25.359460,35.142296,\
    42.082633,46.140080,47.788611,47.887148,47.470293,47.507139,48.682920,51.252751,55.000000,59.307883,63.326941,66.199152,\
    67.286436,66.350198,63.639610,59.866963,56.074110,53.418823,52.928233,55.274323,60.621778,68.582540,78.287693,88.561119,\
    98.156743,106.007043,111.428222,114.237346,114.756729,113.706291,112.009861,110.560634,110.000000,110.560634,112.009861,\
    113.706291,114.756729,114.237346,111.428222,106.007043,98.156743,88.561119,78.287693,68.582540,60.621778,55.274323,\
    52.928233,53.418823,56.074110,59.866963,63.639610,66.350198,67.286436,66.199152,63.326941,59.307883,55.000000,51.252751,\
    48.682920,47.507139,47.470293,47.887148,47.788611,46.140080,42.082633,35.142296,25.359460,13.308217,0.000000,-13.308217,\
    -25.359460,-35.142296,-42.082633,-46.140080,-47.788611,-47.887148,-47.470293,-47.507139,-48.682920,-51.252751,-55.000000,\
    -59.307883,-63.326941,-66.199152,-67.286436,-66.350198,-63.639610,-59.866963,-56.074110,-53.418823,-52.928233,-55.274323,\
    -60.621778,-68.582540,-78.287693,-88.561119,-98.156743,-106.007043,-111.428222,-114.237346,-114.756729,-113.706291,\
    -112.009861,-110.560634,-110.000000,-110.560634,-112.009861,-113.706291,-114.756729,-114.237346,-111.428222,-106.007043,\
    -98.156743,-88.561119,-78.287693,-68.582540,-60.621778,-55.274323,-52.928233,-53.418823,-56.074110,-59.866963,-63.639610,\
    -66.350198,-67.286436,-66.199152,-63.326941,-59.307883,-55.000000,-51.252751,-48.682920,-47.507139,-47.470293,-47.887148,\
    -47.788611,-46.140080,-42.082633,-35.142296,-25.359460,-13.308217,-0.000000,13.308217,25.359460,35.142296,42.082633,\
    46.140080,47.788611,47.887148,47.470293,47.507139,48.682920,51.252751,55.000000,59.307883,63.326941,66.199152,67.286436,\
    66.350198,63.639610,59.866963,56.074110,53.418823,52.928233,55.274323,60.621778,68.582540,78.287693,88.561119,98.156743,\
    106.007043,111.428222,114.237346,114.756729,113.706291,112.009861,110.560634,110.000000,110.560634,112.009861,113.706291,\
    114.756729,114.237346,111.428222,106.007043,98.156743,88.561119,78.287693,68.582540,60.621778,55.274323,52.928233,53.418823,\
    56.074110,59.866963,63.639610,66.350198,67.286436,66.199152,63.326941,59.307883,55.000000,51.252751,48.682920,47.507139,\
    47.470293,47.887148,47.788611,46.140080,42.082633,35.142296,25.359460,13.308217,0.000000,-13.308217,-25.359460,-35.142296,\
    -42.082633,-46.140080,-47.788611,-47.887148,-47.470293,-47.507139,-48.682920,-51.252751,-55.000000,-59.307883,-63.326941,\
    -66.199152,-67.286436,-66.350198,-63.639610,-59.866963,-56.074110,-53.418823,-52.928233,-55.274323,-60.621778,-68.582540,\
    -78.287693,-88.561119,-98.156743,-106.007043,-111.428222,-114.237346,-114.756729,-113.706291,-112.009861,-110.560634,\
    -110.000000,-110.560634,-112.009861,-113.706291,-114.756729,-114.237346,-111.428222,-106.007043,-98.156743,-88.561119,\
    -78.287693,-68.582540,-60.621778,-55.274323,-52.928233,-53.418823,-56.074110,-59.866963,-63.639610,-66.350198,-67.286436,\
    -66.199152,-63.326941,-59.307883,-55.000000,-51.252751,-48.682920,-47.507139,-47.470293,-47.887148,-47.788611,-46.140080,\
    -42.082633,-35.142296,-25.359460,-13.308217,-0.000000,13.308217,25.359460,35.142296,42.082633,46.140080,47.788611,47.887148,\
    47.470293,47.507139,48.682920,51.252751,55.000000,59.307883,63.326941,66.199152,67.286436,66.350198,63.639610,59.866963,\
    56.074110,53.418823,52.928233,55.274323,60.621778,68.582540,78.287693,88.561119,98.156743,106.007043,111.428222,114.237346,\
    114.756729,113.706291,112.009861,110.560634,110.000000,110.560634,112.009861,113.706291,114.756729,114.237346,111.428222,\
    106.007043,98.156743,88.561119,78.287693,68.582540,60.621778,55.274323,52.928233,53.418823,56.074110,59.866963,63.639610,\
    66.350198,67.286436,66.199152,63.326941,59.307883,55.000000,51.252751,48.682920,47.507139,47.470293,47.887148,47.788611,\
    46.140080,42.082633,35.142296,25.359460,13.308217,0.000000,-13.308217,-25.359460,-35.142296,-42.082633,-46.140080,-47.788611,\
    -47.887148,-47.470293,-47.507139,-48.682920,-51.252751,-55.000000,-59.307883,-63.326941,-66.199152,-67.286436,-66.350198,\
    -63.639610,-59.866963,-56.074110,-53.418823,-52.928233,-55.274323,-60.621778,-68.582540,-78.287693,-88.561119,-98.156743,\
    -106.007043,-111.428222,-114.237346,-114.756729,-113.706291,-112.009861,-110.560634,-110.000000,-110.560634,-112.009861,\
    -113.706291,-114.756729,-114.237346,-111.428222,-106.007043,-98.156743,-88.561119,-78.287693,-68.582540,-60.621778,-55.274323,\
    -52.928233,-53.418823,-56.074110,-59.866963,-63.639610,-66.350198,-67.286436,-66.199152,-63.326941,-59.307883,-55.000000,\
    -51.252751,-48.682920,-47.507139,-47.470293,-47.887148,-47.788611,-46.140080,-42.082633,-35.142296,-25.359460,-13.308217,\
    -0.000000,13.308217,25.359460,35.142296,42.082633,46.140080,47.788611,47.887148,47.470293,47.507139,48.682920,51.252751,\
    55.000000,59.307883,63.326941,66.199152,67.286436,66.350198,63.639610,59.866963,56.074110,53.418823,52.928233,55.274323,\
    60.621778,68.582540,78.287693,88.561119,98.156743,106.007043,111.428222,114.237346,114.756729,113.706291,112.009861,\
    110.560634,110.000000,110.560634,112.009861,113.706291,114.756729,114.237346,111.428222,106.007043,98.156743,88.561119,\
    78.287693,68.582540,60.621778,55.274323,52.928233,53.418823,56.074110,59.866963,63.639610,66.350198,67.286436,66.199152,\
    63.326941,59.307883,55.000000,51.252751,48.682920,47.507139,47.470293,47.887148,47.788611,46.140080,42.082633,35.142296,\
    25.359460,13.308217,0.000000,-13.308217,-25.359460,-35.142296,-42.082633,-46.140080,-47.788611,-47.887148,-47.470293,\
    -47.507139,-48.682920,-51.252751,-55.000000,-59.307883,-63.326941,-66.199152,-67.286436,-66.350198,-63.639610,-59.866963,\
    -56.074110,-53.418823,-52.928233,-55.274323,-60.621778,-68.582540,-78.287693,-88.561119,-98.156743,-106.007043,-111.428222,\
    -114.237346,-114.756729,-113.706291,-112.009861,-110.560634,-110.000000,-110.560634,-112.009861,-113.706291,-114.756729,\
    -114.237346,-111.428222,-106.007043,-98.156743,-88.561119,-78.287693,-68.582540,-60.621778,-55.274323,-52.928233,-53.418823,\
    -56.074110,-59.866963,-63.639610,-66.350198,-67.286436,-66.199152,-63.326941,-59.307883,-55.000000,-51.252751,-48.682920,\
    -47.507139,-47.470293,-47.887148,-47.788611,-46.140080,-42.082633,-35.142296,-25.359460,-13.308217,-0.000000,13.308217,\
    25.359460,35.142296,42.082633,46.140080,47.788611,47.887148,47.470293,47.507139,48.682920,51.252751,55.000000,59.307883,\
    63.326941,66.199152,67.286436,66.350198,63.639610,59.866963,56.074110,53.418823,52.928233,55.274323,60.621778,68.582540,\
    78.287693,88.561119,98.156743,106.007043,111.428222,114.237346,114.756729,113.706291,112.009861,110.560634,110.000000,\
    110.560634,112.009861,113.706291,114.756729,114.237346,111.428222,106.007043,98.156743,88.561119,78.287693,68.582540,\
    60.621778,55.274323,52.928233,53.418823,56.074110,59.866963,63.639610,66.350198,67.286436,66.199152,63.326941,59.307883,\
    55.000000,51.252751,48.682920,47.507139,47.470293,47.887148,47.788611,46.140080,42.082633,35.142296,25.359460,13.308217,\
    0.000000,-13.308217,-25.359460,-35.142296,-42.082633,-46.140080,-47.788611,-47.887148,-47.470293,-47.507139,-48.682920,\
    -51.252751,-55.000000,-59.307883,-63.326941,-66.199152,-67.286436,-66.350198,-63.639610,-59.866963,-56.074110,-53.418823,\
    -52.928233,-55.274323,-60.621778,-68.582540,-78.287693,-88.561119,-98.156743,-106.007043,-111.428222,-114.237346,-114.756729,\
    -113.706291,-112.009861,-110.560634,-110.000000,-110.560634,-112.009861,-113.706291,-114.756729,-114.237346,-111.428222,\
    -106.007043,-98.156743,-88.561119,-78.287693,-68.582540,-60.621778,-55.274323,-52.928233,-53.418823,-56.074110,-59.866963,\
    -63.639610,-66.350198,-67.286436,-66.199152,-63.326941,-59.307883,-55.000000,-51.252751,-48.682920,-47.507139,-47.470293,\
    -47.887148,-47.788611,-46.140080,-42.082633,-35.142296,-25.359460,-13.308217,-0.000000,13.308217,25.359460,35.142296,\
    42.082633,46.140080,47.788611,47.887148,47.470293,47.507139,48.682920,51.252751,55.000000,59.307883,63.326941,66.199152,\
    67.286436,66.350198,63.639610,59.866963,56.074110,53.418823,52.928233,55.274323,60.621778,68.582540,78.287693,88.561119,\
    98.156743,106.007043,111.428222,114.237346,114.756729,113.706291,112.009861,110.560634,110.000000,110.560634,112.009861,\
    113.706291,114.756729,114.237346,111.428222,106.007043,98.156743,88.561119,78.287693,68.582540,60.621778,55.274323,52.928233,\
    53.418823,56.074110,59.866963,63.639610,66.350198,67.286436,66.199152,63.326941,59.307883,55.000000,51.252751,48.682920,\
    47.507139,47.470293,47.887148,47.788611,46.140080,42.082633,35.142296,25.359460,13.308217,0.000000,-13.308217,-25.359460,\
    -35.142296,-42.082633,-46.140080,-47.788611,-47.887148,-47.470293,-47.507139,-48.682920,-51.252751,-55.000000,-59.307883,\
    -63.326941,-66.199152,-67.286436,-66.350198,-63.639610,-59.866963,-56.074110,-53.418823,-52.928233,-55.274323,-60.621778,\
    -68.582540,-78.287693,-88.561119,-98.156743,-106.007043,-111.428222,-114.237346,-114.756729,-113.706291,-112.009861,\
    -110.560634,-110.000000,-110.560634,-112.009861,-113.706291,-114.756729,-114.237346,-111.428222,-106.007043,-98.156743,\
    -88.561119,-78.287693,-68.582540,-60.621778,-55.274323,-52.928233,-53.418823,-56.074110,-59.866963,-63.639610,-66.350198,\
    -67.286436,-66.199152,-63.326941,-59.307883,-55.000000,-51.252751,-48.682920,-47.507139,-47.470293,-47.887148,-47.788611,\
    -46.140080,-42.082633,-35.142296,-25.359460,-13.308217,
};

//fft算法来自开源 https://github.com/xiahouzuoxin/fft

const float sin_tb[] =    // 精度(PI PI/2 PI/4 PI/8 PI/16 ... PI/(2^k))
{
        0.000000, 1.000000, 0.707107, 0.382683, 0.195090, 0.098017,
        0.049068, 0.024541, 0.012272, 0.006136, 0.003068, 0.001534,
        0.000767, 0.000383, 0.000192, 0.000096, 0.000048, 0.000024,
        0.000012, 0.000006, 0.000003
        };

const float cos_tb[] =    // 精度(PI PI/2 PI/4 PI/8 PI/16 ... PI/(2^k))
{
        -1.000000, 0.000000, 0.707107, 0.923880, 0.980785, 0.995185,
        0.998795, 0.999699, 0.999925, 0.999981, 0.999995, 0.999999,
        1.000000, 1.000000, 1.000000, 1.000000, 1.000000, 1.000000,
        1.000000, 1.000000, 1.000000
        };


int ones_32(uint32_t n)
{
    unsigned int c = 0 ;
    for(c = 0; n; ++c)
    {
        n &= (n - 1) ;
    }

    return c ;
}


uint32_t floor_log2_32(uint32_t fft_buff)
{
    fft_buff |= (fft_buff >> 1);
    fft_buff |= (fft_buff >> 2);
    fft_buff |= (fft_buff >> 4);
    fft_buff |= (fft_buff >> 8);
    fft_buff |= (fft_buff >> 16);

    return (ones_32(fft_buff >> 1));
}


/*
 * FFT Algorithm
 * === Inputs ===
 * fft_buff : complex numbers
 * N : nodes of FFT. @N should be power of 2, that is 2^(*)
 * === Output ===
 * the @fft_buff contains the result of FFT algorithm, so the original data
 * in @fft_buff is destroyed, please store them before using FFT.
 */
int fft(TYPE_FFT *fft_buff, uint32_t N)
{
    int i, j, l, k, ip;
    static uint32_t M = 0;
    static int le, le2;
    static TYPE_FFT_E sR, sI, tR, tI, uR, uI;

    M = floor_log2_32(N);

    /*
     * bit reversal sorting
     */

    l = N >> 1;
    j = l;
    ip = N - 2;
    for(i = 1; i <= ip; i++)
    {
        if(i < j)
        {
            tR = fft_buff[j].real;
            tI = fft_buff[j].imag;
            fft_buff[j].real = fft_buff[i].real;
            fft_buff[j].imag = fft_buff[i].imag;
            fft_buff[i].real = tR;
            fft_buff[i].imag = tI;
        }

        k = l;

        while(k <= j)
        {
            j = j - k;
            k = k >> 1;
        }

        j = j + k;
    }

    /*
     * For Loops
     */

    for(l = 1; l <= M; l++)  /* loop for ceil{log2(N)} */
    {
        le  = (int)(1 << l);
        le2 = (int)(le >> 1);
        uR = 1;
        uI = 0;

        k = floor_log2_32(le2);
        sR = cos_tb[k];
        sI = -sin_tb[k];

        for(j = 1; j <= le2; j++)  /* loop for each sub DFT */
        {
            for(i = j - 1; i < N; i += le) /* loop for each butterfly */
            {
                ip = i + le2;
                tR = fft_buff[ip].real * uR - fft_buff[ip].imag * uI;
                tI = fft_buff[ip].real * uI + fft_buff[ip].imag * uR;
                fft_buff[ip].real = fft_buff[i].real - tR;
                fft_buff[ip].imag = fft_buff[i].imag - tI;
                fft_buff[i].real += tR;
                fft_buff[i].imag += tI;
            }  /* Next i */

            tR = uR;
            uR = tR * sR - uI * sI;
            uI = tR * sI + uI * sR;

        } /* Next j */

    } /* Next l */

    return 0;
}


/*
 * Inverse FFT Algorithm
 * === Inputs ===
 * fft_buff : complex numbers
 * N : nodes of FFT. @N should be power of 2, that is 2^(*)
 * === Output ===
 * the @fft_buff contains the result of FFT algorithm, so the original data
 * in @fft_buff is destroyed, please store them before using FFT.
 */

int ifft(TYPE_FFT *fft_buff, uint32_t N)
{
    int k = 0;

    for(k = 0; k <= N - 1; k++)
    {
        fft_buff[k].imag = -fft_buff[k].imag;
    }

    fft(fft_buff, N);    /* using FFT */

    for(k = 0; k <= N - 1; k++)
    {
        fft_buff[k].real = fft_buff[k].real / N;
        fft_buff[k].imag = -fft_buff[k].imag / N;
    }

    return 0;
}

static void import_data(void)
{
    int i;

    for(i = 0; i < SAMPLE_NODES; i++)
    {
        fft_buff[i].real = my_sin_wave_table[i];//取前1024个数进行fft变换
        fft_buff[i].imag  = 0.0f;
    }
}


int main(int argc, char *argv[])
{
    int i;
    int f;//频率
    float a;//幅度
    int fd;
    float t;

    printf("FFT\r\n");

    import_data();

    fft(fft_buff, SAMPLE_NODES);

    //fft后的结果在fft_buff
    //将其实部与虚部处理,输出频点与幅度值,导入Excel看效果
    //因为是周期性质,取前半部分即可

    //数据采样频率是360*40,均分到SAMPLE_NODES,则对应频点间隔是 360*40/SAMPLE_NODES
    fd=360*40/SAMPLE_NODES;

    for(i = 0; i < SAMPLE_NODES / 2; i++)
    {
        f = i *fd;
        a = (double)sqrt(fft_buff[i].real * fft_buff[i].real + fft_buff[i].imag * fft_buff[i].imag)/ (SAMPLE_NODES / 2);//转换幅度
        //printf("%d,%f\n", f, a);//>>导入excel查看幅频图效果
    }

    //过滤高频部分
    //将幅度小于某个值的,以30为例过滤
    for(i = 0; i < SAMPLE_NODES; i++)
    {
        a = (double)sqrt(fft_buff[i].real * fft_buff[i].real + fft_buff[i].imag * fft_buff[i].imag)/ (SAMPLE_NODES / 2);
        if(a<30)
        {
            fft_buff[i].real = 0;
            fft_buff[i].imag = 0;
        }
    }

    //再进行逆运算还原
    ifft(fft_buff, SAMPLE_NODES);
    for(i = 0; i < SAMPLE_NODES; i++)
    {
        t=1.0/(360*40)*i;//结合采样频率步进,方便查看波形效果
        //printf("%f,%f\n", t,fft_buff[i].real);//>>导入excel查看还原后的sin效果
    }

    return 0;
}

3.3 幅频图

混合波形使用FFT后的幅频图:

看图可知对应的频点和幅度,与原始混合波形的3种正弦波吻合。

栅栏效应:其中频点与采样频率和FFT转换的采样数有关,需要结合实际定义,因为数字离散信号,频点并不是连续的,如上大概相差14,因为数据点数N为1024点,原始信号的采样频率为fs=14400Hz,则频率间隔fs/N约等于14Hz。这样看频谱图只有fs/N倍数的,例如没有410Hz,只能看到相邻频率406Hz和420Hz处频谱泄漏的数据,相当于透过栅栏观赏风景,只能看到频谱的一部分,而其它频率点看不见,因此很可能部分有用的频率成分被漏掉,这种现象被称为栅栏效应。

3.4 逆变换还原

FFT后的幅频数据过滤其中2个频点的小幅度波形,还原后的效果如下图:

过滤还原后的波形,与原始图,前面python生成的黑色虚线吻合。

傅里叶变换在数字信号领域中应用广泛,在嵌入式系统中,一般用来分析ADC高速采集的数据、滤波处理,或者对音频进行简单分析。但傅里叶变换算法内存消耗有点大,低端设备可能无法运行。


版权声明:与非网经原作者授权转载,版权属于原作者。文章观点仅代表作者本人,不代表与非网立场。文章及其配图仅供工程师学习之用,如有侵权或者其他问题,请联系本站作侵删。 侵权投诉


About Joyk


Aggregate valuable and interesting links.
Joyk means Joy of geeK