Python’s “batteries included” nature makes it easy to interact with just about anything… except speakers and a microphone! As of this moment, there still are not standard libraries which which allow cross-platform interfacing with audio devices. There are some pretty convenient third-party modules, but I hope in the future a standard solution will be distributed with python. I appreciate the differences of Linux architectures such as ALSA and OSS, but toss in Windows and MacOS in the mix and it gets to be a huge mess. For Linux, would I even need anything fancy? I can run “cat file.wav > /dev/dsp” from a command prompt to play audio. There are some standard libraries for operating system specific sound (i.e., winsound), but I want something more versatile. The official audio wiki page on the subject lists a small collection of third-party platform-independent libraries. After excluding those which don’t support microphone access (the ultimate goal of all my poking around in this subject), I dove a little deeper into sounddevice and PyAudio. Both of these I installed with pip (i.e., pip install pyaudio)

For a more modern, cleaner, and more complete GUI-based viewer of realtime audio data (and the FFT frequency data), check out my Python Real-time Audio Frequency Monitor project.

I really like the structure and documentation of sounddevice, but I decided to keep developing with PyAudio for now. Sounddevice seemed to take more system resources than PyAudio (in my limited test conditions: Windows 10 with very fast and modern hardware, Python 3), and would audibly “glitch” music as it was being played every time it attached or detached from the microphone stream. I tried streaming, but after about an hour I couldn’t get clean live access to the microphone without glitching audio playback. Furthermore, every few times I ran this script it crashed my python kernel! I very rarely see this happening. iPython complained: “It seems the kernel died unexpectedly. Use ‘Restart kernel’ to continue using this console” and I eventually moved back to PyAudio. For a less “realtime” application, sounddevice might be a great solution. Here’s the minimal case sounddevice script I tested with (that crashed sometimes). If you have a better one to do live high-speed audio capture, let me know!

import sounddevice #pip install sounddevice

for i in range(30): #30 updates in 1 second
    rec = sounddevice.rec(44100/30)
    sounddevice.wait()
    print(rec.shape)

Here’s a simple demo to show how I get realtime microphone audio into numpy arrays using PyAudio. This isn’t really that special. It’s a good starting point though. Note that rather than have the user define a microphone source in the python script (I had a fancy menu system handling this for a while), I allow PyAudio to just look at the operating system’s default input device. This seems like a realistic expectation, and saves time as long as you don’t expect your user to be recording from two different devices at the same time. This script gets some audio from the microphone and shows the values in the console (ten times).

import pyaudio
import numpy as np

CHUNK = 4096 # number of data points to read at a time
RATE = 44100 # time resolution of the recording device (Hz)

p=pyaudio.PyAudio() # start the PyAudio class
stream=p.open(format=pyaudio.paInt16,channels=1,rate=RATE,input=True,
              frames_per_buffer=CHUNK) #uses default input device

# create a numpy array holding a single read of audio data
for i in range(10): #to it a few times just to see
    data = np.fromstring(stream.read(CHUNK),dtype=np.int16)
    print(data)

# close the stream gracefully
stream.stop_stream()
stream.close()
p.terminate()

01

I tried to push the limit a little bit and see how much useful data I could get from this console window. It turns out that it’s pretty responsive! Here’s a slight modification of the code, made to turn the console window into an impromptu VU meter.

import pyaudio
import numpy as np

CHUNK = 2**11
RATE = 44100

p=pyaudio.PyAudio()
stream=p.open(format=pyaudio.paInt16,channels=1,rate=RATE,input=True,
              frames_per_buffer=CHUNK)

for i in range(int(10*44100/1024)): #go for a few seconds
    data = np.fromstring(stream.read(CHUNK),dtype=np.int16)
    peak=np.average(np.abs(data))*2
    bars="#"*int(50*peak/2**16)
    print("%04d %05d %s"%(i,peak,bars))

stream.stop_stream()
stream.close()
p.terminate()

The results are pretty good! The advantage here is that no libraries are required except PyAudio. For people interested in doing simple math (peak detection, frequency detection, etc.) this is a perfect starting point. Here’s a quick cellphone video:

I’ve made realtime audio visualization (realtime FFT) scripts with Python before, but 80% of that code was creating a GUI. I want to see data in real time while I’m developing this code, but I really don’t want to mess with GUI programming. I then had a crazy idea. Everyone has a web browser, which is a pretty good GUI… with a Python script to analyze audio and save graphs (a lot of them, quickly) and some JavaScript running in a browser to keep refreshing those graphs, I could get an idea of what the audio stream is doing in something kind of like real time. It was intended to be a hack, but I never expected it to work so well! Check this out…

Here’s the python script to listen to the microphone and generate graphs:

import pyaudio
import numpy as np
import pylab
import time

RATE = 44100
CHUNK = int(RATE/20) # RATE / number of updates per second

def soundplot(stream):
    t1=time.time()
    data = np.fromstring(stream.read(CHUNK),dtype=np.int16)
    pylab.plot(data)
    pylab.title(i)
    pylab.grid()
    pylab.axis([0,len(data),-2**16/2,2**16/2])
    pylab.savefig("03.png",dpi=50)
    pylab.close('all')
    print("took %.02f ms"%((time.time()-t1)*1000))

if __name__=="__main__":
    p=pyaudio.PyAudio()
    stream=p.open(format=pyaudio.paInt16,channels=1,rate=RATE,input=True,
                  frames_per_buffer=CHUNK)
    for i in range(int(20*RATE/CHUNK)): #do this for 10 seconds
        soundplot(stream)
    stream.stop_stream()
    stream.close()
    p.terminate()

Here’s the HTML file with JavaScript to keep reloading the image… 

<html>
<script language="javascript">
function RefreshImage(){
document.pic0.src="03.png?a=" + String(Math.random()*99999999);
setTimeout('RefreshImage()',50);
}
</script>
<body onload="RefreshImage()">
<img name="pic0" src="03.png">
</body>
</html>

Here’s the result! I couldn’t believe my eyes. It’s not elegant, but it’s kind of functional!

Why stop there? I went ahead and wrote a microphone listening and processing class which makes this stuff easier. My ultimate goal hasn’t been revealed yet, but I’m sure it’ll be clear in a few weeks. Let’s just say there’s a lot of use in me visualizing streams of continuous data. Anyway, this class is the truly terrible attempt at a word pun by merging the words “SWH”, “ear”, and “Hear”, into the official title “SWHear” which seems to be unique on Google. This class is minimal case, but can be easily modified to implement threaded recording (which won’t cause the rest of the functions to hang) as well as mathematical manipulation of data, such as FFT. With the same HTML file as used above, here’s the new python script and some video of the output:

import pyaudio
import time
import pylab
import numpy as np

class SWHear(object):
    """
    The SWHear class is made to provide access to continuously recorded
    (and mathematically processed) microphone data.
    """

    def __init__(self,device=None,startStreaming=True):
        """fire up the SWHear class."""
        print(" -- initializing SWHear")

        self.chunk = 4096 # number of data points to read at a time
        self.rate = 44100 # time resolution of the recording device (Hz)

        # for tape recording (continuous "tape" of recent audio)
        self.tapeLength=2 #seconds
        self.tape=np.empty(self.rate*self.tapeLength)*np.nan

        self.p=pyaudio.PyAudio() # start the PyAudio class
        if startStreaming:
            self.stream_start()

    ### LOWEST LEVEL AUDIO ACCESS
    # pure access to microphone and stream operations
    # keep math, plotting, FFT, etc out of here.

    def stream_read(self):
        """return values for a single chunk"""
        data = np.fromstring(self.stream.read(self.chunk),dtype=np.int16)
        #print(data)
        return data

    def stream_start(self):
        """connect to the audio device and start a stream"""
        print(" -- stream started")
        self.stream=self.p.open(format=pyaudio.paInt16,channels=1,
                                rate=self.rate,input=True,
                                frames_per_buffer=self.chunk)

    def stream_stop(self):
        """close the stream but keep the PyAudio instance alive."""
        if 'stream' in locals():
            self.stream.stop_stream()
            self.stream.close()
        print(" -- stream CLOSED")

    def close(self):
        """gently detach from things."""
        self.stream_stop()
        self.p.terminate()

    ### TAPE METHODS
    # tape is like a circular magnetic ribbon of tape that's continously
    # recorded and recorded over in a loop. self.tape contains this data.
    # the newest data is always at the end. Don't modify data on the type,
    # but rather do math on it (like FFT) as you read from it.

    def tape_add(self):
        """add a single chunk to the tape."""
        self.tape[:-self.chunk]=self.tape[self.chunk:]
        self.tape[-self.chunk:]=self.stream_read()

    def tape_flush(self):
        """completely fill tape with new data."""
        readsInTape=int(self.rate*self.tapeLength/self.chunk)
        print(" -- flushing %d s tape with %dx%.2f ms reads"%\
                  (self.tapeLength,readsInTape,self.chunk/self.rate))
        for i in range(readsInTape):
            self.tape_add()

    def tape_forever(self,plotSec=.25):
        t1=0
        try:
            while True:
                self.tape_add()
                if (time.time()-t1)>plotSec:
                    t1=time.time()
                    self.tape_plot()
        except:
            print(" ~~ exception (keyboard?)")
            return

    def tape_plot(self,saveAs="03.png"):
        """plot what's in the tape."""
        pylab.plot(np.arange(len(self.tape))/self.rate,self.tape)
        pylab.axis([0,self.tapeLength,-2**16/2,2**16/2])
        if saveAs:
            t1=time.time()
            pylab.savefig(saveAs,dpi=50)
            print("plotting saving took %.02f ms"%((time.time()-t1)*1000))
        else:
            pylab.show()
            print() #good for IPython
        pylab.close('all')

if __name__=="__main__":
    ear=SWHear()
    ear.tape_forever()
    ear.close()
    print("DONE")

I don’t really intend anyone to actually do this, but it’s a cool alternative to recording a small portion of audio, plotting it in a pop-up matplotlib window, and waiting for the user to close it to record a new fraction. I had a lot more text in here demonstrating real-time FFT, but I’d rather consolidate everything FFT related into a single post. For now, I’m happy pursuing microphone-related python projects with PyAudio.

UPDATE: Displaying a single frequency

Use Numpy’s FFT() and FFTFREQ() to turn the linear data into frequency. Set that target and grab the FFT value corresponding to that frequency. I haven’t tested this to be sure it’s working, but it should at least be close…

import pyaudio
import numpy as np
np.set_printoptions(suppress=True) # don't use scientific notation

CHUNK = 4096 # number of data points to read at a time
RATE = 44100 # time resolution of the recording device (Hz)
TARGET = 2100 # show only this one frequency

p=pyaudio.PyAudio() # start the PyAudio class
stream=p.open(format=pyaudio.paInt16,channels=1,rate=RATE,input=True,
              frames_per_buffer=CHUNK) #uses default input device

