//
// Programmer:    Craig Stuart Sapp <craig@ccrma.stanford.edu>
// Creation Date: Wed Jan 10 12:32:52 PST 2001
// Last Modified: Mon Jan 29 10:39:20 PST 2001
// Filename:      ...sig/examples/all/rootspace.cpp
// Web Address:   http://sig.sapp.org/examples/museinfo/humdrum/rootspace.cpp
// Syntax:        C++; museinfo
//
// Description:   Measure the chord information for a given space.
// 

#include "humdrum.h"

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


// function declarations
void   checkOptions             (Options& opts, int argc, char* argv[]);
void   example                  (void);
void   usage                    (const char* command);
int    findMinimum              (Array<double>& values);
void   fillDataArrays           (HumdrumFile& infile, Array<int>& pitch, 
                                 Array<int>& testset, Array<double>& duration, 
                                 Array<double>& level, int& target, 
                                 int& targetindex);
int    generateAnalysis         (Array<int>& pitch, Array<double>& duration, 
                                 Array<double>& level, Array<int>& testset, 
                                 double sx, double e1, double e2,
                                 Array<double>& coef);
int    intcompare               (const void* a, const void* b);
double secondbest               (int bestindex, Array<double>& coef);

// pitch space variables:
int vectorxsq[40] = {
            0, 49,      49, 1000000, 16, 16,      16,  16, 16, 1000000, 
           64,  1,       1,      64, 64, 25,      25,  25, 25,      25, 
      1000000, 81,       4,       4,  4, 81, 1000000,  36, 36,      36, 
           36, 36, 1000000,       9,  9,  9,       9, 100, 49,      49
};

int vectorysq[40] = {
            0,   4,      16, 1000000,  16,   4,       0,   4,  16, 1000000, 
            9,   1,       1,       9,  25,  25,       9,   1,   1,       9, 
      1000000,  16,       4,       0,   4,  16, 1000000,  16,   4,       0,
            4,  16, 1000000,       9,   1,   1,       9,  25,  16,       4
};

// user interface variables
Options      options;            // database for command-line arguments
int          quietQ = 0;         // display only data table or not
int          verboseQ = 1;       // display only data table or not
int          compoundQ  = 1;     // used with the -c option
double       deltamin;           // used with the --deltamin option
double       deltamax;           // used with the --deltamax option
double       deltastep;          // used with the --dstep option
double       lambdamin;          // used with the --lambdamin option
double       lambdamax;          // used with the --lambdamax option
double       lambdastep;         // used with the --lstep option
double       sxmin;              // used with the --sxmin option
double       sxmax;              // used with the --sxmax option
double       xstep;              // used with the --xstep option
double       sy;                 // used with the --sy option


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

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

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

   // figure out the number of input files to process
   int numinputs = options.getArgCount();

   Array<double> tonicscores;
   Array<int> pitch;
   Array<int> testset;
   Array<double> duration;
   Array<double> level;
   Array<double> coef;
   int target = -1;           // the correct root
   int targetindex = -1;      // the correct root's index in the pitches array

   int i, root;
   double ewt, rwt;
   double a, b, c;
   double az, bz, cz;
   az = xstep/1000;
   bz = deltastep/1000;
   cz = lambdastep/1000;
   for (i=0; i<numinputs || i==0; i++) {
      infile.clear();

      // if no command-line arguments read data file from standard input
      if (numinputs < 1) {
         infile.read(cin);
      } else {
         infile.read(options.getArg(i+1));
      }

      fillDataArrays(infile, pitch, testset, duration, level, target, targetindex);
      cout<<"**alpha\t**delta\t**lambd\t**eval"
            "\t**eroot\t**ewht\t**root\t**rwht\t**diff\n";
      for (a=sxmin; a<=sxmax+az; a+=xstep) {
         for (b=deltamin; b<=deltamax+bz; b+=deltastep) {
            for (c=lambdamin; c<=lambdamax+cz; c+=lambdastep) {
               root =generateAnalysis(pitch, duration, level, testset, a, b, c,
                     coef);
               cout << a << "\t" << b << "\t" << c << "\t";
               if (testset[root] == target) {
                  cout << 0;
               } else {
                  cout << 1;
               }
               if (verboseQ) {
                  ewt = coef[root]/pitch.getSize();
                  rwt = coef[targetindex]/pitch.getSize();
                  cout << "\t" << testset[root] 
                       << "\t" << ewt;
                  cout << "\t" << testset[targetindex] 
                       << "\t" << rwt;
                  if (testset[root] == target) {
                     cout << "\t" << secondbest(targetindex, coef);
                  } else {
                     cout << "\t" << rwt - ewt;
                  }
               }
               cout << "\n";
            }
         }
      }
      cout << "*-\t*-\t*-\t*-\t*-\t*-\t*-\t*-\t*-\n";
   }

   return 0;
}


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


