Arduino Forum

Topics => Science and Measurement => Topic started by: akellyirl on Oct 23, 2013, 10:20 pm

Title: Reliable Frequency Measurment Using Autocorrelation
Post by: akellyirl on Oct 23, 2013, 10:20 pm
Accurate Frequency Detection is important for many projects such as Guitar/Piano Tuners, Vibration Analyzers, Heartrate Monitors, MEMs Sensor Analysis and Laboratory Instruments.

There have been many fine examples of projects that try to solve this problem, for example: Arduino Frequency Detection by amandaghassaei on Instructables and Arduino Frequency Counter Library, but they all use Time Domain techniques; analyzing the signal for features such as : Zero-Crossings, Peak Detection, Slope Detection etc..

Waveforms with very strong harmonic content makes the fundamental frequency undiscernable and signals like that will not be identifiable using Time Domain techniques.

I've written an instructable on how to use Autocorrelation to reliably Identify the frequency of signals in the presence of noise and strong harmonic content.

Actually is surprisingly easy.

Here it is.
http://www.instructables.com/id/Reliable-Frequency-Detection-Using-DSP-Techniques/
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: robtillaart on Oct 28, 2013, 09:10 am
Hi Akellyirl,

Nice piece of work, I took the freedom to refactor and optimize your code.

+ created function:  float autoCorrelateFrequency(char * sample, int len, float sampleFreq)
+ made datatype signed char, to remove the -128 from inner loop formula
   (preprocess the data array once)
+ removed the division by 256 from the inner loop
+ added endstate!=3 condition in the outer loop to stop processing asap
+ refactor

The above changes improved the performance from 8613 -> 139 millis to determine the 259.41 Hz in the sample signal.
Given that the UNO can make 8000 10bit  samples / second it should be possible to determine frequencies up to 2Khz
(4000 samples) in less than a second. 

(sample.h file is attached)
Code: [Select]
//    FILE: frequencyDSP.ino
//  AUTHOR: rob tillaart (org version akellyirl)
// VERSION: 0.1.03
// PURPOSE: frequency analysis using DSP
//    DATE: 2013-10-27
//     URL: http://www.instructables.com/id/Reliable-Frequency-Detection-Using-DSP-Techniques

#include "sample.h"

void setup()
{
  Serial.begin(115200);
  Serial.println("FrequencyDSP");
  Serial.println("------------");

  uint32_t start = millis();
  float fr = autoCorrelateFrequency(rawData, sizeof(rawData), 22050);
  uint32_t stop = millis();

  Serial.print("freq:\t");
  Serial.println(fr);
  Serial.print("time:\t");
  Serial.println(stop - start);
}

void loop()
{
}

float autoCorrelateFrequency(char * sample, int len, float sampleFreq)
{
  long sum = 0;
  long sum_old = 0;
  int thresh = 0;
  byte pd_state = 0;
  int period = 0;  // 0 results in inf

  // Autocorrelation
  for (int i=0; i < len && (pd_state != 3); i++)
  {
    sum_old = sum;
    sum = 0;

    for (int k=0; k <len-i; k++)
    {
      sum += (rawData[k]) * (rawData[k+i]);
    }
    sum /= 256;

    // Peak Detect State Machine
    // 0: initial, set threshold
    // 1: detect positive slope
    // 2: detect negative slope
    // 3: peak found
    if (pd_state == 0)
    {
      thresh = sum / 2;
      pd_state = 1;
    }
    else if (pd_state == 1 && (sum > thresh) && (sum - sum_old) > 0)
    {
      pd_state = 2;
    }
    else if (pd_state == 2 && (sum - sum_old) <= 0)
    {
      period = i;
      pd_state = 3;
    }
  }

  return sampleFreq/period;
}


Please give it a try
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: curro92 on Sep 22, 2014, 06:59 am
hi,
in the example as input file is used sample.h raw. How would you read samples from an analog input of Arduino?

Thanks
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: MarkT on Sep 22, 2014, 04:36 pm
Read some samples on a regular basis, shift them right two bits and store
in the char * that autoCorrelateFrequency () expects?
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: jremington on Sep 23, 2014, 01:30 am
A problem with the code as posted is that the autocorrelation sum depends on the number of samples correlated, which become smaller as the index (i) increases. That is, the autocorrelation sum is not properly scaled.