# create a numpy array holding a single read of audio data
for i in range(10): #to it a few times just to see
    data = np.fromstring(stream.read(CHUNK),dtype=np.int16)
    fft = abs(np.fft.fft(data).real)
    fft = fft[:int(len(fft)/2)] # keep only first half
    freq = np.fft.fftfreq(CHUNK,1.0/RATE)
    freq = freq[:int(len(freq)/2)] # keep only first half
    assert freq[-1]>TARGET, "ERROR: increase chunk size"
    val = fft[np.where(freq>TARGET)[0][0]]
    print(val)

# close the stream gracefully
stream.stop_stream()
stream.close()
p.terminate()

UPDATE: Display peak frequency

If your goal is to determine which frequency is producing the loudest tone, use this function. I also added a few lines to graph the output in case you want to observe how it operates. I recommend testing this script with a tone generator, or a YouTube video containing tones of a range of frequencies like this one.

import pyaudio
import numpy as np
import matplotlib.pyplot as plt

np.set_printoptions(suppress=True) # don't use scientific notation

CHUNK = 4096 # number of data points to read at a time
RATE = 44100 # time resolution of the recording device (Hz)

p=pyaudio.PyAudio() # start the PyAudio class
stream=p.open(format=pyaudio.paInt16,channels=1,rate=RATE,input=True,
              frames_per_buffer=CHUNK) #uses default input device

# create a numpy array holding a single read of audio data
for i in range(10): #to it a few times just to see
    data = np.fromstring(stream.read(CHUNK),dtype=np.int16)
    data = data * np.hanning(len(data)) # smooth the FFT by windowing data
    fft = abs(np.fft.fft(data).real)
    fft = fft[:int(len(fft)/2)] # keep only first half
    freq = np.fft.fftfreq(CHUNK,1.0/RATE)
    freq = freq[:int(len(freq)/2)] # keep only first half
    freqPeak = freq[np.where(fft==np.max(fft))[0][0]]+1
    print("peak frequency: %d Hz"%freqPeak)

    # uncomment this if you want to see what the freq vs FFT looks like
    #plt.plot(freq,fft)
    #plt.axis([0,4000,None,None])
    #plt.show()
    #plt.close()

# close the stream gracefully
stream.stop_stream()
stream.close()
p.terminate()

 

This program shows left vs right audio level:

import pyaudio
import numpy as np

maxValue = 2**16
p=pyaudio.PyAudio()
stream=p.open(format=pyaudio.paInt16,channels=2,rate=44100,
              input=True, frames_per_buffer=1024)
while True:
    data = np.fromstring(stream.read(1024),dtype=np.int16)
    dataL = data[0::2]
    dataR = data[1::2]
    peakL = np.abs(np.max(dataL)-np.min(dataL))/maxValue
    peakR = np.abs(np.max(dataR)-np.min(dataR))/maxValue
    print("L:%00.02f R:%00.02f"%(peakL*100, peakR*100))

Output

L:47.26 R:45.17
L:47.55 R:45.63
L:49.44 R:45.98
L:45.27 R:49.80
L:44.39 R:45.75
L:47.50 R:46.96
L:41.49 R:42.64
L:42.95 R:41.39
L:49.56 R:49.62
L:48.29 R:48.80
L:45.03 R:47.62
L:47.99 R:49.35
L:41.58 R:49.21

Or with a tweak…

import pyaudio
import numpy as np

maxValue = 2**16
bars = 35
p=pyaudio.PyAudio()
stream=p.open(format=pyaudio.paInt16,channels=2,rate=44100,
              input=True, frames_per_buffer=1024)
while True:
    data = np.fromstring(stream.read(1024),dtype=np.int16)
    dataL = data[0::2]
    dataR = data[1::2]
    peakL = np.abs(np.max(dataL)-np.min(dataL))/maxValue
    peakR = np.abs(np.max(dataR)-np.min(dataR))/maxValue
    lString = "#"*int(peakL*bars)+"-"*int(bars-peakL*bars)
    rString = "#"*int(peakR*bars)+"-"*int(bars-peakR*bars)
    print("L=[%s]\tR=[%s]"%(lString, rString))

graphical output:





Additional Resources

After using the AVR-ISP mkII for years (actually the cheap eBay knock-offs) to program ATMEL AVR microcontrollers, today I gave the Bus Pirate a shot. Far more than just a microcontroller programmer, this little board is basically a serial interface to basic microcontroller peripherals. In a nutshell, you plug it in via USB and it looks like a serial port which has a command-line interface that lets you do things like turn pins on and off, perform voltage measurements, and it naively supports bidirectional use of common protocols like I2C, SPI, UART, and even HD44780 series LCDs. Note that although you could directly interface with the Bus Pirate using HyperTerminal, I recommend using TeraTerm. It can supply voltages (3.3V and 5V) to power small circuits, and if current draw is too high (indicating something is hooked-up wrong) it automatically turns the supply off. So clever! At <$30, it’s a cool tool to have around. In addition, it’s naively supported as an AVR programmer by AVRDUDE. Although I could write assembly to perform tasks, I almost always write in C for the convenience. For my reference (and that of anyone who may want to do something similar), I’m posting the simplest-case method I use to program AVR microcontrollers with the Bus Pirate on Windows (noting that Linux would be nearly identical). I also wrote a Python script to connect with the Bus Pirate and run simple commands (which turns the power supply on and report the voltage of the VCC line immediately after programming completes).  Yes, there are fancy packages that allow you to interact with Bus Pirate from Python, but the advantage of my method is that it runs from native Python libraries! To get this all up and running for yourself, just install WinAVR (which supplies AVRDUDE and AVR-GCC) and Python 3. I assume this code will work just as well on Python 2, but haven’t tried.

IMG_7092 (1)
the Bus Pirate programming an ATTiny85 microcontroller

 

To ensure my Bus Pirate is working properly, I start off by running the Bus Pirate’s built-in test routine. For full details read the guide. It just involves connecting two pairs of pins together as shown in the picture here, connecting to the Bus Pirate with the serial terminal, and running the command “~”. It will output all sorts of useful information. Once I know my hardware is up and running, I’m good to continue.

Bpv3v2go-pinout

Here’s the code which runs on the microcontroller to twiddle all the pins (saved as main.c). Note that my MCU is an ATTiny85. I’m using standard clock settings (internal RC clock, 8MHz), but if I wanted to modify fuses to do things like use an external clock source or crystal, I’d calculate them with engbedded’s handy dandy fuse calculator (which also shows AVRdude arguments needed to make the change!).

#define	F_CPU (8000000UL)
#include <avr/io.h>
#include <util/delay.h>

int main (void)
{
    DDRB = 255; 
    while(1) 
    {
        PORTB ^= 255;
        _delay_ms(500);
    }
}

To compile the code and program the MCU with it, I always have a bash script in the same folder that I can double-click on to delete old compiled files (so we don’t accidentally re-program our MCU with old code), compile main.c, and load it onto the MCU using the Bus Pirate. You may have to change COM3 to reflect the com port of your Bus Pirate. Note that it is required that you disconnect other terminals from the Bus Pirate before doing this, otherwise you’ll get an “access denied” error.

@echo off
del *.elf
del *.hex
avr-gcc -mmcu=attiny85 -Wall -Os -o main.elf main.c
avr-objcopy -j .text -j .data -O ihex main.elf main.hex
avrdude -c buspirate -p attiny85 -P com3 -e -U flash:w:main.hex
python up.py

Although the programmer briefly supplies my MCU with power from the +5V pin, it’s cut after programming completes. Rather than manually re-opening my terminal program, re-connecting with the bus pirate, re-setting the mode (command “m”) to something random (DIO, command “9”), and re-enableing voltage output (command “W”) just to see my LED blink, I want all that to be automated. Thanks python for making this easy. The last line calls “up.py”. This fancy script even outputs the voltage of the VCC line after it’s turned on!

"""python3 control of buspirate (SWHarden.com)"""

import serial

BUSPIRATE_PORT = 'com3' #customize this! Find it in device manager.

def send(ser,cmd):
    """send the command and listen to the response."""
    ser.write(str(cmd+'\n').encode('ascii')) # send our command
    for line in ser.readlines(): # while there's a response
        print(line.decode('utf-8').strip()) # show it

ser=serial.Serial(BUSPIRATE_PORT, 115200, timeout=1) # is com free?
assert ser.isOpen() #throw an exception if we aren't connected
send(ser,'#') # reset bus pirate (slow, maybe not needed)
send(ser,'m') # change mode (goal is to get away from HiZ)
send(ser,'9') # mode 9 is DIO
send(ser,'W') # turn power supply to ON. Lowercase w for OFF.
send(ser,'v') # show current voltages
ser.close() # disconnect so we can access it from another app
print("disconnected!") # let the user know we're done.

When “burn.cmd” is run, the code is compiled and loaded, the power supply is turned on (and killed if too much current is drawn!), and the voltage on VCC is reported. The output is:

C:\Users\scott\Documents\important\AVR\2016-07-13 ATTiny85 LEDblink>burn.cmd

Detecting BusPirate...
**
**  Bus Pirate v3a
**  Firmware v5.10 (r559)  Bootloader v4.4
**  DEVID:0x0447 REVID:0x3046 (24FJ64GA002 B8)
**  http://dangerousprototypes.com
**
BusPirate: using BINARY mode
avrdude: AVR device initialized and ready to accept instructions

Reading | ################################################## | 100% 0.12s

avrdude: Device signature = 0x1e930b
avrdude: erasing chip
avrdude: reading input file "main.hex"
avrdude: input file main.hex auto detected as Intel Hex
avrdude: writing flash (84 bytes):

Writing | ################################################## | 100% 3.12s

avrdude: 84 bytes of flash written
avrdude: verifying flash memory against main.hex:
avrdude: load data flash data from input file main.hex:
avrdude: input file main.hex auto detected as Intel Hex
avrdude: input file main.hex contains 84 bytes
avrdude: reading on-chip flash data:

Reading | ################################################## | 100% 2.72s

avrdude: verifying ...
avrdude: 84 bytes of flash verified

avrdude: safemode: Fuses OK

avrdude done.  Thank you.

#
RESET

Bus Pirate v3a
Firmware v5.10 (r559)  Bootloader v4.4
DEVID:0x0447 REVID:0x3046 (24FJ64GA002 B8)
http://dangerousprototypes.com
HiZ>
m
1. HiZ
2. 1-WIRE
3. UART
4. I2C
5. SPI
6. 2WIRE
7. 3WIRE
8. LCD
9. DIO
x. exit(without change)

(1)>
9
Ready
DIO>
W
Power supplies ON
DIO>
v
Pinstates:
1.(BR)  2.(RD)  3.(OR)  4.(YW)  5.(GN)  6.(BL)  7.(PU)  8.(GR)  9.(WT)  0.(Blk)
GND     3.3V    5.0V    ADC     VPU     AUX     CLK     MOSI    CS      MISO
P       P       P       I       I       I       I       I       I       I
GND     3.17V   5.00V   0.00V   0.00V   L       L       L       H       L
DIO>
disconnected!


