Arduino pitch detection in realtime (using autocorrelation and peak detection)

Projects General Guidance
43

I've had a bit more time to play with this. There are two significant changes in this version. First the data capture uses the ADC in the continuous mode with an interrupt service routine to move the data into the buffer. Sample rate is fixed at 9615.4 samples per second. The second change is a correction to the interpolation routine which previously had the peak correlation value as one of the endpoints rather than the center of the three points.

I've changed the peak detection algorithm implementation to use "case" statements rather than a cascaded "if" statements. It should be functionally the same, but I think peak detection could be made more robust and this was a first step toward a revision. Eventually I'd like peak detection to evaluate multiple correlation peaks rather than stopping at the first it finds over the threshold.

Finally I opened up the frequency window test to 100-2000 Hz, largely because the phone app doesn't put out much audio power below about 200 Hz, so it gives me more range to test over. In limited testing this implementation seems to mostly put out accurate values if it puts out anything at all.

// Arduino Pitch Detection on A0 with autocorrelation and peak detection
// Original author(s): akellyirl, revised by robtillaart, MrMark, barthulsen
// http://forum.arduino.cc/index.php?topic=540969.15
// Continuous ADC ala http://www.instructables.com/id/Arduino-Audio-Input/ "Step 7"

#define sampleFrequency 9615.4
#define bufferSize 1024

volatile byte  rawData[bufferSize] ;  // Buffer for ADC capture
volatile int sampleCnt ;                    // Pointer to ADC capture buffer
long currentSum, previousSum, twoPreviousSum ;
int threshold = 0;
float frequency = 0;
byte pdState = 0;

void setup() {

  Serial.begin(115200) ;
  pinMode(LED_BUILTIN, OUTPUT);
  cli();//disable interrupts

  //set up continuous sampling of analog pin 0
  //clear ADCSRA and ADCSRB registers
  ADCSRA = 0 ;
  ADCSRB = 0 ;

  ADMUX |= (1 << REFS0) ; //set reference voltage
  ADMUX |= (1 << ADLAR) ; //left align the ADC value- so we can read highest 8 bits from ADCH register only

  ADCSRA |= (1 << ADPS2) | (1 << ADPS1) | (1 << ADPS0); // ADC clock 128 prescaler- 16mHz/128=125kHz->9615 sps
  ADCSRA |= (1 << ADATE); //enable auto trigger
  ADCSRA |= (1 << ADIE) ; //enable interrupts when measurement complete
  ADCSRA |= (1 << ADEN) ; //enable ADC
  ADCSRA |= (1 << ADSC) ; //start ADC measurements
}

ISR(ADC_vect) {     // When ADC sample ready, put in buffer if not full
  if (sampleCnt < bufferSize)
  {
    rawData[sampleCnt] = ADCH ;
    sampleCnt++ ;
  }
}

void readData() {
  sampleCnt = 0 ;
  sei() ;         // Enable interrupts, samples placed in buffer by ISR
  while (sampleCnt < bufferSize) ;  // Spin until buffer is full
  cli() ;         // Disable interrupts
}

void findFrequency() {
  // Calculate mean to remove DC offset
  long meanSum = 0 ;
  for (int k = 0; k < bufferSize; k++) {
    meanSum += rawData[k] ;
  }
  char mean = meanSum / bufferSize ;
  // Remove mean
  for (int k = 0; k < bufferSize; k++) {
    rawData[k] -= mean ;
  }

  // Autocorrelation
  currentSum = 0 ;
  pdState = 0 ;
  for (int i = 0; i < bufferSize && (pdState != 3); i++) {
    // Autocorrelation
    float period = 0 ;
    twoPreviousSum = previousSum ;
    previousSum = currentSum ;
    currentSum = 0 ;
    for (int k = 0; k < bufferSize - i; k++) {
      currentSum += char(rawData[k]) * char(rawData[k + i]) ;
    }
    // Peak detection
    switch (pdState) {
      case 0:   // Set threshold based on zero lag autocorrelation
        threshold = currentSum / 2 ;
        pdState = 1 ;
        break ;
      case 1:   // Look for over threshold and increasing
        if ((currentSum > threshold) && (currentSum - previousSum) > 0) pdState = 2 ;
        break ;
      case 2:   // Look for decreasing (past peak over threshold)
        if ((currentSum - previousSum) <= 0) {
          // quadratic interpolation
          float interpolationValue = 0.5 * (currentSum - twoPreviousSum) / (2 * previousSum - twoPreviousSum - currentSum) ;
          period = i - 1 + interpolationValue ;
          pdState = 3 ;
        }
        break ;
      default:
        pdState = 3 ;
        break ;
    }

    // Frequency identified in Hz
    if (threshold > 100) {
      frequency = sampleFrequency / period;
      if (frequency < 2000) {
        Serial.println(frequency);
      }
    }
  }
}

void loop() {
  digitalWrite(LED_BUILTIN, HIGH);
  readData();
  digitalWrite(LED_BUILTIN, LOW);
  findFrequency();
}