企业官网建设流程全解析

用Python可视化高斯函数：从FWHM到标准差的几何直觉

在信号处理和数据分析中，高斯函数就像空气一样无处不在——从图像处理的高斯模糊到激光雷达波形分解，再到金融时间序列的噪声建模。但大多数教材和教程都停留在公式推导和死记硬背的层面，让学习者错失了理解其几何本质的机会。今天我们将用Python的NumPy和Matplotlib，通过动态可视化的方式，让高斯函数的三个关键参数——半高宽(FWHM)、拐点和标准差(σ)——在你眼前"活"起来。

1. 搭建高斯函数可视化实验室

1.1 创建基础高斯曲线

让我们先构建一个标准高斯函数的"画布"。高斯函数的数学表达式为：

import numpy as np import matplotlib.pyplot as plt def gaussian(x, mu=0, sigma=1): """标准高斯函数""" return np.exp(-((x - mu)**2) / (2 * sigma**2)) / (sigma * np.sqrt(2 * np.pi)) # 生成x轴数据点 x = np.linspace(-5, 5, 1000) y = gaussian(x) # 绘制基础曲线 plt.figure(figsize=(10, 6)) plt.plot(x, y, label='标准高斯分布') plt.title('标准高斯函数曲线') plt.xlabel('x') plt.ylabel('概率密度') plt.grid(True)

这段代码会生成一个μ=0、σ=1的标准正态分布曲线。但静态图像还不足以展示关键特征，我们需要交互式标记这些几何特征。

1.2 动态标记峰值和半高宽

半高宽(FWHM)是峰值一半高度处的全宽度。让我们用代码自动计算并标记它：

peak_height = max(y) half_max = peak_height / 2 # 找到半高宽的两个交点 crossings = np.where(np.diff(np.sign(y - half_max)))[0] fwhm = x[crossings[1]] - x[crossings[0]] # 可视化标记 plt.hlines(half_max, x[crossings[0]], x[crossings[1]], colors='r', linestyles='dashed') plt.vlines(x[crossings[0]], 0, half_max, colors='r', linestyles='dashed') plt.vlines(x[crossings[1]], 0, half_max, colors='r', linestyles='dashed') plt.text(0, half_max, f' FWHM = {fwhm:.2f}', va='bottom', ha='center', color='r')

运行后会看到红色虚线标出了FWHM的宽度。有趣的是，对于σ=1的标准高斯函数，FWHM≈2.355，这正是公式FWHM = 2√(2ln2)σ的计算结果。

2. 探索拐点的几何意义

2.1 计算并标记拐点

拐点是曲线凹凸性改变的点，对应高斯函数二阶导数为零的位置。数学上，高斯函数的拐点出现在x=±σ处：

# 计算二阶导数 def gaussian_second_derivative(x, sigma=1): return ((x**2 - sigma**2) / sigma**4) * gaussian(x, sigma=sigma) y_second_deriv = gaussian_second_derivative(x) # 找到拐点 inflection_points = np.where(np.diff(np.sign(y_second_deriv)))[0] x_infl = x[inflection_points] y_infl = y[inflection_points] # 标记拐点 plt.scatter(x_infl, y_infl, color='g', s=100, label='拐点') for xi, yi in zip(x_infl, y_infl): plt.text(xi, yi, f' ({xi:.1f}, {yi:.2f})', ha='left', va='bottom') # 计算拐点横坐标差值的一半 infl_half_diff = (x_infl[1] - x_infl[0]) / 2 plt.title(f'拐点差值一半 = {infl_half_diff:.1f} (等于σ)')

你会看到绿色点标记出的两个拐点正好位于x=-1和x=1处，验证了"拐点横坐标差值的一半等于σ"这一性质。

2.2 FWHM与σ的几何关系可视化

从图中可以直观看到，FWHM的一半(≈1.177)大于σ(=1)。这是因为：

1. 半高宽测量的是峰值一半高度处的宽度
2. 拐点位于曲线开始快速下降的位置
3. 在峰值到拐点之间，曲线下降速度逐渐加快

这种几何关系解释了为什么在激光雷达波形分解中，直接使用FWHM估计σ时需要除以2√(2ln2)的校正因子。

3. 参数变化对曲线形态的影响

3.1 交互式σ调节观察器

让我们创建一个交互式可视化，动态观察σ变化如何影响曲线特征：