I fixed that, and also eliminated the condition that the program stop as soon as it found a large peak.  A plot of the results shows a number of large peaks in the autocorrelation function, so the result is less clean than the OP claims. (The plot becomes noisier toward the right end, because the number of samples included in the sum is small and the error increases).

The inner loop now looks like this (should probably divide by (k-1), but that is a minor detail)
Code: [Select]
   for (k=0; k <len-i; k++)
   {
     sum += (rawData[k]) * (rawData[k+i]);
   }
   sum /= k ;
   fprintf(fp,"%ld\n",sum);


Note: normally the average value of the data points is subtracted from each term, but here, the average is not large (-4).
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: robtillaart on Sep 23, 2014, 09:12 pm
@jremington
can you post your complete version?
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: jremington on Sep 24, 2014, 05:50 am
See below. This simplified version is for Code::Blocks on a PC, because I wanted to write out a file of the results for plotting and inspection.

Note: another problem with the original code is that the peak detect algorithm detects only the first peak (not necessarily the largest). I did not bother to fix that.

I know this post is old, and several thousand people have looked at it without commenting, so there probably isn't much interest in the technique.

Code: [Select]
#include <stdio.h>
#include <stdlib.h>
#include "sample.h"

float autoCorrelateFrequency(char * sample, int len, float sampleFreq)
{
  long sum = 0;
  long sum_old = 0;
  int thresh = 0,i,k;
  char pd_state = 0;
  int period = 0;  // 0 results in inf
    FILE *fp;
    fp=fopen("test.txt", "w");
  // Autocorrelation
  for (i=0; i < len; i++)
  {
    fprintf(fp,"%d,%d,",i,rawData[i]);
    sum_old = sum;
    sum = 0;

    for (k=0; k <len-i; k++)
    {
      sum += (rawData[k]) * (rawData[k+i]);
    }
    sum /= k ;  //this should be k-1, but if so the last point will result in divide by zero
    fprintf(fp,"%ld\n",sum);
    // Peak Detect State Machine
    // 0: initial, set threshold
    // 1: detect positive slope
    // 2: detect negative slope
    // 3: peak found
    if (pd_state == 0)
    {
      thresh = sum / 2;
      pd_state = 1;
    }
    else if (pd_state == 1 && (sum > thresh) && (sum - sum_old) > 0)
    {
      pd_state = 2;
    }
    else if (pd_state == 2 && (sum - sum_old) <= 0)
    {
      period = i;
      pd_state = 3;
    }
  }
    fclose(fp);
  return sampleFreq/period;
}

int main()
{
  float fr = autoCorrelateFrequency(rawData, sizeof(rawData), 22050);

  printf("freq: %f\n",fr);

    return 0;
}
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: curro92 on Sep 29, 2014, 08:30 am
hi,
something is wrong:
I have compiled the code of robtillaart, I used sample.h of instructables where I added a global variable
const float sample_freq = 22050.0;
(so every sample.h can have different sample rate)

I get the freq = 5512.5Hz result, but the frequency of instructables sample.h  is 260Hz  (with their code I get 259.4Hz)

With the source code of jremington I get the same result.

What am I doing wrong?
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: robtillaart on Sep 29, 2014, 10:43 am
please post the output of the sketch you attached.

which board are you using?
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: curro92 on Sep 29, 2014, 11:16 am
The output is:
Quote
samples: 970, sample_freq: 22050.000000
freq: 5512.500000


I've compiled on XP as
gcc  frequencyDSP_robtillaart.c  -o frequencyDSP_robtillaart.exe

And I have compiled in Linux Mint 17 as

gcc  frequencyDSP_robtillaart.c   -o frequencyDSP_robtillaart

The output is the same in both cases.

Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: curro92 on Sep 29, 2014, 11:29 am

With this other instructables code the output is ok:

len: 970 sample_rate: 22050.000000
freq: 259.411774
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: robtillaart on Sep 29, 2014, 12:05 pm
could be related to size of datatypes?  32 bit - 16 bit
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: KenF on Sep 29, 2014, 12:19 pm
Guys.  This is great work.
I would suggest that the next step would be..
Having found the dominant frequency, subtract it from the raw data.  Then go through the whole process again to find the next most dominant frequency. 
This can be repeated until the rawdata minus all of the frequencies found = 0. That would be a really interesting toy.  Possibly even able to identify someone by their voice.
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: jremington on Sep 29, 2014, 07:51 pm
Just for completeness, I calculated the Fourier transform (as an FFT of length 1024, zero padded) on the original data, and the output is shown on the attached plot.

As you can see, the autocorrelation function does quite a bit better than the FFT to distinguish the frequency of interest, but there are other peaks in the autocorrelation function nearly as large, so it would be difficult to have much confidence in the result.
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: charliesixpack on Oct 08, 2014, 02:12 am
I think you may have a problem with your code.  You have defined rawData to be an array of type char.  When you perform the multiplication in the step which multiplies the two rawData values you are overflowing the type char.  The multiplication is done on type char before it is assigned to variable sum.
Code: [Select]
   for (k=0; k <len-i; k++)
   {
     sum += (rawData[k]) * (rawData[k+i]);
   }

I have been bitten by this before.   Try casting rawData to long before doing the multiplication and see if you get the same answer.

I analyzed your data in Mathematica and got a much different result for the autocorrelation.  Also, you imply that autocorrelation is a frequency domain technique which it is not.  It is done strictly in the time domain.  Taking the FFT of the autocorrelation gives the power spectrum of the signal putting it in the frequency domain.

The plot below is the periodogram of the data.  Mathematica has a function to find the peaks which are found at the following indices of the data (shown as red dots on the plot)...
{13, 24, 35, 47, 59, 70, 82, 93, 105, 116, 128, 139, 151, 162, 174, \
185, 197, 208, 220, 231, 243, 254, 266, 277, 289, 300, 312, 323, 335, \
346, 358, 369, 381, 392, 404, 415, 427, 438, 450, 461, 472, 485, 487, \
500, 511, 522, 534, 545, 557, 568, 580, 591, 603, 614, 626, 637, 649, \
660, 672, 683, 695, 706, 718, 729, 741, 752, 764, 775, 787, 798, 810, \
821, 833, 844, 856, 867, 879, 890, 902, 913, 925, 937, 948, 959}

The fundamental frequency (259.09 Hz)and overtones really stand out when displayed this way.
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: jremington on Oct 09, 2014, 09:28 pm
Not my code, but you have a point. It is certainly easy to overlook those traps in someone else's!

Curiously, using Code::Blocks/gcc on the PC as I was, it does not make any difference to the overall result to cast rawData[] as long, i.e.
Code: [Select]
sum += ( (long) rawData[k]) * ( (long) rawData[k+i]);
It appears that on the PC, Code::Blocks promotes the type char to a larger integer for the purpose of the multiplication. I haven't yet figured out how to view the machine code listing to check this supposition, so I printed out the first few values of rawData[]*rawData[] using the following statement:
Code: [Select]
   fprintf(fp,"%d,%d,%d\n",i,rawData[i],rawData[i]*rawData[i]);and got the following values:
Quote

0,115,13225
1,126,15876
2,62,3844
3,57,3249
4,-80,6400
5,77,5929

Clearly, the multiplication is not overflowing the char data type. So, I don't see why Matlab would give a "much different result" for the autocorrelation function (why not show us what it is?).


Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: jremington on Oct 10, 2014, 07:44 pm
I checked the operations performed by Matlab and found that the built-in function autocorr() calculates the same unweighted function that the OP originally described (except that the mean value of the data is properly subtracted). This is fine for short lags but reduces the weight of peaks associated with large lags. From the autocorr help message "the form adopted here follows that of Box, Jenkins, and Reinsel  ... Time Series Analysis: Forecasting and Control. 3rd ed. Upper Saddle River, NJ: Prentice-Hall, 1994". I've attached a sample plot obtained by using the following Matlab script:
Code: [Select]

