FFTW3.3.9快速傅里叶变换(FFT)库详解

需积分: 47 131 浏览量更新于2024-07-09 收藏 645KB PDF 举报

"FFTW是版本3.3.9的开源傅里叶变换运算库，由Matteo Frigo和Steven G. Johnson编写的文档。该库提供了对任意大小和维度的实值及复值数组进行快速傅里叶变换(FFT)的功能，执行效率高，时间复杂度为O(nlogn)。文档包括对不同类型的FFT操作的教程介绍，如复杂一维、多维DFT，以及针对实数数据的DFT处理。" FFTW是一个高效的开源库，专门用于执行快速傅里叶变换，它是自由软件，遵循特定的许可协议，允许复制和分发原始或修改后的文档，只要保留版权和许可通知。该库的最新版本3.3.9在2020年12月10日发布。 FFTW的教程部分涵盖了多种傅里叶变换类型： 1. 复杂一维DFTs：基础的傅里叶变换，将一维复数序列转换为其频域表示，这是所有FFT的基础。 2. 复杂多维DFTs：扩展到多维空间，适用于处理图像和其他多维数据，将每个维度独立进行一维变换，然后组合结果。 3. 一维DFTs of Real Data：对于实际数据的特殊优化，只计算非对称频谱的一半，从而减少计算量。 4. 多维DFTs of Real Data：在多个维度上应用这些优化，进一步节省计算资源。 5. More DFTs of Real Data：包括更复杂的情况，如半复共轭形式，适用于处理物理问题中的对称性。 FFTW的高效性来源于其自动代码生成技术，它能够为特定平台生成最优的FFT算法，同时支持多种编译器和硬件平台。此外，FFTW还提供各种接口，包括C、C++和Fortran，方便不同编程语言的开发者使用。在实际应用中，FFTW被广泛用于信号处理、图像分析、科学计算和工程领域，因其出色的性能和易用性而备受赞誉。例如，它可以用于滤波、频谱分析、图像处理中的卷积操作，以及在量子力学、流体动力学和统计建模等领域的计算任务。FFTW的灵活性和可扩展性使其成为科研和工程计算的首选工具之一。

10 FFTW 3.3.9

2.5 More DFTs of Real Data

FFTW supports several other transform types via a uniﬁed r2r (real-to-real) interface, so

called because it takes a real (double) array and outputs a real array of the same size.

These r2r transforms currently fall into three categories: DFTs of real input and complex-

Hermitian output in halfcomplex format, DFTs of real input with even/odd symmetry

(a.k.a. discrete cosine/sine transforms, DCTs/DSTs), and discrete Hartley transforms

(DHTs), all described in more detail by the following sections.

The r2r transforms follow the by now familiar interface of creating an fftw_plan, exe-

cuting it with fftw_execute(plan), and destroying it with fftw_destroy_plan(plan).

Furthermore, all r2r transforms share the same planner interface:

fftw_plan fftw_plan_r2r_1d(int n, double *in, double *out,

fftw_r2r_kind kind, unsigned flags);

fftw_plan fftw_plan_r2r_2d(int n0, int n1, double *in, double *out,

fftw_r2r_kind kind0, fftw_r2r_kind kind1,

unsigned flags);

fftw_plan fftw_plan_r2r_3d(int n0, int n1, int n2,

double *in, double *out,

fftw_r2r_kind kind0,

fftw_r2r_kind kind1,

fftw_r2r_kind kind2,

unsigned flags);

fftw_plan fftw_plan_r2r(int rank, const int *n, double *in, double *out,

const fftw_r2r_kind *kind, unsigned flags);

Just as for the complex DFT, these plan 1d/2d/3d/multi-dimensional transforms for con-

tiguous arrays in row-major order, transforming (real) input to output of the same size,

where n speciﬁes the physical dimensions of the arrays. All positive n are supported (with

the exception of n=1 for the FFTW_REDFT00 kind, noted in the real-even subsection below);