//////////////////////////////
//
// secondbest -- find the second best score.  Not necessarily the
//    second beat, but close enough.
//

double secondbest(int bestindex, Array<double>& coef) {
   int secondbest = 0;
   if (coef.getSize() < 2) {
      return 0.0;
   }
   if (bestindex == 0) {
      secondbest = 1;
   } else {
      secondbest = 0;
   }

   int i;
   for (i=0; i<coef.getSize(); i++) {
      if (i != bestindex && coef[i] < coef[secondbest]) {
         secondbest = i;
      }
   }

   return (coef[bestindex] - coef[secondbest])/coef.getSize();
}



//////////////////////////////
//
// generateAnalysis -- only test notes which are present in the
//      chord.
//

int generateAnalysis(Array<int>& pitch, Array<double>& duration, 
      Array<double>& level, Array<int>& testset, double sx, double e1, double e2,
      Array<double>& coef) {

   double suma;
   double sumb;
   coef.setSize(testset.getSize());
   int rootindex = 0;
   double normvalue;
   double sxsq = sx * sx;
   int min = 0;
   int i, j;

   for (i=0; i<testset.getSize(); i++) {
      suma = 0.0;
      sumb = 0.0;
      for (j=0; j<pitch.getSize(); j++) {
         rootindex = ((int)pitch[j] - testset[i] + 40) % 40;
         normvalue = sqrt(vectorxsq[rootindex] * sxsq + vectorysq[rootindex]);
         suma += normvalue * (duration[j] + e1);
         sumb += normvalue * (level[j] + e2);
      }
      coef[i] = sqrt(suma * suma + sumb * sumb);
      if (coef[i] < coef[min]) {
         min = i;
      }
   }
 
   return min;          
}



//////////////////////////////
//
// fillDataArrays -- 
//

