/***************************************************************************
                          Segmenter.cpp  -  description
                             -------------------
    begin                : Don Feb 10 2005
    copyright            : (C) 2005 by klaus frieler
    email                : kgf@omniversum.de
 ***************************************************************************/

/***************************************************************************
 *                                                                         *
 *   This program is free software; you can redistribute it and/or modify  *
 *   it under the terms of the GNU General Public License as published by  *
 *   the Free Software Foundation; either version 2 of the License, or     *
 *   (at your option) any later version.                                   *
 *                                                                         *
 ***************************************************************************/
#include "Segmenter.h"


///////////////////////////////////////
// LBDMBase 
///////////////////////////////////////

dData LBDMBase::getIOIs() 
{
   dData ret;
   if(!m_pMel){
      m_sErrText = "Invalid melody object";
      return ret;
   }
   ret = m_pMel->getIOIs();
   return ret;
}

dData LBDMBase::getOOIs() 
{
   dData ret;
   if(!m_pMel){
      m_sErrText = "Invalid melody object";
      return ret;
   }
   ret = m_pMel->getOOIs();
   return ret;
}

dData LBDMBase::getPitchIntervals()
{
   dData ret;
   if(!m_pMel){
      m_sErrText = "Invalid melody object";
      return ret;
   }
   iData tmp = m_pMel->getPitches();
   if(tmp.empty()){
      m_sErrText = "Empty melody";
      return ret;
   }
   for(uint i = 0; i<tmp.size()-1; ++i){
      double interval = abs((int)tmp[i+1]- (int)tmp[i]);
      ret.push_back(interval);
   }
   return ret;
}

Vector LBDMBase::getBoundaries(const Vector& bounds) 
{
   Vector ret;
   double val = 0;
   if(bounds.size()<3){
      m_sErrText = "LBDMBase: Bounds empty";    
      return ret;
   }
   if(bounds[0]>bounds[1]){
      val = 1.0;
   }
   if(bounds[0]==bounds[1] ){
      val = 0.0;
   }
   ret.push_back(val);

   for(uint i = 1; i<bounds.size()-1; ++i){
      val  = 0;
      double left = bounds[i-1], mid = bounds[i], right = bounds[i+1];
      
      //cout<<setprecision(2)<<"L: ="<<left<<" M: "<<mid<<" R: = "<<right<<endl;
      if(left <= mid && right<= mid){
         val = 1.0;
         if(left==mid || right==mid){
            val = 0.5;
         }
         if(left==mid && right==mid){
            val = 0.0;
         }
      }
      ret.push_back(val);
   }
   uint sz = (uint)bounds.size()-1;
   val = 0;
   if(bounds[sz]>bounds[sz-1]){
      val = 1.0;
   }
   if(bounds[sz]==bounds[sz-1] ){
      val = 0.0;
   }
   ret.push_back(val);

   return ret;
}

Vector LBDMBase::getSecondBoundaries(const Vector& bounds)
{
   Vector ret;
   iData pos1, pos2, pos3;
   Vector locMax     = bounds.getLocalMaxima(pos1);
   Vector loclocMax  = locMax.getLocalMaxima(pos2);

   for(uint i = 0; i<pos2.size();++i){
      pos2[i] = pos1[pos2[i]];
   }
   pos1.clear();
   for(uint i = 0; i<pos2.size();++i){
      for(uint j = 0; i+j+1<pos2.size() && pos2[i+j] == pos2[i+j+1]; ++j){++i;}
      pos1.push_back(pos2[i]);
   }

   for(uint i = 0; i<bounds.size(); ++i){
      ret.push_back(0.0);
   }
   for(uint i = 0; i<bounds.size(); ++i){
      for(uint j = 0; j < pos1.size(); ++j){
         if(i == pos1[j]){ 
            ret[i] = 1.0;
            break;
         }
      }
   }
   return ret;
}

///////////////////////////////////////
//PlainLBDM 
///////////////////////////////////////

int PlainLBDM::getIntervalRelation(const double d1, const double d2) const
{
   return d1>=d2 ? (d1>d2 ? 1 : 0) : -1;
}

