昇腾CANN ops-fft：在 NPU 上做离散傅里叶变换

信号处理、频谱分析、图像滤波——大量科学计算和工程应用都绕不开 FFT。ops-fft 把离散傅里叶变换（DFT）的计算搬到昇腾 NPU 上，用 Vector 单元的并行能力甩开 CPU 几个量级。

但 FFT 在 NPU 上的坑不比数学算子少。基数选择不当、维度顺序混了、频域到时域的归一化遗漏——每一个都能让结果面目全非。

FFT 的硬件前提

昇腾 NPU 的 Vector 单元一次能处理 256 个 float/fp16 的并行运算。FFT 的核心运算是蝶形运算（Butterfly Operation），每一级需要对 256 个数据点做交叉乘加。如果数据长度刚好是 256 的倍数，Vector 单元能跑满。如果不是，需要 padding 或者分段。

ops-fft 内部用两种基数策略：

• Radix-2：每级蝶形 2 点输入 → 适合 2^n 长度的序列，如 256、512、1024、2048
• Mixed-Radix：混合基数（2/4/8），适合任意长度，如 312、625、1250

Mixed-Radix 的计算效率比 Radix-2 低 10-15%，但覆盖了所有输入长度。

四种算子接口

ops-fft 面向四种高频场景：

// 一维 FFTops_fft_1d(input,length,direction=FFT_FORWARD)// 例：频谱分析// 二维 FFT（逐列做一维）ops_fft_2d(input,height,width,direction=FFT_FORWARD)// 例：2D 图像频域滤波// 任意维度 FFT（通过参数指定维度）ops_fftn(input,shape,axes,direction=FFT_FORWARD)// 例：多维信号联合分析// 卷积定理：乘积 → FFT → 逐点乘 → IFFTops_fft_convolution(signal_x,signal_y)// 例：信号卷积、图像滤波

踩坑一：频率分辨率和采样率的关系

FFT 的输出频率分辨率 Δf 由采样率 Fs 和点数 N 决定：Δf = Fs / N。如果这两个参数配错了，频谱的横轴直接歪了。

错误写法：

importtorch_npuimportnumpyasnp# 采样率 44100 Hz，采集 1024 个点Fs=44100N=1024# 生成一个 440 Hz 的正弦波t=np.arange(N)/Fs signal=np.sin(2*np.pi*440*t).astype(np.float32)# 直接做 FFT（忽略了频率分辨率）spectrum=torch_npu.ops_fft_1d(signal.npu(),N,direction='forward')# 错：频谱横轴的单位是“ bins”，不是 Hz# 峰值应该在 bin 440 / (44100/1024) ≈ bin 10 处# 但代码里直接 print(spectrum)，频率轴对不上实际物理意义

正确写法：显示计算并对齐频率轴。

importtorch_npuimportnumpyasnp# 采样率和点数Fs=44100N=1024# 生成信号t=np.arange(N)/Fs signal=np.sin(2*np.pi*440*t).astype(np.float32)# 做 FFTspectrum=torch_npu.ops_fft_1d(signal.npu(),N,direction='forward')# 计算频率轴（Nyquist 频率之前的正值部分）freq_axis=np.fft.fftfreq(N,d=1.0/Fs)[:N//2]# 频谱只取正值部分（负频率对称，不需要重复）spectrum_positive=spectrum[:,:N//2]# 在 freq_axis 上定位峰值peak_idx=np.argmax(np.abs(spectrum_positive))peak_freq=freq_axis[peak_idx]print(f"实际峰值频率:{peak_freq:.1f}Hz")# 输出 ≈ 440 Hz

C++ 侧原理：ops-fft 输出的频率轴信息以元数据形式存储，不是张量的一部分。需要用np.fft.fftfreq(N, d=1.0/Fs)显式计算横轴。

踩坑二：2D FFT 的维度顺序

二维 FFT 有两种理解方式：先对行做一维 FFT，再对列做一维 FFT（Row-First），或者反过来（Col-First）。 ops-fft 默认是 Row-First——但很多图像处理的代码习惯 Col-First，导致结果差了转置。

错误写��：

importtorch_npu# 图像：height × width = 512 × 512# 想做 2D FFT 做频域滤波image=torch.randn(512,512).npu().half()# 按 ops-fft 默认的 Row-First 做 2D FFT# 先对每一行做 FFT，再对每一列做 FFTfreq_image=torch_npu.ops_fft_2d(image,512,512,direction='forward')# 错误：如果代码后面对 ImageJ 的滤波器做卷积# 预期：滤波是在 frequency domain 做的# 实际：这个 frequency domain 的排列顺序取决于 filter 的生成方式# 不一致时，filter 和 image 的频率轴对不上，过滤波段就错了

正确写法：显式控制维度。

importtorch_npuimporttorch.fft# 方法一：用 PyTorch 的 fft2d（底层调用 ops-fft）freq_image=torch.fft.fft2(image)# 默认：Row-First# 方法二：显式指定 axes，确保一致freq_image=torch.fft.fft2(image,dim=(0,1))# 先 dim=0 再 dim=1# 方法三（如果习惯 Col-First）：手动 transposefreq_image=torch.transpose(torch.fft.fft2(image,dim=0),# 先对列 dim=0 做 fft0,1# 再对行 fft，实际变成了 Column-First)# 生成 filter 时确保维度一致filter=torch.ones_like(freq_image)filter[256-10:256+10,256-10:256+10]=0# 去掉低频（中心）filtered_freq=freq_image*filter# 逐点乘filtered_image=torch.fft.ifft2(filtered_freq)# 逆变换回去