void fillDataArrays(HumdrumFile& infile, Array<int>& pitches, Array<int>& testset,
      Array<double>& durations, Array<double>& levels, int& target, 
      int& targetindex) {
   pitches.setSize(0);
   durations.setSize(0);
   levels.setSize(0);
   pitches.allowGrowth(1);
   durations.allowGrowth(1);
   levels.allowGrowth(1);

   Array<int> scorelevels;
   infile.analyzeMetricLevel(scorelevels);
   
   int firsttime = 1;
   int i, j, k;
   int ii, jj;
   int ccount;
   static char buffer[1024] = {0};
   int pitch;
   double duration;
   double level;
   for (i=0; i<infile.getNumLines(); i++) {
      if (infile[i].getType() != E_humrec_data) {
         // ignore non-note data lines
         continue;
      }
      ii = i;
      for (j=0; j<infile[i].getFieldCount(); j++) {
         if (infile[i].getExInterpNum(j) != E_KERN_EXINT) {
            // ignore non-kern data spines
            continue;
         }
         if (firsttime && strcmp(infile[i][j], ".") == 0) {
            // extract the held over note from a previous point in the infile
            ii = infile[i].getDotLine(j);
            jj = infile[i].getDotSpine(j);
         } else {
            ii = i;
            jj = j;
         }
 
         if (strcmp(infile[ii][jj], ".") != 0) {
            // extract all notes in the region of interest, ignoring
            // tied notes.
            ccount = infile[ii].getTokenCount(jj);
            for (k=0; k<ccount; k++) {
               infile[ii].getToken(buffer, jj, k, 128);
               if (strchr(buffer, 'r') != NULL) {
                  // skip over rests
                  continue;
               }
               if (strchr(buffer, '_') != NULL) {
                  // skip over doubly tied notes
                  continue;                

               }
               if (!firsttime && strchr(buffer, ']') != NULL) {
                  // skip over tied notes at the ends of ties.
                  continue;
               }
               // have a note so now extract the metric level and the duration
               pitch = Convert::kernToBase40(buffer);
               if (pitch < 0) {
                  // ignore rests
                  continue;
               }
               // pitch = ((int)pitch - 2 + 40) % 40;
               duration = infile.getTiedDuration(ii, jj, k);
               if (duration == 0.0) {
                  // ignore grace notes and other zero-dur ornaments
                  continue;
               }
               level = scorelevels[ii];
               durations.append(duration);
               levels.append(level);
               pitches.append(pitch);
 
            } // end of a chord
         }
      }  // end of a line
      firsttime = 0;
 
   } // end of the music selection   


   pitches.allowGrowth(0);
   durations.allowGrowth(0);
   levels.allowGrowth(0);

   double temp;
   if (!quietQ) {
      cout << "!! INPUT DATA: \n";
      cout << "!! **p\t**dur\t**level\t**b40pc\t**log2d\t**log2l\n";
      for (i=0; i<pitches.getSize(); i++) {
         temp = -log(durations[i])/log(2.0);
         if (fabs(temp) == 0.0) {
            temp = 0.0;
         }
         cout << "!! " << Convert::base40ToKern(buffer, pitches[i]) << "\t"
              << durations[i] << "\t"
              << 1.0/pow(2.0, levels[i]) << "\t"
              << (pitches[i] - 2 + 40) % 40 << "\t"
              << temp << "\t"
              << levels[i] << "\n";
      }
      cout << "!! *-\t*-\t*-\t*-\t*-\t*-\n";
   }

   for (i=0; i<pitches.getSize(); i++) {
      pitches[i] = (pitches[i] - 2 + 40) % 40;
      durations[i] = -log(durations[i])/log(2.0);
   }

   // make a list of possible roots (assuming only pitches
   // present in the chord can be possible roots):
   testset = pitches;
   qsort(testset.getBase(), testset.getSize(), sizeof(int), intcompare);
   int uniqc = 0;
   Array<int> nodoubles(testset.getSize());
   nodoubles.setSize(0);
   nodoubles.append(testset[0]);
   for (i=1; i<testset.getSize(); i++) {
      if (testset[i] == testset[uniqc]) {
         testset[i] = 1000;
      } else {
         uniqc = i;
         nodoubles.append(testset[i]);
      }
   }
   testset = nodoubles;

//   if (!quietQ) {
//      for (i=0; i<testset.getSize(); i++) {
//         cout << "!! Uniq Pitch: " << testset[i] << endl;
//      }
//   }

   // extract the target value from the **root spine
   for (i=0; i<infile.getNumLines(); i++) {
      if (infile[i].getType() != E_humrec_data) { 
         continue;
      }
      for (j=0; j<infile[i].getFieldCount(); j++) {
         if (strcmp("**root", infile[i].getExInterp(j)) == 0) {
            if (strcmp(".", infile[i][j]) != 0 && strchr(infile[i][j], 'r') == 0) {
               target = Convert::kernToBase40(infile[i][j]);
               target = (target - 2 + 40) % 40;
	       if (!quietQ) {
                  cout << "!!!TARGET-ROOT: " << infile[i][j] << "\n";
               }
               break;
            }
         }
      }
   }

   // figure out which pitch index the target is found at
   targetindex = -1;
   for (i=0; i<testset.getSize(); i++) {
      if (testset[i] == target) {
         targetindex = i;
         break;
      } 
   }
   if (targetindex < 0) {
      cout << "Error: target root is not in the test set of pitches" << endl;
      exit(1);
   }

   if (options.getBoolean("dataonly")) {
      exit(0);
   }

}