This is a minimal-case scenario, but can be obviously expanded to perform some complicated tasks! For example, all commands could be run from a single python program. Considering the Bus Pirate’s ability to communicate with so many different protocols (I2C, 2-write, etc.), being able to naively control it from Python without having to install special additional libraries will certainly prove to be convenient.

PS: I noted there is a surprising delay when initializing programming the AVR with the bus pirate. The process hangs for about 10 seconds after the bus pirate introduces itself with the welcome message, then seems to resume at full speed writing to the flash of the microchip. After a bit of Googling, I believe the delay is due to the Bus Pirate slowly bit-banging SPI to initialize the programming sequence. The AVR has rich SPI functionality, some of which involves its own programming. Satisfied with this answer for now, I’m not going to try to speed it up. It’s a little annoying, but not too bad that I won’t use this to program my AVRs.





Additional Resources

I was recently presented with the need to rename a folder of images based on a timestamp. This way, I can keep saving new files in that folder with overlapping filenames (i.e., 01.jpg, 02.jpg, 03.jpg, etc.), and every time I run this script all images are prepended with a timestamp. I still want the files to be sorted alphabetically, which is why an alphabetical timestamp (rather than a random hash) is preferred.

  • At first I considered a long date such as 2014-04-19-01.jpg, but that adds so much text!
    …also, it doesn’t include time of day.
  • If I include time of day, it becomes 2014-04-19-09-16-23-01.jpg
  • If I eliminate dashes to shorten it, it becomes hard to read, but might work 140419091623-01.jpg
  • If I use Unix Epoch time, it becomes 1397912944-01.jpg

The result I came up with uses base conversion and a string table of numbers and letters (in alphabetical order) to create a second-respecting timestamp hash using an arbitrary number of characters. For simplicity, I used 36 characters: 0-9, and a-z. I then wrote two functions to perform arbitrary base conversion, pulling characters from the hash. Although I could have nearly doubled my available characters by including the full ASCII table, respecting capitalization, I decided to keep it simple. The scheme goes like this:

  • Determine the date / time: 19-Apr-2014 13:08:55
  • Create an integer of Unix Epoch time (seconds past Jan 1, 1970):  1397912935
  • Do a base conversion from a character list: n4a4iv
  • My file name now becomes n4a4iv-01.jpg – I can accept this!
    and when I sort the folder alphabetically, they’re in order by the timestamp

I can now represent any modern time, down to the second, with 6 characters. Here’s some example output:

19-Apr-2014 13:08:55 <-> 1397912935 <-> n4a4iv
19-Apr-2014 13:08:56 <-> 1397912936 <-> n4a4iw
19-Apr-2014 13:08:57 <-> 1397912937 <-> n4a4ix
19-Apr-2014 13:08:58 <-> 1397912938 <-> n4a4iy
19-Apr-2014 13:08:59 <-> 1397912939 <-> n4a4iz
19-Apr-2014 13:09:00 <-> 1397912940 <-> n4a4j0
19-Apr-2014 13:09:01 <-> 1397912941 <-> n4a4j1
19-Apr-2014 13:09:02 <-> 1397912942 <-> n4a4j2
19-Apr-2014 13:09:03 <-> 1397912943 <-> n4a4j3
19-Apr-2014 13:09:04 <-> 1397912944 <-> n4a4j4

Interestingly, if I change my hash characters away from the list of 36 alphanumerics and replace it with just 0 and 1, I can encode/decode the date in binary:

19-Apr-2014 13:27:28 <-> 1397914048 <-> 1010011010100100111100111000000
19-Apr-2014 13:27:29 <-> 1397914049 <-> 1010011010100100111100111000001
19-Apr-2014 13:27:30 <-> 1397914050 <-> 1010011010100100111100111000010
19-Apr-2014 13:27:31 <-> 1397914051 <-> 1010011010100100111100111000011
19-Apr-2014 13:27:32 <-> 1397914052 <-> 1010011010100100111100111000100
19-Apr-2014 13:27:33 <-> 1397914053 <-> 1010011010100100111100111000101
19-Apr-2014 13:27:34 <-> 1397914054 <-> 1010011010100100111100111000110
19-Apr-2014 13:27:35 <-> 1397914055 <-> 1010011010100100111100111000111
19-Apr-2014 13:27:36 <-> 1397914056 <-> 1010011010100100111100111001000
19-Apr-2014 13:27:37 <-> 1397914057 <-> 1010011010100100111100111001001

Here’s the code to generate / decode Unix epoch timestamps in Python:

hashchars='0123456789abcdefghijklmnopqrstuvwxyz'
#hashchars='01' #for binary

def epochToHash(n):
  hash=''
  while n>0:
    hash = hashchars[int(n % len(hashchars))] + hash
    n = int(n / len(hashchars))
  return hash

def epochFromHash(s):
  s=s[::-1]
  epoch=0
  for pos in range(len(s)):
    epoch+=hashchars.find(s[pos])*(len(hashchars)**pos)
  return epoch

import time
t=int(time.time())
for i in range(10):
  t=t+1
  print(time.strftime("%d-%b-%Y %H:%M:%S", time.gmtime(t)),
              "<->", t,"<->",epochToHash(t))




Additional Resources

To maintain high frequency stability, RF oscillator circuits are sometimes “ovenized” where their temperature is raised slightly above ambient room temperature and held precisely at one temperature. Sometimes just the crystal is heated (with a “crystal oven”), and other times the entire oscillator circuit is heated. The advantage of heating the circuit is that other components (especially metal core instructors) are temperature sensitive. Googling for the phrase “crystal oven”, you’ll find no shortage of recommended circuits. Although a more complicated PID (proportional-integral-derivative) controller may seem enticing for these situations, the fact that the enclosure is so well insulated and drifts so little over vast periods of time suggests that it might not be the best application of a PID controller. One of my favorite write-ups is from M0AYF’s site which describes how to build a crystal oven for QRSS purposes. He demonstrates the MK1 and then the next design the MK2 crystal oven controller.  Here are his circuits:

Briefly, desired temperature is set with a potentiometer. An operational amplifier (op-amp) compares the target temperature with measured temperature (using a thermistor – a resistor which varies resistance by tempearture). If the measured temperature is below the target, the op-amp output goes high, and current flows through heating resistors. There are a few differences between the two circuits, but one of the things that struck me as different was the use of negative feedback with the operational amplifier. This means that rather than being on or off (like the air conditioning in your house), it can be on a little bit. I wondered if this would greatly affect frequency stability. In the original circuit, he mentions

The oven then cycles on and off roughly every thirty or forty seconds and hovers around 40 degrees-C thereafter to within better than one degree-C.

I wondered how much this on/off heater cycle affected temperature. Is it negligible, or could it affect frequency of an oscillator circuit? Indeed his application heats an entire enclosure so small variations get averaged-out by the large thermal mass. However in crystal oven designs where only the crystal is heated, such as described by Bill (W4HBK), I’ll bet the effect is much greater. Compare the thermal mass of these two concepts.

How does the amount of thermal mass relate to how well it can be controlled? How important is negative feedback for partial-on heater operation? Can simple ON/OFF heater regulation adequately stabalize a crystal or enclosure? I’d like to design my own heater, pulling the best elements from the rest I see on the internet. My goals are:

  1. use inexpensive thermistors instead of linear temperature sensors (like LM335)
  2. use inexpensive quarter-watt resistors as heaters instead of power resistors
  3. be able to set temperature with a knob
  4. be able to monitor temperature of the heater
  5. be able to monitor power delivered to the heater
  6. maximum long-term temperature stability

Right off the bat, I realized that this requires a PC interface. Even if it’s not used to adjust temperature (an ultimate goal), it will be used to log temperature and power for analysis. I won’t go into the details about how I did it, other than to say that I’m using an ATMEL ATMega8 AVR microcontroller and ten times I second I sample voltage on each of it’s six 10-bit ADC pins (PC0-PC5), and send that data to the computer with USART using an eBay special serial/USB adapter based on FTDI. They’re <$7 (shipped) and come with the USB cable. Obviously in a consumer application I’d etch boards and use the SMT-only FTDI chips, but for messing around at home I a few a few of these little adapters. They’re convenient as heck because I can just add a heater to my prototype boards and it even supplies power and ground. Convenient, right? Power is messier than it could be because it’s being supplied by the PC, but for now it gets the job done. On the software side, Python with PySerial listens to the serial port and copies data to a large numpy array, saving it every once and a while. Occasionally a bit is sent wrong and a number is received incorrectly (maybe one an hour), but the error is recognized and eliminated by the checksum (just the sum of all transmitted numbers). Plotting is done with numpy and matpltolib. Code for all of that is at the bottom of this post.

That’s the data logger circuit I came up with. Reading six channels ten times a second, it’s more than sufficient for voltage measurement. I went ahead and added an op-amp to the board too, since I knew I’d be using one. I dedicated one of the channels to serve as ambient temperature measurement. See the little red thermistor by the blue resistor? I also dedicated another channel to the output of the op-amp. This way I can measure drive to whatever temperature controller circuity I choose to use down the road. For my first test, I’m using a small thermal mass like one would in a crystal oven. Here’s how I made that:

I then build the temperature controller part of the circuit. It’s pretty similar to that previously published. it uses a thermistor in a voltage divider configuration to sense temperature. It uses a trimmer potentiometer to set temperature. An LED indicator light gives some indication of on/off, but keep in mind that a fraction of a volt will turn the Darlington transistor (TIP122) on slightly although it doesn’t reach a level high enough to drive the LED. The amplifier by default is set to high gain (55x), but can be greatly lowered (negative gain actually) with a jumper. This lets me test how important gain is for the circuitry.

controller

When using a crystal oven configuration, I concluded high high gain (cycling the heater on/off) is a BAD idea. While average temperature is held around the same, the crystal oscillates. This is what is occurring above when M0AYF indicates his MK1 heater turns on and off every 40 seconds. While you might be able to get away with it while heating a chassis or something, I think it’s easy to see it’s not a good option for crystal heaters. Instead, look at the low gain (negative gain) configuration. It reaches temperature surprisingly quickly and locks to it steadily. Excellent.

high gain
high gain configuration tends to oscillate every 30 seconds
low gain / negative gain configuration is extremely stable
low gain / negative gain configuration is extremely stable (fairly high temperature)
Here's a similar experiment with a lower target temperature. Noise is due to unregulated USB power supply / voltage reference. Undeniably, this circuit does not oscillate much if any.
Here’s a similar experiment with a lower target temperature. Noise is due to unregulated USB power supply / voltage reference. Undeniably, this circuit does not oscillate much if any.