根因：ops-fft 的 Row-First 是在每个ops_fft_2d调用内部固定的。改变维度顺序只能在 Python 侧通过 transpose 控制，或者调用ops_fftn并显式指定 axes 参数。

踩坑三：频域卷积的归一化遗漏

频域卷积的核心是：时域的卷积等于频域的乘积。用 FFT 做卷积的标准流程是：

FFTx(A) × FFTx(B) → 逐点乘 → IFFTx → 除以 N（归一化）

最后一步「除以 N」经常漏掉——因为数学公式里 IFFT 本身带 1/N，但 FFT 的正向和逆向的 scaling 分配有多种标准，如果不统一，结果会差 N 倍。

错误写法：

importtorch_npu# 两个有限脉冲响应：各 128 点h1=torch.randn(128).npu().half()h2=torch.randn(128).npu().half()N=128# FFT 后逐点乘H1=torch_npu.ops_fft_1d(h1,N,direction='forward')H2=torch_npu.ops_fft_1d(h2,N,direction='forward')product=H1*H2# 错误：直接 IFFT，没有归一化# 结果的值比理论值大了 N=128 倍result_no_norm=torch_npu.ops_fft_1d(product,N,direction='inverse')

正确写法：显式归一化。

importtorch_npuimporttorch# 两个信号h1=torch.randn(128).npu().half()h2=torch.randn(128).npu().half()N=128# FFT 后逐点乘H1=torch_npu.ops_fft_1d(h1,N,direction='forward')H2=torch_npu.ops_fft_1d(h2,N,direction='forward')product=H1*H2# IFFT 后，除以 N 归一化result=torch_npu.ops_fft_1d(product,N,direction='inverse')result_normalized=result/float(N)# 验证：和不使用 FFT 的直接卷积结果对比（参考实现）expected=torch.nn.functional.conv1d(h1.unsqueeze(0).unsqueeze(0),# [1, 1, 128]h2.unsqueeze(0).unsqueeze(0),padding='same').squeeze(0)# 误差应该 < 1e-3asserttorch.allclose(result_normalized,expected,atol=1e-3)

ops-fft 的内部 scaling：默认策略是「Forward 1/sqrt(N)，Inverse sqrt(N)」（Parseval 归一化），但这个默认值在某些信号处理教材里不常见。更通用的约定是「Forward 1，Inverse 1/N」——ops-fft 的 scaling 策略可以通过参数选择。

踩坑四：复数输出的视图混淆

FFT 的输出是复数。但在 Python 侧，复数的实部和虚部访问方式有两种：C++ 风格的 (real, imag）交替存储，和 Fortran 风格的连续实部、连续虚部。ops-fft 默认用 C++ style（即 interleaved），但 NumPy 的 fft 输出�� Fortran style（separate real array then imaginary array）。

错误写法：

importtorch_npuimportnumpyasnp signal=torch.randn(256).npu().half()spectrum=torch_npu.ops_fft_1d(signal,256)# 错误：用访问 NumPy 复数数组的方式访问 ops-fft 的输出# NumPy：spectrum.real 和 spectrum.imag 是分开的 viewnp_spectrum=np.fft.fft(signal.cpu().numpy())real_part_np=np_spectrum.real# ← view，直接可用# ops-fft：spectrum 是一个 interleave 存储的 tensor# 直接取 spectrum.real 会出错（没有这个属性）# 应该用 view_as_complex 或分量的 split

正确写法：确保前后端的复数表示一致。

importtorch_npuimporttorch signal=torch.randn(256).npu().half()spectrum=torch_npu.ops_fft_1d(signal,256)# 方法一：转为 PyTorch 的复数格式（推荐）spectrum_complex=torch.view_as_complex(spectrum.reshape(256,2)# [N, 2] → [N])# 现在 spectrum_complex 是 torch.complex64 类型# 幅度谱magnitude=torch.abs(spectrum_complex)# 方法二：拆分成实部和虚部real_part=spectrum[...,0]# 实部交错存储的第一个imag_part=spectrum[...,1]# 虚部# 如果要和 NumPy 结果对比，也需要转 viewnp_spectrum=np.fft.fft(signal.cpu().numpy())torch_np_complex=torch.from_numpy(np_spectrum)# 比对 magnitudeasserttorch.allclose(torch.abs(spectrum_complex.float()),torch.abs(torch_np_complex.float()),atol=1e-2)

根本原因：ops-fft 为了和 Ascend C 的 Vector 单元数据排布兼容，默认用 interleaved 格式（复数的实部和虚部交错存储）。这个设计对 Vector 单元的 align 最友好，但和 Python 生态的 NumPy/PyTorch 不兼容，需要一次 view 转换。

使用性能

在 Ascend 910 上跑 FFT 的吞吐量：

随着点数增加，Vector 单元的并行度优势越大，加速比从 15× 提升到 63×。瓶颈不在计算（在 Vector 单元上蝶形计算本身很快），而在 HBM 的读写带宽——ops-fft 内部尽可能复用 L1 缓存，数据只过一次 HBM。

FFT 的坑集中在「计算外的辅助信息」上——频率轴的物理单位、2D 变换的维度顺序、IFFT 的归一化因子、复数的存储格式。这四个点每一个都是「算子本身没错，但参数配置和环境适配导致结果偏离」的问题。ops-fft 提供的是算子的底层实现，上层的调用姿势需要调用方自己理清楚——这和其他 CANN 算子库的坑的模式是一致的。