Software: Python, Instrument Control & Data Analysis

← Instruments

Practical notes from my Python measurement and analysis workflows: getting a notebook running, keeping scans organized, plotting every repeat, and processing transport and optical data. Examples use generic paths and explicit units so that they can be adapted to another measurement computer.

1. Start Jupyter in the correct environment

On Windows, open Anaconda Prompt and activate the environment used by the measurement notebook. Replace measurement with the actual environment name shown by conda env list. Start the server from the directory containing the notebook.

conda env list
conda activate measurement
cd /d C:\path\to\notebooks
python -m jupyter notebook --port=8888

Open the local URL printed by Jupyter and select the notebook. A localhost address refers to the computer running the browser; it will not reach a different lab computer by itself. If “jupyter is not recognized” appears, first check environment activation and python -m jupyter --version. Install missing packages into the intended environment, rather than switching interpreters accidentally. In a notebook, check the active kernel with:

import sys
from pathlib import Path
print(sys.executable)
print(Path.cwd())

2. Organize a measurement before plotting

• Keep raw files unchanged. Save figures, fitted parameters, and processed arrays in separate output folders.
• Record axes, units, start time, sweep direction, source settings, averaging, integration time or IF bandwidth, and the script version alongside the data.
• Use pathlib for paths and print the resolved directory when a file seems missing. Relative paths are relative to the current working directory.
• Sort power and repeat indices numerically: rep_10 must follow rep_9, not rep_1. Keep missing repetitions as NaN rather than collapsing the time axis.
• Prefer the frequency axis saved in each file over a reconstructed axis inferred from its filename. Reject inconsistent frequency grids before stacking.

from pathlib import Path
from scipy.io import loadmat
import numpy as np

folder = Path("data/run_001")
print(folder.resolve())
mat = loadmat(folder / "trace.mat")
print([k for k in mat if not k.startswith("__")])
f_Hz = np.asarray(mat["frequency"]).ravel()
s21_dB = np.asarray(mat["S21_dB"]).ravel()
assert f_Hz.shape == s21_dB.shape
assert np.all(np.isfinite(f_Hz)) and np.all(np.diff(f_Hz) > 0)

My repeated VNA files use the fields frequency, S21_dB, and S21_phase, with power and repeat indices in the filename. Confirm the phase unit from the acquisition script before converting it. MATLAB v7.3 files require an HDF5-compatible reader; scipy.io.loadmat is for earlier MAT formats.

3. Every repeat as a time–frequency map

To see drift, switching, or fluctuations, retain every sweep at a fixed power. Use an array with shape (repeat, frequency), then transpose it to put time horizontally and frequency vertically. For the fixed-cadence example below, one repeat takes 5 s. Use recorded timestamps when there are pauses or variable sweep durations.

import matplotlib.pyplot as plt

# Inputs from the loader:
# repeat_indices: complete integer grid, including missing repeats
# map_dB: shape (number_of_repeats, number_of_frequencies)
# f_Hz: frequency axis stored in the files
t_s = 5.0 * (repeat_indices - repeat_indices[0])
assert map_dB.shape == (t_s.size, f_Hz.size)
assert t_s.size >= 2 and f_Hz.size >= 2

cmap = plt.get_cmap("RdBu_r").copy()
cmap.set_bad("lightgray")
fig, ax = plt.subplots(figsize=(12, 3), layout="constrained")
im = ax.pcolormesh(
    t_s, f_Hz / 1e9, np.ma.masked_invalid(map_dB.T),
    shading="nearest", cmap=cmap, rasterized=True
)
ax.set(xlabel="Elapsed time (s)", ylabel="Frequency (GHz)")
fig.colorbar(im, ax=ax, label="S21 magnitude (dB)")
fig.savefig("S21_time_map.png", dpi=300)
plt.show()

The same orientation applies to a phase map. Plot phase in its documented unit; phase wraps are not automatically physical jumps. A diverging color scale is especially useful for changes about a reference. Keep common color limits when comparing datasets. If averaging is needed later, specify whether it averages complex S21, linear magnitude, power, or dB; these operations are different.

