vna calibration

.github/workflows/python_publish.yml

@@ -54,7 +54,13 @@ jobs:
       with:
         name: python-package-distributions
         path: dist/
+    - name: Set up Python
+      uses: actions/setup-python@v5
+      with:
+        python-version: "3.x"
     - name: Publish distribution 📦 to PyPI
       uses: pypa/gh-action-pypi-publish@release/v1
       with:
-        password: ${{ secrets.PYPI_API_TOKEN }}
+        password: ${{ secrets.PYPI_API_TOKEN }}
+        verbose: true
+        print-hash: true

charon_vna/vna.py

@@ -1,10 +1,14 @@
 # %% imports
+import copy
 import time
+from pathlib import Path
 from typing import Optional
 
 import adi
+import iio
 import matplotlib as mpl
 import numpy as np
+import skrf as rf
 import xarray as xr
 from matplotlib import pyplot as plt
 from matplotlib.gridspec import GridSpec
@@ -13,6 +17,11 @@ from matplotlib.ticker import EngFormatter
 from numpy import typing as npt
 from scipy import signal
 
+dir_ = Path(__file__).parent
+
+
+# https://wiki.analog.com/resources/tools-software/linux-drivers/iio-transceiver/ad9361-customization
+
 
 # %% helper functions
 def get_config(sdr: adi.ad9361):
@@ -67,6 +76,25 @@ def generate_tone(f: float, N: int = 1024, fs: Optional[float] = None):
 # %% connection
 sdr = adi.ad9361(uri="ip:192.168.3.1")
 
+# %% verify device configuration
+mode_2r2t = bool(sdr._get_iio_debug_attr("adi,2rx-2tx-mode-enable"))
+if not mode_2r2t:
+    raise ValueError("'adi,2rx-2tx-mode-enable' is not set in pluto. See README.md for instructions for changing this")
+# TODO: it might be possible to change this on the fly. I think we'll actually just fail in __init__ for sdr
+
+# %% switch control outputs
+# NOTE: this doesn't appear to work
+sdr._set_iio_debug_attr_str("adi,gpo-manual-mode-enable", "1")
+sdr._get_iio_debug_attr_str("adi,gpo-manual-mode-enable-mask")
+# but direct register access does
+# https://ez.analog.com/linux-software-drivers/f/q-a/120853/control-fmcomms3-s-gpo-with-python
+ctx = iio.Context("ip:192.168.3.1")
+ctrl = ctx.find_device("ad9361-phy")
+# https://www.analog.com/media/cn/technical-documentation/user-guides/ad9364_register_map_reference_manual_ug-672.pdf
+ctrl.reg_write(0x26, 0x90)  # bit 7: AuxDAC Manual, bit 4: GPO Manual
+ctrl.reg_write(0x27, 0x10)  # bits 7-4: GPO3-0
+
+
 # %% initialization
 sdr.rx_lo = int(2.0e9)
 sdr.tx_lo = int(2.0e9)
@@ -80,8 +108,8 @@ 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.rx_hardwaregain_chan0 = 40
+sdr.rx_hardwaregain_chan1 = 40
 sdr.tx_hardwaregain_chan0 = -10
 
 config = get_config(sdr)
@@ -185,12 +213,11 @@ 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"
+reference_sparams = None
+reference_sparams = dir_ / "RBP-135+_Plus25degC.s2p"
 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])
-    )
+    ref = rf.Network(reference_sparams)
+    rbp135 = xr.DataArray(ref.s, dims=["frequency", "m", "n"], coords=dict(frequency=ref.f, 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")
@@ -198,3 +225,108 @@ if reference_sparams is not None:
     axs[1].legend()
 
 plt.show()
+
+
+# %%
+def s2net(s: xr.DataArray) -> rf.Network:
+    net = rf.Network(frequency=s.frequency)
+    net.s = s.data
+    return net
+
+
+# %% SOL calibration
+cal_frequency = np.linspace(70e6, 600e6, 2001)
+ideal_cal_frequency = rf.Frequency(np.min(cal_frequency), np.max(cal_frequency), len(cal_frequency))
+input("Connect SHORT and press ENTER...")
+short = vna_capture(frequency=cal_frequency)
+input("Connect OPEN and press ENTER...")
+open = vna_capture(frequency=cal_frequency)
+input("Connect LOAD and press ENTER...")
+load = vna_capture(frequency=cal_frequency)
+
+short_net = s2net(short)
+open_net = s2net(open)
+load_net = s2net(load)
+
+cal_ideal = rf.media.DefinedGammaZ0(frequency=ideal_cal_frequency)
+calibration = rf.calibration.OnePort(
+    [short_net, open_net, load_net],
+    [cal_ideal.short(), cal_ideal.open(), cal_ideal.load(0)],
+)
+
+
+# %%
+s = vna_capture(frequency=cal_frequency)
+
+# %%
+ham_bands = [
+    [135.7e3, 137.8e3],
+    [472e3, 479e3],
+    [1.8e6, 2e6],
+    [3.5e6, 4e6],
+    [5332e3, 5405e3],
+    [7e6, 7.3e6],
+    [10.1e6, 10.15e6],
+    [14e6, 14.35e6],
+    [18.068e6, 18.168e6],
+    [21e6, 21.45e6],
+    [24.89e6, 24.99e6],
+    [28e6, 29.7e6],
+    [50e6, 54e6],
+    [144e6, 148e6],
+    [219e6, 220e6],
+    [222e6, 225e6],
+    [420e6, 450e6],
+    [902e6, 928e6],
+    [1240e6, 1300e6],
+    [2300e6, 2310e6],
+    [2390e6, 2450e6],
+    [3400e6, 3450e6],
+    [5650e6, 5925e6],
+    [10e9, 10.5e9],
+    [24e9, 24.25e9],
+    [47e9, 47.2e9],
+    [76e9, 81e9],
+    [122.25e9, 123e9],
+    [134e9, 141e9],
+    [241e9, 250e9],
+    [275e9, np.inf],
+]
+
+# %%
+s_calibrated = calibration.apply_cal(s2net(s))
+
+plt.figure()
+s_calibrated.plot_s_smith()
+# ref.plot_s_smith(m=1, n=1)
+plt.show()
+
+plt.figure()
+for start, stop in ham_bands:
+    plt.axvspan(start, stop, alpha=0.1, color="k")
+s_calibrated.plot_s_db()
+# ref.plot_s_db(m=1, n=1)
+plt.gca().xaxis.set_major_formatter(EngFormatter())
+plt.grid(True)
+plt.xlim(s_calibrated.f[0], s_calibrated.f[-1])
+plt.show()
+
+plt.figure()
+for start, stop in ham_bands:
+    plt.axvspan(start, stop, alpha=0.1, color="k")
+# s_calibrated.plot_s_vswr()
+# drop invalid points
+vswr = copy.deepcopy(s_calibrated.s_vswr[:, 0, 0])
+vswr[vswr < 1] = np.nan
+plt.plot(s_calibrated.f, vswr)
+plt.axhline(1, color="k", linestyle="--")
+plt.ylabel("VSWR")
+plt.xlabel("Frequency [Hz]")
+# ref.plot_s_vswr(m=1, n=1)
+plt.gca().xaxis.set_major_formatter(EngFormatter())
+plt.grid(True)
+plt.ylim(0, 10)
+plt.xlim(s_calibrated.f[0], s_calibrated.f[-1])
+plt.show()
+
+# %%

pyproject.toml

@@ -50,3 +50,6 @@ exclude = '''
   | dist
 )/
 '''
+
+[tool.isort]
+profile = "black"