products of small factors are most eﬃcient (factorizing n-1 and n+1 for FFTW_REDFT00 and

FFTW_RODFT00 kinds, described below), but an O(n log n) algorithm is used even for prime

sizes.

Each dimension has a kind parameter, of type fftw_r2r_kind, specifying the kind of r2r

transform to be used for that dimension. (In the case of fftw_plan_r2r, this is an array

kind[rank] where kind[i] is the transform kind for the dimension n[i].) The kind can

be one of a set of predeﬁned constants, deﬁned in the following subsections.

In other words, FFTW computes the separable product of the speciﬁed r2r transforms

over each dimension, which can be used e.g. for partial diﬀerential equations with mixed

boundary conditions. (For some r2r kinds, notably the halfcomplex DFT and the DHT,

such a separable product is somewhat problematic in more than one dimension, however,

as is described below.)

In the current version of FFTW, all r2r transforms except for the halfcomplex type are com-

puted via pre- or post-processing of halfcomplex transforms, and they are therefore not as

fast as they could be. Since most other general DCT/DST codes employ a similar algorithm,

however, FFTW’s implementation should provide at least competitive performance.

Chapter 2: Tutorial 11

2.5.1 The Halfcomplex-format DFT

An r2r kind of FFTW_R2HC (r2hc) corresponds to an r2c DFT (see Section 2.3

[One-Dimensional DFTs of Real Data], page 6) but with “halfcomplex” format output,

and may sometimes be faster and/or more convenient than the latter. The inverse hc2r

transform is of kind FFTW_HC2R. This consists of the non-redundant half of the complex

output for a 1d real-input DFT of size n, stored as a sequence of n real numbers (double)

in the format:

, r

, . . . , r

n/2

, i

(n+1)/2−1

, . . . , i

, i

Here, r

is the real part of the kth output, and i

is the imaginary part. (Division by 2

is rounded down.) For a halfcomplex array hc[n], the kth component thus has its real

part in hc[k] and its imaginary part in hc[n-k], with the exception of k == 0 or n/2 (the

latter only if n is even)—in these two cases, the imaginary part is zero due to symmetries

of the real-input DFT, and is not stored. Thus, the r2hc transform of n real values is a

halfcomplex array of length n, and vice versa for hc2r.

Aside from the diﬀering format, the output of FFTW_R2HC/FFTW_HC2R is otherwise exactly

the same as for the corresponding 1d r2c/c2r transform (i.e. FFTW_FORWARD/FFTW_BACKWARD

transforms, respectively). Recall that these transforms are unnormalized, so r2hc followed

by hc2r will result in the original data multiplied by n. Furthermore, like the c2r transform,

an out-of-place hc2r transform will destroy its input array.

Although these halfcomplex transforms can be used with the multi-dimensional r2r interface,

the interpretation of such a separable product of transforms along each dimension is prob-

lematic. For example, consider a two-dimensional n0 by n1, r2hc by r2hc transform planned

by fftw_plan_r2r_2d(n0, n1, in, out, FFTW_R2HC, FFTW_R2HC, FFTW_MEASURE). Con-

ceptually, FFTW ﬁrst transforms the rows (of size n1) to produce halfcomplex rows, and

then transforms the columns (of size n0). Half of these column transforms, however, are

of imaginary parts, and should therefore be multiplied by i and combined with the r2hc

transforms of the real columns to produce the 2d DFT amplitudes; FFTW’s r2r transform

does not perform this combination for you. Thus, if a multi-dimensional real-input/output

DFT is required, we recommend using the ordinary r2c/c2r interface (see Section 2.4 [Multi-

Dimensional DFTs of Real Data], page 7).

2.5.2 Real even/odd DFTs (cosine/sine transforms)

The Fourier transform of a real-even function f(−x) = f(x) is real-even, and i times

the Fourier transform of a real-odd function f(−x) = −f(x) is real-odd. Similar results