4. Lock-in voltage to resistance

Subtract a measured offset in X and Y before taking a magnitude. The current must be the actual AC current through the device, with consistent RMS or peak conventions. With a phase-calibrated, predominantly resistive signal, X/I gives the signed in-phase response; the magnitude discards sign and includes the quadrature component.

import numpy as np

def lockin_response(X_V, Y_V, I_A, X0_V=0.0, Y0_V=0.0):
    if I_A <= 0:
        raise ValueError("Use a positive calibrated AC current amplitude")
    Xc = np.asarray(X_V, float) - X0_V
    Yc = np.asarray(Y_V, float) - Y0_V
    return Xc / I_A, Yc / I_A, np.hypot(Xc, Yc) / I_A

# Supply X_V, Y_V, and offsets measured for this run.
# At 10 nA, 30 microvolts corresponds to 3000 ohms.

Offsets are run-specific calibration values, not universal constants. Preserve both channels and label magnitude as |Z| when a reactive response is possible. Separate warming and cooling branches before averaging. Group data by actual stable temperature plateaus, and retain the spread and number of samples in each group. In a two-terminal measurement, contact and lead contributions remain in the measured response.

5. Temperature maps and resonance tracking

• Align acquisition timestamps with the thermometer log. Establish any time offset from the recorded sequence; do not reuse an old offset automatically.
• Keep warming and cooling sweeps separate. For nonuniform temperature spacing, use pcolormesh with the measured coordinates. Mask missing coverage rather than presenting interpolation as measurement.
• Track each resonance in a local window and inspect linecuts. A minimum can jump to a neighboring mode when resonances overlap or broaden.
• Save frequency, linewidth, fit residuals, uncertainties, bounds, and flags. A vanished or unresolved resonance should not become an apparently precise fitted point.
• Distinguish a fitted effective temperature scale from a material assignment. Hybrid modes and temperature-dependent loss require more information than a frequency shift alone.

6. Optical spectra, background subtraction, and array shapes

My optical workflows use NumPy/SciPy and Matplotlib for reflectance, MCD-related spectra, linecuts, and maps. Keep explicit axis names, for example (magnetic_field, gate_voltage, wavelength). Record the exact processing order: reference conversion, background fitting, reference-spectrum subtraction, normalization, and any smoothing.

• Check each axis length against the corresponding array dimension before indexing. A mismatch between the number of linecuts and field values should be fixed at extraction, not hidden by truncating arrays.
• Report the actual nearest reference coordinate selected by the code. A comment describing an old reference point can become stale when settings change.
• Reject a normalization denominator close to zero. Signed normalization can reverse an entire dataset when the reference is negative.
• Show raw and smoothed linecuts together. Reflection padding reduces FFT edge discontinuities; it does not prove a feature is real.
• When converting wavelength to photon energy, E(eV) ≈ 1239.84/λ(nm), reorder the data columns with the energy axis. Do not reverse labels alone.

7. Instrument-control structure

The control code discussed in my workflow combines PyVISA/QCoDeS, NI-DAQmx, and Windows COM/LabVIEW bridges where needed. A useful separation is connection and configuration, read/set methods, scan scheduling, and plotting. Keep hardware-specific addresses and VI paths in local configuration.

• A two-dimensional snake scan has a slow outer axis and a fast inner axis. Reverse alternate acquired rows when storing them on a common increasing fast-axis grid.
• Record requested and returned setpoints separately. Include settling time and a stop condition; a displayed command is not evidence that motion or a sweep has completed.
• For VISA timeout errors, check instrument state, termination, acquisition duration, and stale responses before merely increasing the timeout. PyVISA timeouts are in milliseconds.
• For a Qt plot panel, remove the old colorbar before replacing the image; otherwise repeated updates can leave duplicate axes. Keep lengthy acquisition work out of the GUI event loop.
• Use simulated inputs to check scan ordering and saving before connecting to hardware.

References

• Jupyter: running the notebook
• Matplotlib: pcolormesh
• PyVISA: resources and timeouts

Updated October 2026. Adapted from my measurement and analysis notes; examples require the input arrays and calibration described above.