diff --git a/pyrelacs/app.py b/pyrelacs/app.py index fa2f849..ba772d4 100644 --- a/pyrelacs/app.py +++ b/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() -