Archive for October, 2010

The Goertzel Algorithm

By Kevin Banks
(08/28/02, 09:02:26 AM EDT)

The Goertzel algorithm can perform tone detection using much less CPU horsepower than the Fast Fourier Transform, but many engineers have never heard of it. This article attempts to change that.

Most engineers are familiar with the Fast Fourier Transform (FFT) and would have little trouble using a “canned” FFT routine to detect one or more tones in an audio signal. What many don’t know, however, is that if you only need to detect a few frequencies, a much faster method is available. It’s called the Goertzel algorithm.

Tone detection

Many applications require tone detection, such as:

     

  • DTMF (touch tone) decoding
  • Call progress (dial tone, busy, and so on) decoding
  • Frequency response measurements (sending a tone while simultaneously reading back the result)-if you do this for a range of frequencies, the resulting frequency response curve can be informative. For example, the frequency response curve of a telephone line tells you if any load coils (inductors) are present on that line.
  •  

Although dedicated ICs exist for the applications above, implementing these functions in software costs less. Unfortunately, many embedded systems don’t have the horsepower to perform continuous real-time FFTs. That’s where the Goertzel algorithm comes in.

In this article, I describe what I call a basic Goertzel and an optimized Goertzel.

The basic Goertzel gives you real and imaginary frequency components as a regular Discrete Fourier Transform (DFT) or FFT would. If you need them, magnitude and phase can then be computed from the real/imaginary pair.

The optimized Goertzel is even faster (and simpler) than the basic Goertzel, but doesn’t give you the real and imaginary frequency components. Instead, it gives you the relative magnitude squared. You can take the square root of this result to get the relative magnitude (if needed), but there’s no way to obtain the phase.

In this short article, I won’t try to explain the theoretical background of the algorithm. I do give some links at the end where you can find more detailed explanations. I can tell you that the algorithm works well, having used it in all of the tone detection applications previously listed (and others).

A basic Goertzel

First a quick overview of the algorithm: some intermediate processing is done in every sample. The actual tone detection occurs every Nth sample. (I’ll talk more about N in a minute.)

As with the FFT, you work with blocks of samples. However, that doesn’t mean you have to process the data in blocks. The numerical processing is short enough to be done in the very interrupt service routine (ISR) that is gathering the samples (if you’re getting an interrupt per sample). Or, if you’re getting buffers of samples, you can go ahead and process them a batch at a time.

Before you can do the actual Goertzel, you must do some preliminary calculations:

     

  1. Decide on the sampling rate.
  2. Choose the block size, N.
  3. Precompute one cosine and one sine term.
  4. Precompute one coefficient.
  5.  

These can all be precomputed once and then hardcoded in your program, saving RAM and ROM space; or you can compute them on-the-fly.

Sampling rate

Your sampling rate may already be determined by the application. For example, in telecom applications, it’s common to use a sampling rate of 8kHz (8,000 samples per second). Alternatively, your analog-to-digital converter (or CODEC) may be running from an external clock or crystal over which you have no control.

If you can choose the sampling rate, the usual Nyquist rules apply: the sampling rate will have to be at least twice your highest frequency of interest. I say “at least” because if you are detecting multiple frequencies, it’s possible that an even higher sampling frequency will give better results. What you really want is for every frequency of interest to be an integer factor of the sampling rate.

Block size

Goertzel block size N is like the number of points in an equivalent FFT. It controls the frequency resolution (also called bin width). For example, if your sampling rate is 8kHz and N is 100 samples, then your bin width is 80Hz.

This would steer you towards making N as high as possible, to get the highest frequency resolution. The catch is that the higher N gets, the longer it takes to detect each tone, simply because you have to wait longer for all the samples to come in. For example, at 8kHz sampling, it will take 100ms for 800 samples to be accumulated. If you’re trying to detect tones of short duration, you will have to use compatible values of N.

The third factor influencing your choice of N is the relationship between the sampling rate and the target frequencies. Ideally you want the frequencies to be centered in their respective bins. In other words, you want the target frequencies to be integer multiples of sample_rate/N.

The good news is that, unlike the FFT, N doesn’t have to be a power of two.

Precomputed constants

Once you’ve selected your sampling rate and block size, it’s a simple five-step process to compute the constants you’ll need during processing:

w = (2*π/N)*k
cosine = cos w
sine = sin w
coeff = 2 * cosine

For the per-sample processing you’re going to need three variables. Let’s call them Q0, Q1, and Q2.

Q1 is just the value of Q0 last time. Q2 is just the value of Q0 two times ago (or Q1 last time).

Q1 and Q2 must be initialized to zero at the beginning of each block of samples. For every sample, you need to run the following three equations:

