//
// Copyright:     Copyright 2003 Craig Stuart Sapp
// Programmer:    Craig Stuart Sapp <craig@ccrma.stanford.edu>
// Creation Date: Wed Feb  5 05:29:45 PST 2003
// Last Modified: Wed Feb  5 07:48:29 PST 2003
// Filename:      key.c
// Web Address:   http://peabody.sapp.org/class/dmp2/lab/key/key.c
// Syntax:        C; Max4/MSP2 External Object; CodeWarrior 6.0
// OS:	          Mac OS 9; PPC
//
// Description:   Krumhansl & Schmuckler key-finding algorithm adapted
//                for use in Max for realtime analysis.  Only note
//                attacks are used to calculate the key -- duration
//                of notes are ignored.
//

#include "ext.h"
#include <cmath>              /* included for sqrt() function */

#define BUFFERSZ 1000         /* maximum note count that can be remembered */

typedef struct {              // data structure for note memory
   long ts;                   // timestamp for note arrival
   long note;                 // MIDI note pitch 
} NoteInfo;

typedef struct {
   t_object max_data;         // Max/MSP data, MUST come first in struct
   NoteInfo memory[BUFFERSZ]; // storage buffer for incoming notes
   long     writeindex;       // write index of memory buffer
   long     windowsize;       // time in milliseconds for analysis window
   long     tonic;            // tonic note from analysis
   long     certainty;        // certainty of the answer
   long     mode;             // mode of key (0=major, 1=minor)
   void*    outputTonic;      // outlet for tonic pitch-class (c=0...b=11)
   void*    outputCertainty;  // outlet for certainty (1=un...127=certain)
   void*    outputMode;       // outlet for mode (0=major, 1=minor)
} MyObject;

void* object_data = NULL;

// function declarations:
void*  create_object       (long windowsizeinit);
void   InputBang           (MyObject* mo);
void   InputNote           (MyObject* mo, long value);
void   InputWindowSize     (MyObject* mo, long value);
void   MessageClear        (MyObject* mo);
long   midilimit           (long value);
void   doAnalysis          (MyObject* mo);
void   storeNoteInMemory   (MyObject* mo, long note, long timestamp);
void   makeHistogram       (double histogram[12], NoteInfo memory[BUFFERSZ],
                            long writeindex, long curtime, long windowsize);


/////////////////////////////////////////////////////////////////////////
//
// Initialization functions
//


//////////////////////////////
//
// main -- called once when the object is created in a patcher window.
//

void main(void) {
   setup((t_messlist**)&object_data, (method)create_object,
         NULL, sizeof(MyObject), NULL, A_DEFLONG, A_NOTHING);
   addbang((method)InputBang);            // inlet 1
   addint ((method)InputNote);            // inlet 1
   addinx ((method)InputWindowSize, 1);   // inlet 2
}



//////////////////////////////
//
// create_object -- create the data storage for the object and
//     and setup input 1.  Default value in object is the duration
//     of the analysis window in seconds.
//

void* create_object(long windowsizeinit) {
   MyObject* mo = (MyObject*)newobject(object_data);
   
   if (windowsizeinit == 0) {
      InputWindowSize(mo, 7000);
   } else {
      InputWindowSize(mo, windowsizeinit);
   }
   MessageClear(mo);       // initialize the memory buffer and other data

   mo->outputMode        = intout(mo);   // outlet 3
   mo->outputCertainty   = intout(mo);   // outlet 2
   mo->outputTonic       = intout(mo);   // outlet 1
 
   intin(mo, 1);                         // inlet 2
   return mo;
}


/////////////////////////////////////////////////////////////////////////
//
// Behavior functions
//


//////////////////////////////
//
// InputBang -- behavior of the object when a "bang" message is received.
//

void InputBang(MyObject* mo) {
   doAnalysis(mo);

   outlet_int(mo->outputMode,      mo->mode);         // outlet 3
   outlet_int(mo->outputCertainty, mo->certainty);    // outlet 2
   outlet_int(mo->outputTonic,     mo->tonic);        // outlet 1
}



//////////////////////////////
//
// InputNote -- behavior of the object when a new number comes
//    in on inlet 0.  The input is a MIDI note which will be stored
//    in memory.
//

void InputNote(MyObject* mo, long value) {
   storeNoteInMemory(mo, value, gettime());
}



//////////////////////////////
//
// InputWindowSize -- behavior of the object when a new number comes
//    in on inlet 1.  The input is the window size in milliseconds
//    which will be used to analyze the key.
//

void InputWindowSize(MyObject* mo, long value) {
   if (value < 10) {               // window must be at least 10 ms long
      value = 10;                  // 
   } else if (value > 600000) {    // but not more than 10 minutes long
      value = 600000;
   }

   mo->windowsize = value;
}



//////////////////////////////
//
// MessageClear -- remove all notes from memory buffer.
//

void MessageClear(MyObject* mo) {
   int i;
   for (i=0; i<BUFFERSZ; i++) {
      mo->memory[i].ts   = -1;
      mo->memory[i].note = -1;
   }
   mo->writeindex = 0;
   mo->tonic      = 0;
   mo->certainty  = 0;
   mo->mode       = 0;
}


/////////////////////////////////////////////////////////////////////////
//
// Non-interface functions:
//


//////////////////////////////
//
// midilimit -- limit a number to the range from 0 to 127.
//    if the input is less than 0, return 0.  
//    if the input is greater than 127, return 127.  
//

