0
# Fast Fourier Transform
1

2
GPU-accelerated FFT operations supporting 1D, 2D, and N-dimensional transforms for both complex and real data with comprehensive frequency domain processing capabilities. CuPy provides GPU-optimized FFT implementations using cuFFT library.
3

4
## Capabilities
5

6
### One-dimensional FFT
7

8
Basic 1D FFT operations for complex and real data.
9

10
```python { .api }
11
def fft(a, n=None, axis=-1, norm=None):
12
    """Compute 1D discrete Fourier Transform.
13
    
14
    Args:
15
        a: Input array
16
        n: Length of transformed axis
17
        axis: Axis over which to compute FFT
18
        norm: Normalization mode ('backward', 'ortho', 'forward')
19
        
20
    Returns:
21
        cupy.ndarray: Complex FFT result
22
    """
23

24
def ifft(a, n=None, axis=-1, norm=None):
25
    """Compute 1D inverse discrete Fourier Transform."""
26

27
def rfft(a, n=None, axis=-1, norm=None):
28
    """Compute 1D FFT of real input.
29
    
30
    Returns:
31
        cupy.ndarray: Complex FFT with conjugate symmetry
32
    """
33

34
def irfft(a, n=None, axis=-1, norm=None):
35
    """Compute inverse of rfft.
36
    
37
    Returns:
38
        cupy.ndarray: Real inverse FFT result
39
    """
40

41
def hfft(a, n=None, axis=-1, norm=None):
42
    """Compute FFT of signal with Hermitian symmetry."""
43

44
def ihfft(a, n=None, axis=-1, norm=None):
45
    """Compute inverse of hfft."""
46
```
47

48
### Two-dimensional FFT
49

50
2D FFT operations for image processing and 2D signal analysis.
51

52
```python { .api }
53
def fft2(a, s=None, axes=(-2, -1), norm=None):
54
    """Compute 2D discrete Fourier Transform.
55
    
56
    Args:
57
        a: Input array
58
        s: Shape of transformed axes
59
        axes: Axes over which to compute FFT
60
        norm: Normalization mode
61
        
62
    Returns:
63
        cupy.ndarray: 2D FFT result
64
    """
65

66
def ifft2(a, s=None, axes=(-2, -1), norm=None):
67
    """Compute 2D inverse discrete Fourier Transform."""
68

69
def rfft2(a, s=None, axes=(-2, -1), norm=None):
70
    """Compute 2D FFT of real input."""
71

72
def irfft2(a, s=None, axes=(-2, -1), norm=None):
73
    """Compute 2D inverse of rfft2."""
74
```
75

76
### N-dimensional FFT
77

78
General N-dimensional FFT for multi-dimensional data.
79

80
```python { .api }
81
def fftn(a, s=None, axes=None, norm=None):
82
    """Compute N-dimensional discrete Fourier Transform.
83
    
84
    Args:
85
        a: Input array
86
        s: Shape of transformed axes
87
        axes: Axes over which to compute FFT
88
        norm: Normalization mode
89
        
90
    Returns:
91
        cupy.ndarray: N-D FFT result
92
    """
93

94
def ifftn(a, s=None, axes=None, norm=None):
95
    """Compute N-dimensional inverse discrete Fourier Transform."""
96

97
def rfftn(a, s=None, axes=None, norm=None):
98
    """Compute N-dimensional FFT of real input."""
99

100
def irfftn(a, s=None, axes=None, norm=None):
101
    """Compute N-dimensional inverse of rfftn."""
102
```
103

104
### Frequency Arrays
105

106
Utility functions for creating frequency arrays.
107

108
```python { .api }
109
def fftfreq(n, d=1.0):
110
    """Return discrete Fourier Transform sample frequencies.
111
    
112
    Args:
113
        n: Window length
114
        d: Sample spacing
115
        
116
    Returns:
117
        cupy.ndarray: Array of frequencies
118
    """
119

120
def rfftfreq(n, d=1.0):
121
    """Return sample frequencies for rfft."""
122

123
def fftshift(x, axes=None):
124
    """Shift zero-frequency component to center.
125
    
126
    Args:
127
        x: Input array
128
        axes: Axes to shift
129
        
130
    Returns:
131
        cupy.ndarray: Shifted array
132
    """
133

134
def ifftshift(x, axes=None):
135
    """Inverse of fftshift."""
136
```
137

138
## Usage Examples
139

140
### Basic 1D FFT Operations
141