Q0 = coeff * Q1 – Q2 + sample
Q2 = Q1
Q1 = Q0

After running the per-sample equations N times, it’s time to see if the tone is present or not.

real = (Q1 – Q2 * cosine)
imag = (Q2 * sine)
magnitude2 = real2 + imag2

A simple threshold test of the magnitude will tell you if the tone was present or not. Reset Q2 and Q1 to zero and start the next block.

An optimized Goertzel

The optimized Goertzel requires less computation than the basic one, at the expense of phase information.

The per-sample processing is the same, but the end of block processing is different. Instead of computing real and imaginary components, and then converting those into relative magnitude squared, you directly compute the following:

magnitude2 = Q12 + Q22-Q1*Q2*coeff

This is the form of Goertzel I’ve used most often, and it was the first one I learned about.

Pulling it all together

Listing 1 shows a short demo program I wrote to enable you to test-drive the algorithm. The code was written and tested using the freely available DJGPP C/C++ compiler. You can modify the #defines near the top of the file to try out different values of N, sampling_rate, and target_frequency.

The program does two demonstrations. In the first one, both forms of the Goertzel algorithm are used to compute relative magnitude squared and relative magnitude for three different synthesized signals: one below the target_frequency, one equal to the target_frequency, and one above the target_frequency.

You’ll notice that the results are nearly identical, regardless of which form of the Goertzel algorithm is used. You’ll also notice that the detector values peak near the target frequency.

In the second demonstration, a simulated frequency sweep is run, and the results of just the basic Goertzel are shown. Again, you’ll notice a clear peak in the detector output near the target frequency. Figure 1 shows the output of the code in Listing 1.

Figure 1: Demo output

For SAMPLING_RATE = 8000.000000 N = 205 and FREQUENCY = 941.000000, k = 24 and coeff = 1.482867


For test frequency 691.000000:
real = -360.392059 imag = -45.871609
Relative magnitude squared = 131986.640625
Relative magnitude = 363.299652
Relative magnitude squared = 131986.640625
Relative magnitude = 363.299652

For test frequency 941.000000:
real = -3727.528076 imag = -9286.238281
Relative magnitude squared = 100128688.000000
Relative magnitude = 10006.432617
Relative magnitude squared = 100128680.000000
Relative magnitude = 10006.431641

For test frequency 1191.000000:
real = 424.038116 imag = -346.308716
Relative magnitude squared = 299738.062500
Relative magnitude = 547.483398
Relative magnitude squared = 299738.062500
Relative magnitude = 547.483398
Freq= 641.0 rel mag^2= 146697.87500 rel mag= 383.01160
Freq= 656.0 rel mag^2= 63684.62109 rel mag= 252.35812
Freq= 671.0 rel mag^2= 96753.92188 rel mag= 311.05292
Freq= 686.0 rel mag^2= 166669.90625 rel mag= 408.25226
Freq= 701.0 rel mag^2= 5414.02002 rel mag= 73.58002
Freq= 716.0 rel mag^2= 258318.37500 rel mag= 508.25031
Freq= 731.0 rel mag^2= 178329.68750 rel mag= 422.29099
Freq= 746.0 rel mag^2= 71271.88281 rel mag= 266.96796
Freq= 761.0 rel mag^2= 437814.90625 rel mag= 661.67584
Freq= 776.0 rel mag^2= 81901.81250 rel mag= 286.18494
Freq= 791.0 rel mag^2= 468060.31250 rel mag= 684.14935
Freq= 806.0 rel mag^2= 623345.56250 rel mag= 789.52234
Freq= 821.0 rel mag^2= 18701.58984 rel mag= 136.75375
Freq= 836.0 rel mag^2= 1434181.62500 rel mag= 1197.57324
Freq= 851.0 rel mag^2= 694211.75000 rel mag= 833.19373
Freq= 866.0 rel mag^2= 1120359.50000 rel mag= 1058.47034
Freq= 881.0 rel mag^2= 4626623.00000 rel mag= 2150.95874
Freq= 896.0 rel mag^2= 160420.43750 rel mag= 400.52521
Freq= 911.0 rel mag^2= 19374364.00000 rel mag= 4401.63184
Freq= 926.0 rel mag^2= 81229848.00000 rel mag= 9012.76074
Freq= 941.0 rel mag^2= 100128688.00000 rel mag= 10006.43262
Freq= 956.0 rel mag^2= 43694608.00000 rel mag= 6610.18994
Freq= 971.0 rel mag^2= 1793435.75000 rel mag= 1339.19226
Freq= 986.0 rel mag^2= 3519388.50000 rel mag= 1876.00330
Freq= 1001.0 rel mag^2= 3318844.50000 rel mag= 1821.76965
Freq= 1016.0 rel mag^2= 27707.98828 rel mag= 166.45717
Freq= 1031.0 rel mag^2= 1749922.62500 rel mag= 1322.84644
Freq= 1046.0 rel mag^2= 478859.28125 rel mag= 691.99658
Freq= 1061.0 rel mag^2= 284255.81250 rel mag= 533.15643
Freq= 1076.0 rel mag^2= 898392.93750 rel mag= 947.83594
Freq= 1091.0 rel mag^2= 11303.36035 rel mag= 106.31726
Freq= 1106.0 rel mag^2= 420975.65625 rel mag= 648.82635
Freq= 1121.0 rel mag^2= 325753.78125 rel mag= 570.74841
Freq= 1136.0 rel mag^2= 36595.78906 rel mag= 191.30026
Freq= 1151.0 rel mag^2= 410926.06250 rel mag= 641.03516
Freq= 1166.0 rel mag^2= 45246.58594 rel mag= 212.71245
Freq= 1181.0 rel mag^2= 119967.59375 rel mag= 346.36337
Freq= 1196.0 rel mag^2= 250361.39062 rel mag= 500.36127
Freq= 1211.0 rel mag^2= 1758.44263 rel mag= 41.93379
Freq= 1226.0 rel mag^2= 190195.57812 rel mag= 436.11417
Freq= 1241.0 rel mag^2= 74192.23438 rel mag= 272.38251

