//
//  frequencer.cpp
//  oscSenderExample
//
//  Created by Robert Tubb on 27/11/2012.
//
//

#include "frequencer.h"
Frequencer frequencer;

Frequencer::Frequencer(){
    
    N = 16;

    for(int i = 0; i<N;i++){
        timePoints.push_back(0.0);
        freqMag.push_back(0.0);
        freqPhase.push_back(0.0);
    }
    
}
/////////////////////////////////////////////////////
/////////////////////////////////////////////////////
// No time domain!
void Frequencer::timeEdit(int n, double value){

    // update internal
    timePoints[n] = value;
    
    // calculate transform
    dft();

}

void Frequencer::timeEdit(vector<double> allValues){
    
    // update internal
    vector<double>::iterator iter;
    vector<double>::iterator iter2;
    iter = allValues.begin();
    iter2 = timePoints.begin();
    for (iter = allValues.begin(); iter < allValues.end(); iter++,iter2++){
        *iter = *iter2;
    }
    
    
    // calculate transform
    dft();
    
    //TODO  output freq domain - to PD
    
}

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


vector<double> Frequencer::freqMagEdit(vector<int> indexes, vector<double> values){
    if(indexes.size() != values.size()){
        cout << "freqMagEdit error: indexes dont match values\n";
        return timePoints;
    }
    // put into local store,
    for(int i=0;i<indexes.size();i++){
        freqMag[indexes[i]] = values[i];
        freqMag[N-indexes[i]-1] = freqMag[indexes[i]];
    }
   
    // calculate transform
    idft();

    return timePoints;
}


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

void Frequencer::freqPhaseEdit(int k, double phase){
    N = 16;
    // put into local store
    freqPhase[k] = phase;
    
    // making sure to conjugate
    
    freqPhase[N-k-1] = - freqPhase[k];
    
    // calculate transform
    idft();
    // TODO output time domain
   
}

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

void Frequencer::dft(){
    
    // we have 16 time points. work out freqs and phases
    int n = 0;
    int k = 0;
    double xre = 0.;
    double xim = 0.;
    
    for(k=0; k < N; k++){
        // work out xn
        xre = 0.;
        xim = 0.;
        for(n=0;n < N;n++){
            xre = xre + timePoints[n]*cos( (2.*PI/(float)N)*n*k);
            xim = xim + timePoints[n]*sin( (2.*PI/(float)N)*n*k);
            
        }
        freqMag[k] = sqrt(pow(xre,2) + pow(xim,2.)); // MAGNITUDE
        freqPhase[k] = -atan2(xim,xre);
    }
    
}



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

void Frequencer::idft(){
    
    // we have some freqs and phases. work out real time points (conj)
    int n = 0;
    int k = 0;
    double xre = 0.;
    double xim = 0.;
    
    for(n=0; n < N; n++){
        // work out xn
        xre = 0.;
        xim = 0.;
        for(k=0;k < N;k++){
            xre = xre + freqMag[k]*cos( (2.*PI/(float)N)*n*k + freqPhase[k]);
            
        }
        timePoints[n] = xre/(float)N; // divided by N ???
    }
    
}
/////////////////////////////////////////////////////
/////////////////////////////////////////////////////