Python

The Red Pitaya boards offer an SCPI server over an TCP/IP Socket connection. The makers describe how to use it. But instead of using plain pyVISA, they provide their own SCPI class.

That’s fine, because that class also provides handy functions to set the various in-built applications (signal generator and the likes).

But it is unnecessary complicated for a blinky example. And in my case, where I only needed some scriptable DIOs, it was quite cumbersome.

So, here is the blinky re-written in plain pyVISA:

import pyvisa as visa
from time import sleep

rm = visa.ResourceManager()
rp = rm.open_resource("TCPIP::169.254.XXX.XXX::5000::SOCKET",
                 read_termination="\r\n",
                 write_termination="\r\n"
                 )

print(rp.query("*IDN?"))

while True:
    rp.write("DIG:PIN LED0,1")
    sleep(.5)
    rp.write("DIG:PIN LED0,0")
    sleep(.5)

The magic lies in the read and write terminations. They have to be set to '\r\n'(in that order), or else the communication simply won’t work and time out.

Make sure you install a reasonably recent pyVISA and pyVISA-py (from pip) or libvisa (from your distro’s repository) before you start. For me (Ubuntu) this works as follows:

pip install -U pyvisa pyvisa-py
sudo apt install libvisa

This integrates nicely with existing instrument command structures and allows for quick testing.

SymPy allows quick and (relatively 😜) easy manipulation of algebraic expressions. But it can do so much more! Using the sympy.physics.control toolbox, you can very inspect linear time-invariant systems in a very comfortable way.

I was looking at a paper (Janiszowski, 1993) which features three transfer functions of different properties and wanted to look at their pole-zero diagrams as well as their Bode plots in real frequency. The transfer functions where the following:

(1) $\begin{align*} G_1(s) &= \frac{ 2 + 42s }{ (1 + 2s)(1 + 40s) } \\ G_2(s) &= \frac{ 5 - 60s }{ (1 + 4s)(1 + 40s) } \\ G_3(s) &= \frac{ 4(1 + s) }{ 1 + 4s + 8s^2 } \end{align*}$

Drawing the pole-zero diagrams from this representation is easy:

Zeros: Set numerator to zero
Poles: Set denominator (or its individual terms) to zero

Drawing the Bode diagram however would’ve involved some basic programming. But neither is necessary.

import sympy as sy
from sympy.physics.control.lti import TransferFunction as TF
from sympy.physics.control import bode_plot, pole_zero_plot

s = sy.symbols("s", complex=True)

G1 = (2 + 42*s)/((1 + 2*s)*(1 + 40*s))
tf = TF(*G.as_numer_denom(), s)
bode_plot(tf)
pole_zero_plot(tf)

Bode plot of the transfer function G₁(s)

Pole-Zero diagram of the transfer function G₁(s)

Well, that was easy. TransferFunction (abbreviated with TF in my examples) requires numerator and denominator to be passed as separate arguments. sympy_expr.as_numer_denom() is convenient, as it returns a (numerator, denominator) tuple. Using the asterisk expands this.

Now, what else can we do with this? Look at RLC resonators, for example:

A parallel resonator of resistance R, inductance L and capacitance C

The reactance of the inductor and capacitor are given by

(2) $\begin{equation*} X_C = \frac{1}{sC} \quad \text{and} \quad X_L = sL \end{equation*}$

where s is basically the complex frequency (we’re looking at this with Laplace eyes). Now, off to SymPy we go:

R, C, L = sy.symbols("R, C, L", real=True, positive=True)
Xc, Xl = 1/(s*C), s*L

def par(a, b):
    # Shorthand for a parallel circuit
    return a*b/(a+b)

G = par(Xc, par(R, Xl)).ratsimp()
G = G.subs({R: 1e3, L: 1e-6, C: 1e-6})
tf = TF(*G.as_numer_denom(), s)
bode_plot(tf, 5, 7)
pole_zero_plot(tf)

For visualisation we choose R = 1 kΩ, L = 1 µH and C = 1 µF.

Bode plot of the parallel RLC circuit with R = 1 kΩ, L = 1 µH and C = 1 µF, and hence a resonance frequency of 1 MHz

And its pole-zero plot just for the lolz.

Why is this useful? Well, because it is quick and robust and once you have your transfer function typed out, you get the impulse and step responses for free:

import sympy as sy
from sympy.physics.control import impulse_response_plot
from sympy.physics.control import step_response_plot

s = sy.symbols("s", complex=True)
G1 = (2+42*s)/((1+2*s)*(1+40*s))
G2 = (5-60*s)/((1+4*s)*(1+40*s))
G3 = 4*(1+s)/(1+4*s+8*s**2)

for G in [G1, G2, G3]:
    tf = TF(*G.as_numer_denom(), s)
    bode_plot(tf)
    pole_zero_plot(tf)
    step_response_plot(tf, upper_limit=60)
    impulse_response_plot(tf, upper_limit=60)

Note: SymPy’s adaptive plotting is not particularly good at plotting oscillations, so the impulse response of the RLC circuit above will look ugly.