From 2851be4709a999b5a4a6b8a5232a656df7faa767 Mon Sep 17 00:00:00 2001
From: Sebastien Bourdeauducq <sb@m-labs.hk>
Date: Mon, 10 Aug 2020 14:29:55 +0800
Subject: [PATCH] simplify/optimize

---
 dmi.py     | 12 ++++++------
 noptica.py | 21 +++++----------------
 2 files changed, 11 insertions(+), 22 deletions(-)

diff --git a/dmi.py b/dmi.py
index edc8f18..e6a99c3 100644
--- a/dmi.py
+++ b/dmi.py
@@ -9,14 +9,14 @@ def main():
     # For this PoC code, a small memory leak is less harmful than random overflows.
     gc.disable()
 
-    freq_sample = 5e6
-    freq_base = 1086e6
+    freq_sample = 1e6
+    freq_base = 1088230e3
     bufsize = 4096
 
     induction = InductionHeater("/dev/ttyUSB0", 350e3, 445e3)
     induction.start()
     try:
-        stabilizer = Stabilizer(freq_sample, 0.4, 1088.1e6 - freq_base, 50e-6, induction)
+        stabilizer = Stabilizer(bufsize, 80.0, 1088.1e6 - freq_base, 10e-6, induction)
         position_tracker = PositionTracker(int(0.1*freq_sample/bufsize))
 
         sdr = SoapySDR.Device()
@@ -35,12 +35,12 @@ def main():
 
                     # We can't update faster than the MHS5200A serial interface
                     stabilizer_throttle += 1
-                    if stabilizer_throttle == 40:
+                    if stabilizer_throttle == 8:
                         stabilizer.input(samples_ref)
                         stabilizer_throttle = 0
 
-                    position, leakage = position_tracker.input(samples_ref, samples_meas)
-                    print(np.sum(position)/len(position), leakage)
+                    #position, leakage = position_tracker.input(samples_ref, samples_meas)
+                    #print(np.sum(position)/len(position), leakage)
                 finally:
                     buf_sdr.dispose(buffers)
         finally:
diff --git a/noptica.py b/noptica.py
index 55d7e8b..7e26c67 100644
--- a/noptica.py
+++ b/noptica.py
@@ -104,30 +104,19 @@ class InductionHeater:
         self.thread.join()
         self.serial.close()
 
-
-# https://gist.github.com/endolith/255291
-def parabolic(f, x):
-    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)
-
-
 class Stabilizer:
-    def __init__(self, freq_sample, amp_threshold, freq_target, k, tuner):
-        self.freq_sample = freq_sample
+    def __init__(self, fft_size, amp_threshold, freq_target, k, tuner):
+        self.freqs = np.fft.fftfreq(fft_size)
         self.amp_threshold = amp_threshold
         self.freq_target = freq_target
         self.k = k
         self.tuner = tuner
 
     def input(self, samples):
-        spectrum = np.abs(np.fft.fft(samples*blackmanharris(len(samples)))[0:len(samples)//2])
-        for i in range(len(spectrum)//100):
-            spectrum[i] = 0
-            spectrum[-i] = 0
+        spectrum = np.abs(np.fft.fft(samples*blackmanharris(len(samples))))
         i = np.argmax(spectrum)
-        true_i, amplitude = parabolic(spectrum, i)
-        freq = 0.5 * self.freq_sample * true_i / len(spectrum)
+        freq = self.freqs[i]
+        amplitude = spectrum[i]
 
         if amplitude > self.amp_threshold:
             tuning = (freq - self.freq_target)*self.k