FFT Operations

def fft.fft(a, n=None, axis=-1, norm=None):
    """
    One-dimensional discrete Fourier Transform.
    
    Parameters:
    - a: array-like, input array
    - n: int, length of transformed axis
    - axis: int, axis over which to compute FFT
    - norm: str, normalization mode (None, 'ortho', 'backward', 'forward')
    
    Returns:
    cupy.ndarray, complex-valued FFT result
    """

def fft.ifft(a, n=None, axis=-1, norm=None):
    """
    One-dimensional inverse discrete Fourier Transform.
    
    Parameters:
    - a: array-like, input array  
    - n: int, length of transformed axis
    - axis: int, axis over which to compute IFFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, complex-valued IFFT result
    """

def fft.rfft(a, n=None, axis=-1, norm=None):
    """
    Real-input one-dimensional discrete Fourier Transform.
    
    Parameters:
    - a: array-like, real-valued input array
    - n: int, length of transformed axis
    - axis: int, axis over which to compute FFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, complex-valued FFT result (half spectrum)
    """

def fft.irfft(a, n=None, axis=-1, norm=None):
    """
    Inverse of rfft.
    
    Parameters:
    - a: array-like, complex-valued input array
    - n: int, length of output along transformed axis
    - axis: int, axis over which to compute IFFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, real-valued IFFT result
    """

def fft.hfft(a, n=None, axis=-1, norm=None):
    """
    FFT of Hermitian (conjugate-symmetric) input.
    
    Parameters:
    - a: array-like, Hermitian input array
    - n: int, length of output along transformed axis
    - axis: int, axis over which to compute FFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, real-valued FFT result
    """

def fft.ihfft(a, n=None, axis=-1, norm=None):
    """
    Inverse of hfft.
    
    Parameters:
    - a: array-like, real-valued input array
    - n: int, length of output along transformed axis
    - axis: int, axis over which to compute IFFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, complex-valued Hermitian result
    """

def fft.fft2(a, s=None, axes=(-2, -1), norm=None):
    """
    Two-dimensional discrete Fourier Transform.
    
    Parameters:
    - a: array-like, input array
    - s: tuple of ints, shape along transformed axes
    - axes: tuple of ints, axes over which to compute FFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, complex-valued 2D FFT result
    """

def fft.ifft2(a, s=None, axes=(-2, -1), norm=None):
    """
    Two-dimensional inverse discrete Fourier Transform.
    
    Parameters:
    - a: array-like, input array
    - s: tuple of ints, shape along transformed axes
    - axes: tuple of ints, axes over which to compute IFFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, complex-valued 2D IFFT result
    """

def fft.rfft2(a, s=None, axes=(-2, -1), norm=None):
    """
    Real-input two-dimensional discrete Fourier Transform.
    
    Parameters:
    - a: array-like, real-valued input array
    - s: tuple of ints, shape along transformed axes
    - axes: tuple of ints, axes over which to compute FFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, complex-valued 2D FFT result
    """

def fft.irfft2(a, s=None, axes=(-2, -1), norm=None):
    """
    Inverse of rfft2.
    
    Parameters:
    - a: array-like, complex-valued input array
    - s: tuple of ints, shape along transformed axes
    - axes: tuple of ints, axes over which to compute IFFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, real-valued 2D IFFT result
    """

def fft.fftn(a, s=None, axes=None, norm=None):
    """
    N-dimensional discrete Fourier Transform.
    
    Parameters:
    - a: array-like, input array
    - s: tuple of ints, shape along transformed axes
    - axes: tuple of ints, axes over which to compute FFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, complex-valued N-D FFT result
    """

def fft.ifftn(a, s=None, axes=None, norm=None):
    """
    N-dimensional inverse discrete Fourier Transform.
    
    Parameters:
    - a: array-like, input array
    - s: tuple of ints, shape along transformed axes
    - axes: tuple of ints, axes over which to compute IFFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, complex-valued N-D IFFT result
    """

