方波的Gibbs現(xiàn)象簡(jiǎn)析

今日正文

(1)占空比為50%的方波的傅里葉級(jí)數(shù)展開(kāi)

假設(shè)有一個(gè)方波，周期是T，占空比τ為0.5，幅度為1，如下圖所示。

周期信號(hào)的傅里葉級(jí)數(shù)為：

因此，可以求得，方波的各個(gè)頻率分量所對(duì)應(yīng)的傅里葉系數(shù)Ck，分別為：

也就是說(shuō)，方波用傅里葉級(jí)數(shù)暫開(kāi)后，可以得到：

也就是說(shuō)，上述圖示的方波是由直流和一系列不同幅度的余弦函數(shù)構(gòu)成的。

(2)演示一下，用各個(gè)分量，慢慢疊加，形成方波。

當(dāng)只有直流時(shí)，圖形為：

疊加上頻率為w0的余弦信號(hào)，圖形為：

疊加上頻率為3w0的余弦信號(hào)，圖形為：

疊加上頻率為5w0的余弦信號(hào)，圖形為：

......

疊加上頻率為101w0的余弦信號(hào)，圖形為：

......

疊加上頻率為1001w0的余弦信號(hào)后，圖形為：

(3)吉布斯現(xiàn)象(Gibbs phenomenon)

由上面的疊加圖形可以看到，當(dāng)用余弦波疊加去逼近方波信號(hào)時(shí)，所用的諧波次數(shù)N即使增加到1001后，在不連續(xù)點(diǎn)的附近，仍然會(huì)出現(xiàn)過(guò)沖。

N越大，過(guò)沖的最大值越接近不連續(xù)點(diǎn)，但其峰值并不下降，而是大約等于原函數(shù)在不連續(xù)點(diǎn)處跳變值的9%。

(4) 上面的圖形的Python程序

import numpy as np
import matplotlib.pyplot as plt


def square_wave(T, tau, num_periods):
    """Generate a square wave."""
    t = np.linspace(-T * num_periods / 2, T * num_periods / 2, 100000)
    duty_cycle = tau / T
    waveform = np.zeros_like(t)
    waveform[((t+0.25*T) % T) < (duty_cycle * T )] = 1
    waveform1=np.ones_like(t)*0.5
    return t, waveform,waveform1


def harmonic_component(T, n, amplitude,num_periods):
    """Generate a harmonic component."""
    t = np.linspace(-T * num_periods / 2, T * num_periods / 2, 100000)
    frequency = n / T
    component = amplitude * np.cos(2* np.pi * frequency*t)
    return t, component


def main():
    """Main function."""
    T = 1  # Period
    tau = 0.5 * T  # Pulse width
    num_periods = 5  # Number of periods to plot
    num_harmonics = 1001  # Number of harmonics to include


    # Generate the fundamental square wave
    t, waveform,waveform1 = square_wave(T, tau, num_periods)


    # Plot the fundamental square wave
    plt.figure(figsize=(10, 6))
    # plt.plot(t, waveform, label='Fundamental')


    
    # Generate and add harmonic components
    for n in range(1, num_harmonics + 1):
        _, component = harmonic_component(T, n, np.sinc(n/2),num_periods)
        waveform1 += component
        
    plt.plot(t, waveform1)
    plt.xlabel('Time')
    plt.ylabel('Amplitude')
    plt.title('Square Wave Reconstruction with Harmonic Components')
    plt.legend()
    plt.grid(True)
    plt.show()


if __name__ == "__main__":
    main()