To avoid false detections in production code, you will probably want to qualify the raw detector readings by using a debouncing technique, such as requiring multiple detections in a row before reporting a tone’s presence to the user.

As you can see, the Goertzel algorithm deserves to be added to your signal processing toolbox.

Kevin Banks has been developing embedded systems for 19 years, as a consultant and as an employee at companies including SCI, TxPORT, DiscoveryCom, and Nokia. Currently he is back in consulting mode. He can be reached at kbanks@hiwaay.net.

References

Here are some links to other resources you may find useful:

www.numerix-dsp.com/goertzel.html

ptolemy.eecs.berkeley.edu/papers/96/dtmf_ict/www/node3.html

www.analogdevices.com/library/dspManuals/Using_ADSP-2100_Vol1_books.html

www.analogdevices.com/library/dspManuals/pdf/Volume1/Chapter_14.pdf

www.analogdevices.com/library/dspManuals/Using_ADSP-2100_Vol2_books.html

www.analogdevices.com/library/dspManuals/pdf/2100Chapter_8.pdf

www.delorie.com/djgpp/

Return to the September 2002 Table of Contents

Read the original at:

http://web.archive.org/web/20070926222355/http://www.embedded.com/story/OEG20020819S0057

Advertisements

1 Comment

Step by step to create/modify ramdisk.img

Benno’s blog has an article to change the ramdisk image, I tried it and give more details here.

ramdisk.img is included in the google android sdk, it exists in folder $SDK_ROOT/tools/lib/images/ramdisk.img. The ramdisk.img included in the google android sdk is a gzipped ramdisk.cpio file.

Here is the steps:

  1. Upload the ramdisk.img to your linux machine
  2. Change the ramdisk.img name to ramdisk.cpio.gz, and extract it by: # gzip -d ramdisk.cpio.gz
  3. Create a temporary folder, say tmp, copy ramdisk.cpio to tmp folder
  4. Extract the ramdisk.cpio in the tmp folder with command: # cpio -i -F ramdisk.cpio
  5. Remove the ramdisk.cpio in the tmp folder, and make any changes you want to the extracted ramdisk.cpio in tmp folder
  6. Recreate the ramdisk.cpio with command: # cpio -i -t -F ../ramdisk.cpio | cpio -o -H newc -O ../ramdisk_new.cpio

Some notes:

  1. I change ramdisk.img to ramdisk.cpio.gz, and unzip it. It is because I find the -z parameter is not supported with my cpio. I tried the latest cpio (2.9), it doesn’t work too.
  2. Check cpio version by # cpio –version. I’m using cpio version 2.4.
  3. Find the latest cpio (v 2.9) on site: gnu cpio
  4. Notice that in step 6, the command includes two O’s. First o is lower-case, second is up-case.
  5. Notice in step 6, please remain ramdisk.cpio in up folder of tmp folder. The command need it there.

Read the original at :

http://discuz-android.blogspot.com/2008/01/step-by-step-to-createmodify-ramdiskimg.html

2 Comments

Android on iMx51EVK with 3G modem

You can find the original article here http://sites.google.com/site/fslnovi/imx-demos/imx51-babbage-demos/android-on-bbg25

Android Installation  

1. Install JDK