def fft.rfftn(a, s=None, axes=None, norm=None):
    """
    Real-input N-dimensional discrete Fourier Transform.
    
    Parameters:
    - a: array-like, real-valued input array
    - s: tuple of ints, shape along transformed axes
    - axes: tuple of ints, axes over which to compute FFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, complex-valued N-D FFT result
    """

def fft.irfftn(a, s=None, axes=None, norm=None):
    """
    Inverse of rfftn.
    
    Parameters:
    - a: array-like, complex-valued input array
    - s: tuple of ints, shape along transformed axes
    - axes: tuple of ints, axes over which to compute IFFT
    - norm: str, normalization mode
    
    Returns:
    cupy.ndarray, real-valued N-D IFFT result
    """

def fft.fftfreq(n, d=1.0):
    """
    Return discrete Fourier Transform sample frequencies.
    
    Parameters:
    - n: int, window length
    - d: scalar, sample spacing
    
    Returns:
    cupy.ndarray, array of sample frequencies
    """

def fft.rfftfreq(n, d=1.0):
    """
    Return sample frequencies for rfft.
    
    Parameters:
    - n: int, window length
    - d: scalar, sample spacing
    
    Returns:
    cupy.ndarray, array of sample frequencies for real FFT
    """

def fft.fftshift(x, axes=None):
    """
    Shift zero-frequency component to center of spectrum.
    
    Parameters:
    - x: array-like, input array
    - axes: int or tuple of ints, axes over which to shift
    
    Returns:
    cupy.ndarray, shifted array
    """

def fft.ifftshift(x, axes=None):
    """
    Inverse of fftshift.
    
    Parameters:
    - x: array-like, input array
    - axes: int or tuple of ints, axes over which to shift
    
    Returns:
    cupy.ndarray, shifted array
    """

import cupy.fft.config

# FFT planning and caching configuration
fft.config.enable_nd_planning = True    # Enable N-D FFT planning
fft.config.use_multi_gpus = False       # Multi-GPU FFT support
fft.config.set_cufft_gpus(gpus)         # Set GPUs for cuFFT operations

import cupy as cp
import numpy as np

# Create a test signal (sum of sinusoids)
t = cp.linspace(0, 1, 1000, endpoint=False)
freq1, freq2 = 50, 120  # Hz
signal = cp.sin(2 * cp.pi * freq1 * t) + 0.5 * cp.sin(2 * cp.pi * freq2 * t)

# Add some noise
signal += 0.2 * cp.random.normal(size=t.shape)

# Compute FFT
fft_result = cp.fft.fft(signal)
frequencies = cp.fft.fftfreq(len(signal), d=1/1000)  # Sample rate = 1000 Hz

# Get magnitude spectrum
magnitude = cp.abs(fft_result)
phase = cp.angle(fft_result)

# Find dominant frequencies
dominant_freqs = frequencies[cp.argpartition(magnitude, -10)[-10:]]
print("Dominant frequencies:", cp.asnumpy(dominant_freqs[dominant_freqs >= 0]))

import cupy as cp

# Create a 2D test pattern
x = cp.linspace(-5, 5, 256)
y = cp.linspace(-5, 5, 256)
X, Y = cp.meshgrid(x, y)

# Create image with some periodic pattern
image = cp.exp(-(X**2 + Y**2)/2) * cp.cos(2*cp.pi*X) * cp.sin(2*cp.pi*Y)

# Add noise
image += 0.1 * cp.random.normal(size=image.shape)

# Compute 2D FFT
fft_image = cp.fft.fft2(image)
fft_shifted = cp.fft.fftshift(fft_image)

# Compute power spectrum
power_spectrum = cp.abs(fft_shifted)**2

# Apply low-pass filter in frequency domain
center = cp.array(image.shape) // 2
Y_freq, X_freq = cp.ogrid[:image.shape[0], :image.shape[1]]
distance = cp.sqrt((X_freq - center[1])**2 + (Y_freq - center[0])**2)

