FFT Operations

def fft(a, n=None, axis=-1, norm=None):
    """
    Compute 1-D discrete Fourier Transform.

    Parameters:
    - a: input array
    - n: length of transformed axis
    - axis: axis over which to compute FFT
    - norm: normalization mode ('backward', 'ortho', 'forward')

    Returns:
    cupy.ndarray: transformed array
    """

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

    Parameters:
    - a: input array
    - n: length of transformed axis
    - axis: axis over which to compute IFFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: inverse transformed array
    """

def fft2(a, s=None, axes=(-2, -1), norm=None):
    """
    Compute 2-D discrete Fourier Transform.

    Parameters:
    - a: input array
    - s: shape of transformed axes
    - axes: axes over which to compute FFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: 2-D transformed array
    """

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

    Parameters:
    - a: input array
    - s: shape of transformed axes
    - axes: axes over which to compute IFFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: 2-D inverse transformed array
    """

def fftn(a, s=None, axes=None, norm=None):
    """
    Compute N-D discrete Fourier Transform.

    Parameters:
    - a: input array
    - s: shape of transformed axes
    - axes: axes over which to compute FFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: N-D transformed array
    """

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

    Parameters:
    - a: input array  
    - s: shape of transformed axes
    - axes: axes over which to compute IFFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: N-D inverse transformed array
    """

def rfft(a, n=None, axis=-1, norm=None):
    """
    Compute 1-D discrete Fourier Transform for real input.

    Parameters:
    - a: input array (real)
    - n: length of transformed axis
    - axis: axis over which to compute FFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: transformed array (complex)
    """

def irfft(a, n=None, axis=-1, norm=None):
    """
    Compute 1-D inverse discrete Fourier Transform for real output.

    Parameters:
    - a: input array (complex)
    - n: length of transformed axis
    - axis: axis over which to compute IFFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: inverse transformed array (real)
    """

def rfft2(a, s=None, axes=(-2, -1), norm=None):
    """
    Compute 2-D discrete Fourier Transform for real input.

    Parameters:
    - a: input array (real)
    - s: shape of transformed axes
    - axes: axes over which to compute FFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: 2-D transformed array (complex)
    """

def irfft2(a, s=None, axes=(-2, -1), norm=None):
    """
    Compute 2-D inverse discrete Fourier Transform for real output.

    Parameters:
    - a: input array (complex)
    - s: shape of transformed axes
    - axes: axes over which to compute IFFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: 2-D inverse transformed array (real)
    """

def rfftn(a, s=None, axes=None, norm=None):
    """
    Compute N-D discrete Fourier Transform for real input.

    Parameters:
    - a: input array (real)
    - s: shape of transformed axes
    - axes: axes over which to compute FFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: N-D transformed array (complex)
    """

def irfftn(a, s=None, axes=None, norm=None):
    """
    Compute N-D inverse discrete Fourier Transform for real output.

    Parameters:
    - a: input array (complex)
    - s: shape of transformed axes  
    - axes: axes over which to compute IFFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: N-D inverse transformed array (real)
    """

def hfft(a, n=None, axis=-1, norm=None):
    """
    Compute 1-D discrete Fourier Transform of Hermitian input.

    Parameters:
    - a: input array (Hermitian symmetric)
    - n: length of transformed axis
    - axis: axis over which to compute FFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: transformed array (real)
    """

def ihfft(a, n=None, axis=-1, norm=None):
    """
    Compute 1-D inverse discrete Fourier Transform to produce Hermitian output.

    Parameters:
    - a: input array (real)
    - n: length of transformed axis
    - axis: axis over which to compute IFFT
    - norm: normalization mode

    Returns:
    cupy.ndarray: inverse transformed array (Hermitian)
    """

def fftfreq(n, d=1.0):
    """
    Return discrete Fourier Transform sample frequencies.

    Parameters:
    - n: window length
    - d: sample spacing

    Returns:
    cupy.ndarray: sample frequencies
    """

def rfftfreq(n, d=1.0):
    """
    Return sample frequencies for real FFT.

    Parameters:
    - n: window length
    - d: sample spacing

    Returns:
    cupy.ndarray: sample frequencies for real FFT
    """

def fftshift(x, axes=None):
    """
    Shift zero-frequency component to center of spectrum.

    Parameters:
    - x: input array
    - axes: axes over which to shift

    Returns:
    cupy.ndarray: shifted array
    """

def ifftshift(x, axes=None):
    """
    Inverse of fftshift.

    Parameters:
    - x: input array
    - axes: axes over which to shift

    Returns:
    cupy.ndarray: shifted array
    """

# FFT configuration module
import cupy.fft.config

# Planning and caching controls available through config module

import cupy as cp
import numpy as np

# Create sample signal
t = cp.linspace(0, 1, 1000)
freq1, freq2 = 5, 20  # Hz
signal = cp.sin(2*cp.pi*freq1*t) + 0.5*cp.sin(2*cp.pi*freq2*t)

# Add noise
noise = 0.1 * cp.random.randn(len(t))
noisy_signal = signal + noise

# Compute FFT
fft_result = cp.fft.fft(noisy_signal)
frequencies = cp.fft.fftfreq(len(t), t[1] - t[0])

# Magnitude spectrum
magnitude = cp.abs(fft_result)

# Find peaks (transfer to CPU for analysis)
magnitude_cpu = cp.asnumpy(magnitude)
freq_cpu = cp.asnumpy(frequencies)

import cupy as cp

# Create or load 2D image data
image = cp.random.random((256, 256))

# Compute 2-D FFT
fft2d = cp.fft.fft2(image)

# Shift zero frequency to center
fft2d_shifted = cp.fft.fftshift(fft2d)

# Magnitude and phase
magnitude = cp.abs(fft2d_shifted)
phase = cp.angle(fft2d_shifted)

# Apply frequency domain filter (low-pass)
rows, cols = image.shape
crow, ccol = rows//2, cols//2
mask = cp.zeros((rows, cols))
mask[crow-30:crow+30, ccol-30:ccol+30] = 1

# Apply mask
filtered_fft = fft2d_shifted * mask

# Inverse transform
filtered_fft_ishifted = cp.fft.ifftshift(filtered_fft)
filtered_image = cp.fft.ifft2(filtered_fft_ishifted)
filtered_image = cp.real(filtered_image)  # Take real part

import cupy as cp

# For real-valued signals, use rfft for efficiency
real_signal = cp.random.randn(1024)

# Real FFT (more efficient than complex FFT)
rfft_result = cp.fft.rfft(real_signal)

# Only positive frequencies are computed
n_pos_freq = len(rfft_result)
print(f"Original length: {len(real_signal)}, FFT length: {n_pos_freq}")

# Recover original signal
recovered_signal = cp.fft.irfft(rfft_result)

# Check reconstruction accuracy
reconstruction_error = cp.linalg.norm(real_signal - recovered_signal)
print(f"Reconstruction error: {reconstruction_error}")

import cupy as cp

# 3-D data (e.g., volumetric data)
volume_data = cp.random.random((64, 64, 64))

# 3-D FFT
fft3d = cp.fft.fftn(volume_data)

# Compute magnitude spectrum
magnitude_spectrum = cp.abs(fft3d)

# FFT along specific axes only
fft_xy = cp.fft.fft2(volume_data, axes=(0, 1))  # FFT in X-Y plane for each Z

# Real 3-D FFT for real input
real_volume = cp.random.randn(32, 32, 32)
rfft3d = cp.fft.rfftn(real_volume)

import cupy as cp

# Create signals
signal = cp.random.randn(1000)
kernel = cp.array([1, -1, 1, -1, 1]) / 5  # Simple kernel

# Pad for full convolution
n = len(signal) + len(kernel) - 1
signal_padded = cp.zeros(n)
signal_padded[:len(signal)] = signal
kernel_padded = cp.zeros(n)
kernel_padded[:len(kernel)] = kernel

# FFT-based convolution
signal_fft = cp.fft.fft(signal_padded)
kernel_fft = cp.fft.fft(kernel_padded)
convolved_fft = signal_fft * kernel_fft
convolved = cp.fft.ifft(convolved_fft)
convolved = cp.real(convolved)  # Take real part

# Compare with direct convolution (on CPU)
signal_cpu = cp.asnumpy(signal)
kernel_cpu = cp.asnumpy(kernel)
direct_conv = np.convolve(signal_cpu, kernel_cpu, mode='full')

import cupy as cp

# Generate test signal with multiple frequencies
fs = 1000  # Sampling frequency
t = cp.arange(0, 1, 1/fs)
f1, f2, f3 = 10, 50, 120  # Frequencies in Hz

signal = (cp.sin(2*cp.pi*f1*t) + 
          0.5*cp.sin(2*cp.pi*f2*t) + 
          0.2*cp.sin(2*cp.pi*f3*t))

# Add noise
signal += 0.1 * cp.random.randn(len(t))

# Windowing (Hann window)
window = cp.hanning(len(signal))
windowed_signal = signal * window

# Compute FFT
fft_result = cp.fft.fft(windowed_signal)
freqs = cp.fft.fftfreq(len(signal), 1/fs)

# Power spectral density
psd = cp.abs(fft_result)**2

# Find positive frequencies only
positive_freq_mask = freqs >= 0
positive_freqs = freqs[positive_freq_mask]
positive_psd = psd[positive_freq_mask]

import cupy as cp

# Create noisy signal
t = cp.linspace(0, 1, 1000)
clean_signal = cp.sin(2*cp.pi*10*t)  # 10 Hz signal
noise = 0.5 * cp.sin(2*cp.pi*60*t)  # 60 Hz noise
noisy_signal = clean_signal + noise

# FFT
signal_fft = cp.fft.fft(noisy_signal)
freqs = cp.fft.fftfreq(len(t), t[1] - t[0])

# Design low-pass filter (cutoff at 20 Hz)
cutoff_freq = 20
filter_mask = cp.abs(freqs) <= cutoff_freq

# Apply filter in frequency domain
filtered_fft = signal_fft * filter_mask

# Inverse FFT to get filtered signal
filtered_signal = cp.fft.ifft(filtered_fft)
filtered_signal = cp.real(filtered_signal)

# Compare signals
print(f"Original SNR: {cp.var(clean_signal) / cp.var(noise)}")
residual_noise = filtered_signal - clean_signal
print(f"Filtered SNR: {cp.var(clean_signal) / cp.var(residual_noise)}")