//
// Programmer:    Craig Stuart Sapp <craig@ccrma.stanford.edu>
// Creation Date: Sat Jul 27 14:29:44 PDT 2002
// Last Modified: Sat Jul 27 14:29:47 PDT 2002
// Filename:      ...sig/examples/all/iwsearch.cpp
// Web Address:   http://sig.sapp.org/examples/museinfo/humdrum/iwsearch.cpp
// Syntax:        C++; museinfo
//
// Description:   Search for optimal root-interval weights for chord root
//                identification by Monte Carlo methods.
// 

#include "humdrum.h"

#include <string.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>

typedef Array<double> ArrayDouble;
typedef Array<int>    ArrayInt;

///////////////////////////////////////////////////////////////////////////

// function declarations
void   checkOptions             (Options& opts, int argc, char* argv[]);
void   example                  (void);
void   usage                    (const char* command);
void   getStartingWeights       (Array<double>& weights, 
                                 HumdrumFile& weightfile);
void   getChordInformation      (Array<ArrayInt>& pitches, Array<int>& root, 
                                 Array<double>& count, HumdrumFile& datafile);
void   runSearch                (Array<double>& initialweights, 
                                 Array<ArrayInt>& pitches, Array<int>& root, 
                                 Array<double>& count, int trialsperstep,
                                 double initrange, double decay, int stoplimit);
double getErrors                (Array<double>& weights, 
                                 Array<ArrayInt>& pitches, Array<int>& root, 
                                 Array<double>& count);
void   makeRandomWeights        (Array<ArrayDouble>& stepweights, double range);
void   printWeights             (Array<double> weights);
double getDistance              (Array<double>& a, Array<double>& b);


// global variables
Options      options;            // database for command-line arguments
int          debugQ     = 0;     // used with --debug option
double       rangeStart = 1.0;   // starting point of random variation
double       decay      = 0.95;  // decay of range between each step
int          trialCount = 100;   // number of trials for each step
int          stopLimit  = 10;    // number of trials which can have same error
int          seed       = 0;     // random number seed

///////////////////////////////////////////////////////////////////////////

int main(int argc, char* argv[]) {
   HumdrumFile datafile;
   HumdrumFile weightfile;

   checkOptions(options, argc, argv); // process the command-line options

   datafile.read(options.getArg(1));
   weightfile.read(options.getArg(2));

   Array<double> count;          // the frequency of the chords occurence
   Array<int> root;              // the root of the chord
   Array<ArrayInt> pitches;      // the pitch class set of the chord
   Array<double> initialweights; // starting weights of the search

   getStartingWeights(initialweights, weightfile);
   getChordInformation(pitches, root, count, datafile);
     
   runSearch(initialweights, pitches, root, count, trialCount, rangeStart,
         decay, stopLimit);

   return 0;
}


///////////////////////////////////////////////////////////////////////////


//////////////////////////////
//
// getErrors -- try all chords and count how many root errors occured.
//

double getErrors(Array<double>& weights, Array<ArrayInt>& pitches,
      Array<int>& root, Array<double>& count) {

   Array<double> rootscores;
   rootscores.setSize(40);
   rootscores.allowGrowth(0);
   int i, j, m;
   double errors = 0;
   int min;
   for (m=0; m<root.getSize(); m++) {
      rootscores.setAll(0.0);
      for (i=0; i<rootscores.getSize(); i++) {
         for (j=0; j<pitches[m].getSize(); j++) {
            rootscores[i] += weights[(pitches[m][j]-i+400)%40];
         }
      }
      min = 0;
      for (i=0; i<rootscores.getSize(); i++) {
         if (rootscores[min] > rootscores[i]) {  
            min = i;
         }
      }
      if (root[m] != min+2) {
         if (debugQ) {
            cout << "Error: root=" << root[m] 
                 << "\tbut measured: " << min << endl;
         }
         errors += count[m];
      }
   }

   return errors;
}



//////////////////////////////
//
// runSearch --
//

