-
Notifications
You must be signed in to change notification settings - Fork 1
/
mysdr.py
62 lines (56 loc) · 2.01 KB
/
mysdr.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
#!/usr/bin/env python
# coding: utf-8
import numpy as np
import scipy.signal
from math import floor
from rtlsdr import RtlSdr
# SDR
def derive(fsps, Tmax):
dt = 1.0/fsps # Time step size between samples
nyquist = fsps / 2.0 # Nyquist frequency
N = round(fsps*Tmax) # Total number of samples to collect
return dt, nyquist, N
def downsample(audio, fsps, faudiosps):
dsf = round(fsps/faudiosps) # Downsampling factor
audio = audio[0:floor(len(audio)/dsf)*dsf]
downsampled_audio = np.mean(audio.reshape((-1, dsf)), axis=1).astype('float32')
return downsampled_audio
def get_bpm(fcutoff, nyquist, N):
fcuttoff_nyq = fcutoff / nyquist
midwidth = round(fcuttoff_nyq*N)
zerowidth = int((N-midwidth)/2)
bpm = np.concatenate((np.zeros(zerowidth),np.ones(midwidth),np.zeros(zerowidth)))
return bpm
def process_samples(samples, fsps, fc, Tmax, fcutoff):
### Initialize
dt, nyquist, N = derive(fsps, Tmax)
bpm = get_bpm(fcutoff, nyquist, N)
f_offset = int(250e3) #NEW
shift = np.exp(1.0j * 2.0 * np.pi * f_offset / fsps * np.arange(N))
### Process
shifted_samples = samples * shift
spectrum = np.fft.fftshift(np.fft.fft(shifted_samples))
filteredspectrum = spectrum * bpm
filteredsignal = np.fft.ifft(np.fft.fftshift(filteredspectrum))
# Get angles
theta_og = np.arctan2(filteredsignal.imag,filteredsignal.real)
# Squelch
abssignal = np.abs(filteredsignal)
mean = np.mean(abssignal)
theta = np.where(abssignal < (mean/3), 0, theta_og)
# Get instantaneous frequency
derivtheta = np.convolve([1,-1],theta,'same')
# Clean wrapping artifacts via median filter
audio = scipy.ndimage.median_filter(derivtheta, size=10)
return audio
# Initialize constants
faudiosps = 40960
fsps = 2*256*256*16 # Sampling frequency
fc = 433.5e6 # Center frequency
Tmax = 5 # Sampling time duration
fcutoff = 25e3 # Cutoff for walkie talkie
f_offset = int(250e3)
sdr = RtlSdr()
sdr.sample_rate = fsps
sdr.center_freq = fc + f_offset
sdr.gain = 42.0