long midilimit(long value) {
   if (value < 0)     return   0;
   if (value > 127)   return 127;
   return value;
}


//////////////////////////////
//
// doAnalysis -- calculate the key of the music in memory.
//

void doAnalysis(MyObject* mo) {
   double majorKey[12] = 
      {6.35, 2.23, 3.48, 2.33, 4.38, 4.09, 2.52, 5.19, 2.39, 3.66, 2.29, 2.88};
   double minorKey[12] = 
      {6.33, 2.68, 3.52, 5.38, 2.60, 3.53, 2.54, 4.75, 3.98, 2.69, 3.34, 3.17};
   double mean;
   double major_mean = 3.482500;
   double minor_mean = 3.709167;
   double r_major[12];
   double r_minor[12];
   double maj_numerator;
   double min_numerator;
   double maj_denominator;
   double min_denominator;
   double maj_temp;
   double min_temp;
   double temp;
   int    subscript;
   int    i, j;
   double best_key;
   long   mode = 0;
   double second_best_key;
   int    sec_pitch_class;
   int    pitch_class;
   double confidence;
   double total = 0;
   double histogram[12];
   
   
   makeHistogram(histogram, mo->memory, mo->writeindex, gettime(), mo->windowsize);
   for (i=0; i<12; i++) {
      total += histogram[i];
   }
   mean = total/12.0;

   for (i=0; i<12; i++) {
      maj_numerator   = min_numerator   = 0;
      maj_denominator = min_denominator = 0;
      maj_temp = min_temp = temp = 0;

      // Examine all pitches for each key.
      for (j=0; j<12; j++) {
         subscript = (i+j) % 12;
         temp += (histogram[subscript] - mean) *
                 (histogram[subscript] - mean);

         // For the major keys.
         maj_numerator += (majorKey[j]-major_mean) * 
                          (histogram[subscript]-mean);
         maj_temp += (majorKey[j] - major_mean) * 
                     (majorKey[j] - major_mean);
         
         // For the minor keys.
         min_numerator += (minorKey[j]-minor_mean) * 
                          (histogram[subscript]-mean);
         min_temp += (minorKey[j] - minor_mean) *
                     (minorKey[j] - minor_mean);

      }

      maj_denominator = sqrt(maj_temp * temp);
      min_denominator = sqrt(min_temp * temp);
      
      if (maj_denominator == 0 || min_denominator == 0) {
         mo->tonic      = 0;
         mo->certainty  = 0;
         mo->mode       = 0;
         return;
      }

      r_major[i] = maj_numerator / maj_denominator;
      r_minor[i] = min_numerator / min_denominator;
   }

   // Now determine which correlation is the greatest.
   // Start off with the assumption that C major is the best key.
   best_key = 0.0; 
   pitch_class = 0;            // tonic note of best key
   
   // Compare all the remaining key correlations.
   for (i=0; i<12; i++) {
      if (r_major[i] > best_key) {
         best_key = r_major[i]; 
         mode = 0;                 // major key
         pitch_class = i;
      }
      if (r_minor[i] > best_key) {
         best_key = r_minor[i]; 
         mode = 1;                 // minor key
         pitch_class = i;
      }
   }

   // A confidence measure can be determined by taking the difference
   // between the correlation for the "best key" and subtracting the
   // correlation for the "second best key".  The maximum confidence
   // score is 100; the minimum is zero.
   // First, having found the "best key", find the "second best key."
   second_best_key = 0;
   sec_pitch_class = 0;
   for (i=0; i<12; i++) {
      if (r_major[i] != best_key && r_major[i] > second_best_key) {
         second_best_key = r_major[i];
         sec_pitch_class = i;
      }
      if (r_minor[i] != best_key && r_minor[i] > second_best_key) {
         second_best_key = r_minor[i];
         sec_pitch_class = i;
      }
   }

   // The value 3.0 below is an empirical scaling factor.
   confidence = (best_key - second_best_key) * 100 * 3.0;
   if (confidence > 100.0) confidence = 100.0;


   // set all of the analysis output values
   mo->tonic     = pitch_class;
   mo->certainty = (int)(127.0 * confidence/100.0);
   mo->mode      = mode;
}



///////////////////////////////
//
// storeNoteInMemory --
//

void storeNoteInMemory(MyObject* mo, long note, long timestamp) {
   mo->memory[mo->writeindex].note = midilimit(note);
   mo->memory[mo->writeindex].ts   = timestamp;
   mo->writeindex += 1;
   if (mo->writeindex >= BUFFERSZ) {
      mo->writeindex = 0;
   }
}



//////////////////////////////
//
// makeHistogram -- count the pitch classes in memory
//

void makeHistogram(double histogram[12], NoteInfo memory[BUFFERSZ],
      long writeindex, long curtime, long windowsize) {
   int i;
   int index;  // read index in memory (start from write index going backwards)
   
   // clear the contents of the histogram
   for (i=0; i<12; i++) {
      histogram[i] = 0.0;
   }
   
   for (i=0; i<BUFFERSZ; i++) {
      index = writeindex - 1 - i;
      if (index < 0) {
         index = (index + 2 * BUFFERSZ) % BUFFERSZ;
      }
      if (curtime - memory[index].ts > windowsize) {
         break;
      }
      if (memory[index].note < 0) {
         break;
      }
      histogram[memory[index].note % 12] += 1.0;
   }
}