Real-Time Digital Signal Processing Using pyaudio _ helper and the ipywidgets

The focus of this paper is on teaching real-time digital signal processing to electrical and computer engineers using the Jupyter notebook and the code module pyaudio_helper, which is a component of the package scikitdsp-comm. Specifically, we show how easy it is to design, prototype, and test using PC-based instrumentation, real-time DSP algorithms for processing analog signal inputs and returning analog signal outputs, all within the Jupyter notebook. A key feature is that real-time algorithm prototyping is simplified by configuring a few attributes of a DSP_io_stream object from the pyaudio_helper module, leaving the developer to focus on the real-time DSP code contained in a callback function, using a template notebook cell. Real-time control of running code is provided by ipywidgets. The PC-based instrumentation aspect allows measurement of the analog input/output (I/O) to be captured, stored in text files, and then read back into the notebook to compare with the original design expectations via matplotlib plots. In a typical application slider widgets are used to change variables in the callback. One and two channel audio applications as well as algorithms for complex signal (in-phase/quadrature) waveforms, as found in software-defined radio, can also be developed. The analog I/O devices that can be interfaced are both internal and via USB external sound interfaces. The sampling rate, and hence the bandwidth of the signal that can be processed, is limited by the operating system audio subsystem capabilities, but is at least 48 KHz and often 96 kHz.


Introduction
As the power of personal computer has increased, the dream of rapid prototyping of real-time signal processing, without the need to use dedicated DSP-microprocessors or digital signal processing (DSP) enhanced microcontrollers, such as the ARM Cortex-M4 [cortexM4], can be set aside.Students can focus on the powerful capability of numpy, scipy, and matplotlib, along with packages such as scipy.signal[Scipysignal] and scikit-dsp-comm [DSPComm], to explore real-time signals and systems computing.
The focus of this paper is on teaching real-time DSP to electrical and computer engineers using the Jupyter notebook and the code module pyaudio_helper, which is a component of the package scikit-dsp-comm.To be clear, pyaudio_helper is built upon the well known package [pyaudio], which has its roots in Port Audio [portaudio].Specifically, we show how easy it is to design, prototype, and test using PC-based instrumentation, real-time DSP algorithms for processing analog signal inputs and returning analog signal outputs, all within the Jupyter notebook.Real-time algorithm prototyping is simplified by configuring a DSP_io_stream object from the pyaudio_helper module, allowing the developer to quickly focus on writing a DSP callback function using a template notebook cell.The developer is free to take advantage of scipy.signalfilter functions, write custom classes, and as needed utilize global variables to allow the algorithm to maintain state between callbacks pushed by the underlying PyAudio framework.The PC-based instrumentation aspect allows measurement of the analog input/output (I/O) to be captured, stored in text files, and then read back into the notebook to compare with the original design expectations via matplotlib plots.Real-time control of running code is provided by ipywidgets.In a typical application slider widgets are used to change variables in the callback during I/O streaming.The analog I/O devices that can be interfaced are both internal and via USB external sound interfaces.The sampling rate, and hence the bandwidth of the signal that can be processed, is limited by the operating system audio subsystem capabilities, but is at least 48 KHz and often 96 kHz.
We will ultimately see that to set up an audio stream requires: (1) create and instance of the DSP_io_stream class by assigning valid input and output device ports to it, (2) define a callback function to process the input signal sample frames into output sample frames with a user defined algorithm, and (3) call the method interactive_stream() to start streaming.

Analog Input/Output Using DSP Algorithms
A classic text to learn the theory of digital signal processing is [Opp2010].This book is heavy on the underlying theoretical concepts of DSP, including the mathematical modeling of analog I/O systems as shown in Figure 1.This block diagram is a mathematical abstraction of what will be implemented using [pyaudio] and a PC audio subsystem.An analog or continuous-time signal x(t) enters the system on the left and is converted to the discretetime signal x[n] by the analog to digital block.In practice this block is known as the analog-to-digital converter (ADC).The sampling rate f s , which is the inverse of the sampling period, T , leads to x[n] = x(nT ).To be clear, x[n], denotes a sequence of samples corresponding to the original analog input x(t).The use of brackets versus parentheses differentiates the two signal types as discrete-time and continuous-time respectively.The sampling theorem [Opp2010] tells us that the sampling rate f s must be greater than twice the highest frequency we wish to represent in the discrete-time domain.Violating this condition results in aliasing, which means a signal centered on frequency f 0 > f s /2 will land inside the band of frequencies [0, f s /2].Fortunately, most audio ADCs limit the signal bandwidth of x(t) in such a way that signals with frequency content greater than f s /2 are eliminated from passing through the ADC.Also note in practice, x[n] is a scaled and finite precision version of x(t).In real-time DSP environments the ADC maps the analog signal samples to signed integers, most likely int16.As we shall see in pyaudio, this is indeed the case.