Vector PlainLBDM::assignBoundaryStrengths(const dData& intervals)
{
   Vector    bounds;
   double   sum = 0;
   int      rel1 = 0, rel2  = 0;
   if(intervals.empty()){
      m_sErrText = "PlainLBDM: Internal error";    
      return bounds;
   }
   for(uint i = 0; i<intervals.size(); ++i){
      sum   = 0.0;
      rel1  = 0;
      rel2  = 0;
      if(i>0){
         rel1 = getIntervalRelation(intervals[i-1], intervals[i]); 
         //cout<<"int[i-1] = "<<intervals[i-1]<<" int[i] = "<<intervals[i];
      }
      if(i<intervals.size()-1){
         rel2 = getIntervalRelation(intervals[i], intervals[i+1]); 
         //cout<<" int[i+1] = "<< intervals[i+1]<<endl;
      }
      //cout<<" R1 = "<<rel1<<" R2 = "<<rel2<<endl;
      if(rel2 != 0 && m_param.bUseBonus){
         sum += m_param.bonus;
      }
      
      if(rel1 != 0){
         sum += m_param.diff;
         if(rel1<0) sum += m_param.bigger;
       }
      if(rel2 != 0){
         sum += m_param.diff;
         if(rel2>0) sum += m_param.bigger;
      }
      //cout<<"Sum = "<<sum<<endl;
      bounds.push_back(sum);
   }

   return bounds;
}

Vector PlainLBDM::analyze()
{
   Vector bounds;
   dData intervals;
   switch(m_interval)
   {
      case iois:
         intervals = getIOIs();
         m_param.bUseBonus = true;
         break;
      case pitches:
         intervals = getPitchIntervals();
         m_param.bUseBonus = false;
      break;
      case rests:
         intervals = getOOIs();
         m_param.bUseBonus = false;
      break;
      default:
         m_sErrText = "PlainLBDM: Internal error";    
      break;
   }
   if(intervals.empty()){
      m_sErrText = "PlainLBDM: Internal error";    
      return bounds;
   }
   bounds  = assignBoundaryStrengths(intervals);
   return bounds;
}

///////////////////////////////////////
//RefinedLBDM 
///////////////////////////////////////

Vector RefinedLBDM::analyze()
{
   Vector bounds;
   dData intervals;
   switch(m_interval)
   {
      case iois:
         intervals = getIOIs();
         break;
      case pitches:
         intervals = getPitchIntervals();
      break;
      case rests:
         intervals = getOOIs();
      break;
      default:
         m_sErrText = "RefinedLBDM: Internal error";    
      break;
   }
   if(intervals.empty()){
      m_sErrText = "RefinedLBDM: Internal error";    
      return bounds;
   }
   bounds  = assignBoundaryStrengths(intervals);
   double max = bounds.getMax();
   if(max >0){
     bounds *= 1/max;
   }
   //bounds.dump();
   return bounds;
}

Vector RefinedLBDM::assignBoundaryStrengths(const dData& intervals)
{
   Vector    bounds;
   double   alpha1 = 0, alpha2 = 0, sum = 0;
   
   if(intervals.empty()){
      m_sErrText = "RefinedLBDM: Internal error";    
      return bounds;
   }
   Vector dumper(intervals);
   for(uint i = 0; i<intervals.size(); ++i){
      alpha1 = 0; alpha2 = 0;
      if(i>0){
         alpha1 = intervals[i] * alpha(intervals[i-1], intervals[i]);
         //cout<<"intervals[i] = "<< intervals[i]<<" intervals[i-1] = "<<intervals[i-1]<<" alpha1 = "<<alpha1<<endl;
      }
      if(i<intervals.size()-1){
         alpha2 = intervals[i] * alpha(intervals[i], intervals[i+1]);
         //cout<<"intervals[i] = "<< intervals[i]<<" intervals[i+1] = "<<intervals[i+1]<<" alpha2 = "<<alpha2<<endl;
      }
      sum = alpha1+ alpha2;
      //cout<<"sum = "<<sum<<endl;
      bounds.push_back(sum);
   }
   return bounds;
}

double RefinedLBDM::alpha(double x, double y)const
{
   double diff = x-y, sum = x+y;
   if(diff == 0 || sum == 0) return 0;
   diff = fabs(diff/(x+y));
   if(diff<1e-6) diff = 0.0;
   return diff;
}