Clearly low (or negative) gain is best for crystal heaters. What about chassis / enclosure heaters? Let’s give that a shot. I made an enclosure heater with the same 2 resistors. Again, I’m staying away from expensive components, and that includes power resistors. I used epoxy (gorilla glue) to cement them to the wall of one side of the enclosure.

I put a “heater sensor” thermistor near the resistors on the case so I could get an idea of the heat of the resistors, and a “case sensor” on the opposite side of the case. This will let me know how long it takes the case to reach temperature, and let me compare differences between using near vs. far sensors (with respect to the heating element) to control temperature. I ran the same experiments and this is what I came up with!

heater temperature (blue) and enclosure temperature (green) with low gain (first 20 minutes), then high gain (after) operation. High gain sensor/feedback loop is sufficient to induce oscillation, even with the large thermal mass of the enclosure
CLOSE SENSOR CONTROL, LOW/HIGH GAIN: TOP: heater temperature (blue) and enclosure temperature (green) with low gain (first 20 minutes), then high gain (after) operation. High gain sensor/feedback loop is sufficient to induce oscillation, even with the large thermal mass of the enclosure. BOTTOM: power to the heater (voltage off the op-amp output going into the base of the Darlington transistor). Although I didn’t give the low-gain configuration time to equilibrate, I doubt it would have oscillated on a time scale I am patient enough to see. Future, days-long experimentation will be required to determine if it oscillates significantly.
Even with the far sensor (opposite side of the enclosure as the heater) driving the operational amplifier in high gain mode, oscillations occur. Due to the larger thermal mass and increased distance the heat must travel to be sensed they take much longer to occur, leading them to be slower and larger than oscillations seen earlier when the heater was very close to the sensor.
FAR SENSOR CONTROL, HIGH GAIN: Even with the far sensor (opposite side of the enclosure as the heater) driving the operational amplifier in high gain mode, oscillations occur. Blue is the far sensor temperature. Green is the sensor near the heater temperature. Due to the larger thermal mass and increased distance the heat must travel to be sensed they take much longer to occur, leading them to be slower and larger than oscillations seen earlier when the heater was very close to the sensor.

Right off the bat, we observe that even with the increased thermal mass of the entire enclosure (being heated with two dinky 100 ohm 1/4 watt resistors) the system is prone to temperature oscillation if gain is set too high. For me, this is the final nail in the coffin – I will never use a comparator-type high gain sensor/regulation loop to control heater current. With that out, the only thing to compare is which is better: placing the sensor near the heating element, or far from it. In reality, with a well-insulated device like I seem to have, it seems like it doesn’t make much of a difference! The idea is that by placing it near the heater, it can stabilize quickly. However, placing it far from the heater will give it maximum sensation of “load” temperature. Anywhere in-between should be fine. As long as it’s somewhat thermally coupled to the enclosure, enclosure temperature will pull it slightly away from heater temperature regardless of location. Therefore, I conclude it’s not that critical where the sensor is placed, as long as it has good contact with the enclosure. Perhaps with long-term study (on the order of hours to days) slow oscillations may emerge, but I’ll have to build it in a more permanent configuration to test it out. Lucky, that’s exactly what I plan to do, so check back a few days from now!

Since the data speaks for itself, I’ll be concise with my conclusions:

  • two 1/4 watt 100 Ohm resistors in parallel (50 ohms) are suitable to heat an insulated enclosure with 12V
  • two 1/4 watt 100 Ohm resistors in parallel (50 ohms) are suitable to heat a crystal with 5V
  • low gain or negative gain is preferred to prevent oscillating tempeartures
  • Sensor location on an enclosure is not critical as long as it’s well-coupled to the enclosure and the entire enclosure is well-insulated.

I feel satisfied with today’s work. Next step is to build this device on a larger scale and fix it in a more permanent configuration, then leave it to run for a few weeks and see how it does. On to making the oscillator! If you have any questions or comments, feel free to email me. If you recreate this project, email me! I’d love to hear about it.

Here’s the code that went on the ATMega8 AVR (it continuously transmits voltage measurements on 6 channels).

#define F_CPU 8000000UL
#include <avr/io.h>
#include <util/delay.h>
#include <avr/interrupt.h>

/*
8MHZ: 300,600,1200,2400,4800,9600,14400,19200,38400
1MHZ: 300,600,1200,2400,4800
*/
#define USART_BAUDRATE 38400
#define BAUD_PRESCALE (((F_CPU / (USART_BAUDRATE * 16UL))) - 1)

/*
ISR(ADC_vect)
{
    PORTD^=255;
}
*/

void USART_Init(void){
	UBRRL = BAUD_PRESCALE;
	UBRRH = (BAUD_PRESCALE >> 8);
	UCSRB = (1<<TXEN);
	UCSRC = (1<<URSEL)|(1<<UCSZ1)|(1<<UCSZ0); // 9N1
}

void USART_Transmit( unsigned char data ){
	while ( !( UCSRA & (1<<UDRE)) );
	UDR = data;
}

void sendNum(long unsigned int byte){
	if (byte==0){
		USART_Transmit(48);
	}
	while (byte){
		USART_Transmit(byte%10+48);
		byte-=byte%10;
		byte/=10;
	}
}

int readADC(char adcn){
	ADMUX = 0b0100000+adcn;
	ADCSRA |= (1<<ADSC); // reset value
	while (ADCSRA & (1<<ADSC)) {}; // wait for measurement
	return ADC>>6;
}

int sendADC(char adcn){
	int val;
	val=readADC(adcn);
	sendNum(val);
	USART_Transmit(',');
	return val;
}

int main(void){
	ADCSRA = (1<<ADEN)  | 0b111;
	DDRB=255;
	USART_Init();
	int checksum;

	for(;;){
		PORTB=255;
		checksum=0;
		checksum+=sendADC(0);
		checksum+=sendADC(1);
		checksum+=sendADC(2);
		checksum+=sendADC(3);
		checksum+=sendADC(4);
		checksum+=sendADC(5);
		sendNum(checksum);
		USART_Transmit('n');
		PORTB=0;
		_delay_ms(200);
	}
}

Here’s the command I used to compile the code, set the AVR fuse bits, and load it to the AVR.

del *.elf
del *.hex
avr-gcc -mmcu=atmega8 -Wall -Os -o main.elf main.c -w
pause
cls
avr-objcopy -j .text -j .data -O ihex main.elf main.hex
avrdude -c usbtiny -p m8 -F -U flash:w:"main.hex":a -U lfuse:w:0xe4:m -U hfuse:w:0xd9:m

Here’s the code that runs on the PC to listen to the microchip, match the data to the checksum, and log it occasionally. 

import serial, time
import numpy
ser = serial.Serial("COM16", 38400, timeout=100)

line=ser.readline()[:-1]
t1=time.time()
lines=0

data=[]

def adc2R(adc):
    Vo=adc*5.0/1024.0
    Vi=5.0
    R2=10000.0
    R1=R2*(Vi-Vo)/Vo
    return R1

while True:
    line=ser.readline()[:-1]
    lines+=1
    if "," in line:
        line=line.split(",")
        for i in range(len(line)):
            line[i]=int(line[i][::-1])

    if line[-1]==sum(line[:-1]):
        line=[time.time()]+line[:-1]
        print lines, line
        data.append(line)
    else:
        print  lines, line, "<-- FAIL"

    if lines%50==49:
        numpy.save("data.npy",data)
        print "nSAVINGn%d lines in %.02f sec (%.02f vals/sec)n"%(lines,
            time.time()-t1,lines/(time.time()-t1))

Here’s the code that runs on the PC to graph data.

import matplotlib
matplotlib.use('TkAgg') # <-- THIS MAKES IT FAST!
import numpy
import pylab
import datetime
import time

def adc2F(adc):
    Vo=adc*5.0/1024.0
    K=Vo*100
    C=K-273
    F=C*(9.0/5)+32
    return F

def adc2R(adc):
    Vo=adc*5.0/1024.0
    Vi=5.0
    R2=10000.0
    R1=R2*(Vi-Vo)/Vo
    return R1

def adc2V(adc):
    Vo=adc*5.0/1024.0
    return Vo

if True:
    print "LOADING DATA"
    data=numpy.load("data.npy")
    data=data
    print "LOADED"

    fig=pylab.figure()
    xs=data[:,0]
    tempAmbient=data[:,1]
    tempPower=data[:,2]
    tempHeater=data[:,3]
    tempCase=data[:,4]
    dates=(xs-xs[0])/60.0
    #dates=[]
    #for dt in xs: dates.append(datetime.datetime.fromtimestamp(dt))

    ax1=pylab.subplot(211)
    pylab.title("Temperature Controller - Low Gain")
    pylab.ylabel('Heater (ADC)')
    pylab.plot(dates,tempHeater,'b-')
    pylab.plot(dates,tempCase,'g-')
    #pylab.axhline(115.5,color="k",ls=":")

    #ax2=pylab.subplot(312,sharex=ax1)
    #pylab.ylabel('Case (ADC)')
    #pylab.plot(dates,tempCase,'r-')
    #pylab.plot(dates,tempAmbient,'g-')
    #pylab.axhline(0,color="k",ls=":")

    ax2=pylab.subplot(212,sharex=ax1)
    pylab.ylabel('Heater Power')
    pylab.plot(dates,tempPower)

    #fig.autofmt_xdate()
    pylab.xlabel('Elapsed Time (min)')

    pylab.show()

print "DONE"




Additional Resources

In an effort to resume previous work [A, B, C, D] on developing a crystal oven for radio frequency transmitter / receiver stabilization purposes, the first step for me was to create a device to accurately measure and log temperature. I did this with common, cheap components, and the output is saved to the computer (over 1,000 readings a second). Briefly, I use a LM335 precision temperature sensor ($0.70 on mouser) which outputs voltage with respect to temperature. It acts like a Zener diode where the breakdown voltage relates to temperature. 2.95V is 295K (Kelvin), which is 22ºC / 71ºF. Note that Kelvin is just ºC + 273.15 (the difference between freezing and absolute zero). My goal was to use the ADC of a microcontroller to measure the output. The problem is that my ADC (one of 6 built into the ATMEL ATMega8 microcontroller) has 10-bit resolution, reporting steps from 0-5V as values from 0-1024. Thus, each step represents 0.0049V (0.49ºC / 0.882ºF). While ~1ºF resolution might be acceptable for some temperature measurement or control applications, I want to see fractions of a degree because radio frequency crystal temperature stabilization is critical. Here’s a video overview.

This is the circuit came up with. My goal was to make it cheaply and what I had on hand. It could certainly be better (more stable, more precise, etc.) but this seems to be working nicely. The idea is that you set the gain (the ratio of R2/R1) to increase your desired resolution (so your 5V of ADC recording spans over just several ºF you’re interested in), then set your “base offset” temperature that will produce 0V. In my design, I adjusted so 0V was room temperature, and 5V (maximum) was body temperature. This way when I touched the sensor, I’d watch temperature rise and fall when I let go.  Component values are very non-critical. LM324 is powered 0V GND and +5V Vcc. I chose to keep things simple and use a single rail power supply. It is worth noting that I ended-up using a 3.5V Zener diode for the positive end of the potentiometer rather than 5V.  If your power supply is well regulated 5V will be no problem, but as I was powering this with USB I decided to go for some extra stability by using a Zener reference.

