JP2002117016A - Fast fourier transform method and fast inverse fourier transform method - Google Patents

Abstract

PROBLEM TO BE SOLVED: To attain N=M×2n (M is an odd number)-point FFT and IFFT for input data. SOLUTION: In step S1, (i) is set to 0 first and in step S2, M pieces of data are taken out. In step S3, M-point DFT of the data taken out in step S2 is performed. In step S4, the product of obtained y(k) and a twist coefficient w(i,k) is calculated and the result is stored in y(k). In step S5, the value of y(k) is overwritten to the array (x) where the original data were put. The mentioned processes are repeated N/M times for all pieces of data through steps S6 and S7. In step S9, N/M(=2n)-point FFT is carried out for 0<=k<M. Lastly, rearrangement is carried out in step S12.

Description

DETAILED DESCRIPTION OF THE INVENTION

[0001]

BACKGROUND OF THE INVENTION 1. Field of the Invention The present invention generally relates to a fast Fourier transform method in which the number of points is 2 ⁿ (power of 2), and it is impossible to perform FFT with other number of points. It relates to an inverse Fourier transform method.

[0002]

2. Description of the Related Art Conventional Fast Fourier Transform (Fast Fourier Transform)
Transform (FFT) has a restriction that the number of points N must satisfy N = ^2n. Otherwise, the discrete Fourier transform (Digital) is performed using an existing device.
Fourier Transform (DFT) could not be performed at high speed.

[0003]

SUMMARY OF THE INVENTION N = ^{M.times.2.sup.n} (M is an odd number) points of FFT and IFFT on input data, which were impossible with the conventional fast Fourier transform and fast inverse Fourier transform, are now possible. It is an object of the present invention to provide a fast Fourier transform method and a fast inverse Fourier transform method.

[0004]

According to the fast Fourier transform method of the present invention, in order to solve the above-mentioned problem, M is an odd number, n is an integer, and M × 2 ^n- point complex number data is used as input data. It is characterized by performing fast Fourier transform on input data and outputting M × 2 ^n- point complex number data.

[0005] This fast Fourier transform method performs a predetermined pre-process on the input data, and then performs a 2 ^n- point fast Fourier transform on the M divided areas to obtain M × 2 ^n- points. Outputs the result of discrete Fourier transform of.

In the fast Fourier transform method, in the predetermined preprocessing, data division of the input data, M-point discrete Fourier transform on the divided data, and multiplication of the obtained M-point discrete Fourier transform coefficient Are performed in order.

In the fast Fourier transform method, when the divided data is subjected to an M-point discrete Fourier transform to obtain a discrete Fourier transform coefficient, the amount of calculation is reduced by using the symmetry of the trigonometric function.

[0008] The fast inverse Fourier transform method according to the present invention comprises:
In order to solve the above problem, M is an odd number, n is an integer, and M × 2 ^n- point complex number data is input data, and the input data is subjected to a fast inverse Fourier transform to obtain M × 2
It is characterized by outputting 2 ^n- point complex number data.

In this fast inverse Fourier transform method, a predetermined pre-process is performed on the input data, and then 2 ^n- point fast inverse Fourier transform is performed on the M divided areas, thereby obtaining an M × 2 ^{An n-} point discrete inverse Fourier transform result is output.

In the above-mentioned predetermined preprocessing of the fast inverse Fourier transform method, the data division of the input data, the M-point discrete inverse Fourier transform into the divided data, and the obtained M-point discrete inverse Fourier transform coefficient Are sequentially performed.

In this fast inverse Fourier transform method,
When performing M-point discrete inverse Fourier transform on the divided data to obtain a discrete inverse Fourier transform coefficient, the amount of calculation is reduced by using the symmetry of the trigonometric function.

[0012]

DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS First, the principle of FFT as a premise of the fast Fourier transform method and the fast inverse Fourier transform method of the present invention will be described.

The discrete Fourier transform (DFT) is based on given N data x (0), x (1),..., X (N-1).
This is a transform for obtaining the Fourier coefficient X (k) as shown in the following equation (1).

[0014]

(Equation 1)

In this DFT, the total number of data N is N =
When it can be factored into N ₁ N ₂ , n,
k can be expressed.

[0016]

(Equation 2)

Here, (n ₁ , k ₂ = 0..N ₁ −
1; n ₂ , k ₁ = 0. . . A N ₂ -1). At this time, the above equation (1) can be rewritten as the following equations (2) and (3).

[0018]

(Equation 3)

Further, if the following equations (4) and (5) are used, equation (6) is obtained.

[0020]

(Equation 4)

[0021]

(Equation 5)

