cleanup, fix SDR buffer overflows

stabilize.py

@@ -1,7 +1,8 @@
 import SoapySDR
 import numpy as np
 from scipy.signal import blackmanharris
-import serial
+
+from stabilizer import InductionHeater
 
 fs = 5e6
 base_freq = 1086e6
@@ -9,84 +10,68 @@ threshold = 0.4
 target_freq = 1088.3e6
 k = 200e-6
 
-induction_min = 430e3
-induction_max = 445e3
-
-generator = serial.Serial("/dev/ttyUSB0", 57600, writeTimeout=0)
-
-def set_induction_freq(freq):
-    command = ":s1f{:010d}\n".format(int(freq*1e2))
-    generator.write(command.encode())
-    generator.readline()
-
-def set_induction(amount):
-    amount = max(min(amount, 0.5), -0.5)
-    freq = (induction_min + induction_max)/2 + amount*(induction_max - induction_min)
-    set_induction_freq(freq)
-
-
 def parabolic(f, x):
-    """Quadratic interpolation for estimating the true position of an
-    inter-sample maximum when nearby samples are known.
-    f is a vector and x is an index for that vector.
-    Returns (vx, vy), the coordinates of the vertex of a parabola that goes
-    through point x and its two neighbors.
-    Example:
-    Defining a vector f with a local maximum at index 3 (= 6), find local
-    maximum if points 2, 3, and 4 actually defined a parabola.
-    In [3]: f = [2, 3, 1, 6, 4, 2, 3, 1]
-    In [4]: parabolic(f, argmax(f))
-    Out[4]: (3.2142857142857144, 6.1607142857142856)
-    """
     xv = 1/2. * (f[x-1] - f[x+1]) / (f[x-1] - 2 * f[x] + f[x+1]) + x
     yv = f[x] - 1/4. * (f[x-1] - f[x+1]) * (xv - x)
     return (xv, yv)
 
+
 def fmt_range(x, m):
     pos = int(round(x*80/m))
     return " "*pos + "*"
 
-sdr = SoapySDR.Device()
-samples = np.array([0]*4096, np.complex64)
 
-sdr.setSampleRate(SoapySDR.SOAPY_SDR_RX, 0, fs)
+def main():
-sdr.setFrequency(SoapySDR.SOAPY_SDR_RX, 0, base_freq)
+    induction = InductionHeater("/dev/ttyUSB0", 430e3, 445e3)
+    induction.start()
+    try:
+        sdr = SoapySDR.Device()
+        samples = np.array([0]*4096, np.complex64)
 
-rxStream = sdr.setupStream(SoapySDR.SOAPY_SDR_RX, SoapySDR.SOAPY_SDR_CF32)
+        sdr.setSampleRate(SoapySDR.SOAPY_SDR_RX, 0, fs)
-sdr.activateStream(rxStream)
+        sdr.setFrequency(SoapySDR.SOAPY_SDR_RX, 0, base_freq)
 
-decimation = 0
+        rxStream = sdr.setupStream(SoapySDR.SOAPY_SDR_RX, SoapySDR.SOAPY_SDR_CF32)
+        sdr.activateStream(rxStream)
 
-try:
+        decimation = 0
-    while True:
-        sr = sdr.readStream(rxStream, [samples], len(samples))
-        if sr.ret != len(samples):
-            print("!")
-            continue
 
-        freqs = np.abs(np.fft.fft(samples*blackmanharris(len(samples)))[0:len(samples)//2])
+        try:
-        for i in range(len(freqs)//100):
+            while True:
-            freqs[i] = 0
+                sr = sdr.readStream(rxStream, [samples], len(samples))
-            freqs[-i] = 0
+                if sr.ret != len(samples):
+                    print("!")
+                    continue
 
-        # https://gist.github.com/endolith/255291
+                freqs = np.abs(np.fft.fft(samples*blackmanharris(len(samples)))[0:len(samples)//2])
-        i = np.argmax(freqs)
+                for i in range(len(freqs)//100):
-        true_i = parabolic(freqs, i)[0]
+                    freqs[i] = 0
-        freq = 0.5 * fs * true_i / len(freqs)
+                    freqs[-i] = 0
-        amplitude = freqs[i]
 
-        if amplitude > threshold:
+                # https://gist.github.com/endolith/255291
-            actual_freq = base_freq + freq
+                i = np.argmax(freqs)
-            print("{:8.3f} {:6.2f} {}".format(actual_freq/1e6, amplitude, fmt_range(freq, fs/2)))
+                true_i = parabolic(freqs, i)[0]
-            induction = (actual_freq - target_freq)*k
+                freq = 0.5 * fs * true_i / len(freqs)
-        else:
+                amplitude = freqs[i]
-            print("<no signal>", amplitude)
-            induction = 0.0
 
-        decimation += 1
+                if amplitude > threshold:
-        if decimation == 15:
+                    actual_freq = base_freq + freq
-            set_induction(induction)
+                    print("{:8.3f} {:6.2f} {}".format(actual_freq/1e6, amplitude, fmt_range(freq, fs/2)))
-            decimation = 0
+                    induction_amt = (actual_freq - target_freq)*k
-finally:
+                else:
-    sdr.deactivateStream(rxStream)
+                    print("<no signal>", amplitude)
-    sdr.closeStream(rxStream)
+                    induction_amt = 0.0
+
+                decimation += 1
+                if decimation == 30:
+                    induction.set(induction_amt)
+                    decimation = 0
+        finally:
+            sdr.deactivateStream(rxStream)
+            sdr.closeStream(rxStream)
+    finally:
+        induction.stop()
+
+
+if __name__ == "__main__":
+    main()

stabilizer.py

@@ -0,0 +1,39 @@
+import serial
+import queue
+import threading
+
+
+class InductionHeater:
+    """Interface to the MHS5200A function generator driving the LC tank"""
+    def __init__(self, port, induction_min, induction_max):
+        self.port = port
+        self.induction_min = induction_min
+        self.induction_max = induction_max
+        self.queue = queue.Queue(1)
+
+    def start(self):
+        self.serial = serial.Serial(self.port, 57600)
+        self.thread = threading.Thread(target=self.thread_target)
+        self.thread.start()
+
+    def thread_target(self):
+        while True:
+            amount = self.queue.get()
+            if amount is None:
+                break
+
+            amount = max(min(amount, 0.5), -0.5)
+            freq = ((self.induction_min + self.induction_max)/2
+                    + amount*(self.induction_max - self.induction_min))
+
+            command = ":s1f{:010d}\n".format(int(freq*1e2))
+            self.serial.write(command.encode())
+            self.serial.readline()
+
+    def set(self, amount):
+        self.queue.put(amount, block=False)
+
+    def stop(self):
+        self.queue.put(None, block=True)
+        self.thread.join()
+        self.serial.close()