precision thermometer LM335 LM324 microcontroller

 

On the microcontroller side, analog-to-digital measurement is summed-up pretty well in the datasheet. There is a lot of good documentation on the internet about how to get reliable, stable measurements. Decoupling capacitors, reference voltages, etc etc. That’s outside the scope of today’s topic. In my case, the output of the ADC went into the ATMega8 ADC5 (PC5, pin 28). Decoupling capacitors were placed at ARef and AVcc, according to the datasheet. Microcontroller code is at the bottom of this post.

To get the values to the computer, I used the USART capability of my microcontroller and sent ADC readings (at a rate over 1,000 a second) over a USB adapter based on an FTDI FT232 chip. I got e-bay knock-off FTDI evaluation boards which come with a USB cable too (they’re about $6, free shipping). Yeah, I could have done it cheaper, but this works effortlessly. I don’t use a crystal. I set fuse settings so the MCU runs at 8MHz, and thanks to the nifty online baud rate calculator determined I can use a variety of transfer speeds (up to 38400). At 1MHz (if DIV8 fuse bit is enabled) I’m limited to 4800 baud. Here’s the result, it’s me touching the sensor with my finger (heating it), then letting go.

finger touch
Touching the temperature sensor with my finger, voltage rose exponentially. When removed, it decayed exponentially – a temperature RC circuit, with capacitance being the specific heat capacity of the sensor itself. Small amounts of jitter are expected because I’m powering the MCU from unregulated USB +5V.

I spent a while considering fancy ways to send the data (checksums, frame headers, error correction, etc.) but ended-up just sending it old fashioned ASCII characters. I used to care more about speed, but even sending ASCII it can send over a thousand ADC readings a second, which is plenty for me. I ended-up throttling down the output to 10/second because it was just too much to log comfortable for long recordings (like 24 hours). In retrospect, it would have made sense to catch all those numbers and do averaging on the on the PC side.

I keep my house around 70F at night when I'm there, and you can see the air conditioning kick on and off. In the morning the AC was turned off for the day, temperature rose, and when I got back home I turned the AC on and it started to drop again.
I keep my house around 70F at night when I’m there, and you can see the air conditioning kick on and off. In the morning the AC was turned off for the day, temperature rose, and when I got back home I turned the AC on and it started to drop again.

On the receive side, I have nifty Python with PySerial ready to catch data coming from the microcontroller. It’s decoded, turned to values, and every 1000 receives saves a numpy array as a NPY binary file. I run the project out of my google drive folder, so while I’m at work I can run the plotting program and it loads the NPY file and shows it – today it allowed me to realize that my roomate turned off the air conditioning after I left, because I saw the temperature rising mid-day. The above graph is temperature in my house for the last ~24 hours. That’s about it! Here’s some of the technical stuff.

AVR ATMega8 microcontroller code:

#define F_CPU 8000000UL
#include <avr/io.h>
#include <util/delay.h>
#include <avr/interrupt.h>

/*
8MHZ: 300,600,1200,2400,4800,9600,14400,19200,38400
1MHZ: 300,600,1200,2400,4800
*/
#define USART_BAUDRATE 38400
#define BAUD_PRESCALE (((F_CPU / (USART_BAUDRATE * 16UL))) - 1)

/*
ISR(ADC_vect)
{
    PORTD^=255;
}
*/

void USART_Init(void){
	UBRRL = BAUD_PRESCALE;
	UBRRH = (BAUD_PRESCALE >> 8);
	UCSRB = (1<<TXEN);
	UCSRC = (1<<URSEL)|(1<<UCSZ1)|(1<<UCSZ0); // 9N1
}

void USART_Transmit( unsigned char data ){
	while ( !( UCSRA & (1<<UDRE)) );
	UDR = data;
}

void sendNum(long unsigned int byte){
	if (byte==0){
		USART_Transmit(48);
	}
	while (byte){
		USART_Transmit(byte%10+48);
		byte-=byte%10;
		byte/=10;
	}

}

unsigned int readADC(char adcn){
	ADMUX = 0b0100000+adcn;
	ADCSRA |= (1<<ADSC); // reset value
	while (ADCSRA & (1<<ADSC)) {}; // wait for measurement
	return ADC>>6;
}

void ADC_Init(){
	// ADC Enable, Prescaler 128
	ADCSRA = (1<<ADEN)  | 0b111;
}

int main(void){
	//DDRD=255;
	USART_Init();
	ADC_Init();
	for(;;){
		sendNum(readADC(5));
		USART_Transmit('n');
		_delay_ms(100);
	}
}

Here is the Python code to receive the data and log it to disk:

import serial, time
import numpy
ser = serial.Serial("COM15", 38400, timeout=100)

line=ser.readline()[:-1]
t1=time.time()
lines=0

data=[]

while True:
    line=ser.readline()[:-1]

    if "," in line:
        line=line.split(",")
        for i in range(len(line)):
            line[i]=line[i][::-1]
    else:
        line=[line[::-1]]
    temp=int(line[0])
    lines+=1
    data.append(temp)
    print "#",
    if lines%1000==999:
        numpy.save("DATA.npy",data)
        print
        print line
        print "%d lines in %.02f sec (%.02f vals/sec)"%(lines,
				time.time()-t1,lines/(time.time()-t1))

Here is the Python code to plot the data that has been saved:

import numpy
import pylab

data=numpy.load("DATA.npy")
print data
data=data*.008 #convert to F
xs=numpy.arange(len(data))/9.95  #vals/sec
xs=xs/60.0# minutes
xs=xs/60.0# hours

pylab.plot(xs,data)
pylab.grid(alpha=.5)
pylab.axis([None,None,0*.008,1024*.008])
pylab.ylabel(r'$Delta$ Fahrenheit')
pylab.xlabel("hours")
pylab.show()

If you recreate this project, or have any questions, feel free to email me!





Additional Resources

Here I demonstrate a dirt-cheap method of transmitting data from any microchip to any PC using $3.21 in parts.  I’ve had this idea for a while, but finally got it working tonight. On the transmit side, I’m having a an ATMEL AVR microcontroller (ATMega48) transmit data (every number from 0 to 200 over and over) wirelessly using 433mhz wireless modules. The PC receives the data through the microphone port of a sound card, and a cross-platform Python script I wrote decodes the data from the audio and graphs it on the screen. I did something similar back in 2011, but it wasn’t wireless, and the software wasn’t nearly as robust as it is now.

This is a proof-of-concept demonstration, and part of a larger project. I think there’s a need for this type of thing though! It’s unnecessarily hard to transfer data from a MCU to a PC as it is. There’s USB (For AVR V-USB is a nightmare and requires a precise, specific clock speed, DIP chips don’t have native USB, and some PIC DIP chips do but then you have to go through driver hell), USART RS-232 over serial port works (but who has serial ports these days?), or USART over USB RS-232 interface chips (like FTDI FT-232, but surface mount only), but both also require precise, specific clock speeds. Pretend I want to just measure temperature once a minute. Do I really want to etch circuit boards and solder SMT components? Well, kinda, but I don’t like feeling forced to. Some times you just want a no-nonsense way to get some numbers from your microchip to your computer. This project is a funky out-of-the-box alternative to traditional methods, and one that I hope will raise a few eyebrows.

 

Ultimately, I designed this project to eventually allow multiple “bursting” data transmitters to transmit on the same frequency routinely, thanks to syncing and forced-sync-loss (read on). It’s part of what I’m tongue-in-cheek calling the Scott Harden RF Protocol (SH-RFP). In my goal application, I wish to have about 5 wireless temperature sensors all transmitting data to my PC.  The receive side has some error checking in that it makes sure pulse sizes are intelligent and symmetrical (unlike random noise), and since each number is sent twice (with the second time being in reverse), there’s another layer of error-detection.  This is *NOT* a robust and accurate method to send critical data. It’s a cheap way to send data. It is very range limited, and only is intended to work over a distance of ten or twenty feet. First, let’s see it in action!

The RF modules are pretty simple. At 1.56 on ebay (with free shipping), they’re cheap too! I won’t go into detail documenting the ins and out of these things (that’s done well elsewhere). Briefly, you give them +5V (VCC), 0V (GND), and flip their data pin (ATAD) on and off on the transmitter module, and the receiver module’s DATA pin reflects the same state. The receiver uses a gain circuit which continuously increases gain until signal is detected, so if you’re not transmitting it WILL decode noise and start flipping its output pin. Note that persistent high or low states are prone to noise too, so any protocol you use these things for should have rapid state transitions. It’s also suggested that you maintain an average 50% duty cycle. These modules utilize amplitude shift keying (ASK) to transmit data wirelessly. The graphic below shows what that looks like at the RF level. Transmit and receive is improved by adding a quarter-wavelength vertical antenna to the “ANT” solder pad. At 433MHz, that is about 17cm, so I’m using a 17cm copper wire as an antenna.

Transmitting from the microcontroller is easy as pie! It’s just a matter of copying-in a few lines of C.  It doesn’t rely on USART, SPI, I2C, or any other protocol. Part of why I developed this method is because I often use ATTiny44A which doesn’t have USART for serial interfacing. The “SH-RFP” is easy to implement just by adding a few lines of code. I can handle that.  How does it work? I can define it simply by a few rules:

 

SHRFP (Scott Harden RF Protocol)

Pulses can be one of 3 lengths: A (0), B (1), or C (break).

Each pulse represents high, then low of that length.

Step 1: prime synchronization by sending ten ABCs

Step 2: indicate we’re starting data by sending C.

Step 3: for each number you want to send:

A: send your number bit by bit (A=0, B=1)

B: send your number bit by bit (A=1, B=0)

C: indicate number end by sending C.

 Step 4: tell PC to release the signal by sending ten Cs.

Decoding is the same thing in reverse. I use an eBay sound card at $1.29 (with free shipping) to get the signal into the PC. Syncronization is required to allow the PC to know that real data (not noise) is starting. Sending the same number twice (once with reversed bit polarity) is a proofchecking mechanisms that lets us throw-out data that isn’t accurate.

From a software side, I’m using PyAudio to collect data from the sound card, and the PythonXY distribution to handle analysis with numpy, scipy, and plotting with QwtPlot, and general GUI functionality with PyQt. I think that’s about everything.

 The demonstration interface is pretty self-explanatory. The top-right shows a sample piece of data. The top left is a histogram of the number of samples of each pulse width. A clean signal should have 3 pulses (A=0, B=1, C=break). Note that you’re supposed to look at the peaks to determine the best lengths to tell the software to use to distinguish A, B, and C. This was intentionally not hard-coded because I want to rapidly switch from one microcontroller platform to another which may be operating at a different clock speed, and if all the sudden it’s running 3 times slower it will be no problem to decide on the PC side. Slick, huh? The bottom-left shows data values coming in. The bottom-right graphs those values. Rate reporting lets us know that I’m receiving over 700 good data points a second. That’s pretty cool, especially considering I’m recording at 44,100 Hz. 