In the above equation (6), if k ₁ is regarded as a constant, N 1
_{It is a two-} point DFT. That is, by performing the calculations of the above equations (4) and (5), the N-point DFT is calculated as N ₂
It is that it can be divided into _one point DFTN.
Then, if further factorized N _2, the smaller the number of points by applying repeatedly the above method DF
T can be divided. Thus, in the FFT, the DFT operation is efficiently performed by the division method.

Since the definition of the inverse transformation is the following equation (7),
The division calculation can be performed in the same manner as in the forward conversion.

[0024]

(Equation 6)

Next, an FFT apparatus having the number of points M × ²ⁿ , which is an embodiment of the present invention, will be described. First, the forward conversion will be described.

F of the number of points N = M × 2 ⁿ (M is an odd number)
FT is the M to _{N 1} described above, the ^{_{N 2 N / M (= 2}} n)
Is applied to the existing point number N = 2 ⁿ FF
This can be realized using T. FIG. 1 is a flowchart thereof. In the following description and FIG. 1, it is assumed that the array is a complex number, and that the real part and the imaginary part are substituted in the substitution between the elements of the array. For example, if x '(i) = x (i), x
The real and imaginary parts of (i) are substituted into the real and imaginary parts of x '(i), respectively. Also, N input data are x (0), x
(1),..., X (N−1) are stored in the array x.

The FFT method with the number of points M × ²ⁿ will be described below with reference to FIG. First, in step S1, i =
Then, in step S2, the M pieces of data are converted to x ′ (j).
(0 ≦ j <M). Next, in step S3, the data obtained in step S2 is subjected to the calculation of the above equation (4), that is, the M-point DFT. That is, the following equation (8) is calculated.

[0028]

(Equation 7)

The calculation of the equation (8) can be performed at high speed by performing the procedure described later. Next, step S4
Then, the product of the obtained y (k) and the twist coefficient w (i, k) is calculated (calculation of the above equation (5)), and the result is stored in y (k). The twist coefficient w (i, k) is a value defined by the following equation (9).

[0030]

(Equation 8)

Since the twist coefficient w (i, k) becomes 1 when k = 0, this operation may be performed in the range of 1 ≦ k ≦ M−1.

In step S5, the value of y (k) is overwritten on the array x containing the original data. Instead of overwriting, another array may be prepared and stored there.

FIG. 2 shows steps S1 to S in FIG.
5 is a schematic representation of the operation of FIG. From array x to N /
M pieces of data are taken out at every M, an operation of the product of the DFT and the twist coefficient is performed, and the result is returned to the array based on the result. This calculation is performed by N / M through steps S6 and S7.
Repeat for all input data.

Next, the operation of the above equation (6) is performed by a normal FFT
It should be done in. Then, after step S8, step S8
9, the FFT of N / M (= 2 ⁿ ) points is defined as 0 ≦ k <M
About going. If the result of this FFT is returned to the original array, the value to be obtained for the array x is obtained. The relationship between the position where data is present and the index of the actual array is as shown in FIG. Therefore, after steps S10 and S11 in FIG. 1, rearrangement is performed in step S12 as follows.

[0035]

(Equation 9)

In this array X (k) (0 ≦ k <N), discrete Fourier coefficients of x (i) are obtained.

Next, the inverse transformation will be described. In the inverse transform, the DFT of 2 is an IDFT shown in the following equation in FIG.

[0038]

(Equation 10)

In step S4 of FIG. 1, w (i, k) is set as shown in the following equation.

[0040]

[Equation 11]

Then, the FFT of step S9 in FIG.
By using the FFT, the calculation can be performed in the same procedure as the forward transform.

Next, the high-speed M-point DFT will be described.

[0043]

(Equation 12)

M (M is an odd number) point D shown in the above equation
The FT can be calculated at high speed by using the symmetry of the trigonometric function. Complex number x (n) (0 ≦ n ≦ M−1)
Is given by the following equation, y (k) becomes the following equation (10).

[0045]

(Equation 13)

[0046]

[Equation 14]

Here, the terms of the above equation (10) are expressed by the following equations (11), (12), (13) and (14).

[0048]

(Equation 15)

[0049]

(Equation 16)

[0050]

[Equation 17]

[0051]

(Equation 18)

Then, the above equation (10) becomes the following equation (15)
It becomes an expression.

[0053]

[Equation 19]

Further, the following equation holds when 1 ≦ n ≦ M−1.

[0055]

(Equation 20)

Therefore, the above equations (11) to (14) are calculated as in the following equations (16), (17), (18) and (19), so that the amount of calculation can be reduced. .

[0057]

(Equation 21)

[0058]

(Equation 22)

[0059]

(Equation 23)

[0060]

(Equation 24)

On the other hand, for k of 1 ≦ k ≦ (M−1) / 2, the following equation (20) holds.

