From a4c78ff0cc4c85d75580803ddfc0f0b080370c47 Mon Sep 17 00:00:00 2001
From: Maximilian Gruber <maximilian.gruber@ptb.de>
Date: Fri, 26 Apr 2024 15:09:43 +0200
Subject: [PATCH] wip: enhance output plot

---
 app/cocal_methods.py | 45 ++++++++++++++++++++++++++++++--------------
 app/testing.py       | 26 ++++++++++++++++---------
 2 files changed, 48 insertions(+), 23 deletions(-)

diff --git a/app/cocal_methods.py b/app/cocal_methods.py
index 2c5820e..c0fe483 100644
--- a/app/cocal_methods.py
+++ b/app/cocal_methods.py
@@ -378,18 +378,30 @@ class CocalMethods:
         return tree
 
     def generate_spectral_domain_visualization(self):
+        
+        # provide shortcuts for plotting
+        b = self.transfer_behavior["IIR"]["b"]
+        a = self.transfer_behavior["IIR"]["a"]
+        ba_cov = self.transfer_behavior["IIR"]["Uba"]
+
+        # further relevant parameters
+        time_delta = np.mean(np.diff(self.ref_times[0]))
+        frame_rate = 1.0 / time_delta
+        
+        Nb = len(b)
+        ba = np.concatenate((b, a[1:]))
+        mask = np.logical_and(self.ref_frequency < 5000.0, self.ref_frequency >= 0.0)
 
         # visualize coefficient uncertainties by plotting the spectrum uncertainty via MC
-        draw_samples = lambda size: np.random.multivariate_normal(
-            ba, ba_cov, size=size
-        )
+        def draw_samples(size):
+            return np.random.multivariate_normal(ba, ba_cov, size=size)
 
         w = np.linspace(0, np.pi, 50)
         f = w / (2 * np.pi) * frame_rate
         w_exp = np.exp(-1j * w)  # ???
-        evaluate = lambda sample: complex_2_real_imag(
-            self.freqz_core(sample, Nb, w_exp)
-        )
+
+        def evaluate(sample):
+            return complex_2_real_imag(self.freqz_core(sample, Nb, w_exp))
 
         umc_kwargs = {
             "draw_samples": draw_samples,
@@ -406,32 +418,37 @@ class CocalMethods:
         h_unc = real_imag_2_complex(np.sqrt(np.diag(h_cov)))
 
         # visualize result
-        fig, ax = plt.subplots(2, 1)
-        ax[0].plot(f, np.abs(h))
+        fig, ax = plt.subplots(2, 1, sharex=True)
+        ax[0].plot(f, np.abs(h), label="fitted TF")
         ax[0].fill_between(
             f,
             np.abs(h) - np.abs(h_unc),
             np.abs(h) + np.abs(h_unc),
             alpha=0.3,
-            label="unc",
+            label="unc of fitted TF",
         )
         ax[0].scatter(
             self.ref_frequency[mask],
             np.abs(self.dut_spectrum[mask] / self.ref_spectrum[mask]),
-            label="ref",
+            label="empirical TF",
+            s=1,
         )
+        ax[0].plot(f, np.ones_like(f), "--r", label="ideal")
+        ax[0].plot(f, np.abs(1.0/h), label="compensated empirical")
         ax[0].legend()
+        ax[0].set_xscale("log")
+        ax[0].set_yscale("log")
         ax[1].plot(f, h_unc)
-        plt.savefig(self.result_image_path)
-        #plt.show()
-            
+        #plt.savefig(self.result_image_path)
+        plt.show()
+
 
     def perform_dummy_computations(self):
         self.start_date = datetime.date.today().isoformat()
 
         self.transfer_behavior = {
             "IIR": {
-                "a": np.ones((1,)), 
+                "a": np.ones((2,)), 
                 "b": np.ones((1,)), 
                 "Uba": np.eye(2),
                 }
diff --git a/app/testing.py b/app/testing.py
index fdf737f..985b25e 100644
--- a/app/testing.py
+++ b/app/testing.py
@@ -1,30 +1,38 @@
 # call as `python -m app/testing`
 
 import sys
+
 sys.path.append(".")
 
 from pathlib import Path
 from app import cocal_methods
 
+
 class session:
-    hash = "abc"
+    hash = "ad554ee7119047559a2250e7e9678ee4"
+
 
 cocal_session = session()
 
 cocal = cocal_methods.CocalMethods(cocal_session, generate_plots=True)
 cocal.testsignal_path = Path("./app/audio/Rosa_Rauschen.mp3")
 
-#cocal.record_and_save_reference()
+# cocal.record_and_save_reference()
 
+# cocal.record_and_save_reference(35)
+cocal.ref_paths = [
+    "sessions/ad554ee7119047559a2250e7e9678ee4/reference_1.wav",
+    "sessions/ad554ee7119047559a2250e7e9678ee4/reference_2.wav",
+]
+cocal.dut_path = "sessions/ad554ee7119047559a2250e7e9678ee4/upload_2.wav"
 
-#cocal.record_and_save_reference(35)
-#cocal.ref_paths = ["sessions\\abc\\reference_demo_22.wav", "sessions\\abc\\reference_demo_34.wav"]
-cocal.ref_paths = ["sessions/ad18bf881c964efba09bca831fa7bd45/reference_1.wav", "sessions/ad18bf881c964efba09bca831fa7bd45/reference_3.wav"]
-cocal.dut_path = "sessions/ad18bf881c964efba09bca831fa7bd45/upload.wav"
-
-#cocal.perform_computations()
-cocal.perform_dummy_computations()
+# cocal.perform_computations()
+cocal.load_files()
+cocal.perform_computations()
 
 cocal.generate_dcc()
 
+cocal.generate_spectral_domain_visualization()
+
+
 print("Fin.")
-- 
GitLab