//////////////////////////////
//
// intcompare --
//
 
int intcompare(const void* a, const void* b) {
   if (*((int*)a) < *((int*)b)) {
      return -1;
   } else if (*((int*)a) > *((int*)b)) {
      return 1;
   } else {
      return 0;
   }
} 



//////////////////////////////
//
// findMinimum --
//

int findMinimum(Array<double>& values) {
   int i;
   int min = 0;
   for (i=1; i<values.getSize(); i++) {
      if (values[i] < values[min]) {
         min = i;
      }
   }
   return min;
}



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

void checkOptions(Options& opts, int argc, char* argv[]) {
   opts.define("a|append=b",        "append analysis to data in output");   
   opts.define("q|quiet=b",         "suppress extra information");   
   opts.define("v|noverbose=b",     "suppress extra information");   
   opts.define("C|compound=b",      "don't try to use compound meters");   
   opts.define("deltamin|e1min=d:-10.0",   "first rhythm scalar min");   
   opts.define("deltamax|e1max=d:10.0",    "first rhythm scalar max");   
   opts.define("deltastep|e1step=d:0.25",  "first rhythm scalar step");   
   opts.define("lambdamin|e2min=d:-10.0",  "second rhythm scalar min");   
   opts.define("lambdamax|e2max=d:10.0",   "second rhythm scalar min");   
   opts.define("lambdastep|e2step=d:0.25", "second rhythm scalar step");   
   opts.define("sxmin=d:0.5",       "first pitch scalar min");   
   opts.define("sxmax=d:1.0",       "first pitch scalar max");   
   opts.define("xstep=d:0.015625",  "first pitch scalar step");   
   opts.define("sy=d:1",            "second pitch scalar");   
   opts.define("dataonly=b",        "show input data for calculation only");
   opts.define("s|scores=b",        "show scores");   
   opts.define("n|numeric=b",       "show result in numeric root/mode");   

   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, March 2001" << endl;
      exit(0);
   } else if (opts.getBoolean("version")) {
      cout << argv[0] << ", version: 26 Mar 2001" << endl;
      cout << "compiled: " << __DATE__ << endl;
      cout << MUSEINFO_VERSION << endl;
      exit(0);
   } else if (opts.getBoolean("help")) {
      usage(opts.getCommand().c_str());
      exit(0);
   } else if (opts.getBoolean("example")) {
      example();
      exit(0);
   }

   // compound diabled for now
   if (opts.getBoolean("compound")) {
      compoundQ = 0;
   } else {
      compoundQ = 0;
   }

   quietQ   = opts.getBoolean("quiet");
   verboseQ = !opts.getBoolean("noverbose");

   deltamin   = opts.getDouble("deltamin");
   deltamax   = opts.getDouble("deltamax");
   deltastep  = opts.getDouble("deltastep");
   lambdamin  = opts.getDouble("lambdamin");
   lambdamax  = opts.getDouble("lambdamax");
   lambdastep = opts.getDouble("lambdastep");
   sxmin      = opts.getDouble("sxmin");
   sxmax      = opts.getDouble("sxmax");
   xstep      = opts.getDouble("xstep");
   sy         = opts.getDouble("sy");

   if (!quietQ) {
      cout << "!!!dim: " 
           << (sxmax - sxmin)/xstep+1 << ":"
           << (deltamax - deltamin)/deltastep+1 << ":"
           << (lambdamax - lambdamin)/lambdastep+1 << "\n";
      cout << "!!!alpha:" << sxmin << ":" 
           << xstep << ":" << sxmax << "\n";
      cout << "!!!delta:" << deltamin << ":" 
           << deltastep << ":" << deltamax << "\n";
      cout << "!!!lambda:" << lambdamin << ":" 
           << lambdastep << ":" << lambdamax << "\n";
   }

}



//////////////////////////////
//
// 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: bffd718c30435f0ce7373e932643de8e rootspace.cpp [20160320]