[0062]

(Equation 25)

As is clear from comparison between the above equations (15) and (20), y (k) is represented by Y ₁ (k), Y ₂ (k), Y ₃ (k),
By calculating Y ₄ (k) after performing partial calculation, y (M−k) can be obtained only from addition / subtraction operations using the partial calculation. When k = 0, the following equation is established.

[0064]

(Equation 26)

It can be seen from this equation that the calculation can be performed only by addition.

From the above, it can be understood that the high-speed calculation of the M-point DFT can be summarized as shown in the flowchart of FIG.

First, y (0) = x (0) + x (1). . . x (M-
1) Calculate. Next, the following (a) and (b) are expressed as 1 ≦ k
Repeat for k of ≤ (M-1) / 2. Here, (a) represents Y ₁ (k), Y ₂ (k), and Y ₁ (k) in accordance with equations (16) to (19).
Y ₃ (k) and Y ₄ (k) are obtained, and (b) is (a)
Using the value obtained in step (15) and (20), y
(k) and y (Mk).

Similarly, the above high-speed calculation method can be applied to the IFFT expressed by the following equation.

[0069]

[Equation 27]

Then, assuming that y (n) = y _r (n) + ji _i (n) and 1 ≦ k ≦ (M−1) / 2, the following (2)
1), (22), (23), (24), (25), (2)
6)

[0071]

[Equation 28]

[0072]

(Equation 29)

[0073]

[Equation 30]

[0074]

(Equation 31)

[0075]

(Equation 32)

[0076]

[Equation 33]

First, x〜 (0) = y (0) + y
(1). . . Calculate y (M-1). Next, 1 ≦ k ≦ (M-1) / 2
For k of X _1, according to the equations (21) to (24), X ₁
(k), X ₂ (k), X ₃ (k), and X ₄ (k) are obtained, and x〜 (k) and x〜 (Mk) are obtained according to the equations (15) and (20). As described above, the IFFT can be calculated.

Next, the calculation amount will be described. In general, the operation costs of real number addition / subtraction and real number multiplication depend on the architecture of a computer, but it is assumed that all the operation costs are 1 in order to simplify the estimation of the operation amount.
Further, it is assumed that the operation cost of addition and subtraction of a complex number is two real number additions and subtractions (total 2), and the operation cost of complex number multiplication is four times of real number multiplication and two times of real number addition and subtraction (total 6).

Here, from the concept of the above calculation cost,
The calculation amount of the high-speed M-point DFT will be described. The calculation amount C _DFT (M) of the high-speed M-point DFT is
It is determined according to FIG.

First, after setting k = 0 in step S11 in FIG. 4, in step S12, M-1 complex number additions and subtractions are performed. Next, after adding k + 1 in step S13, in step S14, (M−1) / 2 × (4M−6) real number additions and subtractions and (M−1) / 2 × (2M−1) times are performed. Perform real number multiplication. Then, in step S15, (M-1)
Perform real number addition / subtraction twice / four times. The above processing is repeated until k <(M−1) / 2 is not satisfied in step S16.

Therefore, the computation amount C _DFT (M) of the high-speed M-point DFT can be obtained by the following equation.

[0082]

(Equation 34)

Next, the calculation amount of the M × 2 ^n- point FFT will be described. Since the calculation amount does not change between FFT and IFFT, only the calculation amount of FFT will be considered below.

First, a calculation amount C _FFT (N) when the N = 2 ^n- point FFT is performed in the radix-2 is obtained. N = 2 ⁿ
In the point FFT, the product with the twist coefficient is N
/ 2 · log ₂ (N) complex number multiplications and N · log ₂ (N) complex number additions for the butterfly operation, C
_FFT (N) is obtained by the following equation.

[0085]

(Equation 35)

Subsequently, N = M × 2ⁿPoint FFT
Calculation amount C ′ when performed by the method of the present invention_FFT(N)
Confuse. When the amount of calculation is determined along FIG.
In step S3, since DFT is performed N / M times, N / M
M × C _DFTThe calculation amount of (M) is obtained. Next, FIG.
In step S4, N / M × (M−1) complex number multiplications
Therefore, the calculation amount of 6 × N / M × (M−1) is obtained.
You. Then, in step S9 in FIG.
Since the FFT is performed M times, M × C_FFT(N / M)
Find the amount of computation. Then, the following equation is obtained.

[0087]

[Equation 36]

Therefore, in the FFT of the present invention,
It can be seen that if M is sufficiently smaller than N, the calculation amount on the order of N log ₂ (N) is sufficient.