Here’s the MCU code I used. It’s an ATMega48 ATMEL AVR microcontroller. Easy code!

#define F_CPU 8000000UL

#include <avr/io.h>
#include <util/delay.h>

void tick(char ticks){
	while (ticks>0){
		_delay_us(100);
		ticks--;
	}
}

void pulse(char ticks){
	PORTB=255;
	tick(ticks);
	PORTB=0;
	tick(ticks);
}

void send_sync(){
	char i;
	for (i=0;i<10;i++){
		pulse(1);
		pulse(2);
		pulse(3);
	}
	pulse(3);
}

void send_lose(){
	char i;
	for (i=0;i<5;i++){
		pulse(3);
	}
}

void sendByte(int val){
	// TODO - make faster by only sending needed bytes
	char i;
	for (i=0;i<8;i++){
		if ((val>>i)&1){pulse(2);}
		else{pulse(1);}
	}
}

void send(int val){
	sendByte(val);  // regular
	sendByte(~val); // inverted
	pulse(3);
}

int main (void)
{
    DDRB = 255;
	int i;

    while(1) {
		send_sync();
		for (i=0;i<200;i++){
			send(i);
		}
		send_lose();
	}
}

Here’re some relevant snippits of the PC code. Download the full project below if you’re interested.

import matplotlib
matplotlib.use('TkAgg') # -- THIS MAKES IT FAST!
import numpy
import pyaudio
import threading
import pylab
import scipy
import time
import sys

class SwhRecorder:
    """Simple, cross-platform class to record from the microphone.
    This version is optimized for SH-RFP (Scott Harden RF Protocol)
    Pulse data extraction. It's dirty, but it's easy and it works.

    BIG PICTURE:
    continuously record sound in buffers.
    if buffer is detected:

        ### POPULATE DELAYS[] ###
        downsample data
        find Is where data>0
        use ediff1d to get differences between Is
        append >1 values to delays[]
        --if the old audio[] ended high, figure out how many
        --on next run, add that number to the first I

        ### PLUCK DELAYS, POPULATE VALUES ###
        only analyze delays through the last 'break'
        values[] is populated with decoded delays.

    ."""

    def __init__(self):
        """minimal garb is executed when class is loaded."""
        self.RATE=44100
        self.BUFFERSIZE=2**10
        print "BUFFER:",self.BUFFERSIZE
        self.threadsDieNow=False
        self.newAudio=[]
        self.lastAudio=[]
        self.SHRFP=True
        self.dataString=""
        self.LEFTOVER=[]

        self.pulses=[]
        self.pulsesToKeep=1000

        self.data=[]
        self.dataToKeep=1000

        self.SIZE0=5
        self.SIZE1=10
        self.SIZE2=15
        self.SIZEF=3

        self.totalBits=0
        self.totalNumbers=0
        self.totalSamples=0
        self.totalTime=0

        self.nothingNewToShow=True

    def setup(self):
        """initialize sound card."""
        #TODO - windows detection vs. alsa or something for linux
        #TODO - try/except for sound card selection/initiation
        self.p = pyaudio.PyAudio()
        self.inStream = self.p.open(input_device_index=None,
                                    format=pyaudio.paInt16,channels=1,
                                    rate=self.RATE,input=True,
                                    frames_per_buffer=self.BUFFERSIZE)

    def close(self):
        """cleanly back out and release sound card."""
        self.p.close(self.inStream)

    def decodeBit(self,s):
        "given a good string 1001101001 etc, return number or None"
        if len(s)<2:return -2
        s=s[::-1]
        A=s[:len(s)/2] #INVERTED
        A=A.replace("0","z").replace("1","0").replace("z","1")
        B=s[len(s)/2:] #NORMAL

        if A<>B:
            return -1
        else:
            return int(A,2)

    def analyzeDataString(self):
        i=0
        bit=""
        lastB=0
        while i<len(self.dataString):
            if self.dataString[i]=="B":
                self.data.append(self.decodeBit(bit))
                self.totalNumbers+=1
                lastB=i
            if self.dataString[i] in ['B','?']:
                bit=""
            else:
                bit+=self.dataString[i]
            i+=1
        self.dataString=self.dataString[lastB+1:]
        if len(self.data)>self.dataToKeep:
            self.data=self.data[-self.dataToKeep:]

    def continuousAnalysis(self):
        """keep watching newAudio, and process it."""
        while True:
            while len(self.newAudio)< self.BUFFERSIZE:
                time.sleep(.1)

            analysisStart=time.time()

            audio=self.newAudio

            # TODO - insert previous audio sequence here

            # GET Is where data is positive
            Ipositive=numpy.nonzero(audio>0)[0]
            diffs=numpy.ediff1d(Ipositive)
            Idiffs=numpy.where(diffs>1)[0]
            Icross=Ipositive[Idiffs]
            pulses=diffs[Idiffs]

            # remove some of the audio buffer, leaving the overhang

            if len(Icross)>0:
                processedThrough=Icross[-1]+diffs[Idiffs[-1]]
            else:
                processedThrough=len(audio)

            self.lastAudio=self.newAudio[:processedThrough]
            self.newAudio=self.newAudio[processedThrough:]

            if False:
                # chart audio data (use it to check algorythm)
                pylab.plot(audio,'b')
                pylab.axhline(0,color='k',ls=':')

                for i in range(len(Icross)):
                    # plot each below-zero pulse whose length is measured
                    pylab.axvspan(Icross[i],Icross[i]+diffs[Idiffs[i]],
                                  color='b',alpha=.2,lw=0)

                # plot the hangover that will be carried to next chunk
                pylab.axvspan(Icross[i]+diffs[Idiffs[i]],len(audio),
                              color='r',alpha=.2)
                pylab.show()
                return

            # TODO - histogram of this point to assess quality
            s=''
            for pulse in pulses:
                if (self.SIZE0-self.SIZEF)<pulse<(self.SIZE0+self.SIZEF):
                    s+="0"
                elif (self.SIZE1-self.SIZEF)<pulse<(self.SIZE1+self.SIZEF):
                    s+="1"
                elif (self.SIZE2-self.SIZEF)<pulse<(self.SIZE2+self.SIZEF):
                    s+="B"
                else:
                    s+="?"

            self.pulses=pulses
            self.totalBits+=len(pulses)

            print "[%.02f ms took %.02f ms] T: 0=%d 1=%d B=%d ?=%d"%(
                          len(audio)*1000.0/self.RATE,
                          time.time()-analysisStart,
                          s.count('0'),s.count('1'),s.count('B'),s.count('?'))

            self.dataString+=s
            self.analyzeDataString()

            self.totalSamples+=self.BUFFERSIZE
            self.totalTime=self.totalSamples/float(self.RATE)
            self.totalBitRate=self.totalBits/self.totalTime
            self.totalDataRate=self.totalNumbers/self.totalTime

            self.nothingNewToShow=False

    def continuousRecord(self):
        """record forever, adding to self.newAudio[]. Thread this out."""
        while self.threadsDieNow==False:
            maxSecBack=5
            while len(self.newAudio)>(maxSecBack*self.RATE):
                print "DELETING NEW AUDIO!"
                self.newAudio=self.newAudio[self.BUFFERSIZE:]
            audioString=self.inStream.read(self.BUFFERSIZE)
            audio=numpy.fromstring(audioString,dtype=numpy.int16)
            self.newAudio=numpy.hstack((self.newAudio,audio))

    def continuousDataGo(self):
        self.t = threading.Thread(target=self.continuousRecord)
        self.t.start()
        self.t2 = threading.Thread(target=self.continuousAnalysis)
        self.t2.start()

    def continuousEnd(self):
        """shut down continuous recording."""
        self.threadsDieNow=True

if __name__ == "__main__":
    SHR=SwhRecorder()
    SHR.SHRFP_decode=True
    SHR.setup()
    SHR.continuousDataGo()

    #SHR.DataStart()

    print "---DONE---"

Finally, if you’re interested, here’s the full code (and demo audio WAV files):

DOWNLOAD: SCOTT HARDEN RF PROTOCOL DEMO.zip

If you use these concepts, hardware, or ideas in your project, let me know about it! Send me an email showing me your project – I’d love to see it. Good luck!





Additional Resources

WARNING: this project is largely outdated, and some of the modules are no longer supported by modern distributions of Python.

For a more modern, cleaner, and more complete GUI-based viewer of realtime audio data (and the FFT frequency data), check out my Python Real-time Audio Frequency Monitor project.

I’m no stranger to visualizing linear data in the frequency-domain. Between the high definition spectrograph suite I wrote in my first year of dental school (QRSS-VD, which differentiates tones to sub-Hz resolution), to the various scripts over the years (which go into FFT imaginary number theory, linear data signal filtering with python, and real time audio graphing with wckgraph), I’ve tried dozens of combinations of techniques to capture data, analyze it, and display it with Python. Because I’m now branching into making microcontroller devices which measure and transfer analog data to a computer, I need a way to rapidly visualize data obtained in Python. Since my microcontroller device isn’t up and running yet, linear data from a PC microphone will have to do.  Here’s a quick and dirty start-to-finish project anyone can tease apart to figure out how to do some of these not-so-intuitive processes in Python. To my knowledge, this is a cross-platform solution too. For the sound card interaction, it relies on the cross-platform sound card interface library PyAudio. My python distro is 2.7 (python xy), but pythonxy doesn’t [yet] supply PyAudio.

The code behind it is a little jumbled, but it works. For recording, I wrote a class “SwhRecorder” which uses threading to continuously record audio and save it as a numpy array. When the class is loaded and started, your GUI can wait until it sees newAudio become True, then it can grab audio directly, or use fft() to pull the spectral component (which is what I do in the video). Note that my fft() relies on numpy.fft.fft(). The return is a nearly-symmetrical mirror image of the frequency components, which (get ready to cringe mathematicians) I simply split into two arrays, reverse one of them, and add together. To turn this absolute value into dB, I’d take the log10(fft) and multiply it by 20. You know, if you’re into that kind of thing, you should really check out a post I made about FFT theory and analyzing audio data in python.

Here’s the meat of the code. To run it, you should really grab the zip file at the bottom of the page. I’ll start with the recorder class:

import matplotlib
matplotlib.use('TkAgg') # THIS MAKES IT FAST!
import numpy
import scipy
import struct
import pyaudio
import threading
import pylab
import struct

