下面主要看106k、212k、424k、848k速率下基带信号的波形及其频谱分析
import numpy as np
import matplotlib.pyplot as plt
import scipy.signal as signal
import random
fs=13560000
bit_0=[-1]*128
bit_1=[1]*128
start_bit=bit_0
stop_bit=bit_1
CEB_SOF=[1]*512+[-1]*1300+[1]*250
CEB_EOF=[-1]*1300+[1]*1200
ceb_frame_length=16
ceb_frame_bit=[]
ceb_frame_pattern=[]
for i in range(ceb_frame_length):
ceb_frame_bit.append(random.choice([0,1]))
ceb_frame_pattern=ceb_frame_pattern+CEB_SOF
cnt=0
for i in range(ceb_frame_length):
if cnt == 8:
if ceb_frame_bit[i]==0:
ceb_frame_pattern=ceb_frame_pattern+bit_0
else:
ceb_frame_pattern=ceb_frame_pattern+bit_1
ceb_frame_pattern=ceb_frame_pattern+stop_bit
cnt=0
if cnt == 0:
ceb_frame_pattern=ceb_frame_pattern+start_bit
if ceb_frame_bit[i]==0:
ceb_frame_pattern=ceb_frame_pattern+bit_0
else:
ceb_frame_pattern=ceb_frame_pattern+bit_1
else :
if ceb_frame_bit[i]==0:
ceb_frame_pattern=ceb_frame_pattern+bit_0
else:
ceb_frame_pattern=ceb_frame_pattern+bit_1
cnt=cnt+1
ceb_frame_pattern=ceb_frame_pattern+CEB_EOF
#design a low-pass filter
cutoff_freq=5.4e6
filter_order=21
nyquist_freq=fs/2
cutoff_norm=cutoff_freq/nyquist_freq
taps=signal.firwin(filter_order,cutoff=cutoff_norm,window='hamming')
#106k analysis
# FFT analysis
fft_result=np.fft.fft(np.array(ceb_frame_pattern))
freqs=np.fft.fftfreq(len(fft_result),1/fs)
length_fft=len(fft_result)//2
# find 95% energy frequency distribution
energy_threshold = 0.95 * np.sum(np.abs(fft_result[1:length_fft])**2)
cumulative_energy = 0
freq_range = []
for i in range(len(fft_result)//2):
if i==0:
continue
if cumulative_energy < energy_threshold:
freq_range.append(freqs[i])
cumulative_energy += np.abs(fft_result[i])**2
plt.figure(figsize=(12,6))
plt.subplot(2,1,1)
plt.plot(ceb_frame_pattern,marker='.')
plt.title("CEB 106k bba signal")
plt.subplot(2,1,2)
plt.plot(freqs[:length_fft],np.abs(fft_result)[:length_fft])
plt.axvline(freq_range[0], color='r', linestyle='--', label='95% Energy Range')
plt.axvline(freq_range[-1], color='r', linestyle='--')
plt.title("Fequency Spectrum")
plt.xlabel("Fequency(Hz))")
plt.tight_layout()
#212k analysis
ceb_frame_pattern_dwon_sample=[]
for i in range(len(ceb_frame_pattern)):
if i%2==0:
ceb_frame_pattern_dwon_sample.append(ceb_frame_pattern[i])
# FFT analysis
fft_result=np.fft.fft(np.array(ceb_frame_pattern_dwon_sample))
freqs=np.fft.fftfreq(len(fft_result),1/fs)
length_fft=len(fft_result)//2
# find 95% energy frequency distribution
energy_threshold = 0.95 * np.sum(np.abs(fft_result[1:length_fft])**2)
cumulative_energy = 0
freq_range = []
for i in range(len(fft_result)//2):
if i==0:
continue
if cumulative_energy < energy_threshold:
freq_range.append(freqs[i])
cumulative_energy += np.abs(fft_result[i])**2
plt.figure(figsize=(12,6))
plt.subplot(2,1,1)
plt.plot(ceb_frame_pattern_dwon_sample,marker='.')
plt.title("CEB 212k bba signal")
plt.subplot(2,1,2)
plt.plot(freqs[:length_fft],np.abs(fft_result)[:length_fft])
plt.axvline(freq_range[0], color='r', linestyle='--', label='95% Energy Range')
plt.axvline(freq_range[-1], color='r', linestyle='--')
plt.title("Fequency Spectrum")
plt.xlabel("Fequency(Hz))")
plt.tight_layout()
#424k analysis
ceb_frame_pattern_dwon_sample=[]
for i in range(len(ceb_frame_pattern)):
if i%4==0:
ceb_frame_pattern_dwon_sample.append(ceb_frame_pattern[i])
# FFT analysis
fft_result=np.fft.fft(np.array(ceb_frame_pattern_dwon_sample))
freqs=np.fft.fftfreq(len(fft_result),1/fs)
length_fft=len(fft_result)//2
# find 95% energy frequency distribution
energy_threshold = 0.95 * np.sum(np.abs(fft_result[1:length_fft])**2)
cumulative_energy = 0
freq_range = []
for i in range(len(fft_result)//2):
if i==0:
continue
if cumulative_energy < energy_threshold:
freq_range.append(freqs[i])
cumulative_energy += np.abs(fft_result[i])**2
plt.figure(figsize=(12,6))
plt.subplot(2,1,1)
plt.plot(ceb_frame_pattern_dwon_sample,marker='.')
plt.title("CEB 424k bba signal")
plt.subplot(2,1,2)
plt.plot(freqs[:length_fft],np.abs(fft_result)[:length_fft])
plt.axvline(freq_range[0], color='r', linestyle='--', label='95% Energy Range')
plt.axvline(freq_range[-1], color='r', linestyle='--')
plt.title("Fequency Spectrum")
plt.xlabel("Fequency(Hz))")
plt.tight_layout()
#848k analysis
ceb_frame_pattern_dwon_sample=[]
for i in range(len(ceb_frame_pattern)):
if i%4==0:
ceb_frame_pattern_dwon_sample.append(ceb_frame_pattern[i])
# FFT analysis
fft_result=np.fft.fft(np.array(ceb_frame_pattern_dwon_sample))
freqs=np.fft.fftfreq(len(fft_result),1/fs)
length_fft=len(fft_result)//2
# find 95% energy frequency distribution
energy_threshold = 0.95 * np.sum(np.abs(fft_result[1:length_fft])**2)
cumulative_energy = 0
freq_range = []
for i in range(len(fft_result)//2):
if i==0:
continue
if cumulative_energy < energy_threshold:
freq_range.append(freqs[i])
cumulative_energy += np.abs(fft_result[i])**2
plt.figure(figsize=(12,6))
plt.subplot(2,1,1)
plt.plot(ceb_frame_pattern_dwon_sample,marker='.')
plt.title("CEB 848k bba signal")
plt.subplot(2,1,2)
plt.plot(freqs[:length_fft],np.abs(fft_result)[:length_fft])
plt.axvline(freq_range[0], color='r', linestyle='--', label='95% Energy Range')
plt.axvline(freq_range[-1], color='r', linestyle='--')
plt.title("Fequency Spectrum")
plt.xlabel("Fequency(Hz))")
plt.tight_layout()
plt.show()