For example, comparing the operation amounts of the FFT and the DFT of the present invention, when N = 5 × 2 ⁹ (= 2560), for example, C _DFT (2560) / C ′ _FFT (2560) = 208, and DFT The calculation can be performed about 200 times faster than in the case of.

Further, in order to compare the amount of operation between the FFT of the present invention and the conventional FFT using the radix of 2, the rate of increase in the amount of operation can be defined as in the following equation.

[0091]

(37)

In practice, the value N which is not a power of 2 is a radix 2
Although the function C _FFT (N) of the _FFT operation amount cannot be applied, this is used in order to obtain a reference for comparison. The following are the values of R _FFT (N, M) for some N and M combinations.

[0093]

[Table 1]

From the above table, it can be seen that if M does not become extremely large, an increase of several tens of percent is sufficient.

[0095]

According to the present invention, N = M, which is impossible with the conventional device.
FFT of × 2n points (M is an odd number) becomes possible. In addition, if there is an FFT device with a high speed, by utilizing the device, the development cost can be reduced as compared with a case where a completely new device is manufactured.

[Brief description of the drawings]

FIG. 1 is an FFT of the number of points N = M × 2 ⁿ (M is an odd number)
It is a flowchart for realizing.

FIG. 2 is a diagram schematically showing an operation of the flowchart shown in FIG. 1;

FIG. 3 is a diagram illustrating a relationship between a position of array data and an appropriate index position.

FIG. 4 is a flowchart of high-speed calculation of DFT.

[Explanation of symbols]

S2: a process of extracting M data from the array x, S3
DFT processing for M data, S4 DFT coefficient multiplication by a twist coefficient, S5 Processing to overwrite the result obtained in S4 on array x

 ────────────────────────────────────────────────── ─── Continued on the front page (72) Inventor Masayuki Nishiguchi 6-35 Kita Shinagawa, Shinagawa-ku, Tokyo Sony Corporation F-term (reference) 5B056 AA02 BB13 CC03 DD04

Claims (10)

[Claims] 1. When M is an odd number and n is an integer, M × 2
The complex number data for ⁿ points as input data, performs a fast Fourier transform on the input data, Fast Fourier transform method and outputting the complex data of the M × 2 ⁿ points. 2. Performing a predetermined pre-process on the input data, and then performing a 2 ^n- point fast Fourier transform on the M-divided region to obtain an M × 2 ^n- point discrete Fourier transform result 2. The fast Fourier transform method according to claim 1, wherein 3. The predetermined preprocessing includes, in order, data division of the input data, M-point discrete Fourier transform into the divided data, and multiplication of the obtained M-point discrete Fourier transform coefficients by a predetermined coefficient. 3. The fast Fourier transform method according to claim 2, wherein each method is performed. 4. The predetermined preprocessing includes extracting data of M points for each N / M obtained by dividing N points constituting an array x by an odd number M, performing a discrete Fourier transform on the data of the M points, and The M-point discrete Fourier transform coefficient is multiplied by a twist coefficient, and the multiplication result is represented by the array x
4. The fast Fourier transform method according to claim 3, wherein 5. The method according to claim 1, wherein the predetermined preprocessing is to reduce an operation amount by using a symmetry of a trigonometric function when performing M-point discrete Fourier transform on the divided data to obtain a discrete Fourier transform coefficient. 4. The fast Fourier transform method according to claim 3, wherein: 6. When M is an odd number and n is an integer, M × 2
A fast inverse Fourier transform method, characterized in that ^n- point complex number data is used as input data, the input data is subjected to fast inverse Fourier transform, and M × 2 ^n- point complex number data is output. 7. After performing predetermined pre-processing on the input data, a discrete inverse Fourier transform of M × 2 ⁿ points is performed by performing a 2 ^n- point fast inverse Fourier transform on the M divided regions. 7. The fast inverse Fourier transform method according to claim 6, wherein a conversion result is output. 8. The predetermined pre-processing includes data division of the input data, M-point discrete inverse Fourier transform to the divided data, and multiplication of the obtained M-point discrete inverse Fourier transform coefficient by a predetermined coefficient. 8. The fast inverse Fourier transform method according to claim 7, wherein 9. The predetermined preprocessing includes extracting M-point data for each N / M obtained by dividing N points constituting an array x by an odd number M, and performing a discrete inverse Fourier transform on the M-point data; 9. The fast inverse Fourier transform method according to claim 8, wherein the obtained M-point inverse discrete Fourier transform coefficient is multiplied by a twist coefficient, and the multiplication result is returned to the array x. 10. The predetermined pre-processing reduces an operation amount by using the symmetry of a trigonometric function when performing M-point discrete inverse Fourier transform on the divided data to obtain a discrete inverse Fourier transform coefficient. 9. The fast inverse Fourier transform method according to claim 8, wherein: