[app] adding pyqtplot for displaying the beat

pyrelacs/app.py

@@ -1,8 +1,8 @@
+import pathlib
 from PyQt6.QtGui import QAction
 import sys
-import pathlib
 
-from PyQt6.QtCore import QProcess, QSize, QThreadPool, Qt, QSettings
+from PyQt6.QtCore import QSize, QThreadPool, QSettings
 from PyQt6.QtWidgets import (
     QApplication,
     QGridLayout,
@@ -12,9 +12,12 @@ from PyQt6.QtWidgets import (
     QMainWindow,
     QPlainTextEdit,
 )
+import pyqtgraph as pg
 import uldaq
 from IPython import embed
+from scipy.signal import welch, find_peaks
 import numpy as np
+import nixio as nix
 
 from pyrelacs.util.logging import config_logging
 import pyrelacs.info as info
@@ -29,6 +32,7 @@ class PyRelacs(QMainWindow):
         super().__init__()
         self.setWindowTitle("PyRelacs")
         self.setMinimumSize(1000, 1000)
+        self.plot_graph = pg.PlotWidget()
 
         self.threadpool = QThreadPool()
         self.repros = Repro()
@@ -41,13 +45,18 @@ class PyRelacs(QMainWindow):
         self.daq_disconnect_button.setCheckable(True)
         self.daq_disconnect_button.clicked.connect(self.disconnect_dac)
 
+        self.plot_calibration_button = QPushButton("Plot Calibration")
+        self.plot_calibration_button.setCheckable(True)
+        self.plot_calibration_button.clicked.connect(self.plot_calibration)
+
         self.text = QPlainTextEdit()
         self.text.setReadOnly(True)
 
         layout = QGridLayout()
-        layout.addWidget(self.daq_connect_button, 0, 0)
+        layout.addWidget(self.plot_calibration_button, 0, 0)
         layout.addWidget(self.daq_disconnect_button, 0, 1)
-        layout.addWidget(self.text, 2, 0, 1, 2)
+        layout.addWidget(self.text, 3, 0, 1, 2)
+        layout.addWidget(self.plot_graph, 2, 0, 1, 2)
 
         self.toolbar = QToolBar("Repros")
         self.addToolBar(self.toolbar)
@@ -58,6 +67,67 @@ class PyRelacs(QMainWindow):
         widget.setLayout(layout)
         self.setCentralWidget(widget)
 
+        self.nix_file = nix.File.open(
+            str(pathlib.Path(__file__).parent / "data"), nix.FileMode.ReadOnly
+        )
+
+    def plot_calibration(self):
+        def decibel(power, ref_power=1.0, min_power=1e-20):
+            """Transform power to decibel relative to ref_power.
+
+            \\[ decibel = 10 \\cdot \\log_{10}(power/ref\\_power) \\]
+            Power values smaller than `min_power` are set to `-np.inf`.
+
+            Parameters
+            ----------
+            power: float or array
+                Power values, for example from a power spectrum or spectrogram.
+            ref_power: float or None or 'peak'
+                Reference power for computing decibel.
+                If set to `None` or 'peak', the maximum power is used.
+            min_power: float
+                Power values smaller than `min_power` are set to `-np.inf`.
+
+            Returns
+            -------
+            decibel_psd: array
+                Power values in decibel relative to `ref_power`.
+            """
+            if np.isscalar(power):
+                tmp_power = np.array([power])
+                decibel_psd = np.array([power])
+            else:
+                tmp_power = power
+                decibel_psd = power.copy()
+            if ref_power is None or ref_power == "peak":
+                ref_power = np.max(decibel_psd)
+            decibel_psd[tmp_power <= min_power] = float("-inf")
+            decibel_psd[tmp_power > min_power] = 10.0 * np.log10(
+                decibel_psd[tmp_power > min_power] / ref_power
+            )
+            if np.isscalar(power):
+                return decibel_psd[0]
+            else:
+                return decibel_psd
+
+        block = self.nix_file.blocks[0]
+        for stim, fish in zip(
+            list(block.data_arrays)[::2], list(block.data_arrays)[1::2]
+        ):
+            beat = stim[:] + fish[:]
+            beat_squared = beat**2
+
+            f, powerspec = welch(beat, fs=40_000.0)
+            powerspec = decibel(powerspec)
+
+            f_sq, powerspec_sq = welch(beat_squared, fs=40_000.0)
+            powerspec_sq = decibel(powerspec_sq)
+            peaks = find_peaks(powerspec_sq, prominence=20)[0]
+            pen = pg.mkPen()
+            self.plot_graph.plot(
+                np.arange(0, len(beat)) / 40_000.0, beat_squared, pen=pen
+            )
+
     def connect_dac(self):
         devices = uldaq.get_daq_device_inventory(uldaq.InterfaceType.USB)
         try:
@@ -93,7 +163,7 @@ class PyRelacs(QMainWindow):
 
     def run_repro(self, n, fn):
         self.text.appendPlainText(f"started Repro {n}, {fn}")
-        worker = Worker(self.repros.run_repro, n, fn)
+        worker = Worker(self.repros.run_repro, self.nix_file, n, fn)
         worker.signals.result.connect(self.print_output)
         worker.signals.finished.connect(self.thread_complete)
         worker.signals.progress.connect(self.progress_fn)
@@ -142,4 +212,3 @@ def main():
 
 if __name__ == "__main__":
     main()
-