Install JDK 5 from here :
http://java.sun.com/javase/downloads/index_jdk5.jsp

NOTE : JDK 6 is not compatible. JDK 5 is essential

Installation notes for JDK 5 are here:
http://java.sun.com/j2se/1.5.0/install-linux.html
Download the .bin file into /usr/java and install

2. Setup other packages for Android based on instructions here :
http://source.android.com/source/download.html

3. Edit ~/.profile and add the following commands :

export ARCH=arm
export CROSS_COMPILE=/opt/gcc-4.1.2-glibc-2.5-nptl-3/arm-none-linux-gnueabi/bin/arm-none-linux-gnueabi-
export JAVA_HOME=/usr/java/jdk1.5.0_22
export ANDROID_JAVA_HOME=$JAVA_HOME
export PATH=$JAVA_HOME/bin:$PATH

4. gcc-4.3 and g++4.3 must be installed to compile the SDK.

However, gcc and g++ installed using sudo apt-get install will pull in gcc4.4 and g++4.4
This will not lead to issues in compilation.

Follow all the steps in User guide.

Errors during build : https://sites.google.com/site/fslnovi/imx-demos/imx51-babbage-demos/android-on-bbg25/android-error

Hardware for USB Modem

Huawei E180
Bought from : http://cgi.ebay.com.au/ws/eBayISAPI.dll?ViewItem&item=220459919467
Reference manual : http://www.meteor.ie/rubylith/files/mbb/Meteor%20Mobile%20Broadband%20Trial%20Quick%20start%20guide.pdf

SIM Card
Tmobile 3G with voice and data plan

How to get Log messages ?

Connect 3G module with SIM inserted and then bootup to Android shell:
$ dmesg
$ ls /dev/ttyUSB*
$ logcat -b radio &
$ logcat pppd:* *:S &

* E180 is identified by kernel as two USB serial ports:
# ls /dev/ttyUSB*
/dev/ttyUSB0
/dev/ttyUSB1

But our default configuration for RIL daemon is to listen on ttyUSB3. So change /init.rc as below:
service ril-daemon .... -d /dev/ttyUSB3 -u /dev/ttyUSB0
->
service ril-daemon --- -d /dev/ttyUSB1 -u /dev/ttyUSB0

init.rc is in the root partition of Android root filesystem. After you build Android, it's under ~/myAndroid/out/target/product/imx51_BBG/

How to change init.rc and integrate into rootfs ?
On host machine
$ cd ~/myAndroid
$ rm out/target/imx51_BBG/ramdisk.img
$ CHANGE out/target/imx51_BBG/root/init.rc per your needs
$ make
ramdisk.img will be re-generated
or if you don't want to make whole Android, you can:
$
out/host/linux-x86/bin/mkbootfs out/target/product/imx51_BBG/root | out/host/linux-x86/bin/minigzip > out/target/product/imx51_BBG/ramdisk.img

If you only have ramdisk.img but don't have corresponding "root/" directory, you can get it by:
$ gunzip ramdisk.img
$ mkdir root
$ cd root
$ sudo cpio -i --make-directories < ../ramdisk.img
$ ls
Now you have your "root/". You can change anything (e.g. init.rc) and the re-create the ramdisk.img with method above.

How to detect that you are connected to the network ?

The entire log is attached. But the following sections indicate connection was successful. 

D/GSM     ( 2057): [GsmDataConnectionTracker] setState: INITING
D/AT      ( 1955): AT< OK
D/AT      ( 1955): AT> AT+CGQREQ=1
D/AT      ( 1955): AT< OK
D/AT      ( 1955): AT> AT+CGQMIN=1
D/AT      ( 1955): AT< OK
D/AT      ( 1955): AT> AT+CGEREP=1,0
D/AT      ( 1955): AT< OK
D/AT      ( 1955): AT> AT+CGACT=0,1
D/AT      ( 1955): AT< OK
D/AT      ( 1955): AT> ATD*99***1#
D/AT      ( 1955): AT< CONNECT 236800
D/RIL     ( 1955): PPPd started

The other obvious check would be to check if we get an IP address. 

Sometimes, if you are booting from SD, the excetution mode of "/etc/init.gprs-pppd" is not set which will cause it to fail while bringing up PPP.
To fix this, you need to change the mode in system partition on the SD card after connecting it to the linux host. 

$ chmod 777 /etc/init.gprs-pppd

Radio Interface Layer Sites
----------------------------
a. http://www.netmite.com/android/mydroid/development/pdk/docs/telephony.html
b. http://i-miss-erin.blogspot.com/2009/11/radio-layer-interface-in-android.html 

Android Developers Site
-----------------------
a. http://developer.android.com/guide/index.html

1 Comment