Analog
To Digital

Digital
To Analog Digital Signal Processing Algorithms The DSP algorithms block can be any operation on samples x[n] that makes sense.Ultimately, once we discuss frame-based processing in the next section, we will see how Python code fulfills this.At this beginning stage, the notion is that the samples flow through the algorithm one at a time, that is one input results in one output sample.The output samples are converted back to analog signal y(t) by placing the samples into a digital-to-analog converter (DAC).The DAC does not simply set y(nT ) = y[n], as a continuous function time t must be output.A reconstruction operation takes place inside the DAC which interpolates the y[n] signal samples over continuous time.In most DACs this is accomplished with a combination of digital and analog filters, the details of which is outside the scope of this paper.The use of In a DSP theory class the algorithm for producing y[n] from x[n] is typically a causal linear time-invariant (LTI) system/filter, implemented via a difference equation, i.e., where a k , k = 1, 2, . . ., N and b m , m = 0, 1, . . ., M are the filter coefficients.The filter coefficients that implement a particular filter design can be obtained using design tools in [DSPComm].Other algorithms of course are possible.We might have a two channel system and perform operations on both signals, say combining them, filtering, and locally generating time varying periodic signals to create audio special effects.When first learning about real-time DSP it is important to start with simple algorithm configurations, so that external measurements can be used to characterize the systems and verify that the intended results are realized.Developing a real-time DSP project follows along the lines of, design, implement, and test using external test equipment.The Jupyter notebook allows all of this to happen in one place, particularly if the test instrumentation is also PC-based, since PC-based instrument results can be exported as csv and then imported in Jupyter notebook using loadtxt.Here we advocate the use of PC-based instruments, so that all parties, student/instructor/tinkerer, can explore real-time DSP from most anywhere at any time.In this paper we use the Analog Discovery 2 [AD2] for signal generation (two function generator channels), signal measurement (two scope channels, with fast Fourier transform (FFT) spectrum analysis included).It is also helpful to have a signal generator cell phone app available, and of course music from a cell phone or PC.All of the cabling is done using 3.5mm stereo patch cables and small pin header adapters [3p5mm] to interface to the AD2.

Frame-based Real-Time DSP Using the DSP_io_stream class
The block diagram of Figure 2 illustrates the essence of this paper.Implementing the structure of this figure relies upon the class DSP_io_stream, which is housed in sk_dsp_comm.pyaudio_helper.py.To make use of this class requires the scipy stack (numpy, scipy, and matplotlib), as well as [DSPComm] and [pyaudio].PyAudio is multi-platform, with the configuration platform dependent.The set-up is documented at [pyaudio] and SPCommTutorial.The classes and functions of pyaudio_helper are detailed in Figure 3.We will make reference to the classes, methods, and functions throughout the remainder of this paper.With DSP_io_stream one or two channel streaming is possible, as shown in Figure 2. The ADCs and DACs can be internal to the PC or external, say using a USB interface.In a modern PC the audio subsystem has a microphone hardwired to the ADCs and the DACs are connected to the speakers and 3.5mm headphone jack.To provide more flexibility in doing real-time DSP, an external USB audio interface is essential.Two worthy options are the Sabrent at less than $10 and the iMic at under $40.You get what you pay for.The iMic is ideal for full two channel audio I/O processing and also has a line-in/mic switch setting, while the Sabrent offers a single channel input and two channel output.Both are very capable for their intended purposes.A photograph of the AD2 with the iMic interface, 3.5mm splitters and the pin header interfaces mentioned earlier, is shown in Figure 4.The 3.5mm audio splitters are optional, but allow headphones to be plugged into the output while leaving the AD2 scope connected, and the ability to input music/function generator from a cellphone while leaving the AD2 input cable connected (pins wires may need to be pulled off the AD2 to avoid interaction between the two devices in parallel).
When a DSP_io_stream is created (top of Figure 3) it needs to know which input and output devices to connect to.If you just want and input or just an out, you still need to supply a valid output or input device, respectively.To list the internal/external devices available on a given PC we use the function available_devices() from Figure 3.If you add or remove devices while the notebook kernel is running, you will need to restart the kernel to get an accurate listing of devices.The code block below was run with the iMic plugged into a USB hub: The output list can be viewed as a look-up table (LUT) for how to patch physical devices into the block diagram of Figure 2.
We now shift the focus to the interior of Figure 2   throughput.On a PC, with its multitasking OS, there is a lot going on.To get reasonable audio sample throughput the PC audio subsystem fills or packs an input buffer with frame_length samples (or two times frame_length), sample for a two channel stream) originating as 16-bit signed integers (i.e., int16), before calling the callback function.The details of the callback function is the subject of the next section.As the callback prepares to exit, an output buffer of 16-bit signed integers is formed, again of length frame_length, and the buffer is absorbed by the PC audio subsystem.In the context of embedded systems programming, the callback can be thought of as an interrupt service routine.To the PC audio community the frame or buffer just described, is also known as a CHUNK.In a two-channel stream the frame holds an interleaving of left and right channels, ...LRLRL... in the buffer formed/absorbed by the PC audio system.Understand that the efficiency of frame-based processing comes with a price.The buffering either side of the callback block of Figure 2 introduces a latency or processing time delay of at least two times the frame_length times the sampling period.
Moving along with this top level discussion, the central block of Figure 2 is labeled Frame-Based DSP Callback, and as we have alluded to already, is where the real-time DSP code resides.Global variables are needed inside the callback, as the input/output signature is fixed by [pyaudio].The globals allow algorithm parameters to be available inside the callback, e.g., filter coefficients, and in the case of a digital filter, the filter state must be maintained from frame-to-frame.We will see in the examples section how scipy.signal.lfilter(),which implements (1), conveniently supports frame-based digital filtering.To allow interactive control of parameters of the DSP algorithm we can use ipywidgets.We will also see later the sliders widgets are particularly suited to this task.

Anatomy of a PyAudio Callback function
Before writing the callback we first need to instantiate a DSP_io_stream object, as shown in the following code block: The constructor for DSP_io_stream of Figure 3 and the code block above confirm that most importantly we need to supply a function callback name, and most likely provide custom input/output device numbers, choose a sampling rate, and optionally choose the length of the capture buffer.
A basic single channel loop through callback function, where the input samples are passed to the output, is shown in the code block below: The frame_length has been set to 1024, and of the four required inputs from [pyaudio], the first, in_data, is the input buffer which we first convert to a int16 ndarray using np.frombuffer, and then as a working array convert to float32.Note to fill the full dynamic range of the fixedpoint signal samples, means that the x[n] sample values can range over [−2 15 , 2 15 − 1].Passing over the comments we set y=x, and finally convert the output array y back to int16 and then in the return line back to a byte-string buffer using .tobytes().
In general when y is converted from float back to int16, clipping/overflow will occur unless the dynamic range mentioned above is observed.Along the way code instrumentation methods from Figure 3 are included to record time spent in the callback (DSP_callback_tic() and DSP_callback_toc()) and store samples for later analysis in the attribute capture_buffer (DSP_capture_add_samples).These features will be examined in an upcoming example.
To start streaming we need to call the method interactive_stream(), which runs the stream in a thread and displays ipywidgets start/stop buttons below the code cell as shown in Figure 5.

Performance Measurements
The loop through example is good place to explore some performance metrics of 2, and take a look at some of the instrumentation that is part of the DSP_io_stream class.The methods DSP_callback_tic() and DSP_callback_toc() store time stamps in attributes of the class.Another attribute stores samples in the attribute data_capture.For the instrumentation to collect operating data we need to set Tcapture greater than zero.We will also set the total run time to 2s: DSP_IO = pah.DSP_io_stream(callback,2,2,fs=48000, Tcapture=2) DSP_IO.interactive_stream(2,1) Running the above in Jupyter notebook cell will capture 2s of data.The method stream_stats() displays the following: Ideal Callback period = 21.33 (ms) Average Callback Period = 21.33 (ms) Average Callback process time = 0.40 (ms) which tells us that as expected for a sampling rate of 48 kHz, and a frame length of 1024 is simply The time spent in the callback should be very small, as very little processing is being done.We can also examine the callback latency by first having the AD2 input a low duty cycle pulse train at a 2 Hz rate, thus having 500 ms between pules.We then use the scope to measure the time difference between the input (scope channel C2) and output (scope channel C1) waveforms.
The resulting plot is shown in Figure 6.We see that PyAudio and and the PC audio subsystem introduces about 70.7ms of latency.
A hybrid iMic ADC and builtin DAC results in 138 ms on macOS.
Moving to Win 10 latency increases to 142 ms, using default USB drivers.
Latency ~ 70.7ms Fig. 6: Callback latency measurement using the AD2 where C2 is the input and C1 is the output, of a 2 Hz pulse train in the loop through app.
The frequency response magnitude of an LTI system can be measured using the fact that [Opp2010] at the output of a system driven by white noise, the measured power output spectrum is a scaled version of the underlying system frequency response magnitude squared, i.e.,

S y,measured
where σ 2 x is the variance of the input white noise signal.Here we use this technique to first estimate the frequency response magnitude of the input path (ADC only) using the attribute DSP_IO.capture_buffer, and secondly take end-to-end (ADC-DAC) measurements using the AD2 spectrum analyzer in dB average mode (500 records).In both cases the white noise input is provided by the AD2 function generator.Finally, the AD2 measurement is saved to a CSV file and imported into the Jupyter notebook, as shown in the code block below.This allows an overlay of the ADC and ADC-DAC measurements, entirely in the Jupyter notebook.The results are compared in Figure 7, where we see a roll-off of about 3 dB at about 14 kHz in both the ADC path and the composite ADC-DAC path.The composite ADC-DAC begins to rise above 17 kHz and flattens to 2 dB down from 18-20 kHz.As a practical matter, humans do not hear sound much above 16 kHz, so the peaking is not much of an issue.Testing of the Sabrent device the composite ADC-DAC 3 dB roll-off occurs at about 17 kHz.The native PC audio output can for example be tested in combination with the iMic or Sabrent ADCs.Fig. 7: Gain flatness of the loop through app of just the ADC path via the DSP_IO.capture_bufferand then the ADC-DAC path using the AD2 spectrum analyzer to average the noise spectrum.

Examples
In this section we consider a collection of applications examples.This first is a simple two channel loop-through with addition of left and right gain sliders.The second is again two channel, but now cross left-right panning is developed.In of these examples the DSP is memoryless, so there is no need to maintain state using Python globals.The third example is an equal-ripple bandpass filter, which utilizes sk_dsp_comm.fir_design_helperto design the filter.The final example develops a three-band audio equalizer using peaking filters to raise and lower the gain over a narrow band of frequencies.

Left and Right Gain Sliders
In this first example the signal processing is again minimal, but now two-channel (stereo) processing is utilized, and left and right channel gain slider using ipywidgets are introduced.Since the audio stream is running in a thread, the ipywidgets can freely run and interactively control parameters inside the callback function.The two slider widgets are created below, followed by the callback, and finally calling the interactive_stream method to run without limit in two channel mode.A 1 kHz sinusoid test signal is input to the left channel and a 5 kHz sinusoid is input to the right channel.While viewing the AD2 scope output in real-time, the gain sliders are adjusted and the signal levels move up and down.A screenshot taken from the Jupyter notebook is combined with a screenshot of the scope output to verify the correlation between the observed signal amplitudes and the slider positions is given in Figure 8.The callback listing, including the set-up of the ipywidgets gain sliders, is given below: Note for this two channel stream, the audio subsystem interleaves left and right samples, so now the class methods get_LR and pack_LR of Figure 3 are utilized to unpack the left and right samples and then repack them, respectively.A screenshot of the gain sliders app, including an AD2 scope capture, with C1 on the left channel and C2 on the right channel, is given in Figure 8.
The ability to control the left and right audio level are as expected, especially when listening.

Cross Left-Right Channel Panning
This example again works with a two channel signal flow.The application is to implement a cross channel panning system.Ordinarily panning moves a single channel of audio from 100% left to 100% right as a slider moves from 0% to 100% of its range.At 50% the single channel should have equal amplitude in both channels.In cross channel panning two input channels are super imposed, but such that at 0% the left and right channels are fully in their own channel.At 50% the left and right outputs are equally mixed.At 100% the input channels are now swapped.Assuming that a represents the panning values on the interval [0, 100], a mathematical model of the cross panning app is where L in and L out are the left channel inputs and outputs respectively, and similarly R in and R out for the right channel.In code we have: For dissimilar left and right audio channels, the action of the slider creates a spinning effect when listening.It is possible to extend this app with an automation, so that a low frequency sinusoid or other waveform changes the panning value at a rate controlled by a slider.

