我正在尝试计算 Bremmstrahlung 的频谱，其中涉及计算傅里叶变换的加速度。我正在求解非线性 ODE 以数值计算时域中的加速度。在使用 Numpy 的 fft 进行傅里叶变换后，得到的频谱看起来非常不平滑和“非物理”。我无法粘贴整个代码，所以我发布了我认为相关的片段。有人可以指出我做错了什么吗？

注意：我的加速度是两个变量（beta 和 b，影响参数）的函数，我想为不同的 b 绘制它，为了排除 beta 项，我只是将所有值相加对于给定的冲击参数，不同 beta 的加速度 b。

我的频谱也是傅里叶变换加速度平方乘以常数因子（拉莫尔公式）

'''
Fourier Transformation of the acceleration
'''
N = 2**8
sampling_frequency = 10
a_w_normalized = [[] for a in range(len(impact_parameter))]
a_w =[[] for a in range(len(impact_parameter))]
intensity_normalized = [[] for a in range(len(impact_parameter))]
acc_summed_over_velocity = []
intensity_summed_over_velocity =[]

for index,b in enumerate(impact_parameter):
    acc_sum =[0 for x in range(N)]
    intensity_sum = [0 for x in range(N)] 
    for j in range(len(velocity_z_component)):
        #window_kaiser = signal.kaiser(N, 15)
        #window_hann = signal.hann(N,sym=True)
        window =1
        fft_input = acc_normalized[index][j]*window
        ft_acc_normalized = np.abs(np.fft.fft(fft_input,norm=None))

        acc_sum =np.add(ft_acc_normalized,acc_sum )
        intensity_list = [power_spectrum_factor * (a ** 2) for a in ft_acc_normalized]
        intensity_sum = np.add(intensity_list,intensity_sum)
        intensity_normalized[index].append(intensity_list)
        a_w_normalized[index].append(ft_acc_normalized)
    acc_summed_over_velocity.append(acc_sum)
    intensity_summed_over_velocity.append(intensity_sum *acceleration_factor**2)

#intensity_summed_over_velocity=  ma.masked_less_equal(intensity_summed_over_velocity,1e-5)

'''
Plotting acceleration for selected values of impact paramters in time and frequency domain
'''

if plot is True:
    plt.figure(figsize=(12, 8))
    for i in range(len(impact_parameter)):
        plt.plot(t,acc_timedomain_summed_over_velocity[i], label='b={:.3f}'.format(impact_parameter[i]), )
    plt.legend()
    plt.ylabel(r'$ a(\tilde t) $', fontsize=14)
    plt.xlabel(r'$ \tilde t $', fontsize=14)
    ax.spines['left'].set_position('center')
    ax.spines['bottom'].set_position('center')
    ax.spines['top'].set_visible(False)
    ax.spines['right'].set_visible(False)
    ax.spines['left'].set_smart_bounds(True)
    ax.spines['bottom'].set_smart_bounds(True)
    ax.xaxis.set_ticks_position('bottom')
    ax.yaxis.set_ticks_position('left')
    plt.xticks(fontsize=15)
    plt.yticks(fontsize=15)
    if screening is False:
        plt.title('a(t) vs time without screening',fontsize=15)
        plt.savefig('./Plots/acc_without_screening.png')
    else:
        plt.title('a(t) vs time with screening',fontsize=15)
        plt.savefig('./Plots/acc_with_screening.png')

    '''
    Fourier transformed intensity for different impact paramater
    '''

    plt.figure(figsize=(12, 8))
    freq_normalized = np.fft.fftfreq(N)*(2*np.pi*sampling_frequency)
    for index,b in enumerate(plasma.impact_parameter):
        plt.plot(np.abs(freq_normalized), intensity_summed_over_velocity[index], label=r'$ \tilde b={:.2f}$'.format(b), )
    plt.xlabel(r'$ \tilde \omega $', fontsize=14)
    plt.ylabel(r'$ I_{\omega} $', fontsize=16)
    plt.legend(loc="lower left")
    plt.xscale('log')
    plt.yscale('log')
    plt.xticks(fontsize=12)
    plt.yticks(fontsize=12)

    if screening is False:
        plt.title('Single particle spectrum without screening')
        plt.savefig('./Plots/spectrum_no_screening.png')
    else:
        plt.title('Single particle spectrum with screening')
        plt.savefig('./Plots/spectrum_debye_screening.png')

    plt.show()