uct logo PfP logo

Pyscellania: a selection of Python programs

Miscellaneous codes in Python.

Oscilloscope and spectrum analyser

Dual oscilloscope in Python

This implements an oscilloscope and spectrum analyser running from the computer’s sound card.

This is an example of how an application with a nontrivial graphical user interface and real-time data acquisition can be implemented in Python.

There are two panels, one displaying the sound card input as a function of time (two traces are possible for a stereo input) and the other displaying the power spectrum of the input.

Possible sample rates depend on the sound card. Most should be able to handle 48kHz.

Requires:

  • PyQt4 and PyQwt for the GUI framework.
  • PyAudio to access the sound card.
  • numpy for numerics.

Tested on Linux and Windows.

Download the code for the scope and its icons.

Non-linear least squares fit

LRC resonance curve fit

Something for the second or third year lab course.

An example of how to do a nonlinear least squares fit, together with estimating the uncertainties and correlations in the fit parameters. Fitting is via the Levenberg-Marquard method, as implemented in scipy.

Experimental data was obtained from a LRC resonance measurement with uncertainties estimated from the oscilloscope settings. The well-known resonance function is fitted to the data to obtain the resonance frequency, the resistance of the inductor and the Q value. The value of R and C are known, so L can be calculated.

Requires:

  • scipy, numpy and pylab.

Download the code for fitting and data.

Solar system in thirty lines of code

Numerical integration of solar system

Simulation of the solar system using numerical integration of the equations of motion by the symplectic Euler method, a simple first order method.

This is sufficient to show that one can get realistic behaviour using the gravitational force and the equations of motion.

Planets, and trails to trace out the orbits, are visualised using VPython. Planet sizes are not to scale; in fact they are all the same. Size is chosen as a compromise with visibility in the simulation. You can add a more realistic depiction by adding a line of data with sizes, and using this to set the initial size of the planets.

Requires:

  • visual.
  • solardat.py: a file of initial data for the solar system. This is a lot longer than 30 lines! (Well, as written ...)

The integration loop is:

for i in range(nbodyused):  # N is number of bodies: loop over all
    a[i]=vector(0,0,0) # initialise acceleration of body i to 0
    for j in range(nbodyused): # loop over bodies exerting force
        if i == j: continue  # body doesn't exert force on self
        R = r[i]-r[j]  # displacement
        a[i]=a[i]-Gm[j]*R/mag(R)**3
# update v and r
v=v+a*dt
r=r+v*dt

Download the code for solar system and solar system initial data

Wavelet analysis (continuous wavelet transform)

Wavelet power

Wavelet analysis gives information on both location and scale of features. It provides an alternative to Fourier analysis.

The Python code illustrated here implements the continuous wavelet transform. A basic wavelet class is defined, that performs the transform upon instantiation, and different wavelet families are defined by subclassing this basic class.

The plot shown uses the Morlet wavelet.

Requires:

  • numpy.
  • pylab for the display of the results.

Download the code for wavelets

Solving the Time Dependent Schroedinger Equation

TDSE

This program solves the TDSE using the pseudo-spectral approach. (e.g. DeVries, A first course in computational physics, Wiley).

A Gaussian wavepacket in 1D scatters from a potential barrier. The solution involves flipping between position and momentum space representations using the FFT. Because of this, the boundary conditions are periodic, and the visualisation highlights this by plotting the position on a circle. It is usual to plot either the real and imaginary parts of the wave function, or the density (square of the wave function). In this visualisation, the (complex) wave function at a point on the circle is plotted on a local complex plane normal to the circle at that point. The real part of the wave function can be obtained by projecting onto a cylindrical extrusion of the circle; the imaginary part by projecting onto the plane of the circle.

The result is more Pythonic than the usual approach.

An example of a more usual approach, using finite differences, and plotting the real and imaginary parts of the wave function, is on the SciPy website. Reflective boundary conditions are used here, so that the system is effectively inside an infinite square well.

Requires:

  • visual.

Download the code for the TDSE