///////////////////////////////////////
//TemperleySegmenter implementation
///////////////////////////////////////

dData TemperleySegmenter::getIOIs() 
{
   dData ret;
   if(!m_pMel){
      m_sErrText = "Invalid melody object";
      return ret;
   }
   ret = m_pMel->getIOIs();
   return ret;
}

dData TemperleySegmenter::getOOIs() 
{
   dData ret;
   if(!m_pMel){
      m_sErrText = "Invalid melody object";
      return ret;
   }
   ret = m_pMel->getOOIs();
   return ret;
}


Vector TemperleySegmenter::analyze()
{
   Vector ret;
   Vector gapScores = assignGapScores();
   if(gapScores.empty()) return ret;
   double *analysis  = new double[gapScores.size()];
   uint   *best      = new uint[gapScores.size()];
   if(!analysis || !best){
      m_sErrText = "Out of memory";
      return ret;
   }
   memset(analysis, 0, sizeof(double)*gapScores.size());
   memset(best, 0, sizeof(uint)*gapScores.size());
   for(uint i = 1; i<gapScores.size(); ++i){
      double bestScore = -1e3;
      //cout<<endl;
      //m_pMel->at(i).dump();
      //cout<<"Gap "<<i<<" "<<gapScores[i]<<endl;
      for(uint j = 1; i>=j && j <m_params.testLength; ++j){
         double lenPenalty = getLengthPenalty(j);
         double lenFactor = sqrt((double)j);
         double phasePenalty = getPhasePenalty(m_pMel->at(i), m_pMel->at(i-j));
         double gap = (gapScores[i-j] + gapScores[i])/2.0;
         //if(gapScores[i-j] == 0 || gapScores[i] == 0) gap *=2;
         double score = lenFactor*( gap - lenPenalty - phasePenalty);
         
         //double score = lenFactor*( 0- lenPenalty );
         //cout<<"Testing i="<<i<<" with "<<j<<endl;
         //cout<<i<<" "<<j<<" LenPenalty = "<<lenPenalty<<" LenFactor = "<<lenFactor<<" Phasepenalty = "<<phasePenalty<<endl;
         //cout<<i<<" "<<j<<" gap[i-j] = "<<gapScores[i-j]<<" gap[i] = "<<gapScores[i]<<" gap = "<<gap<<endl;
         //cout<<i<<" "<<j<<" score = "<<score<<" ";
         
         score +=analysis[i-j];
         
         //cout<<"analysis = "<<analysis[i-j]<<" score = "<<score<<endl;
         //cout<<" total score "<<score<<" bestscore = "<<bestScore<<endl;
         if(score>bestScore){
            best[i] = j;
            analysis[i] = score;
            bestScore   = score;
            
            //cout<<"analysis = "<<analysis[i]<<" score = "<<score<<" Bestscore= "<<bestScore<<" Best i = "<<best[i]<<endl;
         }
      }
   }
   for(uint i = 0; i<gapScores.size(); ++i){
      //cout<<"Old: Analysis["<<i<<"] = "<<analysis[i]<<" best = "<<best[i]<<endl;
   }
   memset(analysis, 0, sizeof(double)*gapScores.size());

   uint k = gapScores.size()-1;
   while(k>0){
      //cout<<"k = "<<k<<" best[k] = "<<best[k]<<endl;
      analysis[k-best[k]] = 1.0;
      k = k-best[k];
   }
   analysis[0] = 0;
   for(uint i = 1; i<gapScores.size(); ++i){
      //cout<<"i = "<<i<<" ret = "<<analysis[i]<<endl;
      ret.push_back(analysis[i]);
   }
   delete [] analysis;
   delete [] best;
   return ret;
}

Vector TemperleySegmenter::assignGapScores() 
{
   dData iois = getIOIs();
   dData oois = getOOIs();
   Vector ret;
   if(!iois.size() || iois.size() != oois.size()) return ret;
   double sumIOI = 0;
   ret.push_back(0);
   for(uint i = 0; i<iois.size(); ++i){
      sumIOI += iois[i];
      double gapScore = (i+1)*m_params.gapWeight*(iois[i] + oois[i])/sumIOI;
      //cout<<"ioi = "<<iois[i]<<" ooi = "<<oois[i]<<" sumIOI = "<<sumIOI/(i+1)<<" gapScore = "<<gapScore<<endl;
      ret.push_back(gapScore);
   }
   return ret;
}

