1

我发现了一篇关于如何使用 Python 和 rtlsdr 库从 RTLSDR 获取广播 FM Radio 的优秀且可能过于详细的帖子。 https://witestlab.poly.edu/blog/capture-and-decode-fm-radio/

通过一些修改,我能够直接播放无线电音频剪辑:

#from pylab import psd, xlabel, ylabel, show
import matplotlib.pyplot as plt
from rtlsdr import RtlSdr
from scipy import signal
from time import sleep
import numpy as np
from scipy.io import wavfile
import sounddevice

sdr = RtlSdr()

# configure device
Freq = 91.7e6 #mhz
Fs = 1140000
F_offset = 250000
Fc = Freq - F_offset
sdr.sample_rate = Fs
sdr.center_freq = Fc
sdr.gain = 'auto'
samples = sdr.read_samples(8192000)
print(samples)
print(samples.shape)
print(np.max(samples))
#continue
#sdr.close()
x1 = np.array(samples).astype("complex64")
fc1 = np.exp(-1.0j*2.0*np.pi* F_offset/Fs*np.arange(len(x1))) 
x2 = x1 * fc1
# plt.specgram(x2, NFFT=2048, Fs=Fs)  
# plt.title("x2")  
# plt.xlabel("Time (s)")  
# plt.ylabel("Frequency (Hz)")  
# plt.ylim(-Fs/2, Fs/2)  
# plt.xlim(0,len(x2)/Fs)  
# plt.ticklabel_format(style='plain', axis='y' )  
# plt.show()
bandwidth = 200000#khz broadcast radio.
n_taps = 64
# Use Remez algorithm to design filter coefficients
lpf = signal.remez(n_taps, [0, bandwidth, bandwidth+(Fs/2-bandwidth)/4, Fs/2], [1,0], Hz=Fs)  
x3 = signal.lfilter(lpf, 1.0, x2)
dec_rate = int(Fs / bandwidth)
x4 = x3[0::dec_rate]
Fs_y = Fs/dec_rate
f_bw = 200000  
dec_rate = int(Fs / f_bw)  
x4 = signal.decimate(x2, dec_rate) 
# Calculate the new sampling rate
# Fs_y = Fs/dec_rate  
# plt.scatter(np.real(x4[0:50000]), np.imag(x4[0:50000]), color="red", alpha=0.05)  
# plt.title("x4")  
# plt.xlabel("Real")  
# plt.xlim(-1.1,1.1)  
# plt.ylabel("Imag")  
# plt.ylim(-1.1,1.1)   
# plt.show() 
y5 = x4[1:] * np.conj(x4[:-1])
x5 = np.angle(y5)
# plt.psd(x5, NFFT=2048, Fs=Fs_y, color="blue")  
# plt.title("x5")  
# plt.axvspan(0,             15000,         color="red", alpha=0.2)  
# plt.axvspan(19000-500,     19000+500,     color="green", alpha=0.4)  
# plt.axvspan(19000*2-15000, 19000*2+15000, color="orange", alpha=0.2)  
# plt.axvspan(19000*3-1500,  19000*3+1500,  color="blue", alpha=0.2)  
# plt.ticklabel_format(style='plain', axis='y' )  
# plt.show()


# The de-emphasis filter
# Given a signal 'x5' (in a numpy array) with sampling rate Fs_y
d = Fs_y * 75e-6   # Calculate the # of samples to hit the -3dB point  
x = np.exp(-1/d)   # Calculate the decay between each sample  
b = [1-x]          # Create the filter coefficients  
a = [1,-x]  
x6 = signal.lfilter(b,a,x5)  
audio_freq = 48100.0  
dec_audio = int(Fs_y/audio_freq)  
Fs_audio = Fs_y / dec_audio
x7 = signal.decimate(x6, dec_audio) 

sounddevice.play(x7,Fs_audio) 
x7 *= 10000 / np.max(np.abs(x7))  

#sounddevice.play(x7,Fs_audio)
x7.astype("int16").tofile("wbfm-mono.raw")  #Raw file.
wavfile.write('wav.wav',int(Fs_audio), x7.astype("int16"))
print('playing...')
sounddevice.play(x7.astype("int16"),Fs_audio,blocking=True)
#print('try2')
#Sounds no good
#sounddevice.play(x7, blocking=True)
# use matplotlib to estimate and plot the PSD
#psd(samples, NFFT=1024, Fs=sdr.sample_rate/1e6, Fc=sdr.center_freq/1e6)
#xlabel('Frequency (MHz)')
#ylabel('Relative power (dB)')
#print("showing")
#show()