142
```python
143
import cupy as cp
144
import matplotlib.pyplot as plt
145

146
# Create a test signal
147
t = cp.linspace(0, 1, 1000, endpoint=False)
148
signal = cp.sin(2 * cp.pi * 50 * t) + cp.sin(2 * cp.pi * 120 * t)
149

150
# Add some noise
151
signal += 0.2 * cp.random.randn(len(t))
152

153
# Compute FFT
154
fft_result = cp.fft.fft(signal)
155
freqs = cp.fft.fftfreq(len(signal), d=t[1]-t[0])
156

157
# Get magnitude spectrum
158
magnitude = cp.abs(fft_result)
159

160
print(f"Signal length: {len(signal)}")
161
print(f"FFT result shape: {fft_result.shape}")
162
print(f"Peak frequencies found at: {freqs[cp.argsort(magnitude)[-3:-1]]}")
163
```
164

165
### Real FFT for Efficiency
166

167
```python
168
# For real-valued signals, use rfft for efficiency
169
real_signal = cp.cos(2 * cp.pi * 10 * t) + cp.cos(2 * cp.pi * 25 * t)
170

171
# Real FFT (more efficient for real inputs)
172
rfft_result = cp.fft.rfft(real_signal)
173
rfreqs = cp.fft.rfftfreq(len(real_signal), d=t[1]-t[0])
174

175
# Reconstruct original signal
176
reconstructed = cp.fft.irfft(rfft_result)
177

178
print(f"Original signal shape: {real_signal.shape}")
179
print(f"RFFT result shape: {rfft_result.shape}")
180
print(f"Reconstruction error: {cp.mean(cp.abs(real_signal - reconstructed)):.2e}")
181
```
182

183
### 2D FFT for Image Processing
184

185
```python
186
# Create a 2D test pattern
187
x = cp.linspace(-5, 5, 256)
188
y = cp.linspace(-5, 5, 256)
189
X, Y = cp.meshgrid(x, y)
190

191
# Create a 2D sinusoidal pattern
192
image = cp.sin(2 * cp.pi * 0.5 * X) * cp.cos(2 * cp.pi * 0.3 * Y)
193

194
# Compute 2D FFT
195
fft2_result = cp.fft.fft2(image)
196

197
# Shift zero frequency to center for visualization
198
fft2_shifted = cp.fft.fftshift(fft2_result)
199

200
# Get magnitude spectrum
201
magnitude_2d = cp.abs(fft2_shifted)
202

203
print(f"Image shape: {image.shape}")
204
print(f"2D FFT shape: {fft2_result.shape}")
205
print(f"Max magnitude: {cp.max(magnitude_2d):.2f}")
206
```
207

208
### Frequency Domain Filtering
209

210
```python
211
# Create noisy signal
212
t = cp.linspace(0, 2, 2000, endpoint=False)
213
clean_signal = cp.sin(2 * cp.pi * 5 * t)  # 5 Hz signal
214
noisy_signal = clean_signal + 0.5 * cp.sin(2 * cp.pi * 50 * t)  # Add 50 Hz noise
215

216
# FFT of noisy signal
217
fft_noisy = cp.fft.fft(noisy_signal)
218
freqs = cp.fft.fftfreq(len(noisy_signal), d=t[1]-t[0])
219

220
# Create low-pass filter (remove frequencies above 10 Hz)
221
filter_mask = cp.abs(freqs) < 10
222
fft_filtered = fft_noisy * filter_mask
223

224
# Inverse FFT to get filtered signal
225
filtered_signal = cp.fft.ifft(fft_filtered).real
226

227
# Compare results
228
original_power = cp.mean(clean_signal**2)
229
noisy_power = cp.mean(noisy_signal**2)
230
filtered_power = cp.mean(filtered_signal**2)
231

232
print(f"Original signal power: {original_power:.3f}")
233
print(f"Noisy signal power: {noisy_power:.3f}")
234
print(f"Filtered signal power: {filtered_power:.3f}")
235
print(f"Filter effectiveness: {cp.mean(cp.abs(clean_signal - filtered_signal)):.4f}")
236
```
237

238
### Convolution using FFT
239

