rusty_psf/fft/
transform.rs

1use ndarray::{Array2, ArrayView2, s};
2use num_complex::Complex64;
3use rustfft::FftPlanner;
4
5/// Perform 2D FFT on a complex-valued array.
6/// 
7/// This function performs a 2D Fast Fourier Transform using row-wise and column-wise
8/// 1D FFTs. The input array is processed in-place to minimize memory usage.
9/// 
10/// # Arguments
11/// * `input` - Complex-valued 2D array to transform
12/// 
13/// # Returns
14/// A new 2D array containing the Fourier transform
15pub fn fft2d(input: &ArrayView2<Complex64>) -> Array2<Complex64> {
16    let (rows, cols) = input.dim();
17    let mut output = input.to_owned();
18    let mut planner = FftPlanner::new();
19    
20    // Row-wise FFT
21    let row_fft = planner.plan_fft_forward(cols);
22    for mut row in output.rows_mut() {
23        // Convert row to Vec to ensure contiguous memory
24        let mut row_vec = row.to_vec();
25        row_fft.process(&mut row_vec);
26        // Copy back to row
27        row.iter_mut().zip(row_vec.iter()).for_each(|(dst, src)| *dst = *src);
28    }
29    
30    // Column-wise FFT
31    let col_fft = planner.plan_fft_forward(rows);
32    for mut col in output.columns_mut() {
33        // Convert column to Vec to ensure contiguous memory
34        let mut col_vec = col.to_vec();
35        col_fft.process(&mut col_vec);
36        // Copy back to column
37        col.iter_mut().zip(col_vec.iter()).for_each(|(dst, src)| *dst = *src);
38    }
39    
40    output
41}
42
43/// Perform 2D inverse FFT on a complex-valued array.
44/// 
45/// This function performs a 2D inverse Fast Fourier Transform using row-wise and column-wise
46/// 1D inverse FFTs. The result is normalized by dividing by the total number of elements.
47/// 
48/// # Arguments
49/// * `input` - Complex-valued 2D array to transform
50/// 
51/// # Returns
52/// A new 2D array containing the inverse Fourier transform
53pub fn ifft2d(input: &ArrayView2<Complex64>) -> Array2<Complex64> {
54    let (rows, cols) = input.dim();
55    let mut output = input.to_owned();
56    let mut planner = FftPlanner::new();
57    
58    // Row-wise IFFT
59    let row_fft = planner.plan_fft_inverse(cols);
60    for mut row in output.rows_mut() {
61        // Convert row to Vec to ensure contiguous memory
62        let mut row_vec = row.to_vec();
63        row_fft.process(&mut row_vec);
64        // Copy back to row
65        row.iter_mut().zip(row_vec.iter()).for_each(|(dst, src)| *dst = *src);
66    }
67    
68    // Column-wise IFFT
69    let col_fft = planner.plan_fft_inverse(rows);
70    for mut col in output.columns_mut() {
71        // Convert column to Vec to ensure contiguous memory
72        let mut col_vec = col.to_vec();
73        col_fft.process(&mut col_vec);
74        // Copy back to column
75        col.iter_mut().zip(col_vec.iter()).for_each(|(dst, src)| *dst = *src);
76    }
77    
78    // Normalize
79    let n = (rows * cols) as f64;
80    output.mapv_inplace(|x| x / n);
81    
82    output
83}
84
85/// Shift zero-frequency component to center of spectrum.
86/// 
87/// This function rearranges the output of an FFT so that the zero-frequency component
88/// is at the center of the spectrum. This is useful for visualization and further processing.
89/// 
90/// # Arguments
91/// * `input` - Complex-valued 2D array to shift
92/// 
93/// # Returns
94/// A new 2D array with the zero-frequency component at the center
95pub fn fftshift(input: &ArrayView2<Complex64>) -> Array2<Complex64> {
96    let (rows, cols) = input.dim();
97    let row_mid = rows / 2;
98    let col_mid = cols / 2;
99    
100    let mut output = Array2::zeros((rows, cols));
101    
102    // Top-left quadrant → Bottom-right
103    output.slice_mut(s![row_mid.., col_mid..])
104        .assign(&input.slice(s![..row_mid, ..col_mid]));
105    
106    // Top-right quadrant → Bottom-left
107    output.slice_mut(s![row_mid.., ..col_mid])
108        .assign(&input.slice(s![..row_mid, col_mid..]));
109    
110    // Bottom-left quadrant → Top-right
111    output.slice_mut(s![..row_mid, col_mid..])
112        .assign(&input.slice(s![row_mid.., ..col_mid]));
113    
114    // Bottom-right quadrant → Top-left
115    output.slice_mut(s![..row_mid, ..col_mid])
116        .assign(&input.slice(s![row_mid.., col_mid..]));
117    
118    output
119}
120
121/// Inverse shift of zero-frequency component from center.
122/// 
123/// This function is the inverse operation of `fftshift`. For even-sized arrays,
124/// this operation is identical to `fftshift`.
125/// 
126/// # Arguments
127/// * `input` - Complex-valued 2D array to shift
128/// 
129/// # Returns
130/// A new 2D array with the zero-frequency component moved from the center
131pub fn ifftshift(input: &ArrayView2<Complex64>) -> Array2<Complex64> {
132    // For even-sized arrays, fftshift and ifftshift are the same operation
133    fftshift(input)
134}