double  TemperleySegmenter::getPhasePenalty(const Tone& tone1, const Tone& tone2) const
{
   BarPosition pos1  = tone1.getBarPosition(), pos2 = tone2.getBarPosition();
   Signature   sig   = pos1.getSignature();
   if(pos1 == invalidBarPosition || pos2 == invalidBarPosition) return 0.0;

   double penalty = m_params.phase4Penalty;
   if(pos1.getTicks() != pos2.getTicks()) return penalty;

   uint beat1 = pos1.getBeat(), beat2 = pos2.getBeat();
   if(beat1 == beat2) return 0.0;
   for(uint num = 2; num<6; ++num){
      if(sig.count()>num && sig.count()%num == 0){
         if((beat1-beat2)%(sig.count()/num)== 0){
            penalty = m_params.phase3Penalty;
            break;
         }
      }
   }
   return penalty;
}

template<class T>
vector<T> getModus(const vector<T> &vec, uint &modFreq)
{
   uint i = 0, j = 0;
   vector<T> modi;
   vector<pair<T, uint> > freq;
   uint modus = 0;
   modFreq = 1;
   if(vec.empty()) return modi;
   for( i = 0; i<vec.size();++i){
      T val = vec[i];
      //val in modi?
      bool bFound = false;
      for(j= 0; j< freq.size(); ++j){
         if(freq[j].first == val){
            //found: increase frequency
            freq[j].second++;
            if(freq[j].second>modFreq) modFreq = freq[j].second;
            bFound = true;
         }
      }
      //if not found add to list
      if(!bFound){
         pair<T, uint> n;
         n.first = val;
         n.second = 1;
         freq.push_back(n);
      }
   }
   //return vector of values with highest modus
   for(i = 0; i<freq.size(); ++i){
      if(freq[i].second == modFreq)
         modi.push_back(freq[i].first);
   }
   return modi;
}
template<class T>
T getMinimum(const vector<T>& vec)
{
   uint dummy;
   return getMinimum<T>(vec, dummy);
}

template<class T>
T getMinimum(const vector<T>& vec, uint& pos)
{
   uint i = 0;
   if(vec.empty()) return T(0);
   T min = vec[0];

   for(i = 0; i < vec.size();++i){
      if(vec[i]<min){
          min = vec[i];
          pos = i;
      }
   }
   return min;
}

Vector SimpleSegmenter::analyze()
{
   Vector  phrases;
   uint  i = 0;
   uint dummy = 0;
   uint size = m_pMel->size();
   iData durTicks    = m_pMel->getDurationTicks();
   dData dur         = m_pMel->getDurations();
   iData modi = getModus<uint>(durTicks, dummy);
   uint modus = getMinimum<uint>(modi);
   iData::iterator durIter = durTicks.begin();
   if(size == 0){
      return phrases;
   }
   phrases.assign( size, 0);
   phrases[size-1] = 1.0;
   if(size<3) return phrases;

   for(i = 1; i< durTicks.size()-2; ++i){
      if(m_params.m_bRunningModus){
         durIter++;
         if(i >= 1){
            iData durTicksSeq;
            durTicksSeq.insert(durTicksSeq.begin(), durTicks.begin(), durIter);  
            //Vector tmp(durTicksSeq);
            //tmp.dump();
            iData tmpModi = getModus<uint>(durTicksSeq, dummy);
            modus = getMinimum<uint>(tmpModi);
         }
      }
      //cout<<"Modus: "<<modus<<endl;
      if(dur[i]>m_params.m_gapTime|| durTicks[i]>=m_params.m_modusFactor*modus){
         phrases[i]   = 1.0;
      }
   }
   //dumpData("Durations", m_durations, out());
   //out()<<"Modus: "<<modus<<endl;
   //dumpData("Phrases  ", phrases, out());
   //out()<<"Lengths: "<<phrases.size()<<" "<<m_durations.size()<<" "<<m_pitches.size()<<endl;
   return phrases;   
}
