图像傅里叶变换信号振幅随中心距离分布的计算方法

你要实现的是射电天文数据处理中非常常规的径向振幅谱统计，按以下步骤操作即可：

• 预处理：对齐信号零频位置
  标准FFT输出的零频分量默认在矩阵角点，不符合射电领域傅里叶平面零频在中心的惯例，首先对原始复矩阵做fftshift操作，将零频平移到矩阵正中心。如果你的输入矩阵已经完成零频居中，跳过这步即可。
  二维M×N矩阵的中心索引为(M//2, N//2)。
• 生成径向距离矩阵
  构建与信号矩阵同尺寸的坐标网格，逐像素计算到中心索引的L2范数（欧氏距离），公式为对任意像素位置(i,j)，距离值r = sqrt( (i - center_i)² + (j - center_j)² )。
  由于距离是连续浮点数，需要先做离散分箱才能做分组统计：常规处理是直接将r值取整为整数，按单个像素为步长划分径向距离区间；如果需要更高的统计精度，也可以按0.2/0.5像素的自定义步长设置分箱边界。
• 分箱统计对应距离的信号振幅
  遍历每个径向距离分箱，提取箱内所有位置的复信号值，按你定义的规则计算振幅：
  • 单像素振幅直接取复值的模即可
  • 同一分箱内包含多个像素时，直接对箱内所有复值构成的向量求L2范数，就是该距离处的总振幅；如果需要平均振幅，将总振幅除以箱内像素数量即可。
    基于Python+Numpy的最简参考实现如下：

  import numpy as np
  
  def get_radial_amplitude(raw_fft_signal, zero_centered=False):
      # 零频居中处理
      sig = np.fft.fftshift(raw_fft_signal) if not zero_centered else raw_fft_signal
      h, w = sig.shape
      center_y, center_x = h // 2, w // 2
      # 计算逐像素到中心的L2距离
      y_grid, x_grid = np.ogrid[:h, :w]
      dist_map = np.sqrt((y_grid - center_y)**2 + (x_grid - center_x)**2)
      dist_map = dist_map.astype(np.int32)  # 按整数像素分箱
      max_dist = dist_map.max()
      radial_amp = np.zeros(max_dist + 1, dtype=np.float64)
      # 逐半径统计振幅
      for r in range(max_dist + 1):
          bin_mask = (dist_map == r)
          radial_amp[r] = np.linalg.norm(sig[bin_mask])
          # 若需计算平均振幅，替换为下一行
          # radial_amp[r] = np.linalg.norm(sig[bin_mask]) / bin_mask.sum()
      return radial_amp, np.arange(max_dist + 1)
  

• 注意事项
  • 如果你需要物理单位的uv距离而非像素距离，将统计得到的像素半径乘以傅里叶平面的单像素分辨率（单位通常为波长）即可
  • 如果原始信号做过填零、加窗预处理，大半径区间的振幅需要做对应的响应修正，否则结果会存在系统偏差
  • 若矩阵尺寸为奇数，中心位置的计算直接向下取整即可，和射电领域常用软件的处理逻辑一致

内容的提问来源于stack exchange，提问作者haricash