class SwhRecorder:
    """Simple, cross-platform class to record from the microphone."""

    def __init__(self):
        """minimal garb is executed when class is loaded."""
        self.RATE=48100
        self.BUFFERSIZE=2**12 #1024 is a good buffer size
        self.secToRecord=.1
        self.threadsDieNow=False
        self.newAudio=False

    def setup(self):
        """initialize sound card."""
        #TODO - windows detection vs. alsa or something for linux
        #TODO - try/except for sound card selection/initiation

        self.buffersToRecord=int(self.RATE*self.secToRecord/self.BUFFERSIZE)
        if self.buffersToRecord==0: self.buffersToRecord=1
        self.samplesToRecord=int(self.BUFFERSIZE*self.buffersToRecord)
        self.chunksToRecord=int(self.samplesToRecord/self.BUFFERSIZE)
        self.secPerPoint=1.0/self.RATE

        self.p = pyaudio.PyAudio()
        self.inStream = self.p.open(format=pyaudio.paInt16,channels=1,
            rate=self.RATE,input=True,frames_per_buffer=self.BUFFERSIZE)
        self.xsBuffer=numpy.arange(self.BUFFERSIZE)*self.secPerPoint
        self.xs=numpy.arange(self.chunksToRecord*self.BUFFERSIZE)*self.secPerPoint
        self.audio=numpy.empty((self.chunksToRecord*self.BUFFERSIZE),dtype=numpy.int16)

    def close(self):
        """cleanly back out and release sound card."""
        self.p.close(self.inStream)

    ### RECORDING AUDIO ###

    def getAudio(self):
        """get a single buffer size worth of audio."""
        audioString=self.inStream.read(self.BUFFERSIZE)
        return numpy.fromstring(audioString,dtype=numpy.int16)

    def record(self,forever=True):
        """record secToRecord seconds of audio."""
        while True:
            if self.threadsDieNow: break
            for i in range(self.chunksToRecord):
                self.audio[i*self.BUFFERSIZE:(i+1)*self.BUFFERSIZE]=self.getAudio()
            self.newAudio=True
            if forever==False: break

    def continuousStart(self):
        """CALL THIS to start running forever."""
        self.t = threading.Thread(target=self.record)
        self.t.start()

    def continuousEnd(self):
        """shut down continuous recording."""
        self.threadsDieNow=True

    ### MATH ###

    def downsample(self,data,mult):
        """Given 1D data, return the binned average."""
        overhang=len(data)%mult
        if overhang: data=data[:-overhang]
        data=numpy.reshape(data,(len(data)/mult,mult))
        data=numpy.average(data,1)
        return data

    def fft(self,data=None,trimBy=10,logScale=False,divBy=100):
        if data==None:
            data=self.audio.flatten()
        left,right=numpy.split(numpy.abs(numpy.fft.fft(data)),2)
        ys=numpy.add(left,right[::-1])
        if logScale:
            ys=numpy.multiply(20,numpy.log10(ys))
        xs=numpy.arange(self.BUFFERSIZE/2,dtype=float)
        if trimBy:
            i=int((self.BUFFERSIZE/2)/trimBy)
            ys=ys[:i]
            xs=xs[:i]*self.RATE/self.BUFFERSIZE
        if divBy:
            ys=ys/float(divBy)
        return xs,ys

    ### VISUALIZATION ###

    def plotAudio(self):
        """open a matplotlib popup window showing audio data."""
        pylab.plot(self.audio.flatten())
        pylab.show()

And now here’s the GUI launcher:

import ui_plot
import sys
import numpy
from PyQt4 import QtCore, QtGui
import PyQt4.Qwt5 as Qwt
from recorder import *

def plotSomething():
    if SR.newAudio==False:
        return
    xs,ys=SR.fft()
    c.setData(xs,ys)
    uiplot.qwtPlot.replot()
    SR.newAudio=False

if __name__ == "__main__":
    app = QtGui.QApplication(sys.argv)

    win_plot = ui_plot.QtGui.QMainWindow()
    uiplot = ui_plot.Ui_win_plot()
    uiplot.setupUi(win_plot)
    uiplot.btnA.clicked.connect(plotSomething)
    #uiplot.btnB.clicked.connect(lambda: uiplot.timer.setInterval(100.0))
    #uiplot.btnC.clicked.connect(lambda: uiplot.timer.setInterval(10.0))
    #uiplot.btnD.clicked.connect(lambda: uiplot.timer.setInterval(1.0))
    c=Qwt.QwtPlotCurve()
    c.attach(uiplot.qwtPlot)

    uiplot.qwtPlot.setAxisScale(uiplot.qwtPlot.yLeft, 0, 1000)

    uiplot.timer = QtCore.QTimer()
    uiplot.timer.start(1.0)

    win_plot.connect(uiplot.timer, QtCore.SIGNAL('timeout()'), plotSomething)

    SR=SwhRecorder()
    SR.setup()
    SR.continuousStart()

    ### DISPLAY WINDOWS
    win_plot.show()
    code=app.exec_()
    SR.close()
    sys.exit(code)

Note that by commenting-out the FFT line and using “c.setData(SR.xs,SR.audio)” you can plot linear PCM data to visualize sound waves like this:

Finally, here’s the zip file. It contains everything you need to run the program on your own computer (including the UI scripts which are not written on this page)

DOWNLOADSWHRecorder.zip

If you make a cool project based on this one, I’d love to hear about it. Good luck!

 





Additional Resources

 

WARNING: this project is largely outdated, and some of the modules are no longer supported by modern distributions of Python.

For a more modern, cleaner, and more complete GUI-based viewer of realtime audio data (and the FFT frequency data), check out my Python Real-time Audio Frequency Monitor project.

I love using python for handing data. Displaying it isn’t always as easy. Python fast to write, and numpy, scipy, and matplotlib are an incredible combination. I love matplotlib for displaying data and use it all the time, but when it comes to realtime data visualization, matplotlib (admittedly) falls behind. Imagine trying to plot sound waves in real time. Matplotlib simply can’t handle it. I’ve recently been making progress toward this end with PyQwt with the Python X,Y distribution. It is a cross-platform solution which should perform identically on Windows, Linux, and MacOS. Here’s an example of what it looks like plotting some dummy data (a sine wave) being transformed with numpy.roll().

How did I do it? Easy. First, I made the GUI with QtDesigner (which comes with Python x,y). I saved the GUI as a .ui file. I then used the pyuic4 command to generate a python script from the .ui file. In reality, I use a little helper script I wrote designed to build .py files from .ui files and start a little “ui.py” file which imports all of the ui classes. It’s overkill for this, but I’ll put it in the ZIP anyway.  Here’s what the GUI looks like in QtDesigner:

 

After that, I tie everything together in a little script which updates the plot in real time. It takes inputs from button click events and tells a clock (QTimer) how often to update/replot the data. Replotting it involves just rolling it with numpy.roll().  Check it out:

import ui_plot #this was generated by pyuic4 command
import sys
import numpy
from PyQt4 import QtCore, QtGui
import PyQt4.Qwt5 as Qwt

numPoints=1000
xs=numpy.arange(numPoints)
ys=numpy.sin(3.14159*xs*10/numPoints) #this is our data

def plotSomething():
    global ys
    ys=numpy.roll(ys,-1)
    c.setData(xs, ys)
    uiplot.qwtPlot.replot()

if __name__ == "__main__":
    app = QtGui.QApplication(sys.argv)
    win_plot = ui_plot.QtGui.QMainWindow()
    uiplot = ui_plot.Ui_win_plot()
    uiplot.setupUi(win_plot)

    # tell buttons what to do when clicked
    uiplot.btnA.clicked.connect(plotSomething)
    uiplot.btnB.clicked.connect(lambda: uiplot.timer.setInterval(100.0))
    uiplot.btnC.clicked.connect(lambda: uiplot.timer.setInterval(10.0))
    uiplot.btnD.clicked.connect(lambda: uiplot.timer.setInterval(1.0))

    # set up the QwtPlot (pay attention!)
    c=Qwt.QwtPlotCurve()  #make a curve
    c.attach(uiplot.qwtPlot) #attach it to the qwtPlot object
    uiplot.timer = QtCore.QTimer() #start a timer (to call replot events)
    uiplot.timer.start(100.0) #set the interval (in ms)
    win_plot.connect(uiplot.timer, QtCore.SIGNAL('timeout()'), plotSomething)

    # show the main window
    win_plot.show()
    sys.exit(app.exec_())

I’ll put all the files in a ZIP to help out anyone interested in giving this a shot. Clicking different buttons updates the graph at different speeds. If you make something cool with this concept, let me know! I’d love to see it.

DOWNLOAD PROJECT: realtime_python_graph.zip





Additional Resources

I recently migrated to Python(x,y) and noticed my matplotlib graphs are resizing unacceptably slowly when I use the pan/zoom button. I’m quite a fan of numpy, scipy, matplotlib, the python imaging library (PIL), and GUI platforms like Tk/TkInter, pyGTK, and pyQT, but getting them all to play nicely is a sometimes pain. I’m considering migrating entirely to Python(x,y) because, as a single distribution, it’s designed to install all these libraries (and many more) in a compatible way out of the box. However, when I did, I noticed matplotlib graphs would resize, rescale, and drag around the axes very slowly. After a lot of digging on the interweb, I figured out what was going wrong. I’ll show you by plotting 20 random data points the slow way (left) then the fast way (right).

THE PROBLEM: See the difference between the two plots? The one on the left (SLOW!) uses the Qt4Agg backend, which renders the matplotlib plot on a QT4 canvas. This is slower than the one on the right, which uses the more traditional TkAgg backend to draw the plot on a Tk canvas with tkinter (FASTER!). Check out matplotlib’s official description of what a backend is and which ones you can use. When you just install Python and matplotlib, Tk is used by default. 

import numpy
import matplotlib
matplotlib.use('TkAgg') # <-- THIS MAKES IT FAST!
import pylab
pylab.plot(numpy.random.random_integers(0,100,20))
pylab.title("USING: "+matplotlib.get_backend())
pylab.show()

THE FIX: Tell matplotlib to stop using QT to draw the plot, and let it plot with Tk. This can be done immediately after importing matplotlib, but must be done before importing pylab using the line matplotlib.use('TkAgg'). Here’s the full example I used to generate the demonstration plots above. Change TkAgg to Qt4Agg (or comment-out the ‘use’ line if you’re using PythonXY) and you will see performance go down the tube. Alternatively, make a change to the matplotlib rc file to customize default behavior when the package is loaded.





Additional Resources

UPDATE: An improved ECG design was posted in August, 2016.
Check out: http://www.swharden.com/wp/2016-08-08-diy-ecg-with-1-op-amp/

Of the hundreds of projects I’ve shared over the years, none has attracted more attention than my DIY ECG machine on the cheap posted almost 4 years ago. This weekend I re-visited the project and made something I’m excited to share!  The original project was immensely popular, my first featured article on Hack-A-Day, and today “ECG” still represents the second most searched term by people who land on my site. My gmail account also has had 194 incoming emails from people asking details about the project. A lot of it was by frustrated students trying to recreate the project running into trouble because it was somewhat poorly documented. Clearly, it’s a project that a wide range of people are interested in, and I’m happy to revisit it bringing new knowledge and insight to the project. I will do my best to document it thoroughly so anyone can recreate it!