m=csvread('test.txt');
ACF=autocorr(m(:,2),969);
plot(ACF);


One can see how the peaks diminish as the lag increases, which I consider to be misleading. The Fourier transform of this particular autocorrelation function will not be the same as the power spectrum plot of the zero-filled FFT that I posted earlier.
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: pantherturm on Nov 14, 2014, 11:53 am
Hello!

I've tested the code of robtillaart but i'm getting some weird characters on the serial monitor. Any idea why is this happening?
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: robtillaart on Nov 14, 2014, 07:35 pm
No, with so little information ...
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: pantherturm on Nov 15, 2014, 08:27 am
I'm using the Uno board. I have uploaded the following code...


Code: [Select]


#include "sample.h"

void setup()
{
  Serial.begin(115200);
  Serial.println("FrequencyDSP");
  Serial.println("------------");

  uint32_t start = millis();
  float fr = autoCorrelateFrequency(rawData, sizeof(rawData), 22050);
  uint32_t stop = millis();

  Serial.print("freq:\t");
  Serial.println(fr);
  Serial.print("time:\t");
  Serial.println(stop - start);
}

void loop()
{
}

float autoCorrelateFrequency(char * sample, int len, float sampleFreq)
{
  long sum = 0;
  long sum_old = 0;
  int thresh = 0;
  byte pd_state = 0;
  int period = 0;  // 0 results in inf

  // Autocorrelation
  for (int i=0; i < len && (pd_state != 3); i++)
  {
    sum_old = sum;
    sum = 0;

    for (int k=0; k <len-i; k++)
    {
      sum += (rawData[k]) * (rawData[k+i]);
    }
    sum /= 256;

    // Peak Detect State Machine
    // 0: initial, set threshold
    // 1: detect positive slope
    // 2: detect negative slope
    // 3: peak found
    if (pd_state == 0)
    {
      thresh = sum / 2;
      pd_state = 1;
    }
    else if (pd_state == 1 && (sum > thresh) && (sum - sum_old) > 0)
    {
      pd_state = 2;
    }
    else if (pd_state == 2 && (sum - sum_old) <= 0)
    {
      period = i;
      pd_state = 3;
    }
  }

  return sampleFreq/period;
}


but instead of a frequency what i'm getting on the serial monitor is this

Quote

üÕ&
K±μ
Yþ)¤

Sorry i'm not an experienced arduino user not sure what else i should provide.
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: wildbill on Nov 15, 2014, 01:10 pm
Check that you have set the serial monitor's baud rate to match that in your sketch (115200).
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: pantherturm on Nov 15, 2014, 01:54 pm
Yep that fixed it! Thank you!
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: michinyon on Dec 01, 2014, 06:58 am
Quote
When you perform the multiplication in the step which multiplies the two rawData values you are overflowing the type char.  The multiplication is done on type char before it is assigned to variable sum.
I am unconvinced by this.   The multiplication of 2 single-byte numbers on the Arduino results in a 16-bit result.
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: michinyon on Dec 01, 2014, 07:38 am
The other problem with the autoCorrelateFrequency( )  function in Reply # 1,  is that this function gets passed a pointer to the data array,  as a function argument,   and then doesn't use it -  it goes straight to the global data.

This problem becomes more obvious if you put each function in a separate file.
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: michinyon on Dec 01, 2014, 08:31 am
Quote
I get the freq = 5512.5Hz result,
So do I.

I tried to run this program using Netbeans using cygwin/gcc on a windows computer,  and had all sorts of grief.

I had to define the data as unsigned char to get it to compile ,    and then I got 5512 Hz

But if I defined the data as char,  I got 259 Hz.

And if I defined the data as signed char,   I also got 259 Hz.