现在上面的工作就像他们设计的那样,但是当我做了一些改变让它读取一个业余无线电中继器(5khz 宽度而不是 200khz)时,我让它工作了,但有几个问题:

  1. 与在 GQRX 中收听相比,音频具有更高的静态音调
  2. 没有静噪 - 它出现的恒定广播 FM 不需要静噪,当没有人说话时,我将如何处理静噪以不产生响亮的静电?

尝试解码业余无线电频率时,我忘记更改此处是否有任何步骤或处理?

#from pylab import psd, xlabel, ylabel, show
import matplotlib.pyplot as plt
from rtlsdr import RtlSdr
from scipy import signal
from time import sleep
import numpy as np
from scipy.io import wavfile
import sounddevice

sdr = RtlSdr()

# configure device
Freq = 440.712e6 #mhz
Fs = 1140000
F_offset = 2500
Fc = Freq - F_offset
sdr.sample_rate = Fs
sdr.center_freq = Fc
sdr.gain = 'auto'
samples = sdr.read_samples(8192000)
print(samples)
print(samples.shape)
print(np.max(samples))
#continue
#sdr.close()
x1 = np.array(samples).astype("complex64")
fc1 = np.exp(-1.0j*2.0*np.pi* F_offset/Fs*np.arange(len(x1))) 
x2 = x1 * fc1
plt.specgram(x2, NFFT=2048, Fs=Fs)  
plt.title("x2")  
plt.xlabel("Time (s)")  
plt.ylabel("Frequency (Hz)")  
plt.ylim(-Fs/2, Fs/2)  
plt.xlim(0,len(x2)/Fs)  
plt.ticklabel_format(style='plain', axis='y' )  
plt.show()
bandwidth = 2500#khz broadcast radio.
n_taps = 64
# Use Remez algorithm to design filter coefficients
lpf = signal.remez(n_taps, [0, bandwidth, bandwidth+(Fs/2-bandwidth)/4, Fs/2], [1,0], Hz=Fs)  
x3 = signal.lfilter(lpf, 1.0, x2)
dec_rate = int(Fs / bandwidth)
x4 = x3[0::dec_rate]
Fs_y = Fs/dec_rate
f_bw = 200000  
dec_rate = int(Fs / f_bw)  
x4 = signal.decimate(x2, dec_rate) 
# Calculate the new sampling rate
Fs_y = Fs/dec_rate  
# plt.scatter(np.real(x4[0:50000]), np.imag(x4[0:50000]), color="red", alpha=0.05)  
# plt.title("x4")  
# plt.xlabel("Real")  
# plt.xlim(-1.1,1.1)  
# plt.ylabel("Imag")  
# plt.ylim(-1.1,1.1)   
# plt.show() 
y5 = x4[1:] * np.conj(x4[:-1])
x5 = np.angle(y5)
# plt.psd(x5, NFFT=2048, Fs=Fs_y, color="blue")  
# plt.title("x5")  
# plt.axvspan(0,             15000,         color="red", alpha=0.2)  
# plt.axvspan(19000-500,     19000+500,     color="green", alpha=0.4)  
# plt.axvspan(19000*2-15000, 19000*2+15000, color="orange", alpha=0.2)  
# plt.axvspan(19000*3-1500,  19000*3+1500,  color="blue", alpha=0.2)  
# plt.ticklabel_format(style='plain', axis='y' )  
# plt.show()


# The de-emphasis filter
# Given a signal 'x5' (in a numpy array) with sampling rate Fs_y
d = Fs_y * 75e-6   # Calculate the # of samples to hit the -3dB point  
x = np.exp(-1/d)   # Calculate the decay between each sample  
b = [1-x]          # Create the filter coefficients  
a = [1,-x]  
x6 = signal.lfilter(b,a,x5)  
audio_freq = 41000.0  
dec_audio = int(Fs_y/audio_freq)  
Fs_audio = Fs_y / dec_audio
x7 = signal.decimate(x6, dec_audio) 

sounddevice.play(x7,Fs_audio) 
x7 *= 10000 / np.max(np.abs(x7))  

#sounddevice.play(x7,Fs_audio)
x7.astype("int16").tofile("wbfm-mono.raw")  #Raw file.
wavfile.write('wav.wav',int(Fs_audio), x7.astype("int16"))
print('playing...')
sounddevice.play(x7.astype("int16"),Fs_audio,blocking=True)
#print('try2')
#Sounds no good
#sounddevice.play(x7, blocking=True)
# use matplotlib to estimate and plot the PSD
#psd(samples, NFFT=1024, Fs=sdr.sample_rate/1e6, Fc=sdr.center_freq/1e6)
#xlabel('Frequency (MHz)')
#ylabel('Relative power (dB)')
#print("showing")
#show()
4

0 回答 0