void runSearch(Array<double>& initialweights, Array<ArrayInt>& pitches, 
      Array<int>& root, Array<double>& count, int trialsperstep,
      double initrange, double decay, int stoplimit) {

   double range = initrange;
   Array<ArrayDouble> stepweights;
   Array<double>      stepscores;
   stepweights.setSize(trialsperstep);
   stepscores.setSize(trialsperstep);
   int i;
   for (i=0; i<stepweights.getSize(); i++) {
      stepweights[i].setSize(40);
      stepweights[i].allowGrowth(0);
   }
   Array<double> bestweight(40);
   bestweight = initialweights;
   double bestscore = getErrors(bestweight, pitches, root, count);

   cout << "!! The initial best score is: " << bestscore << endl;
   cout << "!! Starting weights: \n";
   printWeights(bestweight);
   cout << endl;
   

   int currentlimit = 0;
   int currentbest = 0;
   int step = 0;
   while (currentlimit < stoplimit) {
      step++;
      stepweights[0] = bestweight;
      stepscores[0]  = bestscore;
      
      makeRandomWeights(stepweights, range);

      for (i=1; i<trialsperstep; i++) {
         stepscores[i] = getErrors(stepweights[i], pitches, root, count);
         if (stepscores[i] <= stepscores[currentbest]) { 
            // include equals to wander on same level
            currentbest = i;
         }
      }
      if (stepscores[currentbest] <= bestscore) {
         bestscore    = stepscores[currentbest];
         bestweight   = stepweights[currentbest];
         currentlimit = 0;
      } else {
         currentlimit++;
      }

      cout << "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" << "\n";
      cout << "!! Step:\t"          << step << "\n";
      cout << "!! Trials:\t"        << trialCount << "\n";
      cout << "!! Stop Limit:\t"    << stoplimit << "\n";
      cout << "!! Start Range:\t"   << initrange << "\n";
      cout << "!! Range:\t"         << range << "\n";
      cout << "!! Decay:\t"         << decay << "\n";
      cout << "!! Best Score:\t"    << bestscore << "\n";
      cout << "!! Step Distance:\t" << getDistance(bestweight, stepweights[0]) 
           << "\n";
      cout << "!! Random Seed:\t"   << seed << "\n";
      printWeights(bestweight);
      cout << endl;

      range *= decay;
   }
}



//////////////////////////////
//
// getDistance -- find the Euclidian distance between two vectors.
//

double getDistance(Array<double>& a, Array<double>& b)  {
   Array<double> c(40);
   int i;
   double sum = 0.0;
   for (i=0; i<c.getSize(); i++) {
      if (i==5 || i==11 || i==22 || i==28 || i == 34) {
         continue;
      }
      c[i] = a[i] - b[i];
      sum += c[i] * c[i];
   }
   
   return sqrt(sum);
}



//////////////////////////////
//
// makeRandomWeights --
//

void makeRandomWeights(Array<ArrayDouble>& stepweights, double range) {
   int i, j;
   double delta = 0.0;
   for (i=1; i<stepweights.getSize(); i++) {
      for (j=0; j<stepweights[i].getSize(); j++) {
         delta = drand48() * 2 * range - range;
         stepweights[i][j] = stepweights[0][j] + delta;
      }
   }
}



//////////////////////////////
//
// getStartingWeights --
//

void getStartingWeights(Array<double>& weights, HumdrumFile& weightfile) {
   weights.setSize(40);
   weights.setAll(100000);
   weights.allowGrowth(0);

   int i, j;
   int root;
   double weight;
   for (i=0; i<weightfile.getNumLines(); i++) {
      root = -1;
      weight = 10000;
      if (!weightfile[i].isData()) {
         continue;
      }
      for (j=0; j<weightfile[i].getFieldCount(); j++) {
         if (strcmp("**kern", weightfile[i].getExInterp(j)) == 0) {
            root = Convert::kernToBase40(weightfile[i][j]) % 40;
         } else if (strcmp("**weight", weightfile[i].getExInterp(j)) == 0) {
            weight = strtod(weightfile[i][j], NULL);
         }
      }
      if (root >= 0) {
         weights[root] = weight;
      }
   }

   if (debugQ) {
      printWeights(weights);
   }
}


//////////////////////////////
//
// printWeights --
//

void printWeights(Array<double> weights) {
   char buffer[128] = {0};
   cout << "**kern\t**weight\n";
   int i;
   for (i=0; i<weights.getSize(); i++) {
      if (i==5 || i==11 || i==22 || i==28 || i==34) {
         continue;   // invalid interval
      }
      cout << Convert::base40ToKern(buffer, i+4*40);
      cout << "\t" << weights[i] << "\n";
   }
   cout << "*-\t*-\n";
}



//////////////////////////////
//
// getChordInformation --
//