from ipywidgets import interact, FloatSlider def plot_gaussian_with_sigma(sigma=1): y = gaussian(x, sigma=sigma) plt.figure(figsize=(10, 6)) plt.plot(x, y) # 计算并标记FWHM peak = max(y) half_max = peak / 2 crossings = np.where(np.diff(np.sign(y - half_max)))[0] fwhm = x[crossings[1]] - x[crossings[0]] plt.hlines(half_max, x[crossings[0]], x[crossings[1]], colors='r', linestyles='dashed') # 计算并标记拐点 y_second_deriv = gaussian_second_derivative(x, sigma) infl_points = np.where(np.diff(np.sign(y_second_deriv)))[0] plt.scatter(x[infl_points], y[infl_points], color='g', s=100) plt.title(f'σ={sigma:.1f}, FWHM={fwhm:.2f}, 拐点差值一半={sigma:.1f}') plt.grid(True) plt.ylim(0, 0.5) plt.show() interact(plot_gaussian_with_sigma, sigma=FloatSlider(min=0.5, max=2, step=0.1, value=1))

拖动滑块时，你会看到：

• σ增大 → 曲线变宽，FWHM增大
• 拐点始终位于x=±σ处
• FWHM与σ的比例关系保持不变

3.2 多σ对比分析

为了更清晰地比较不同σ值的影响，我们可以绘制多条曲线：

sigmas = [0.5, 1.0, 1.5] plt.figure(figsize=(10, 6)) for sigma in sigmas: y = gaussian(x, sigma=sigma) plt.plot(x, y, label=f'σ={sigma}') # 标记拐点 infl_points = np.where(np.diff(np.sign(gaussian_second_derivative(x, sigma))))[0] plt.scatter(x[infl_points], y[infl_points], s=50) plt.title('不同σ值的高斯曲线比较') plt.legend() plt.grid(True)

这个对比清晰地展示了σ如何控制曲线的"胖瘦"程度，而拐点始终忠实地位于±σ处。

4. 从理论到实践：激光雷达波形分解案例

4.1 模拟激光雷达回波波形

在实际激光雷达应用中，回波波形常被建模为多个高斯函数的叠加。让我们模拟一个双峰波形：

def multi_gaussian(x, params): """多高斯函数叠加 params: [(mu1, sigma1, amp1), (mu2, sigma2, amp2), ...] """ y = np.zeros_like(x) for mu, sigma, amp in params: y += amp * np.exp(-(x - mu)**2 / (2 * sigma**2)) return y # 模拟双峰波形 x_lidar = np.linspace(0, 100, 500) params = [(30, 5, 1), (60, 8, 0.7)] y_lidar = multi_gaussian(x_lidar, params) # 添加噪声 noise = np.random.normal(0, 0.02, len(x_lidar)) y_lidar += noise plt.figure(figsize=(12, 6)) plt.plot(x_lidar, y_lidar, label='模拟激光雷达回波') plt.title('含噪声的双峰激光雷达波形') plt.xlabel('时间/ns') plt.ylabel('强度') plt.grid(True)

4.2 波形分解中的参数估计

在实际应用中，我们需要从这样的波形中估计各个高斯成分的参数。虽然完整分解算法复杂，但理解FWHM和σ的关系可以帮助我们：

1. 峰值检测：找到局部最大值点
2. FWHM估计：测量每个峰值的半高宽
3. σ转换：使用σ = FWHM / (2√(2ln2))计算标准差
4. 幅度估计：测量峰值高度

from scipy.signal import find_peaks # 峰值检测 peaks, _ = find_peaks(y_lidar, height=0.2, distance=50) print(f"检测到的峰值位置: {x_lidar[peaks]}") # 对每个峰值分析FWHM for peak in peaks: peak_height = y_lidar[peak] half_max = peak_height / 2 # 找到当前峰值附近的交叉点 left = np.where((y_lidar[:peak] < half_max) & (y_lidar[:peak] > half_max - 0.05))[0] right = np.where((y_lidar[peak:] < half_max) & (y_lidar[peak:] > half_max - 0.05))[0] + peak if len(left) > 0 and len(right) > 0: fwhm = x_lidar[right[0]] - x_lidar[left[-1]] sigma_est = fwhm / (2 * np.sqrt(2 * np.log(2))) plt.hlines(half_max, x_lidar[left[-1]], x_lidar[right[0]], colors='r', linestyles='dashed') plt.text(x_lidar[peak], half_max, f' σ≈{sigma_est:.1f}', va='bottom') plt.scatter(x_lidar[peaks], y_lidar[peaks], color='g', s=100, label='检测峰值') plt.legend()

这个简化示例展示了如何从实际波形中估计高斯参数。在真实场景中，还需要考虑基线校正、噪声滤波和更精确的拟合算法。