c@5: //
c@5: //  NoveltyCurve.cpp
c@5: //  Tempogram
c@5: //
c@5: //  Created by Carl Bussey on 10/07/2014.
c@5: //  Copyright (c) 2014 Carl Bussey. All rights reserved.
c@5: //
c@5: 
c@7: #include "NoveltyCurve.h"
c@5: using namespace std;
c@5: 
c@7: NoveltyCurve::NoveltyCurve(float samplingFrequency, int fftLength, int numberOfBlocks, int compressionConstant) :
c@5:     m_samplingFrequency(samplingFrequency),
c@7:     m_fftLength(fftLength),
c@7:     m_blockSize(fftLength/2 + 1),
c@5:     m_numberOfBlocks(numberOfBlocks),
c@5:     m_compressionConstant(compressionConstant),
c@5:     m_numberOfBands(5),
c@5:     m_bandBoundaries(NULL),
c@5:     m_hannLength(65),
c@5:     m_bandSum(NULL)
c@5: {
c@5:     initialise();
c@5: }
c@5: 
c@5: NoveltyCurve::~NoveltyCurve(){
c@5:     cleanup();
c@5: }
c@5: 
c@9: //allocate all space and set variable
c@5: void
c@5: NoveltyCurve::initialise(){
c@5:     data = vector<float>(m_numberOfBlocks);
c@5:     
c@9:     // for bandwise processing, the band is split into 5 bands. m_bandBoundaries contains the upper and lower bin boundaries for each band.
c@9:     m_bandBoundaries = new int[m_numberOfBands+1];
c@5:     m_bandBoundaries[0] = 0;
c@5:     for (int band = 1; band < m_numberOfBands; band++){
c@5:         float lowFreq = 500*pow(2.5, band-1);
c@7:         m_bandBoundaries[band] = m_fftLength*lowFreq/m_samplingFrequency;
c@5:     }
c@7:     m_bandBoundaries[m_numberOfBands] = m_blockSize;
c@5:     
c@7:     m_bandSum = new float [m_numberOfBands];
c@5: }
c@5: 
c@9: //delete space allocated in initialise()
c@5: void
c@5: NoveltyCurve::cleanup(){
c@5:     delete []m_bandBoundaries;
c@5:     m_bandBoundaries = NULL;
c@5:     delete []m_bandSum;
c@5:     m_bandSum = NULL;
c@5: }
c@5: 
c@9: //calculate max of spectrogram
c@7: float NoveltyCurve::calculateMax(vector< vector<float> > &spectrogram){
c@5:     float max = 0;
c@5:     
c@5:     for (int j = 0; j < m_numberOfBlocks; j++){
c@7:         for (int i = 0; i < m_blockSize; i++){
c@7:             max = max > fabs(spectrogram[i][j]) ? max : fabs(spectrogram[i][j]);
c@5:         }
c@5:     }
c@5:     
c@5:     return max;
c@5: }
c@5: 
c@9: //subtract local average of novelty curve
c@9: //uses m_hannWindow as filter
c@7: void NoveltyCurve::subtractLocalAverage(vector<float> &noveltyCurve){
c@5:     vector<float> localAverage(m_numberOfBlocks);
c@5:     
c@9:     float * m_hannWindow = new float[m_hannLength];
c@9:     WindowFunction::hanning(m_hannWindow, m_hannLength, true);
c@9:     
c@5:     FIRFilter *filter = new FIRFilter(m_numberOfBlocks, m_hannLength);
c@5:     filter->process(&noveltyCurve[0], m_hannWindow, &localAverage[0]);
c@5:     delete filter;
c@7:     filter = NULL;
c@5:     
c@7:     assert(noveltyCurve.size() == m_numberOfBlocks);
c@5:     for (int i = 0; i < m_numberOfBlocks; i++){
c@5:         noveltyCurve[i] -= localAverage[i];
c@5:         noveltyCurve[i] = noveltyCurve[i] >= 0 ? noveltyCurve[i] : 0;
c@5:     }
c@9:     
c@9:     delete m_hannWindow;
c@9:     m_hannWindow = NULL;
c@5: }
c@5: 
c@9: //smoothed differentiator filter. Flips upper half of hanning window about y-axis to create coefficients.
c@7: void NoveltyCurve::smoothedDifferentiator(vector< vector<float> > &spectrogram, int smoothLength){
c@7:     
c@7:     float * diffHannWindow = new float [smoothLength];
c@7:     WindowFunction::hanning(diffHannWindow, smoothLength, true);
c@7:     
c@7:     if(smoothLength%2) diffHannWindow[(smoothLength+1)/2 - 1] = 0;
c@7:     for(int i = (smoothLength+1)/2; i < smoothLength; i++){
c@7:         diffHannWindow[i] = -diffHannWindow[i];
c@7:     }
c@7:     
c@7:     FIRFilter *smoothFilter = new FIRFilter(m_numberOfBlocks, smoothLength);
c@7:     
c@7:     for (int i = 0; i < m_blockSize; i++){
c@7:         smoothFilter->process(&spectrogram[i][0], diffHannWindow, &spectrogram[i][0]);
c@7:     }
c@7:     
c@7:     delete smoothFilter;
c@7:     smoothFilter = NULL;
c@7: }
c@7: 
c@9: //half rectification (set negative to zero)
c@9: void NoveltyCurve::halfWaveRectify(vector< vector<float> > &spectrogram){
c@7:     
c@7:     for (int block = 0; block < m_numberOfBlocks; block++){
c@7:         for (int k = 0; k < m_blockSize; k++){
c@7:             if (spectrogram[k][block] < 0.0) spectrogram[k][block] = 0.0;
c@7:         }
c@7:     }
c@7: }
c@7: 
c@9: //process method
c@5: vector<float>
c@7: NoveltyCurve::spectrogramToNoveltyCurve(vector< vector<float> > &spectrogram){
c@7:     
c@7:     assert(spectrogram.size() == m_blockSize);
c@7:     assert(spectrogram[0].size() == m_numberOfBlocks);
c@5:     
c@9:     //normalise and log spectrogram
c@5:     float normaliseScale = calculateMax(spectrogram);
c@5:     for (int block = 0; block < m_numberOfBlocks; block++){
c@7:         for (int k = 0; k < m_blockSize; k++){
c@7:             if(normaliseScale != 0.0) spectrogram[k][block] /= normaliseScale; //normalise
c@7:             spectrogram[k][block] = log(1+m_compressionConstant*spectrogram[k][block]);
c@7:         }
c@7:     }
c@7: 
c@9:     //smooted differentiator
c@7:     smoothedDifferentiator(spectrogram, 5); //make smoothLength a parameter!
c@9:     //halfwave rectification
c@7:     halfWaveRectify(spectrogram);
c@7:     
c@9:     //bandwise processing
c@7:     for (int block = 0; block < m_numberOfBlocks; block++){
c@5:         for (int band = 0; band < m_numberOfBands; band++){
c@7:             int k = m_bandBoundaries[band];
c@7:             int bandEnd = m_bandBoundaries[band+1];
c@7:             m_bandSum[band] = 0;
c@5:             
c@7:             while(k < bandEnd){
c@7:                 m_bandSum[band] += spectrogram[k][block];
c@7:                 k++;
c@5:             }
c@5:         }
c@5:         float total = 0;
c@5:         for(int band = 0; band < m_numberOfBands; band++){
c@7:             total += m_bandSum[band];
c@5:         }
c@7:         data[block] = total/m_numberOfBands;
c@5:     }
c@5:     
c@9:     //subtract local averages
c@7:     subtractLocalAverage(data);
c@5:     
c@5:     return data;
c@7: }