240
```python
241
# FFT-based convolution is efficient for large arrays
242
signal_length = 10000
243
kernel_length = 100
244

245
# Create test signal and kernel
246
signal = cp.random.randn(signal_length)
247
kernel = cp.exp(-cp.linspace(-2, 2, kernel_length)**2)  # Gaussian kernel
248

249
# Pad for proper convolution
250
total_length = signal_length + kernel_length - 1
251
signal_padded = cp.zeros(total_length)
252
signal_padded[:signal_length] = signal
253

254
kernel_padded = cp.zeros(total_length)
255
kernel_padded[:kernel_length] = kernel
256

257
# FFT-based convolution
258
signal_fft = cp.fft.fft(signal_padded)
259
kernel_fft = cp.fft.fft(kernel_padded)
260
convolved_fft = signal_fft * kernel_fft
261
convolved = cp.fft.ifft(convolved_fft).real
262

263
# Compare with direct convolution (for verification on smaller arrays)
264
if signal_length < 1000:
265
    direct_conv = cp.convolve(signal, kernel, mode='full')
266
    fft_conv_trimmed = convolved[:len(direct_conv)]
267
    error = cp.mean(cp.abs(direct_conv - fft_conv_trimmed))
268
    print(f"FFT convolution error: {error:.2e}")
269

270
print(f"Original signal length: {signal_length}")
271
print(f"Convolved signal length: {len(convolved)}")
272
```
273

274
### Spectral Analysis
275

276
```python
277
# Power spectral density estimation
278
fs = 1000  # Sampling frequency
279
t = cp.arange(0, 2, 1/fs)  # 2 seconds of data
280

281
# Create signal with multiple frequency components
282
signal = (cp.sin(2*cp.pi*50*t) + 
283
         0.5*cp.sin(2*cp.pi*120*t) + 
284
         0.2*cp.sin(2*cp.pi*200*t) +
285
         0.1*cp.random.randn(len(t)))
286

287
# Compute power spectral density
288
fft_signal = cp.fft.fft(signal)
289
freqs = cp.fft.fftfreq(len(signal), 1/fs)
290
power_spectrum = cp.abs(fft_signal)**2
291

292
# Find peaks in spectrum
293
positive_freqs = freqs[:len(freqs)//2]
294
positive_power = power_spectrum[:len(power_spectrum)//2]
295

296
# Get indices of largest peaks
297
peak_indices = cp.argsort(positive_power)[-5:]
298
peak_freqs = positive_freqs[peak_indices]
299
peak_powers = positive_power[peak_indices]
300

301
print("Top 5 frequency components:")
302
for freq, power in zip(peak_freqs, peak_powers):
303
    print(f"  {freq:.1f} Hz: {power:.2e}")
304
```
305

306
### Multi-dimensional FFT
307

308
```python
309
# 3D FFT for volumetric data
310
volume_shape = (64, 64, 64)
311
volume_data = cp.random.randn(*volume_shape)
312

313
# Add some structured pattern
314
x, y, z = cp.meshgrid(cp.arange(64), cp.arange(64), cp.arange(64), indexing='ij')
315
pattern = cp.sin(2*cp.pi*x/16) * cp.cos(2*cp.pi*y/16) * cp.sin(2*cp.pi*z/16)
316
volume_data += pattern
317

318
# 3D FFT
319
fft3d = cp.fft.fftn(volume_data)
320

321
# Get magnitude
322
magnitude_3d = cp.abs(fft3d)
323

324
print(f"3D volume shape: {volume_data.shape}")
325
print(f"3D FFT shape: {fft3d.shape}")
326
print(f"Peak magnitude: {cp.max(magnitude_3d):.2f}")
327
print(f"Memory usage: {fft3d.nbytes / 1024**2:.1f} MB")
328
```
329

330
### Performance Optimization
331

332
```python
333
# FFT performance depends on array size
334
# Powers of 2 are typically fastest
335
sizes = [1000, 1024, 2000, 2048, 4000, 4096]
336

337
for size in sizes:
338
    data = cp.random.randn(size)
339
    
340
    # Time the FFT
341
    start = cp.cuda.Event()
342
    end = cp.cuda.Event()
343
    
344
    start.record()
345
    result = cp.fft.fft(data)
346
    end.record()
347
    end.synchronize()
348
    
349
    elapsed = start.elapsed_time(end)
350
    print(f"Size {size:4d}: {elapsed:.2f} ms")
351

352
# For repeated FFTs of the same size, CuPy caches plans
353
data = cp.random.randn(1024)
354
for i in range(5):
355
    start = cp.cuda.Event()
356
    end = cp.cuda.Event()
357
    
358
    start.record()
359
    result = cp.fft.fft(data)
360
    end.record()
361
    end.synchronize()
362
    
363
    elapsed = start.elapsed_time(end)
364
    print(f"Iteration {i+1}: {elapsed:.2f} ms")
365
```