diff --git a/charon_vna/vna.py b/charon_vna/vna.py new file mode 100644 index 0000000..9fd6a91 --- /dev/null +++ b/charon_vna/vna.py @@ -0,0 +1,200 @@ +# %% imports +import time +from typing import Optional + +import adi +import matplotlib as mpl +import numpy as np +import xarray as xr +from matplotlib import pyplot as plt +from matplotlib.gridspec import GridSpec +from matplotlib.patches import Circle +from matplotlib.ticker import EngFormatter +from numpy import typing as npt +from scipy import signal + + +# %% helper functions +def get_config(sdr: adi.ad9361): + config = dict() + config["rx_lo"] = sdr.rx_lo + config["rx_rf_bandwidth"] = sdr.rx_rf_bandwidth + config["rx_enabled_channels"] = sdr.rx_enabled_channels + for chan in config["rx_enabled_channels"]: + config[f"rx_hardwaregain_chan{chan}"] = getattr(sdr, f"rx_hardwaregain_chan{chan}") + config[f"gain_control_mode_chan{chan}"] = getattr(sdr, f"gain_control_mode_chan{chan}") + + config["tx_lo"] = sdr.tx_lo + config["tx_rf_bandwidth"] = sdr.tx_rf_bandwidth + config["tx_cyclic_buffer"] = sdr.tx_cyclic_buffer + config["tx_enabled_channels"] = sdr.tx_enabled_channels + for chan in config["tx_enabled_channels"]: + config[f"tx_hardwaregain_chan{chan}"] = getattr(sdr, f"tx_hardwaregain_chan{chan}") + + config["filter"] = sdr.filter + config["sample_rate"] = sdr.sample_rate + config["loopback"] = sdr.loopback + + return config + + +def db10(p): + return 10 * np.log10(np.abs(p)) + + +def db20(p): + return 20 * np.log10(np.abs(p)) + + +def minmax(x): + return (np.min(x), np.max(x)) + + +def generate_tone(f: float, N: int = 1024, fs: Optional[float] = None): + if fs is None: + fs = sdr.sample_rate + fs = int(fs) + fc = int(f / (fs / N)) * (fs / N) + ts = 1 / float(fs) + t = np.arange(0, N * ts, ts) + i = np.cos(2 * np.pi * t * fc) * 2**14 + q = np.sin(2 * np.pi * t * fc) * 2**14 + iq = i + 1j * q + + return iq + + +# %% connection +sdr = adi.ad9361(uri="ip:192.168.3.1") + +# %% initialization +sdr.rx_lo = int(2.0e9) +sdr.tx_lo = int(2.0e9) +sdr.sample_rate = 30e6 +sdr.rx_rf_bandwidth = int(4e6) +sdr.tx_rf_bandwidth = int(4e6) +sdr.rx_destroy_buffer() +sdr.tx_destroy_buffer() +sdr.rx_enabled_channels = [0, 1] +sdr.tx_enabled_channels = [0] +sdr.loopback = 0 +sdr.gain_control_mode_chan0 = "manual" +sdr.gain_control_mode_chan1 = "manual" +sdr.rx_hardwaregain_chan0 = 10 +sdr.rx_hardwaregain_chan1 = 10 +sdr.tx_hardwaregain_chan0 = -10 + +config = get_config(sdr) +config + +# %% +sdr.tx_destroy_buffer() # must destroy buffer before changing cyclicity +sdr.tx_cyclic_buffer = True +sdr.tx(generate_tone(1e6)) + +# %% +sdr.rx_destroy_buffer() +data = sdr.rx() + +# %% Plot in time +fig, axs = plt.subplots(2, 1, sharex=True, tight_layout=True) +axs[0].plot(np.real(data).T) +axs[1].plot(np.imag(data).T) +axs[0].set_ylabel("Real") +axs[1].set_ylabel("Imag") +axs[-1].set_xlabel("Time") +fig.show() + +# %% Plot in frequency +f, Pxx_den = signal.periodogram(data, sdr.sample_rate, axis=-1, return_onesided=False) +plt.figure() +for cc, chan in enumerate(sdr.rx_enabled_channels): + plt.semilogy(f, Pxx_den[cc], label=f"Channel {chan}") +plt.legend() +plt.ylim([1e-7, 1e2]) +plt.xlabel("frequency [Hz]") +plt.ylabel("PSD [V**2/Hz]") +plt.grid(True) +plt.show() + + +# %% TX helper functions +def set_output_power(power: float): + if power == -5: + # FIXME: this is a hack because I don't want to go through re-calibration + tx_gain = -8 + else: + raise NotImplementedError() + # # TODO: correct over frequency + # tx_gain_idx = np.abs(pout.sel(tx_channel=0) - power).argmin(dim="tx_gain") + # tx_gain = pout.coords["tx_gain"][tx_gain_idx] + sdr.tx_hardwaregain_chan0 = float(tx_gain) + + +def set_output(frequency: float, power: float, offset_frequency: float = 1e6): + sdr.tx_destroy_buffer() + set_output_power(power) + sdr.tx_lo = int(frequency - offset_frequency) + offset_frequency = frequency - sdr.tx_lo + sdr.tx_cyclic_buffer = True + sdr.tx(generate_tone(offset_frequency)) + + +# %% +def vna_capture(frequency: npt.ArrayLike): + s = xr.DataArray( + np.empty(len(frequency), dtype=np.complex128), + dims=["frequency"], + coords=dict( + frequency=frequency, + ), + ) + for freq in s.frequency.data: + set_output(frequency=freq, power=-5) + sdr.rx_destroy_buffer() + sdr.rx_lo = int(freq) + sdr.rx_enabled_channels = [0, 1] + sdr.gain_control_mode_chan0 = "manual" + sdr.gain_control_mode_chan1 = "manual" + sdr.rx_hardwaregain_chan0 = 40 + sdr.rx_hardwaregain_chan1 = 40 + rx = sdr.rx() + s.loc[dict(frequency=freq)] = np.mean(rx[1] / rx[0]) + + return s + + +# %% +s = vna_capture(frequency=np.linspace(70e6, 200e6, 101)) + +# %% Plot Logmag +fig, axs = plt.subplots(2, 1, sharex=True, tight_layout=True) + +axs[0].plot(s.frequency, db20(s), label="Measured") +axs[1].plot(s.frequency, np.rad2deg(np.angle((s))), label="Measured") + +axs[0].grid(True) +axs[1].grid(True) + +axs[0].set_ylim(-80, 0) +axs[1].set_ylim(-200, 200) +axs[1].set_xlim(np.min(s.frequency), np.max(s.frequency)) +axs[1].xaxis.set_major_formatter(EngFormatter(places=1)) +axs[1].set_xlabel("Frequency") + +axs[0].set_ylabel("|S11| [dB]") +axs[1].set_ylabel("∠S11 [deg]") + +reference_sparams = "/home/brendan/Documents/projects/bh_instruments/rbp135.npz" +if reference_sparams is not None: + rbp135 = np.load(reference_sparams) + rbp135 = xr.DataArray( + rbp135["s"], dims=["frequency", "m", "n"], coords=dict(frequency=rbp135["frequency"], m=[1, 2], n=[1, 2]) + ) + + axs[0].plot(rbp135.frequency, db20(rbp135.sel(m=1, n=1)), label="Datasheet") + axs[1].plot(rbp135.frequency, np.rad2deg(np.angle(rbp135.sel(m=2, n=1))), label="Datasheet") + axs[0].legend() + axs[1].legend() + +plt.show()