hold for a discrete Fourier transform, and thus for these symmetries the need for complex

inputs/outputs is entirely eliminated. Moreover, one gains a factor of two in speed/space

from the fact that the data are real, and an additional factor of two from the even/odd

symmetry: only the non-redundant (ﬁrst) half of the array need be stored. The result is

the real-even DFT (REDFT) and the real-odd DFT (RODFT), also known as the discrete

cosine and sine transforms (DCT and DST), respectively.

(In this section, we describe the 1d transforms; multi-dimensional transforms are just a

separable product of these transforms operating along each dimension.)

12 FFTW 3.3.9

Because of the discrete sampling, one has an additional choice: is the data even/odd around

a sampling point, or around the point halfway between two samples? The latter corresponds

to shifting the samples by half an interval, and gives rise to several transform variants

denoted by REDFTab and RODFTab: a and b are 0 or 1, and indicate whether the input

(a) and/or output (b) are shifted by half a sample (1 means it is shifted). These are also

known as types I-IV of the DCT and DST, and all four types are supported by FFTW’s

r2r interface.

The r2r kinds for the various REDFT and RODFT types supported by FFTW, along with

the boundary conditions at both ends of the input array (n real numbers in[j=0..n-1]),

are:

• FFTW_REDFT00 (DCT-I): even around j = 0 and even around j = n − 1.

• FFTW_REDFT10 (DCT-II, “the” DCT): even around j = −0.5 and even around j =

n − 0.5.

• FFTW_REDFT01 (DCT-III, “the” IDCT): even around j = 0 and odd around j = n.

• FFTW_REDFT11 (DCT-IV): even around j = −0.5 and odd around j = n − 0.5.

• FFTW_RODFT00 (DST-I): odd around j = −1 and odd around j = n.

• FFTW_RODFT10 (DST-II): odd around j = −0.5 and odd around j = n − 0.5.

• FFTW_RODFT01 (DST-III): odd around j = −1 and even around j = n − 1.

• FFTW_RODFT11 (DST-IV): odd around j = −0.5 and even around j = n − 0.5.

Note that these symmetries apply to the “logical” array being transformed; there are no

constraints on your physical input data. So, for example, if you specify a size-5 REDFT00

(DCT-I) of the data abcde, it corresponds to the DFT of the logical even array abcdedcb of

size 8. A size-4 REDFT10 (DCT-II) of the data abcd corresponds to the size-8 logical DFT

of the even array abcddcba, shifted by half a sample.

All of these transforms are invertible. The inverse of R*DFT00 is R*DFT00; of R*DFT10 is

R*DFT01 and vice versa (these are often called simply “the” DCT and IDCT, respectively);

and of R*DFT11 is R*DFT11. However, the transforms computed by FFTW are unnor-

malized, exactly like the corresponding real and complex DFTs, so computing a transform

followed by its inverse yields the original array scaled by N , where N is the logical DFT

size. For REDFT00, N = 2(n − 1); for RODFT00, N = 2(n + 1); otherwise, N = 2n.

Note that the boundary conditions of the transform output array are given by the input

boundary conditions of the inverse transform. Thus, the above transforms are all inequiva-

lent in terms of input/output boundary conditions, even neglecting the 0.5 shift diﬀerence.

FFTW is most eﬃcient when N is a product of small factors; note that this diﬀers from the

factorization of the physical size n for REDFT00 and RODFT00! There is another oddity:

n=1 REDFT00 transforms correspond to N = 0, and so are not deﬁned (the planner will

return NULL). Otherwise, any positive n is supported.

For the precise mathematical deﬁnitions of these transforms as used by FFTW, see