I believe there is an "issue" with C, on different platforms,  concerning whether char numbers are signed or unsigned if not explicitly stated.   That is probably where your 5512 Hz comes from.





Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: nenea_dani on Apr 12, 2018, 12:29 am
Hi Akellyirl,
In the version that collects data in A0 why just at low frequencies the result is correct? For 1000hz real they get 891hz and for 2000 they get 1783hz. Is the problem in autocorrelation or in machine state? I used http://www.akellyirl.com/arduino-guitar-tuner/
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: MrMark on Apr 12, 2018, 03:07 am
Hi Akellyirl,
In the version that collects data in A0 why just at low frequencies the result is correct? For 1000hz real they get 891hz and for 2000 they get 1783hz. Is the problem in autocorrelation or in machine state? I used http://www.akellyirl.com/arduino-guitar-tuner/
The original post is several years old, so I'll be surprised if the original poster responds.

The issue is that the measured frequency resolution decreases as the frequency increases.  This is because the autocorrelation peaks are measured at discrete steps.  Resolution is approximately f2/fs where "f" is the frequency being measured and "fs is the sample rate.  So at 2000 Hz with a sample rate of 8919 Hz, the resolution is only about 448 Hz.

A brute force way of getting better resolution is to use a higher sampling rate, but RAM to store the sample then becomes an issue.  A more sophisticated approach would be to interpolate between correlation points by using some sort of curve fit around the discrete peak.
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: nenea_dani on Apr 12, 2018, 10:43 pm
I migrated from a FHT frequency meter for 150-4000Hz that has an error of less than 10 percent in the idea that I can measure it more precisely.  The author said he can get even decimal but not specified as just for two hundred or three hundred Hz. Looks like I'm out of math. The formula " f square / fs"  " is very weird for me and I still want to know where it comes from unless I ask too much.  Do I understand that autocorrelation is convenient to implement in general only for low frequencies ?  Anyway thanks for the information. For now I understand that I am in a direction that surpasses me and I return to FHT that even if I do not understand how it works, at least it gives very good results.
Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: tim77777 on Jul 09, 2018, 02:34 am
I found two major problems with the changes made to akellyirl's original code.

1. Although the changes work OK on the supplied data file, when applied to a changing input such as from an A/D it becomes obvious that often a frequency exactly half of the actual, is being outputted.

2. Changing the raw data from 'byte' to 'signed char', increases the overall processing time for each outputted frequency from 100mS to about 160mS.

The modifications involved a lot of formatting/structural changes from the original code (added function etc), so maybe I messed it up somewhere.  This is how I combined the changes with Akellyirl's original code.

Code: [Select]

#define LENGTH 256

signed char rawData[LENGTH];


// Sample Frequency in kHz
const float sample_freq = 8919;

int len = sizeof(rawData);
int i,k;
long sum, sum_old;
int thresh = 0;
float freq_per = 0;
byte pd_state = 0;
int count;


void setup()
{
  analogReference(DEFAULT);   // Connect to 3.3V
  analogRead(A0);
  Serial.begin(115200);
  count = 0;
}


void loop(){

    if (count < LENGTH)
    {   
      rawData[count] = analogRead(A0)>>2; 
      count++;
    }
 
    else
    {
       sum = 0;
       pd_state = 0;
       int period = 0;
       
       for(i=0; i < len; i++)
       {       
            sum_old = sum;
            sum = 0;
           
            for (k=0; k <len-i; k++)
            {
               sum += (rawData[k]) * (rawData[k+i]);
            }
            sum /= k ;
               
            if (pd_state == 0)
            {
              thresh = sum / 2;
              pd_state = 1;
            }
            else if (pd_state == 1 && (sum > thresh) && (sum - sum_old) > 0)
            {
              pd_state = 2;
            }
            else if (pd_state == 2 && (sum - sum_old) <= 0)
            {
              period = i;
              pd_state = 3;
            }
 
       }
 
        if (thresh >100 && period > 0)
        {
          freq_per = sample_freq/period;
          Serial.println(freq_per);
        }
     
        count = 0;
   
    }
 
}

Title: Re: Reliable Frequency Measurment Using Autocorrelation
Post by: MarkT on Jul 10, 2018, 02:09 pm
The accuracy of frequency samples without interpolation should be fs/N where N is the FFT size.

However in the time domain the time error is 1/fs, which is f/fs wave periods, so the frequency error is about
f * f/fs, or f^2/fs.

Once you have an estimate of the frequency you can use other methods to refine the estimate, perhaps as
simple as crude band-pass filtering and zero-crossing counting.