void getChordInformation(Array<ArrayInt>& pitches, Array<int>& root, 
     Array<double>& count, HumdrumFile& datafile) {

   count.setSize(datafile.getNumLines());
   count.setSize(0);
   root.setSize(datafile.getNumLines());
   root.setSize(0);
   pitches.setSize(datafile.getNumLines());
   pitches.setSize(0);

   char buffer[1024] = {0};
   int i, j, k;
   int troot = -1;
   int tpitch = -1;
   double tcount = 0.0;
   Array<int> tpitches;
   tpitches.setSize(100);
   tpitches.setSize(0);
   for (i=0; i<datafile.getNumLines(); i++) {
      tpitches.setSize(0);
      troot = -1;
      tcount = 0.0;
      for (j=0; j<datafile[i].getFieldCount(); j++) {
         if (strcmp("**root", datafile[i].getExInterp(j)) == 0) {
            troot = Convert::kernToBase40(datafile[i][j]) % 40;
         } else if (strcmp("**count", datafile[i].getExInterp(j)) == 0) {
            tcount = strtod(datafile[i][j], NULL);
         } else if (strcmp("**kern", datafile[i].getExInterp(j)) == 0) {
            int notecount = datafile[i].getTokenCount(j);
            for (k=0; k<notecount; k++) { 
               datafile[i].getToken(buffer, j, k);
               tpitch = Convert::kernToBase40(buffer);
               tpitches.append(tpitch);
            }
         }
      }

      if (troot < 0 || tcount <= 0.0 || tpitches.getSize() == 0) {
         continue;
      }
      root.append(troot);
      count.append(tcount);
      pitches.append(tpitches);
   }


   if (debugQ) {
      cout << "**count\t**root\t**kern\n";
      for (i=0; i<count.getSize(); i++) {
         cout << count[i] << "\t";
         cout << Convert::base40ToKern(buffer, root[i]+3*40) << "\t";
         for (j=0; j<pitches[i].getSize(); j++) {
            cout << Convert::base40ToKern(buffer, pitches[i][j]);
            if (j < pitches[i].getSize() - 1) {
               cout << " ";
            }
         }
         cout << "\n";
      }
      cout << "*-\t*-\t*-\n";
   }

}



//////////////////////////////
//
// checkOptions -- validate and process command-line options.
//

void checkOptions(Options& opts, int argc, char* argv[]) {
   opts.define("r|range=d:1.0",          "starting point of random variation");
   opts.define("d|decay=d:0.99",         "decay of range between each step");
   opts.define("s|seed=i:0",             "random number seed");
   opts.define("c|count|trials=i:10000", "number of trials for each step");
   opts.define("l|limit=i:100",          "number of trials with same errors");

   opts.define("debug=b",         "trace input parsing");   
   opts.define("author=b",        "author of the program");   
   opts.define("version=b",       "compilation information"); 
   opts.define("example=b",       "example usage"); 
   opts.define("h|help=b",        "short description"); 

   opts.process(argc, argv);
   
   // handle basic options:
   if (opts.getBoolean("author")) {
      cout << "Written by Craig Stuart Sapp, "
           << "craig@ccrma.stanford.edu, July 2002" << endl;
      exit(0);
   } else if (opts.getBoolean("version")) {
      cout << argv[0] << ", version: 28 July 2002" << endl;
      cout << "compiled: " << __DATE__ << endl;
      cout << MUSEINFO_VERSION << endl;
      exit(0);
   } else if (opts.getBoolean("help")) {
      usage(opts.getCommand());
      exit(0);
   } else if (opts.getBoolean("example")) {
      example();
      exit(0);
   }

   rangeStart = opts.getDouble("range");
   decay      = opts.getDouble("decay");
   trialCount = opts.getInteger("count");
   stopLimit  = opts.getInteger("limit");
   debugQ     = opts.getBoolean("debug");
   seed       = opts.getInteger("seed");

   if (seed == 0) {
     seed = time(NULL);
   }

   srand48(seed);
}



//////////////////////////////
//
// example -- example usage of the maxent program
//

void example(void) {
   cout <<
   "                                                                        \n"
   << endl;
}



//////////////////////////////
//
// usage -- gives the usage statement for the quality program
//

void usage(const char* command) {
   cout <<
   "                                                                        \n"
   << endl;
}



// md5sum: 8bc0eeac6bd40b6cb75af3b19852c871 iwsearch.cpp [20050403]