Section 4.8 [What FFTW Really Computes], page 42. (For people accustomed to the

There are also type V-VIII transforms, which correspond to a logical DFT of odd size N , independent of

whether the physical size n is odd, but we do not support these variants.

Chapter 2: Tutorial 13

DCT/DST, FFTW’s deﬁnitions have a coeﬃcient of 2 in front of the cos/sin functions so

that they correspond precisely to an even/odd DFT of size N. Some authors also include

additional multiplicative factors of

√

2 for selected inputs and outputs; this makes the

transform orthogonal, but sacriﬁces the direct equivalence to a symmetric DFT.)

Which type do you need?

Since the required ﬂavor of even/odd DFT depends upon your problem, you are the best

judge of this choice, but we can make a few comments on relative eﬃciency to help you

in your selection. In particular, R*DFT01 and R*DFT10 tend to be slightly faster than

R*DFT11 (especially for odd sizes), while the R*DFT00 transforms are sometimes signiﬁ-

cantly slower (especially for even sizes).

Thus, if only the boundary conditions on the transform inputs are speciﬁed, we generally

recommend R*DFT10 over R*DFT00 and R*DFT01 over R*DFT11 (unless the half-sample

shift or the self-inverse property is signiﬁcant for your problem).

If performance is important to you and you are using only small sizes (say n < 200), e.g. for

multi-dimensional transforms, then you might consider generating hard-coded transforms

of those sizes and types that you are interested in (see Section 10.4 [Generating your own

code], page 100).

We are interested in hearing what types of symmetric transforms you ﬁnd most useful.

2.5.3 The Discrete Hartley Transform

If you are planning to use the DHT because you’ve heard that it is “faster” than the DFT

(FFT), stop here. The DHT is not faster than the DFT. That story is an old but enduring

misconception that was debunked in 1987.

The discrete Hartley transform (DHT) is an invertible linear transform closely related to

the DFT. In the DFT, one multiplies each input by cos − i ∗ sin (a complex exponential),

whereas in the DHT each input is multiplied by simply cos+sin. Thus, the DHT transforms

n real numbers to n real numbers, and has the convenient property of being its own inverse.

In FFTW, a DHT (of any positive n) can be speciﬁed by an r2r kind of FFTW_DHT.

Like the DFT, in FFTW the DHT is unnormalized, so computing a DHT of size n followed

by another DHT of the same size will result in the original array multiplied by n.

The DHT was originally proposed as a more eﬃcient alternative to the DFT for real data,

but it was subsequently shown that a specialized DFT (such as FFTW’s r2hc or r2c trans-

forms) could be just as fast. In FFTW, the DHT is actually computed by post-processing

an r2hc transform, so there is ordinarily no reason to prefer it from a performance per-

R*DFT00 is sometimes slower in FFTW because we discovered that the standard algorithm for computing

this by a pre/post-processed real DFT—the algorithm used in FFTPACK, Numerical Recipes, and other

sources for decades now—has serious numerical problems: it already loses several decimal places of accuracy

for 16k sizes. There seem to be only two alternatives in the literature that do not suﬀer similarly: a recursive

decomposition into smaller DCTs, which would require a large set of codelets for eﬃciency and generality,

or sacriﬁcing a factor of ∼ 2 in speed to use a real DFT of twice the size. We currently employ the latter

technique for general n, as well as a limited form of the former method: a split-radix decomposition when

n is odd (N a multiple of 4). For N containing many factors of 2, the split-radix method seems to recover

most of the speed of the standard algorithm without the accuracy tradeoﬀ.

剩余115页未读，继续阅读

TracelessLe

粉丝: 5w+
资源: 466

FFTW3.3.9快速傅里叶变换(FFT)库详解

傅里叶变换

FFTW开源lib库文件及源码

QRealFourier:基于Qt的FFT实现-开源

Mingw下的数据处理与分析：使用开源库进行数学运算

QT+fftw3库进行傅里叶变换

FFTW3实现快速傅里叶变换

傅立叶变换，离散傅立叶变换，快速傅立叶变换的区别。

qt调用fftw库对调制信号傅里叶变换

斯坦福大学傅里叶变换及应用笔记.pdf

C语言处理傅里叶变换的库

最新资源