/*----------------------------------------------------------------------------
fft.c - fast Fourier transform and its inverse (both recursively)
Copyright (C) 2004, Jerome R. Breitenbach. All rights reserved.
The author gives permission to anyone to freely copy, distribute, and use
this file, under the following conditions:
- No changes are made.
- No direct commercial advantage is obtained.
- No liability is attributed to the author for any damages incurred.
----------------------------------------------------------------------------*/
/******************************************************************************
* This file defines a C function fft that, by calling another function *
* fft_rec (also defined), calculates an FFT recursively. Usage: *
* fft(N, x, X); *
* Parameters: *
* N: number of points in FFT (must equal 2^n for some integer n >= 1) *
* x: pointer to N time-domain samples given in rectangular form (Re x, *
* Im x) *
* X: pointer to N frequency-domain samples calculated in rectangular form *
* (Re X, Im X) *
* Similarly, a function ifft with the same parameters is defined that *
* calculates an inverse FFT (IFFT) recursively. Usage: *
* ifft(N, x, X); *
* Here, N and X are given, and x is calculated. *
******************************************************************************/
#include <stdlib.h>
#include <math.h>
/* macros */
#define TWO_PI (6.2831853071795864769252867665590057683943L)
/* function prototypes */
void fft(int N, double (*x)[2], double (*X)[2]);
void fft_rec(int N, int offset, int delta,
double (*x)[2], double (*X)[2], double (*XX)[2]);
void ifft(int N, double (*x)[2], double (*X)[2]);
/* FFT */
void fft(int N, double (*x)[2], double (*X)[2])
{
/* Declare a pointer to scratch space. */
double (*XX)[2] = malloc(2 * N * sizeof(double));
/* Calculate FFT by a recursion. */
fft_rec(N, 0, 1, x, X, XX);
/* Free memory. */
free(XX);
}
/* FFT recursion */
void fft_rec(int N, int offset, int delta,
double (*x)[2], double (*X)[2], double (*XX)[2])
{
int N2 = N/2; /* half the number of points in FFT */
int k; /* generic index */
double cs, sn; /* cosine and sine */
int k00, k01, k10, k11; /* indices for butterflies */
double tmp0, tmp1; /* temporary storage */
if(N != 2) /* Perform recursive step. */
{
/* Calculate two (N/2)-point DFT's. */
fft_rec(N2, offset, 2*delta, x, XX, X);
fft_rec(N2, offset+delta, 2*delta, x, XX, X);
/* Combine the two (N/2)-point DFT's into one N-point DFT. */
for(k=0; k<N2; k++)
{
k00 = offset + k*delta; k01 = k00 + N2*delta;
k10 = offset + 2*k*delta; k11 = k10 + delta;
cs = cos(TWO_PI*k/(double)N); sn = sin(TWO_PI*k/(double)N);
tmp0 = cs * XX[k11][0] + sn * XX[k11][1];
tmp1 = cs * XX[k11][1] - sn * XX[k11][0];
X[k01][0] = XX[k10][0] - tmp0;
X[k01][1] = XX[k10][1] - tmp1;
X[k00][0] = XX[k10][0] + tmp0;
X[k00][1] = XX[k10][1] + tmp1;
}
}
else /* Perform 2-point DFT. */
{
k00 = offset; k01 = k00 + delta;
X[k01][0] = x[k00][0] - x[k01][0];
X[k01][1] = x[k00][1] - x[k01][1];
X[k00][0] = x[k00][0] + x[k01][0];
X[k00][1] = x[k00][1] + x[k01][1];
}
}
/* IFFT */
void ifft(int N, double (*x)[2], double (*X)[2])
{
int N2 = N/2; /* half the number of points in IFFT */
int i; /* generic index */
double tmp0, tmp1; /* temporary storage */
/* Calculate IFFT via reciprocity property of DFT. */
fft(N, X, x);
x[0][0] = x[0][0]/N; x[0][1] = x[0][1]/N;
x[N2][0] = x[N2][0]/N; x[N2][1] = x[N2][1]/N;
for(i=1; i<N2; i++)
{
tmp0 = x[i][0]/N; tmp1 = x[i][1]/N;
x[i][0] = x[N-i][0]/N; x[i][1] = x[N-i][1]/N;
x[N-i][0] = tmp0; x[N-i][1] = tmp1;
}
}
fft.zip_fft_fft linux_linux fft
版权申诉
84 浏览量
2022-09-23
09:55:25
上传
评论
收藏 2KB ZIP 举报
我虽横行却不霸道
- 粉丝: 75
- 资源: 1万+
最新资源
- Python 程序语言设计模式思路-行为型模式:策略模式:将算法封装成独立的类,并使它们可以互相替换及支付模式数据压缩
- main.py
- Last Loaded Test.DBK
- Screenshot_20240520_163011.jpg
- ubuntu-python3-whisper-tornado docker镜像 Dockerfile
- ubuntu-python3-whisper-tornado docker镜像07
- 新录音 8.m4a
- ubuntu-python3-whisper-tornado docker镜像
- ubuntu-python3-whisper-tornado docker镜像
- ubuntu-python3-whisper-tornado docker镜像09
资源上传下载、课程学习等过程中有任何疑问或建议,欢迎提出宝贵意见哦~我们会及时处理!
点击此处反馈