FIR Bandpass Filter
In this example we design a high-order FIR bandpass filter using sk_dsp_comm.fir_design_helperand then implement the design to operate at f s = 48 kHz.Here we choose the bandpass critical frequencies to be 2700, 3200, 4800, and 5300 Hz, with a passband ripple of 0.5 dB and stopband attenuation of 50 dB (see fir_d).Theory is compared with AD2 measurements using, again using noise excitation.When implementing a digital filter using frame-based processing, scipy.signal.lfilterworks nicely.The key is to first create a zero initial condition array zi and hold this in a global variable.Each time lfilter is used in the callback the old initial condition zi is passed in, then the returned zi is held until the next time through the callback.
import sk_dsp_comm.fir_design_helperas fir_d import scipy.signalas signal b = fir_d.fir_remez_bpf(2700,3200,4800,5300, .5,50,48000,18)Following the call to DSP_io.interactive_stream() the start button is clicked and the AD2 spectrum analyzer estimates the power spectrum.The estimate is saved as a CSV file and brought into the Jupyter notebook to overlay the theoretical design.The comparison results are given in Figure 10.
The theory and measured magnitude response plots are in very close agreement, making the end-to-end design, implement, test very satisfying.

Three Band Equalizer
Here we consider the second-order peaking filter, which has infinite impulse response, and place three of them in cascade with a ipywidgets slider used to control the gain of each filter.The peaking filter is used in the design of audio equalizer, where perhaps each filter is centered on octave frequency spacings running from from 10 Hz up to 16 kHz, or so.Each peaking filter can be implemented as a 2nd-order difference equation, i.e., N = 2 in equation ( 1).The design equations for a single peaking filter are given below using z-transform [Opp2010] notation: Fig. 10: An overlay plot of the theoretical frequency response with the measured using an AD2 noise spectrum capture import to the Jupyter notebook.
which has coefficients and f c is the center frequency in Hz relative to sampling rate f s in Hz, and G dB is the peaking filter gain in dB.Conveniently, the function peaking is available in the module sk_dsp_comm.sigsys.The app code is given below starting with the slider creation: Following the call to DSP_io.interactive_stream() the start button is clicked and the FFT spectrum analyzer estimates the power spectrum.The estimate is saved as a CSV file and brought into the Jupyter notebook to overlay the theoretical design.The comparison results are given in Figure 11.
Reasonable agreement is achieved, but listening to music is a more effective way of evaluating the end result.To complete the design more peaking filters should be added.

Conclusions and Future Work
In this paper we have described an approach to implement real-time DSP in the Jupyter notebook.This real-time capability rests on top of PyAudio and the wrapper class DSP_io_stream contained in sk_dsp_comm.pyaudio_helper.The ipywidgets allow for interactivity while real-time DSP code is running.The callback function does the work using frame-based algorithms, which takes some getting used to.By working through examples we have shown that much can be accomplished with little coding.
A limitation of using PyAudio is the input-to-output latency.At a 48 kHz sampling rate a simple loop though app has around  70 ms of delay.For the application discussed in the paper latency is not a show stopper.
In the future we hope to easily develop algorithms that can demodulate software-defined radio (SDR) streams and send the recovered modulation signal out the computer's audio interface via PyAudio.Environments such as GNURadio companion already support this, but being able to do this right in the Jupyter notebook is our desire.

Fig. 2 :
Fig. 2: Two channel analog signal processing implemented using frame-based real-time DSP.

Fig. 5 :
Fig.5: Setting up an interactive stream for the simple y = x loop through, using a run time of 0, which implies run forever.

Fig. 9 :
Fig. 9: Cross left/right panning control: (a) launching the app in the Jupyter notebook and (b)-(d) a sequence of scope screenshots as the panning slider is moved from 0% to 50%, and then to 100%.

Fig. 11 :
Fig. 11: Three band equalizer: (a) launching the app in the Jupyter notebook and (b) an overlay plot of the theoretical log-frequency response with the measured using an AD2 noise spectrum capture import to the Jupyter notebook.
to discuss frame-based DSP and the Frame-Based DSP Callback.When a DSP microcontroller is configured for real-time DSP, it can focus on just this one task very well.Sample-by-sample processing is possible with low I/O latency and overall reasonable audio sample