钟形谱异构频率范围离散化算法优化：简化采样实现

核心思路 #

通过非线性分位数映射实现频率轴的非均匀采样，让峰值频率f_peak附近的采样密度最高，完全匹配钟形谱的能量分布特性。核心仅依赖4个必要参数：频率范围起止边界、峰值频率、总采样点数，无需额外冗余参数。

数学逻辑：

1. 在标准正态分布的分位数区间（覆盖99.7%分布范围）上均匀采样
2. 利用逆高斯累积分布函数（CDF）将均匀采样点映射到概率密度更高的区域（对应频率峰值附近）
3. 将映射结果线性缩放至目标频率范围，确保f_peak对应分布中心

Python 实现 #

import numpy as np
from scipy.special import erfinv

def generate_bell_spectrum(f_start, f_end, f_peak, num_points):
    # 生成覆盖99.7%高斯分布的均匀分位数值
    u = np.linspace(0.0015, 0.9985, num_points)
    # 逆高斯CDF映射，得到高密区集中的采样点
    z = np.sqrt(2) * erfinv(2*u - 1)
    # 线性缩放至目标频率范围，对齐峰值
    z_range = z.max() - z.min()
    f = f_peak + (z / z_range) * (f_end - f_start)
    # 生成高斯型钟形谱（可替换为洛伦兹等其他钟形）
    sigma = (f_end - f_start) / 6
    spectrum = np.exp(-(f - f_peak)**2 / (2*sigma**2))
    return f, spectrum

# 示例调用
f_start, f_end, f_peak = 100, 1000, 500
num_points = 50
freqs, spec_vals = generate_bell_spectrum(f_start, f_end, f_peak, num_points)

MATLAB 实现 #

function [f, spectrum] = generate_bell_spectrum(f_start, f_end, f_peak, num_points)
    % 生成覆盖99.7%高斯分布的均匀分位数值
    u = linspace(0.0015, 0.9985, num_points);
    % 逆高斯CDF映射，聚焦峰值附近采样
    z = sqrt(2) * erfcinv(2*(1 - u));
    % 线性缩放至目标频率范围，对齐峰值
    z_range = max(z) - min(z);
    f = f_peak + (z / z_range) * (f_end - f_start);
    % 生成高斯型钟形谱
    sigma = (f_end - f_start) / 6;
    spectrum = exp(-(f - f_peak).^2 / (2*sigma^2));
end

% 示例调用
f_start = 100; f_end = 1000; f_peak = 500; num_points = 50;
[freqs, spec_vals] = generate_bell_spectrum(f_start, f_end, f_peak, num_points);

补充说明 #

• 若需替换钟形谱类型（如洛伦兹谱），只需修改spectrum的计算表达式：
  洛伦兹谱示例：spectrum = 1 / (1 + ((f - f_peak)/sigma)^2)
• 采样点的密度分布可通过调整分位数区间的上下限微调，但默认的0.0015-0.9985已足够覆盖绝大多数场景，无需额外参数

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