plot.py

import numpy as np
import matplotlib.pyplot as plt

def make_debug_plots(
        times, data, sampling_rate,
        save_fig : str = './debug_output.pdf'):
    """
    Make some debug plots
    """
    blocks = data.reshape(-1, 2000, 2)
    std = np.std(blocks, axis=1)
    abs_spectra = np.abs(np.fft.rfft(blocks, axis=1))
    mean_spectra = np.mean(abs_spectra, axis=0)
    spectrum_quantiles = np.quantile(abs_spectra, [0.05, 0.16, 0.84, 0.95], axis=0)

    freqs = np.fft.rfftfreq(2000, 1/sampling_rate)

    fig, axs = plt.subplots(
        2, 3, sharex='col', layout='constrained', figsize=(12, 6),
        width_ratios=(2, 1, 2)
        )

    for i in range(2): # both polarizations
        axs[i, 0].plot(times, std[:, i], ls='', marker='.', alpha=.2)
        axs[i, 1].hist(std[:, i], bins=max(int(len(std) / 100), 10), orientation='horizontal')
        axs[i, 2].plot(freqs, mean_spectra[:,i])
        for quantile in spectrum_quantiles:
            axs[i, 2].plot(freqs, quantile[:,i], color='grey', lw=.5)

        axs[i, 0].set_ylabel('XY'[i] + ' STD [ADC]')
        axs[i, 2].set_ylabel('XY'[i] + ' Spectrum')

    axs[-1, 0].set_xlabel('Time [ns]')
    axs[-1, 1].set_xlabel('N')
    axs[-1, 2].set_xlabel('Frequency [GHz]')

    plt.savefig(save_fig)
    plt.close()