Mock power dataset

import pandas as pd
import numpy as np

import os
import sys
from radiocalibrationtoolkit import *

[INFO] LFmap: Import successful.

# This ensures Plotly output works in multiple places:
# plotly_mimetype: VS Code notebook UI
# notebook: "Jupyter: Export to HTML" command in VS Code
# See https://plotly.com/python/renderers/#multiple-renderers
import plotly.io as pio
pio.renderers.default = "plotly_mimetype+notebook"

piko = 1e-12

# read HW response
hw_file_path = "./antenna_setup_files/HardwareProfileList_realistic.xml"
# hw_file_path = "./antenna_setup_files/HardwareProfileList_flat.xml"

hw_dict = read_hw_file(hw_file_path, interp_args={"fill_value": "extrapolate"})

hw_reponse_1 = hw_dict["RResponse"]["LNA"]
hw_reponse_2 = hw_dict["RResponse"]["digitizer"]
hw_reponse_3 = hw_dict["RResponse"]["cable_fromLNA2digitizer"]
hw_reponse_4 = hw_dict["RResponse"]["impedance_matching_EW"]

# merge all hw responses to one function
def hw_response_func(x):
    return dB2PowerAmp(
        hw_reponse_1(x) + hw_reponse_2(x) + hw_reponse_3(x) + hw_reponse_4(x)
    )


# impedance function
impedance_func = hw_dict["IImpedance"][
    "antenna_EW"
]

# read sidereal voltage square spectral density
sidereal_voltage2_density_DF = pd.read_csv(
    "./voltage2_density/voltage2_density_Salla_EW_GSM16.csv",
    # "./voltage2_density/voltage2_density_isoAnt_GSM16.csv",
    index_col=0,
)
sidereal_voltage2_density_DF.columns = sidereal_voltage2_density_DF.columns.astype(
    float
)

<?xml version="1.0" encoding="iso-8859-1"?>
<Element HardwareProfileList at 0x7fa19afe2300>

def hw_response_func(x):
    return dB2PowerAmp(
        hw_reponse_1(x) + hw_reponse_2(x) + hw_reponse_3(x) + hw_reponse_4(x)
    )

mock_trace_generator = Mock_trace_generator(
    sidereal_voltage2_density_DF=sidereal_voltage2_density_DF,
    hw_response_func=hw_response_func,
    impedance_func=impedance_func,
    voltage2ADC=2048,
    time_trace_size=2048,
    sampling_frequency_MHz=250,
)
freq_MHz_bins = mock_trace_generator.get_frequency_bins()

additional_noise = 5e-4*piko
number_of_traces = 5000
mock_traces_DF = mock_trace_generator.generate_mock_trace(
    number_of_traces,
    # temp_celsius=30,
    additional_noise=additional_noise,
    # nbi={"67.25": 1},
    # nbi_err=0.3,
    rounding=True
)

100%|███████████████████████████████████████| 5000/5000 [03:11<00:00, 26.14it/s]

# system parameters
sampling_frequency_MHz = 250
N = mock_traces_DF.columns[2:].size
ADC2Volts = 1/2048
trace_time_length_sec = N/(sampling_frequency_MHz*1e+6)

# FFT to spectra
spectra_df = time_trace_df_2_spectra_df(mock_traces_DF, DSC=2, sampling_frequency_MHz=250)

The power \(P_s\) within band from \(f\) to \(f+\delta f\) is calculated using the following formula:

\begin{equation} P_s = 2\frac{1}{T} \sum_{k=f}^{f+\delta f} \frac{|X(k)|^2}{R(f)} \Delta f \end{equation}

where \(\Delta t\) is the sampling time, \(T\) is the time trace length, \(N\) number of samples, \(f_s\) sampling frequency and \(R\) the impedance,

and \(X(k)\) is defined as

\begin{equation} X(k) = X(k)_{DFT} \Delta t = \frac{X(k)_{DFT}}{fs} \end{equation}

meaning that \(X(k)_{DFT}\) is direct output using discrete Fourier transform on the time trace. Furthermore, this relations hold:

\begin{equation} \Delta f = \frac{fs}{N} = \frac{1}{N\Delta t} = \frac{1}{T} \end{equation}

# use the formula, create the integrand first
integrand_df = ((spectra_df* ADC2Volts / (sampling_frequency_MHz*1e+6))**2).divide(
    impedance_func(spectra_df.columns.values)
)

# integrate
mock_power_unbinned_DF = (2 / trace_time_length_sec) *integrate_spectral_density(
    integrand_df,
    # integrated_MHz_bands=np.linspace(0, 125, 126),
    integrated_MHz_bands=np.linspace(30, 81, 52),
    integrating_method='on_discontinuous_function',
)

mock_power_unbinned_DF = pd.concat((mock_traces_DF.iloc[:,:2], mock_power_unbinned_DF), axis=1)

mock_power_DF = bin_df_rows(mock_power_unbinned_DF, binning_column='lst', bins=list(range(25)))
mock_power_DF.index.name = 'lst'
mock_power_DF = mock_power_DF.drop(['temp_c', 'lst'], axis=1)

fig = px.imshow(
    mock_power_DF.T * 1e12, width=600, aspect="cube", color_continuous_scale="jet"
)
fig.update_layout(
    title='<b>Mock power dataset: {} mock traces</b>'.format(number_of_traces),
    xaxis=dict(title="<b>LST</b>", tickprefix="<b>", ticksuffix="</b>", dtick=2),
    yaxis=dict(
        title="<b>frequency [MHz]</b>",
        tickprefix="<b>",
        ticksuffix="</b>",
        range=(30, 80),
        tick0=0,
        dtick=10,
        autorange=False,
    ),
    coloraxis=dict(
        colorbar=dict(
            title=dict(
                text="<b>power [pW]</b>",
                side="right",
            ),
        tickprefix="<b>",
        ticksuffix="</b>",
        ),
    ),
    font=dict(
        # family=font,
        size=20,
        color="black",
    ),
)
fig.update_layout(
    coloraxis=dict(colorbar=dict(title=dict(text="<b>Power [pW]</b>", side="right"))
 ,cmin=0, cmax=20))
fig.show()

# mock_power_DF.to_csv('./mock_power_dataset-{}.csv'.format(number_of_traces))
# mock_power_DF.to_csv('./mock_power_dataset-{}_fixed_temp_30C_norounding_flathwgain_isoant.csv'.format(number_of_traces))
# mock_power_DF.to_csv('./mock_power_dataset-{}_fixed_temp_30C_norounding_flathwgain_isoant_extranoise.csv'.format(number_of_traces))
# mock_power_DF.to_csv('./mock_power_dataset-{}_fixed_temp_30C_rounded_flathwgain_isoant.csv'.format(number_of_traces))

# mock_power_DF.to_csv('./mock_power_dataset-{}_fixed_temp_30C_norounding_realhwgain_realant.csv'.format(number_of_traces))
# mock_power_DF.to_csv('./mock_power_dataset-{}_fixed_temp_30C_rounded_realhwgain_realant.csv'.format(number_of_traces))