# Create circular low-pass filter
cutoff_freq = 30
filter_mask = distance <= cutoff_freq
filtered_fft = fft_shifted * filter_mask

# Inverse FFT to get filtered image
filtered_image = cp.fft.ifft2(cp.fft.ifftshift(filtered_fft))
filtered_image = cp.real(filtered_image)

import cupy as cp

# For real-valued signals, use rfft for efficiency
real_signal = cp.cos(2 * cp.pi * 10 * cp.linspace(0, 1, 1000))

# Real FFT (more memory efficient for real inputs)
rfft_result = cp.fft.rfft(real_signal)
rfreq = cp.fft.rfftfreq(len(real_signal))

# The result contains only positive frequencies
print(f"Original signal length: {len(real_signal)}")
print(f"RFFT result length: {len(rfft_result)}")

# Reconstruct original signal
reconstructed = cp.fft.irfft(rfft_result)
print(f"Reconstruction error: {cp.mean(cp.abs(real_signal - reconstructed)):.2e}")

import cupy as cp

# Cross-correlation using FFT convolution
def cross_correlate_fft(a, b):
    """Compute cross-correlation using FFT."""
    # Pad to avoid circular correlation
    n = len(a) + len(b) - 1
    
    # FFT of both signals
    fft_a = cp.fft.fft(a, n)
    fft_b = cp.fft.fft(b, n)
    
    # Cross-correlation in frequency domain
    cross_corr_fft = fft_a * cp.conj(fft_b)
    
    # Inverse FFT to get cross-correlation
    cross_corr = cp.fft.ifft(cross_corr_fft)
    return cp.real(cross_corr)

# Example signals
signal1 = cp.random.randn(1000)
signal2 = cp.roll(signal1, 100) + 0.1 * cp.random.randn(1000)

# Compute cross-correlation
correlation = cross_correlate_fft(signal1, signal2)
delay = cp.argmax(correlation) - len(signal1) + 1
print(f"Detected delay: {delay} samples")

# Spectral differentiation using FFT
def spectral_derivative(f, dx=1.0):
    """Compute derivative using spectral method."""
    n = len(f)
    k = cp.fft.fftfreq(n, dx) * 2j * cp.pi
    
    fft_f = cp.fft.fft(f)
    fft_df = k * fft_f
    df = cp.fft.ifft(fft_df)
    
    return cp.real(df)

# Test with a known function
x = cp.linspace(0, 2*cp.pi, 100, endpoint=False)
f = cp.sin(x)
df_analytical = cp.cos(x)
df_spectral = spectral_derivative(f, dx=x[1]-x[0])

error = cp.mean(cp.abs(df_analytical - df_spectral))
print(f"Spectral derivative error: {error:.2e}")

import cupy as cp
import time

# Batch processing multiple signals
n_signals = 100
signal_length = 4096
signals = cp.random.randn(n_signals, signal_length)

# Method 1: Process each signal individually
start = time.time()
individual_ffts = []
for i in range(n_signals):
    individual_ffts.append(cp.fft.fft(signals[i]))
individual_time = time.time() - start

# Method 2: Batch processing along axis
start = time.time()
batch_ffts = cp.fft.fft(signals, axis=1)
batch_time = time.time() - start

print(f"Individual processing: {individual_time:.4f} seconds")
print(f"Batch processing: {batch_time:.4f} seconds")
print(f"Speedup: {individual_time/batch_time:.2f}x")

# Use appropriate data types
signal_float32 = cp.random.randn(10000).astype(cp.complex64)
signal_float64 = signal_float32.astype(cp.complex128)

# float32 FFT (faster, less precision)
start = time.time()
fft32 = cp.fft.fft(signal_float32)
time32 = time.time() - start

# float64 FFT (slower, higher precision)
start = time.time()
fft64 = cp.fft.fft(signal_float64)
time64 = time.time() - start

print(f"Complex64 FFT: {time32:.4f} seconds")
print(f"Complex128 FFT: {time64:.4f} seconds")