The goal of this project is to collect heartbeat information on a computer with minimal cost and minimal complexity.  I accomplished this with fewer than a dozen components (all of which can be purchased at RadioShack). It serves both as a light-based heartbeat monitor (similar to a pulse oximeter, though it’s not designed to quantitatively measure blood oxygen saturation), and an electrocardiogram (ECG) to visualize electrical activity generated by heart while it contracts. Let’s jump right to the good part – this is what comes out of the machine:

That’s my actual heartbeat. Cool, right? Before I go into how the circuit works, let’s touch on how we measure heartbeat with ECG vs. light (like a pulse oximeter).  To form a heartbeat, the pacemaker region of the heart (called the SA node, which is near the upper right of the heart) begins to fire and the atria (the two top chambers of the heart) contract. The SA node generates a little electrical shock which stimulated a synchronized contraction. This is exactly what defibrillators do when a heart has stopped beating. When a heart attack is occurring and a patient is undergoing ventricular fibrillation, it means that heart muscle cells are contracting randomly and not in unison, so the heart quivers instead of pumping as an organ. Defibrillators synchronize the heart beat with a sudden rush of current over the heart to reset all of the cells to begin firing at the same time (thanks Ron for requesting a more technical description).  If a current is run over the muscle, the cells (cardiomyocytes) all contract at the same time, and blood moves. The AV node (closer to the center of the heart) in combination with a slow conducting pathway (called the bundle of His) control contraction of the ventricles (the really large chambers at the bottom of the heart), which produce the really large spikes we see on an ECG.  To measure ECG, optimally we’d place electrodes on the surface of the heart. Since that would be painful, we do the best we can by measuring voltage changes (often in the mV range) on the surface of the skin. If we amplify it enough, we can visualize it. Depending on where the pads are placed, we can see different regions of the heart contract by their unique electrophysiological signature. ECG requires sticky pads on your chest and is extremely sensitive to small fluctuations in voltage. Alternatively, a pulse oximeter measures blood oxygenation and can monitor heartbeat by clipping onto a finger tip. It does this by shining light through your finger and measuring how much light is absorbed. This goes up and down as blood is pumped through your finger. If you look at the relationship between absorbency in the red vs. infrared wavelengths, you can infer the oxygenation state of the blood. I’m not doing that today because I’m mostly interested in detecting heart beats.

For operation as a pulse oximeter-type optical heartbeat detector (a photoplethysmograph which produces a photoplethysmogram), I use a bright red LED to shine light through my finger and be detected by a phototransistor (bottom left of the diagram). I talk about how this works in more detail in a previous post. Basically the phototransistor acts like a variable resistor which conducts different amounts of current depending on how much light it sees. This changes the voltage above it in a way that changes with heartbeats. If this small signal is used as the input, this device acts like a pulse oximeter.

For operation as an electrocardiograph (ECG), I attach the (in) directly to a lead on my chest. One of them is grounded (it doesn’t matter which for this circuit – if they’re switched the ECG just looks upside down), and the other is recording. In my original article, I used pennies with wires soldered to them taped to my chest as leads. Today, I’m using fancier sticky pads which are a little more conductive. In either case, one lead goes in the center of your chest, and the other goes to your left side under your arm pit. I like these sticky pads because they stick to my skin better than pennies taped on with electrical tape. I got 100 Nikomed Nikotabs EKG Electrodes 0315 on eBay for $5.51 with free shipping (score!). Just gator clip to them and you’re good to go!

In both cases, I need to build a device to amplify small signals. This is accomplished with the following circuit. The core of the circuit is an LM324 quad operational amplifier.  These chips are everywhere, and extremely cheap. It looks like Thai Shine sells 10 for $2.86 (with free shipping). That’s about a quarter each. Nice!  A lot of ECG projects use instrumentation amplifiers like the AD620 (which I have used with fantastic results), but these are expensive (about $5.00 each). The main difference is that instrumentation amplifiers amplify the difference between two points (which reduces noise and probably makes for a better ECG machine), but for today an operational amplifier will do a good enough job amplifying a small signal with respect to ground. I get around the noise issue by some simple filtering techniques. Let’s take a look at the circuit.

This project utilizes one of the op-amps as a virtual ground. One complaint of using op-amps in simple projects is that they often need + and – voltages. Yeah, this could be done with two 9V batteries to generate +9V and -9V, but I think it’s easier to use a single power source (+ and GND). A way to get around that is to use one of the op-amps as a current source and feed it half of the power supply voltage (VCC), and use the output as a virtual ground (allowing VCC to be your + and 0V GND to be your -). For a good description of how to do this intelligently, read the single supply op amps web page. The caveat is that your signals should remain around VCC/2, which can be done if it is decoupled by feeding it through a series capacitor. The project works at 12V or 5V, but was designed for (and has much better output) at 12V. The remaining 3 op-amps of the LM324 serve three unique functions:

STAGE 1: High gain amplifier. The input signals from either the ECG or pulse oximeter are fed into a chain of 3 opamp stages. The first is a preamplifier. The output is decoupled through a series capacitor to place it near VCC/2, and amplified greatly thanks to the 1.8Mohm negative feedback resistor. Changing this value changes initial gain.

STAGE 2: active low-pass filter. The 10kOhm variable resistor lets you adjust the frequency cutoff. The opamp serves as a unity gain current source / voltage follower that has high input impedance when measuring the output f the low-pass filter and reproduces its voltage with a low impedance output. There’s some more information about active filtering on this page. It’s best to look at the output of this stage and adjust the potentiometer until the 60Hz noise (caused by the AC wiring in the walls) is most reduced while the lower-frequency component of your heartbeat is retained. With the oximeter, virtually no noise gets through. Because the ECG signal is much smaller, this filter has to be less aggressive, and this noise is filtered-out by software (more on this later).

STAGE 3: final amplifier with low-pass filter. It has a gain of ~20 (determined by the ratio of the 1.8kOhm to 100Ohm resistors) and lowpass filtering components are provided by the 22uF capacitor across the negative feedback resistor. If you try to run this circuit at 5V and want more gain (more voltage swing), consider increasing the value of the 1.8kOhm resistor (wit the capacitor removed). Once you have a good gain, add different capacitor values until your signal is left but the noise reduced. For 12V, these values work fine. Let’s see it in action!

Now for the second half – getting it into the computer. The cheapest and easiest way to do this is to simply feed the output into a sound card! A sound card is an analog-to-digital converter (ADC) that everybody has and can sample up to 48 thousand samples a second! (overkill for this application) The first thing you should do is add an output potentiometer to allow you to drop the voltage down if it’s too big for the sound card (in the case of the oximeter) but but also allow full-volume in the case of sensitive measurements (like ECG). Then open-up sound editing software (I like GoldWave for Windows or Audacity for Linux, both of which are free) and record the input. You can do filtering (low-pass filter at 40Hz with a sharp cutoff) to further eliminate any noise that may have sneaked through. Re-sample at 1,000 Hz (1kHz) and save the output as a text file and you’re ready to graph it! Check it out.

Here are the results of some actual data recorded and processed with the method shown in the video. let’s look at the pulse oximeter first.

That looks pretty good, certainly enough for heartbeat detection. There’s obvious room for improvement, but as a proof of concept it’s clearly working. Let’s switch gears and look at the ECG. It’s much more challenging because it’s signal is a couple orders of magnitude smaller than the pulse oximeter, so a lot more noise gets through. Filtering it out offers dramatic improvements!

Here’s the code I used to generate the graphs from the text files that GoldWave saves. It requires Python, Matplotlib (pylab), and Numpy. In my case, I’m using 32-bit 2.6 versions of everything.

# DIY Sound Card ECG/Pulse Oximeter
# by Scott Harden (2013) http://www.SWHarden.com

import pylab
import numpy

f=open("light.txt")
raw=f.readlines()[1:]
f.close()

data = numpy.array(raw,dtype=float)
data = data-min(data) #make all points positive
data = data/max(data)*100.0 #normalize
times = numpy.array(range(len(data)))/1000.0
pylab.figure(figsize=(15,5))
pylab.plot(times,data)
pylab.xlabel("Time Elapsed (seconds)")
pylab.ylabel("Amplitude (% max)")
pylab.title("Pulse Oximeter - filtered")
pylab.subplots_adjust(left=.05,right=.98)
pylab.show()

Future directions involve several projects I hope to work on soon. First, it would be cool to miniaturize everything with surface mount technology (SMT) to bring these things down to the size of a postage stamp. Second, improved finger, toe, or ear clips (or even taped-on sensors) over long duration would provide a pretty interesting way to analyze heart rate variability or modulation in response to stress, sleep apnea, etc. Instead of feeding the signal into a computer, one could send it to a micro-controller for processing. I’ve made some darn-good progress making multi-channel cross-platform USB option for getting physiology data into a computer, but have some work still to do. Alternatively, this data could be graphed on a graphical LCD for an all-in-one little device that doesn’t require a computer. Yep, lots of possible projects can use this as a starting point.

Notes about safety: If you’re worried about electrical shock, or unsure of your ability to make a safe device, don’t attempt to build an ECG machine. For an ECG to work, you have to make good electrical contact with your skin near your heart, and some people feel this is potentially dangerous. Actually, some people like to argue about how dangerous it actually is, as seen on Hack-A-Day comments and my previous post comments. Some people have suggested the danger is negligible and pointed-out that it’s similar to inserting ear-bud headphones into your ears. Others have suggested that it’s dangerous and pointed-out that milliamps can kill a person. Others contest that pulses of current are far more dangerous than a continuous applied current. Realists speculate that virtually no current would be delivered by this circuit if it is wired properly. Rational, cautionary people worried about it reduce risk of accidental current by applying bidirectional diodes at the level of the chest leads, which short any current (above 0.7V) similar to that shown here. Electrically-savvy folks would design an optically decoupled solution. Intelligent folks who abstain from arguing on the internet would probably consult the datasheets regarding ECG input protection. In all cases, don’t attach electrical devices to your body unless you are confident in their safety. As a catch-all, I present the ECG circuit for educational purposes only, and state that it may not be safe and should not be replicated  There, will that cover me in court in case someone tapes wires to their chest and plugs them in the wall socket?

LET ME KNOW WHAT YOU THINK! If you make this, I’m especially interested to see how it came out. Take pictures of your projects and send them my way! If you make improvements, or take this project further, I’d be happy to link to it on this page. I hope this page describes the project well enough that anyone can recreate it, regardless of electronics experience. Finally, I hope that people are inspired by the cool things that can be done with surprisingly simple electronics. Get out there, be creative, and go build something cool!