【经典算法实现 40】一维傅里叶变换 及 逆变换代码实现
本文链接：《【经典算法实现 40】一维傅里叶变换 及 逆变换代码实现》
一、一维傅里叶变换 及 其逆变换
一维傅里叶变换：
F ( u ) = ∑ M = 0 M − 1 f ( x ) e − j u x 2 π M ， u = 0 , 1 , 2 , ⋯ , M − 1 F(u) = \sum_{M=0}^{M-1} f(x) e^{-j u x \frac{2\pi} M}， u=0,1,2,\cdots,M-1 F(u)=M=0∑M−1f(x)e−juxM2π，u=0,1,2,⋯,M−1
一维傅里叶逆变换：
f
(
x
)
=
1
M
∑
u
=
0
M
−
1
F
(
u
)
e
j
u
x
2
π
M
，
x
=
0
,
1
,
2
,
⋯
,
M
−
1
f(x) = {\frac 1 M} \sum_{u=0}^{M-1} F(u) e^{j u x \frac{2\pi} M}， x=0,1,2,\cdots,M-1
f(x)=M1u=0∑M−1F(u)ejuxM2π，x=0,1,2,⋯,M−1
代入欧拉公式:
e
j
x
=
cos
(
x
)
+
j
sin
(
x
)
e^{jx} = \cos(x) + j \sin(x)
ejx=cos(x)+jsin(x)
e
−
j
x
=
cos
(
x
)
−
j
sin
(
x
)
e^{-jx} = \cos(x) - j \sin(x)
e−jx=cos(x)−jsin(x)
e
j
π
+
1
=
0
e^{j\pi} + 1 = 0
ejπ+1=0
sin
(
x
)
=
e
i
x
−
e
−
j
x
2
i
\sin(x) = \frac {e^{ix} - e^{-jx}} {2i}
sin(x)=2ieix−e−jx
cos
(
x
)
=
e
i
x
+
e
−
j
x
2
i
\cos(x) = \frac {e^{ix} + e^{-jx}} {2i}
cos(x)=2ieix+e−jx
我们可以得到如下：
一维傅里叶变换：
F
(
u
)
=
∑
M
=
0
M
−
1
f
(
x
)
(
e
−
j
u
x
2
π
M
)
，
u
=
0
,
1
,
2
,
⋯
,
M
−
1
⇒
F
(
u
)
=
∑
M
=
0
M
−
1
f
(
x
)
(
cos
(
u
x
2
π
M
)
−
j
sin
(
u
x
2
π
M
)
)
，
u
=
0
,
1
,
2
,
⋯
,
M
−
1
F(u) = \sum_{M=0}^{M-1} f(x) \left( e^{-j u x \frac{2\pi} M} \right)， u=0,1,2,\cdots,M-1 \\ \Rightarrow F(u) = \sum_{M=0}^{M-1} f(x) \left( \cos(u x \frac{2\pi} M) - j \sin(u x \frac{2\pi} M) \right)， u=0,1,2,\cdots,M-1
F(u)=M=0∑M−1f(x)(e−juxM2π)，u=0,1,2,⋯,M−1⇒F(u)=M=0∑M−1f(x)(cos(uxM2π)−jsin(uxM2π))，u=0,1,2,⋯,M−1
一维傅里叶逆变换：
f
(
x
)
=
1
M
∑
u
=
0
M
−
1
F
(
u
)
e
j
u
x
2
π
M
，
x
=
0
,
1
,
2
,
⋯
,
M
−
1
⇒
f
(
x
)
=
1
M
∑
u
=
0
M
−
1
F
(
u
)
(
cos
(
u
x
2
π
M
)
+
j
sin
(
u
x
2
π
M
)
)
，
x
=
0
,
1
,
2
,
⋯
,
M
−
1
f(x) = {\frac 1 M} \sum_{u=0}^{M-1} F(u) e^{j u x \frac{2\pi} M}， x=0,1,2,\cdots,M-1 \\ \Rightarrow f(x) = {\frac 1 M} \sum_{u=0}^{M-1} F(u) \left( \cos(u x \frac{2\pi} M) + j \sin(u x \frac{2\pi} M) \right) ， x=0,1,2,\cdots,M-1
f(x)=M1u=0∑M−1F(u)ejuxM2π，x=0,1,2,⋯,M−1⇒f(x)=M1u=0∑M−1F(u)(cos(uxM2π)+jsin(uxM2π))，x=0,1,2,⋯,M−1
在计算每个
F
(
u
)
F(u)
F(u)时，
由于 一维傅里叶逆变换中的
u
u
u 是一个复数，
假设复数
F
(
u
)
=
a
+
b
j
F(u) = a + b j
F(u)=a+bj，将
u
u
u 代入逆变换中，得到如下：
f
(
x
)
=
1
M
∑
u
=
0
M
−
1
F
(
u
)
(
cos
(
u
x
2
π
M
)
+
j
sin
(
u
x
2
π
M
)
)
，
x
=
0
,
1
,
2
,
⋯
,
M
−
1
f(x) = {\frac 1 M} \sum_{u=0}^{M-1} F(u) \left( \cos(u x \frac{2\pi} M) + j \sin(u x \frac{2\pi} M) \right) ， x=0,1,2,\cdots,M-1
f(x)=M1u=0∑M−1F(u)(cos(uxM2π)+jsin(uxM2π))，x=0,1,2,⋯,M−1
对于每个
F
(
u
)
F(u)
F(u)，实部虚部计算结果如下：
f
(
x
)
=
(
a
+
b
j
)
(
cos
(
u
x
2
π
M
)
+
j
sin
(
u
x
2
π
M
)
)
f(x) = ( a + b j) \left( \cos(u x \frac{2\pi} M) + j \sin(u x \frac{2\pi} M) \right)
f(x)=(a+bj)(cos(uxM2π)+jsin(uxM2π))
f
(
x
)
=
(
a
cos
(
u
x
2
π
M
)
+
b
j
cos
(
u
x
2
π
M
)
+
a
j
sin
(
u
x
2
π
M
)
+
b
j
⋅
j
sin
(
u
x
2
π
M
)
)
f(x) =\left( a\cos(u x \frac{2\pi} M) + bj\cos(u x \frac{2\pi} M) + aj \sin(u x \frac{2\pi} M) + bj\cdot j \sin(u x \frac{2\pi} M) \right)
f(x)=(acos(uxM2π)+bjcos(uxM2π)+ajsin(uxM2π)+bj⋅jsin(uxM2π))
分离虚部，实部结果如下：
f
(
x
)
.
r
e
a
l
=
(
a
cos
(
u
x
2
π
M
)
−
b
sin
(
u
x
2
π
M
)
)
f(x).real =\left( a\cos(u x \frac{2\pi} M) - b \sin(u x \frac{2\pi} M) \right)
f(x).real=(acos(uxM2π)−bsin(uxM2π))
f
(
x
)
.
i
m
a
g
=
(
a
sin
(
u
x
2
π
M
)
+
b
cos
(
u
x
2
π
M
)
)
⋅
j
f(x).imag =\left( a \sin(u x \frac{2\pi} M) + b\cos(u x \frac{2\pi} M) \right) \cdot j
f(x).imag=(asin(uxM2π)+bcos(uxM2π))⋅j
二、C代码实现
// 一维傅 里叶变换代码实现
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#define NUM 15
typedef struct Complex_{
double x;
double y;
}Complex_;
float f[NUM];
Complex_ F[NUM];
Complex_ IF[NUM];
// 一维傅里叶变换
void DFT(void)
{
int i=0, j=0;
double tmp;
clock_t start, end;
start = clock();
srand(time(NULL));
// 初始化 F[NUM]
for(i = 0; i<NUM; i++){
//f[i] = sin(M_PI * i/(NUM*2));
f[i] = rand()%99;
F[i].x = 0; // 实数
F[i].y = 0; // 虚数
}
// M 求和
for(j = 0; j<NUM; j++)
{
// 计算每个数的 傅里叶变换结果
for(i = 0; i<NUM; i++)
{
tmp = -2 * M_PI * j * i / NUM;
F[j].x += f[i] * cos(tmp); // 实数
F[j].y += f[i] * sin(tmp); // 虚数
//printf("f[%d]=%f a[%d]=%9f b[%d]=%9f \n",i, f[i], j, F[j].x, j, F[j].y);
}
}
end = clock();
printf("\n\n一维傅里叶变换 完毕, 使用时间:%lf NUM:%d\n", (double)(end-start)/CLOCKS_PER_SEC, NUM);
for(i = 0; i<NUM; i++){
printf("f[%d]=%9f ---> F[%d]=%12f + %12fj\n", i, f[i], i, F[i].x, F[i].y);
}
}
// 一维傅里叶逆变换
void IDFT(void)
{
int i=0, j=0;
double tmp;
clock_t start, end;
start = clock();
// 初始化 IF[NUM]
for(i = 0; i<NUM; i++){
IF[i].x = 0; // 实数
IF[i].y = 0; // 虚数
}
// M 求和
for(j = 0; j<NUM; j++)
{
// 计算每个数的 傅里叶变换结果
for(i = 0; i<NUM; i++)
{
tmp = 2 * M_PI * j * i / NUM;
IF[j].x += F[i].x * cos(tmp) - F[i].y * sin(tmp); // 实数
IF[j].y += F[i].x * sin(tmp) + F[i].y * cos(tmp); // 虚数
//printf("F[%d].x=%9f F[%d].y=%9f ---> IF[%d].x=%9f IF[%d].y=%9f\n", i, F[i].x, i, F[i].y, j, IF[j].x, j, IF[j].y );
}
IF[j].x /= NUM;
IF[j].y /= NUM;
//printf("F[%d].x=%9f F[%d].y=%9f ---> IF[%d].x=%9f IF[%d].y=%9f\n", j, F[j].x, j, F[j].y, j, IF[j].x, j, IF[j].y );
}
end = clock();
printf("\n\n一维傅里叶逆变换 完毕, 使用时间:%lf NUM:%d\n", (double)(end-start)/CLOCKS_PER_SEC, NUM);
for(i = 0; i<NUM; i++){
printf("IF[%d] = %9f + %9fj\n", i, IF[i].x, IF[i].y);
}
}
int main(void)
{
DFT();
printf("\n\n");
IDFT();
printf("\n\n");
return 0;
}
三、奇数位的变换结果
一维傅里叶变换 完毕, 使用时间:0.000000 NUM:15
f[ 0]=61.000000 ---> F[ 0]= 879.000000 + 0.000000j
f[ 1]=82.000000 ---> F[ 1]= 90.893436 + 13.541617j
f[ 2]=65.000000 ---> F[ 2]= -79.685020 + 5.178500j
f[ 3]=48.000000 ---> F[ 3]= 105.981226 + -28.670453j
f[ 4]=49.000000 ---> F[ 4]= -112.202453 + -67.242233j
f[ 5]=96.000000 ---> F[ 5]= -6.000000 + 6.928203j
f[ 6]=48.000000 ---> F[ 6]= -40.481226 + -21.662298j
f[ 7]=26.000000 ---> F[ 7]= 59.494037 + -155.029158j
f[ 8]= 1.000000 ---> F[ 8]= 59.494037 + 155.029158j
f[ 9]=88.000000 ---> F[ 9]= -40.481226 + 21.662298j
f[10]=45.000000 ---> F[10]= -6.000000 + -6.928203j
f[11]=88.000000 ---> F[11]= -112.202453 + 67.242233j
f[12]=44.000000 ---> F[12]= 105.981226 + 28.670453j
f[13]=89.000000 ---> F[13]= -79.685020 + -5.178500j
f[14]=49.000000 ---> F[14]= 90.893436 + -13.541617j
一维傅里叶逆变换 完毕, 使用时间:0.000000 NUM:15
IF[ 0] = 61.000000 + 0.000000j
IF[ 1] = 82.000000 + -0.000000j
IF[ 2] = 65.000000 + -0.000000j
IF[ 3] = 48.000000 + -0.000000j
IF[ 4] = 49.000000 + -0.000000j
IF[ 5] = 96.000000 + 0.000000j
IF[ 6] = 48.000000 + -0.000000j
IF[ 7] = 26.000000 + 0.000000j
IF[ 8] = 1.000000 + 0.000000j
IF[ 9] = 88.000000 + -0.000000j
IF[10] = 45.000000 + -0.000000j
IF[11] = 88.000000 + 0.000000j
IF[12] = 44.000000 + 0.000000j
IF[13] = 89.000000 + -0.000000j
IF[14] = 49.000000 + -0.000000j
四、偶数位的变换结果
一维傅里叶变换 完毕, 使用时间:0.000000 NUM:16
f[ 0]=74.000000 ---> F[ 0]= 907.000000 + 0.000000j
f[ 1]=83.000000 ---> F[ 1]= 74.197457 + -62.069951j
f[ 2]=48.000000 ---> F[ 2]= 62.468037 + 55.887302j
f[ 3]=34.000000 ---> F[ 3]= 81.074579 + 0.703911j
f[ 4]=96.000000 ---> F[ 4]= 35.000000 + 64.000000j
f[ 5]=29.000000 ---> F[ 5]= -53.802501 + -55.656479j
f[ 6]=87.000000 ---> F[ 6]= -70.468037 + -118.112698j
f[ 7]=53.000000 ---> F[ 7]= -21.469535 + 121.569659j
f[ 8]=54.000000 ---> F[ 8]= 63.000000 + -0.000000j
f[ 9]=41.000000 ---> F[ 9]= -21.469535 + -121.569659j
f[10]=21.000000 ---> F[10]= -70.468037 + 118.112698j
f[11]=75.000000 ---> F[11]= -53.802501 + 55.656479j
f[12]=36.000000 ---> F[12]= 35.000000 + -64.000000j
f[13]=26.000000 ---> F[13]= 81.074579 + -0.703911j
f[14]=69.000000 ---> F[14]= 62.468037 + -55.887302j
f[15]=81.000000 ---> F[15]= 74.197457 + 62.069951j
一维傅里叶逆变换 完毕, 使用时间:0.000000 NUM:16
IF[ 0] = 74.000000 + 0.000000j
IF[ 1] = 83.000000 + -0.000000j
IF[ 2] = 48.000000 + -0.000000j
IF[ 3] = 34.000000 + -0.000000j
IF[ 4] = 96.000000 + -0.000000j
IF[ 5] = 29.000000 + -0.000000j
IF[ 6] = 87.000000 + 0.000000j
IF[ 7] = 53.000000 + -0.000000j
IF[ 8] = 54.000000 + -0.000000j
IF[ 9] = 41.000000 + -0.000000j
IF[10] = 21.000000 + -0.000000j
IF[11] = 75.000000 + -0.000000j
IF[12] = 36.000000 + 0.000000j
IF[13] = 26.000000 + 0.000000j
IF[14] = 69.000000 + 0.000000j
IF[15] = 81.000000 + -0.000000j
从上面的奇数 和 偶数位的傅里叶运算结果中，大家有没有发现什么规律呢？
根据这个规律，我们引出了快速傅里叶变换，请看下文。
《傅里叶变换（由公式推导出来的c++源代码）》
《傅里叶变换和正弦函数和欧拉公式》
《傅里叶变换及其实现（MATLAB）》
890
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
