rt300@0: //
rt300@0: //  dsptools.cpp
rt300@0: //  wablet
rt300@0: //
rt300@0: //  Created by Robert Tubb on 21/06/2011.
rt300@0: //  Copyright 2011 __MyCompanyName__. All rights reserved.
rt300@0: //
rt300@0: 
rt300@0: #include "dsptools.h"
rt300@0: #include "ofMain.h"
rt300@0: 
rt300@0: DSPTools::DSPTools(){
rt300@0:     for(int i = 0; i < 3; i++){
rt300@0:         xv[i] = 0.0;
rt300@0:         yv[i] = 0.0;
rt300@0:     }
rt300@0:     getHPCoefficientsButterworth2Pole(ax, by);
rt300@0: }
rt300@0: 
rt300@0: double DSPTools::highpass1(double ax){
rt300@0:     static double xm1 = 0.0;
rt300@0:     double sample;
rt300@0:     sample = 2*(ax - xm1);
rt300@0:     xm1 = ax;
rt300@0:     return sample;    
rt300@0:     
rt300@0: }
rt300@0: double DSPTools::lowpass1(double ax){
rt300@0:     static double xm1 = 0.0;
rt300@0:     double sample;
rt300@0:     sample = 0.5*(ax + xm1);
rt300@0:     xm1 = ax;
rt300@0:     return sample;
rt300@0:     
rt300@0: }
rt300@0: 
rt300@0: void DSPTools::getLPCoefficientsButterworth2Pole(const int samplerate, const double cutoff, double* const ax, double* const by)
rt300@0: {
rt300@0:     double sqrt2 = 1.4142135623730950488;
rt300@0:     
rt300@0:     double QcRaw  = (2 * PI * cutoff) / samplerate; // Find cutoff frequency in [0..PI]
rt300@0:     double QcWarp = tan(QcRaw); // Warp cutoff frequency
rt300@0:     
rt300@0:     double gain = 1 / (1+sqrt2/QcWarp + 2/(QcWarp*QcWarp));
rt300@0:     by[2] = (1 - sqrt2/QcWarp + 2/(QcWarp*QcWarp)) * gain;
rt300@0:     by[1] = (2 - 2 * 2/(QcWarp*QcWarp)) * gain;
rt300@0:     by[0] = 1;
rt300@0:     ax[0] = 1 * gain;
rt300@0:     ax[1] = 2 * gain;
rt300@0:     ax[2] = 1 * gain;
rt300@0: }
rt300@0: void DSPTools::getHPCoefficientsButterworth2Pole(double* const ax, double* const by)
rt300@0: {
rt300@0:     ax[0] = 0.997987115675119;
rt300@0:     ax[1] = -1.995974231350238;
rt300@0:     ax[2] = 0.997987115675119;
rt300@0: 
rt300@0:     by[0] = 1.000000000000000;
rt300@0:     by[1] = -1.995970179642828;
rt300@0:     by[2] = 0.995978283057647;
rt300@0: }
rt300@0: 
rt300@0: double DSPTools::butter(double sample)
rt300@0: {
rt300@0: 
rt300@0:         xv[2] = xv[1]; 
rt300@0:         xv[1] = xv[0];
rt300@0:         xv[0] = sample;
rt300@0:         yv[2] = yv[1]; 
rt300@0:         yv[1] = yv[0];
rt300@0:         
rt300@0:         yv[0] =   (ax[0] * xv[0] + ax[1] * xv[1] + ax[2] * xv[2]
rt300@0:                    - by[1] * yv[1]
rt300@0:                    - by[2] * yv[2]);
rt300@0:         
rt300@0:         return yv[0];
rt300@0: 
rt300@0: }