本文是自己对傅里叶变换的学习记录。

正交

参考[1]，我们可以将定义域在$[a, b]$上的函数$f(x)$表达为一个$\mathbb{R}^{M}$空间的向量$F$：

$$F=[f(a), f(a+\Delta x), \dots, f(a+(M-1)\Delta x)]$$

当$\Delta x\to 0$时，我们可以得到$F\to f$。

假如我们在$\mathbb{R}^{M}$空间中有一组正交向量基，$\{g_0,\dots,g_{M-1}\}$，则$F$可以表达为：

$$F = a_0 g_0 + \dots + a_{M-1}g_{M-1}$$

那么系数$a_k=\frac{\int_a^b F*g_k}{\int_{a}^b g_k*g_k}$。

而三角函数$\sin m\omega x$，$\cos m\omega x$则是两两正交的一组基底，其中$\omega=2*\pi/(b-a)$

傅里叶变换

从上面的正交中我们可以看到，假如我们可以找到一组基底部，那么前面的系数就是我们该基底向量在整个组合中的权重（频率）。这也就对应了常说的，傅里叶变换是时域和频域的变换。

参考[1]，根据上述找到的三角函数基底，我们可以通过欧拉公式，三角函数性质等等的变换，得到有限定义域，无穷情况下的傅里叶变换公式。

$$\hat{f}(\omega) = \int_{-\infty}^{\infty} f(t)\exp^{-i\omega t} dt$$

离散情况下

$$\hat{x}_k = \frac{1}{M-1}\sum_{s=0}^{M-1} x_s \exp^{-i\frac{2\pi}{M}ks}$$

同样也可以将频域转换为时域，即逆傅里叶变换。

傅里叶变换与卷积

给定两个函数，其卷积

$$\int f(x)*g(t-x) dx$$

等价于这两个函数在频域的乘积

$$\hat{f}(x) * \hat{g}(x)$$

证明参见[2]。

快速傅里叶变换

快速傅里叶变换是加快离散形式的傅里叶变换的方法。

当我们把$x_s$看成其在是多项式在$s$幂的系数时，求得的$\hat x_k$恰好是多项式代入$\exp^{-i\frac{2\pi}{M}k}$的结果。

快速傅里叶变换是利用奇偶性，通过递归的方式求解傅里叶变换值，参考[4]。

例如，对于如下多项式

$$f(x) = a_0+a_1x+\dots+a_7x^7$$

按照奇偶性，写成如下形式

$$\begin{aligned} f(x) &= G(x^2) + xH(x^2)\\ G(x) &= a_0+a_2x+a_4x^2+a_6x^3\\ H(x) & = a_1+a_3x +a_5x^2+a_7x^3\\ \end{aligned}$$

同时根据欧拉公式，当$\omega_M^k=\exp^{-i\frac{2\pi}{M}k}=\cos\frac{k}{M}2\pi +i\sin\frac{k}{M}2\pi$，得到$\omega_{M}^k=-\omega_{M}^{k+\frac{M}{2}}$，同时利用$G(x^2)$和$H(x^2)$均为偶函数，可以得到如下关系：

$$\begin{aligned} f(\omega_{M}^k) &=G(\omega_{M/2}^{k}) + \omega_{M}^k \times H(\omega_{M/2}^k),\\ f(\omega_{M}^{k+M/2}) &=G(\omega_{M/2}^{k}) - \omega_{M}^k \times H(\omega_{M/2}^k). \end{aligned}$$

这样我们知道$G(\omega_{M/2}^k)$和$H(\omega_{M/2}^k)$即可。同时根据$G$和$H$函数的定义，它们也是一个多项式（只不过幂数均减为一半），因此也可以调用快速傅里叶求的，形成递归。

python的代码来源为[3]，

def FFT(x):
    """A recursive implementation of the 1D Cooley-Tukey FFT"""
    x = np.asarray(x, dtype=float)
    N = x.shape[0]
    
    if N % 2 > 0:
        raise ValueError("size of x must be a power of 2")
    elif N <= 32:  # this cutoff should be optimized
        return DFT_slow(x)
    else:
        X_even = FFT(x[::2])
        X_odd = FFT(x[1::2])
        factor = np.exp(-2j * np.pi * np.arange(N) / N)
        return np.concatenate([X_even + factor[:N // 2] * X_odd,
                               X_even + factor[N // 2:] * X_odd])

堆叠成矩阵形式

def FFT_vectorized(x):
    """A vectorized, non-recursive version of the Cooley-Tukey FFT"""
    x = np.asarray(x, dtype=float)
    N = x.shape[0]

    if np.log2(N) % 1 > 0:
        raise ValueError("size of x must be a power of 2")

    # N_min here is equivalent to the stopping condition above,
    # and should be a power of 2
    N_min = min(N, 32)
    
    # Perform an O[N^2] DFT on all length-N_min sub-problems at once
    n = np.arange(N_min)
    k = n[:, None]
    M = np.exp(-2j * np.pi * n * k / N_min)
    X = np.dot(M, x.reshape((N_min, -1)))

    # build-up each level of the recursive calculation all at once
    while X.shape[0] < N:
        X_even = X[:, :X.shape[1] // 2]
        X_odd = X[:, X.shape[1] // 2:]
        factor = np.exp(-1j * np.pi * np.arange(X.shape[0])
                        / X.shape[0])[:, None]
        X = np.vstack([X_even + factor * X_odd,
                       X_even - factor * X_odd])

    return X.ravel()

应用

见例题替换字符后匹配，

class Solution:
    def matchReplacement(self, s: str, sub: str, mappings: List[List[str]]) -> bool:
        n, m = len(s), len(sub)
        ok = set((x, y) for x, y in mappings)
        mp1 = defaultdict(lambda: [0] * n)  # 用于记录每种字符出现的索引
        mp2 = defaultdict(lambda: [1] * m)  # 用于记录哪些索引处的字符不能变成key
        for i, c in enumerate(s):
            mp1[c][i] = 1
        for j, c1 in enumerate(sub):
            for c2 in mp1:
                if c1 == c2 or (c1, c2) in ok:
                    mp2[c2][j] = 0

        # 卷积求出所有不合法的起始位置 i-j
        # 将第二个数组 reverse 之后 i+j 变成 i+(m-j-1)
        bad = [False] * (n - m + 1)
        for key in mp1:
            nums1, nums2 = mp1[key], mp2[key][::-1]
            conv = convolve(nums1, nums2)
            for i in range(m - 1, n):
                bad[i - m + 1] |= conv[i] != 0

        return not all(bad)

        
import numpy as np

def convolve(nums1: Any, nums2: Any) -> np.ndarray:
    """fft求卷积"""
    fftLen = 1
    while 2 * fftLen < len(nums1) + len(nums2) - 1:
        fftLen *= 2
    fftLen *= 2
    Fa = np.fft.rfft(nums1, fftLen)
    Fb = np.fft.rfft(nums2, fftLen)
    Fc = Fa * Fb
    res = np.fft.irfft(Fc, fftLen)
    res = np.rint(res).astype(np.int64)
    return res[: len(nums1) + len(nums2) - 1]

参考

• [1] Ai有爱：从数学底层深度理解信号与系统

• [2] DBinary：卷积与傅里叶变换

• [3] Cory Maklin: Fast Fourier Transform

• [4] Oi wiki:快速傅里叶变换