# HG changeset patch
# User Carl Bussey <c.bussey@se10.qmul.ac.uk>
# Date 1408546837 -3600
# Node ID fe23998968b4b85c6a2e108daf138d482f573c96
# Parent  957b83524c06e66419bbb8fcc4c228c1b718e91f
* Added tempogram via autocorrelation feature, using AutocorrelationProcessor
* Moved calculateMax() from NoveltyCurveProcessor to SpectrogramProcessor

diff -r 957b83524c06 -r fe23998968b4 AutocorrelationProcessor.cpp
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/AutocorrelationProcessor.cpp	Wed Aug 20 16:00:37 2014 +0100
@@ -0,0 +1,52 @@
+//
+//  AutocorrelationProcessor.cpp
+//  Tempogram
+//
+//  Created by Carl Bussey on 20/08/2014.
+//  Copyright (c) 2014 Carl Bussey. All rights reserved.
+//
+
+#include "AutocorrelationProcessor.h"
+using namespace std;
+
+AutocorrelationProcessor::AutocorrelationProcessor(const size_t &windowLength, const unsigned int &hopSize, const unsigned int &lagIncrement) :
+    m_windowLength(windowLength),
+    m_hopSize(hopSize),
+    m_lagIncrement(lagIncrement)
+{
+    //Nothing to do here
+}
+
+AutoCorrelation AutocorrelationProcessor::process(float * input, const size_t &inputLength) const
+{
+    int readBlockPointerIndex = 0;
+    AutoCorrelation autocorrelation;
+    
+    while(readBlockPointerIndex <= (int)inputLength) {
+        int readPointer = readBlockPointerIndex - m_windowLength/2;
+        
+        autocorrelation.push_back(processBlock(&input[readPointer], min(inputLength-readPointer, m_windowLength)));
+        readBlockPointerIndex += m_hopSize;
+    }
+    
+    return autocorrelation;
+}
+
+vector<float> AutocorrelationProcessor::processBlock(float * blockInput, const size_t &blockLength) const
+{
+    vector<float> autocorrelation;
+    
+    int N = m_windowLength/m_lagIncrement;
+    
+    for (int lag = 0; lag < N; lag++){
+        float sum = 0;
+        int sampleLag = m_lagIncrement*lag;
+        
+        for (int n = sampleLag; n < (int)blockLength; n++){
+            sum += blockInput[n-sampleLag]*blockInput[n];
+        }
+        autocorrelation.push_back(sum/(2*N + 1 - lag));
+    }
+    
+    return autocorrelation;
+}
\ No newline at end of file
diff -r 957b83524c06 -r fe23998968b4 AutocorrelationProcessor.h
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/AutocorrelationProcessor.h	Wed Aug 20 16:00:37 2014 +0100
@@ -0,0 +1,29 @@
+//
+//  AutocorrelationProcessor.h
+//  Tempogram
+//
+//  Created by Carl Bussey on 20/08/2014.
+//  Copyright (c) 2014 Carl Bussey. All rights reserved.
+//
+
+#ifndef __Tempogram__Autocorrelation__
+#define __Tempogram__Autocorrelation__
+
+#include <iostream>
+#include <vector>
+
+typedef std::vector< std::vector<float> > AutoCorrelation;
+
+class AutocorrelationProcessor{
+public:
+    AutocorrelationProcessor(const size_t &windowLength, const unsigned int &hopSize, const unsigned int &lagIncrement);
+    AutoCorrelation process(float * input, const size_t &inputLength) const;
+    std::vector<float> processBlock(float * input, const size_t &inputLength) const;
+private:
+    size_t m_windowLength;
+    unsigned int m_hopSize;
+    unsigned int m_lagIncrement;
+    
+};
+
+#endif /* defined(__Tempogram__Autocorrelation__) */
diff -r 957b83524c06 -r fe23998968b4 Makefile
--- a/Makefile	Tue Aug 19 18:20:26 2014 +0100
+++ b/Makefile	Wed Aug 20 16:00:37 2014 +0100
@@ -21,11 +21,11 @@
 
 # Edit this to list the .cpp or .c files in your plugin project
 #
-PLUGIN_SOURCES := TempogramPlugin.cpp FIRFilter.cpp WindowFunction.cpp plugins.cpp NoveltyCurveProcessor.cpp SpectrogramProcessor.cpp
+PLUGIN_SOURCES := TempogramPlugin.cpp FIRFilter.cpp WindowFunction.cpp plugins.cpp NoveltyCurveProcessor.cpp SpectrogramProcessor.cpp AutocorrelationProcessor.cpp
 
 # Edit this to list the .h files in your plugin project
 #
-PLUGIN_HEADERS := TempogramPlugin.h FIRFilter.h WindowFunction.h NoveltyCurveProcessor.h SpectrogramProcessor.h
+PLUGIN_HEADERS := TempogramPlugin.h FIRFilter.h WindowFunction.h NoveltyCurveProcessor.h SpectrogramProcessor.h AutocorrelationProcessor.h
 
 # Edit this to the location of the Vamp plugin SDK, relative to your
 # project directory
diff -r 957b83524c06 -r fe23998968b4 NoveltyCurveProcessor.cpp
--- a/NoveltyCurveProcessor.cpp	Tue Aug 19 18:20:26 2014 +0100
+++ b/NoveltyCurveProcessor.cpp	Wed Aug 20 16:00:37 2014 +0100
@@ -52,23 +52,6 @@
     m_pBandSum = 0;
 }
 
-//calculate max of spectrogram
-float NoveltyCurveProcessor::calculateMax(const Spectrogram &spectrogram) const
-{
-    float max = 0;
-    
-    int length = spectrogram.size();
-    int height = spectrogram[0].size();
-    
-    for (int i = 0; i < length; i++){
-        for (int j = 0; j < height; j++){
-            max = max > fabs(spectrogram[i][j]) ? max : fabs(spectrogram[i][j]);
-        }
-    }
-    
-    return max;
-}
-
 //subtract local average of novelty curve
 //uses m_hannWindow as filter
 void NoveltyCurveProcessor::subtractLocalAverage(vector<float> &noveltyCurve, const size_t &smoothLength) const
@@ -94,7 +77,6 @@
 //smoothed differentiator filter. Flips upper half of hanning window about y-axis to create coefficients.
 void NoveltyCurveProcessor::smoothedDifferentiator(SpectrogramTransposed &spectrogramTransposed, const size_t &smoothLength) const
 {
-    
     int numberOfBlocks = spectrogramTransposed[0].size();
     
     float * diffHannWindow = new float [smoothLength];
@@ -116,7 +98,7 @@
 void NoveltyCurveProcessor::halfWaveRectify(Spectrogram &spectrogram) const
 {
     int length = spectrogram.size();
-    int height = spectrogram[0].size();
+    int height = length > 0 ? spectrogram[0].size() : 0;
     
     for (int i = 0; i < length; i++){
         for (int j = 0; j < height; j++){
@@ -131,14 +113,15 @@
 {
     int numberOfBlocks = spectrogram.size();
     std::vector<float> noveltyCurve(numberOfBlocks);
-    SpectrogramTransposed spectrogramTransposed(spectrogram[0].size(), vector<float>(spectrogram.size()));
+    SpectrogramTransposed spectrogramTransposed(m_blockSize, vector<float>(spectrogram.size()));
     
     //normalise and log spectrogram
-    float normaliseScale = calculateMax(spectrogram);
+    float normaliseScale = SpectrogramProcessor::calculateMax(spectrogram);
     for (int block = 0; block < (int)numberOfBlocks; block++){
         for (int k = 0; k < (int)m_blockSize; k++){
-            spectrogramTransposed[k][block] = log(1+m_compressionConstant*spectrogram[block][k]);
-            if(normaliseScale != 0.0) spectrogramTransposed[k][block] /= normaliseScale; //normalise
+            float magnitude = spectrogram[block][k];
+            if(normaliseScale != 0.0) magnitude /= normaliseScale; //normalise
+            spectrogramTransposed[k][block] = log(1+m_compressionConstant*magnitude);
         }
     }
     
diff -r 957b83524c06 -r fe23998968b4 NoveltyCurveProcessor.h
--- a/NoveltyCurveProcessor.h	Tue Aug 19 18:20:26 2014 +0100
+++ b/NoveltyCurveProcessor.h	Wed Aug 20 16:00:37 2014 +0100
@@ -31,7 +31,6 @@
     
     void initialise();
     void cleanup();
-    float calculateMax(const Spectrogram &spectrogram) const;
     void subtractLocalAverage(std::vector<float> &noveltyCurve, const size_t &smoothLength) const;
     void smoothedDifferentiator(SpectrogramTransposed &spectrogram, const size_t &smoothLength) const;
     void halfWaveRectify(SpectrogramTransposed &spectrogram) const;
diff -r 957b83524c06 -r fe23998968b4 SpectrogramProcessor.cpp
--- a/SpectrogramProcessor.cpp	Tue Aug 19 18:20:26 2014 +0100
+++ b/SpectrogramProcessor.cpp	Wed Aug 20 16:00:37 2014 +0100
@@ -56,8 +56,25 @@
     return spectrogramT;
 }
 
+//calculate max of spectrogram
+float SpectrogramProcessor::calculateMax(const Spectrogram &spectrogram)
+{
+    float max = 0;
+    
+    int length = spectrogram.size();
+    int height = length > 0 ? spectrogram[0].size() : 0;
+    
+    for (int i = 0; i < length; i++){
+        for (int j = 0; j < height; j++){
+            max = max > fabs(spectrogram[i][j]) ? max : fabs(spectrogram[i][j]);
+        }
+    }
+    
+    return max;
+}
+
 //process method
-Spectrogram SpectrogramProcessor::process(const float * const pInput, const size_t &inputLength, const float * pWindow, const bool &transposeOutput) const
+Spectrogram SpectrogramProcessor::process(const float * const pInput, const size_t &inputLength, const float * pWindow) const
 {
     Spectrogram spectrogram;
     
@@ -96,6 +113,5 @@
         writeBlockPointer++;
     }
     
-    if(transposeOutput) return transpose(spectrogram);
-    else return spectrogram;
+    return spectrogram;
 }
diff -r 957b83524c06 -r fe23998968b4 SpectrogramProcessor.h
--- a/SpectrogramProcessor.h	Tue Aug 19 18:20:26 2014 +0100
+++ b/SpectrogramProcessor.h	Wed Aug 20 16:00:37 2014 +0100
@@ -31,8 +31,9 @@
     SpectrogramProcessor(const size_t &windowLength, const size_t &fftLength, const size_t &hopSize);
     ~SpectrogramProcessor();
     
-    Spectrogram process(const float * const pInput, const size_t &inputLength, const float * pWindow, const bool &transposeOutput = false) const;
+    Spectrogram process(const float * const pInput, const size_t &inputLength, const float * pWindow) const;
     static SpectrogramTransposed transpose(const Spectrogram &spectrogram);
+    static float calculateMax(const Spectrogram &spectrogram);
 };
 
 #endif /* defined(__Tempogram__Spectrogram__) */
diff -r 957b83524c06 -r fe23998968b4 TempogramPlugin.cpp
--- a/TempogramPlugin.cpp	Tue Aug 19 18:20:26 2014 +0100
+++ b/TempogramPlugin.cpp	Wed Aug 20 16:00:37 2014 +0100
@@ -5,8 +5,7 @@
 
 
 #include "TempogramPlugin.h"
-#include <sstream>
-#include <stdexcept>
+
 
 using Vamp::FFT;
 using Vamp::RealTime;
@@ -323,54 +322,73 @@
     float d_sampleRate;
     float tempogramInputSampleRate = (float)m_inputSampleRate/m_inputStepSize;
     
+    OutputDescriptor d1;
+    d1.identifier = "cyclicTempogram";
+    d1.name = "Cyclic Tempogram";
+    d1.description = "Cyclic Tempogram";
+    d1.unit = "";
+    d1.hasFixedBinCount = true;
+    d1.binCount = m_cyclicTempogramOctaveDivider > 0 && !isnan(m_cyclicTempogramOctaveDivider) ? m_cyclicTempogramOctaveDivider : 0;
+    d1.hasKnownExtents = false;
+    d1.isQuantized = false;
+    d1.sampleType = OutputDescriptor::FixedSampleRate;
+    d_sampleRate = tempogramInputSampleRate/m_tempogramHopSize;
+    d1.sampleRate = d_sampleRate > 0.0 && !isnan(d_sampleRate) ? d_sampleRate : 0;
+    d1.hasDuration = false;
+    list.push_back(d1);
+    
+    OutputDescriptor d2;
+    d2.identifier = "tempogramDFT";
+    d2.name = "Tempogram via DFT";
+    d2.description = "Tempogram via DFT";
+    d2.unit = "BPM";
+    d2.hasFixedBinCount = true;
+    d2.binCount = m_tempogramMaxBin - m_tempogramMinBin + 1;
+    d2.hasKnownExtents = false;
+    d2.isQuantized = false;
+    d2.sampleType = OutputDescriptor::FixedSampleRate;
+    d_sampleRate = tempogramInputSampleRate/m_tempogramHopSize;
+    d2.sampleRate = d_sampleRate > 0.0 && !isnan(d_sampleRate) ? d_sampleRate : 0.0;
+    for(int i = m_tempogramMinBin; i <= (int)m_tempogramMaxBin; i++){
+        float w = ((float)i/m_tempogramFftLength)*(tempogramInputSampleRate);
+        d2.binNames.push_back(floatToString(w*60));
+    }
+    d2.hasDuration = false;
+    list.push_back(d2);
+    
     OutputDescriptor d3;
-    d3.identifier = "cyclicTempogram";
-    d3.name = "Cyclic Tempogram";
-    d3.description = "Cyclic Tempogram";
-    d3.unit = "";
+    d3.identifier = "tempogramACT";
+    d3.name = "Tempogram via ACT";
+    d3.description = "Tempogram via ACT";
+    d3.unit = "BPM";
     d3.hasFixedBinCount = true;
-    d3.binCount = m_cyclicTempogramOctaveDivider > 0 && !isnan(m_cyclicTempogramOctaveDivider) ? m_cyclicTempogramOctaveDivider : 0;
+    d3.binCount = m_tempogramMaxBin - m_tempogramMinBin + 1;
     d3.hasKnownExtents = false;
     d3.isQuantized = false;
     d3.sampleType = OutputDescriptor::FixedSampleRate;
     d_sampleRate = tempogramInputSampleRate/m_tempogramHopSize;
-    d3.sampleRate = d_sampleRate > 0.0 && !isnan(d_sampleRate) ? d_sampleRate : 0;
+    d3.sampleRate = d_sampleRate > 0.0 && !isnan(d_sampleRate) ? d_sampleRate : 0.0;
+    for(int i = m_tempogramMinBin; i <= (int)m_tempogramMaxBin; i++){
+        float w = ((float)i/m_tempogramFftLength)*(tempogramInputSampleRate);
+        d3.binNames.push_back(floatToString(w*60));
+    }
     d3.hasDuration = false;
     list.push_back(d3);
     
-    OutputDescriptor d1;
-    d1.identifier = "tempogram";
-    d1.name = "Tempogram";
-    d1.description = "Tempogram";
-    d1.unit = "BPM";
-    d1.hasFixedBinCount = true;
-    d1.binCount = m_tempogramMaxBin - m_tempogramMinBin + 1;
-    d1.hasKnownExtents = false;
-    d1.isQuantized = false;
-    d1.sampleType = OutputDescriptor::FixedSampleRate;
-    d_sampleRate = tempogramInputSampleRate/m_tempogramHopSize;
-    d1.sampleRate = d_sampleRate > 0.0 && !isnan(d_sampleRate) ? d_sampleRate : 0.0;
-    for(int i = m_tempogramMinBin; i <= (int)m_tempogramMaxBin; i++){
-        float w = ((float)i/m_tempogramFftLength)*(tempogramInputSampleRate);
-        d1.binNames.push_back(floatToString(w*60));
-    }
-    d1.hasDuration = false;
-    list.push_back(d1);
-    
-    OutputDescriptor d2;
-    d2.identifier = "nc";
-    d2.name = "Novelty Curve";
-    d2.description = "Novelty Curve";
-    d2.unit = "";
-    d2.hasFixedBinCount = true;
-    d2.binCount = 1;
-    d2.hasKnownExtents = false;
-    d2.isQuantized = false;
-    d2.sampleType = OutputDescriptor::FixedSampleRate;
+    OutputDescriptor d4;
+    d4.identifier = "nc";
+    d4.name = "Novelty Curve";
+    d4.description = "Novelty Curve";
+    d4.unit = "";
+    d4.hasFixedBinCount = true;
+    d4.binCount = 1;
+    d4.hasKnownExtents = false;
+    d4.isQuantized = false;
+    d4.sampleType = OutputDescriptor::FixedSampleRate;
     d_sampleRate = tempogramInputSampleRate;
-    d2.sampleRate = d_sampleRate > 0 && !isnan(d_sampleRate) ? d_sampleRate : 0;
-    d2.hasDuration = false;
-    list.push_back(d2);
+    d4.sampleRate = d_sampleRate > 0 && !isnan(d_sampleRate) ? d_sampleRate : 0;
+    d4.hasDuration = false;
+    list.push_back(d4);
     
     return list;
 }
@@ -442,7 +460,7 @@
         Feature noveltyCurveFeature;
         noveltyCurveFeature.values.push_back(noveltyCurve[i]);
         noveltyCurveFeature.hasTimestamp = false;
-        featureSet[2].push_back(noveltyCurveFeature);
+        featureSet[3].push_back(noveltyCurveFeature);
         assert(!isnan(noveltyCurveFeature.values.back()));
     }
     
@@ -452,23 +470,30 @@
     //initialise spectrogram processor
     SpectrogramProcessor spectrogramProcessor(m_tempogramWindowLength, m_tempogramFftLength, m_tempogramHopSize);
     //compute spectrogram from novelty curve data (i.e., tempogram)
-    Tempogram tempogram = spectrogramProcessor.process(&noveltyCurve[0], numberOfBlocks, hannWindow);
+    Tempogram tempogramDFT = spectrogramProcessor.process(&noveltyCurve[0], numberOfBlocks, hannWindow);
     delete []hannWindow;
     hannWindow = 0;
     
-    int tempogramLength = tempogram.size();
+    int tempogramLag = 1;
+    AutocorrelationProcessor autocorrelationProcessor(m_tempogramWindowLength, m_tempogramHopSize, tempogramLag);
+    Tempogram tempogramACT = autocorrelationProcessor.process(&noveltyCurve[0], numberOfBlocks);
+    
+    int tempogramLength = tempogramDFT.size();
     
     //push tempogram data to featureset 0 and set timestamps.
     for (int block = 0; block < tempogramLength; block++){
-        Feature tempogramFeature;
+        Feature tempogramDFTFeature;
+        Feature tempogramACTFeature;
         
-        assert(tempogram[block].size() == (m_tempogramFftLength/2 + 1));
+        assert(tempogramDFT[block].size() == (m_tempogramFftLength/2 + 1));
         for(int k = m_tempogramMinBin; k < (int)m_tempogramMaxBin; k++){
-            tempogramFeature.values.push_back(tempogram[block][k]);
-            assert(!isnan(tempogramFeature.values.back()));
+            tempogramDFTFeature.values.push_back(tempogramDFT[block][k]);
+            tempogramACTFeature.values.push_back(tempogramACT[block][k]);
         }
-        tempogramFeature.hasTimestamp = false;
-        featureSet[1].push_back(tempogramFeature);
+        tempogramDFTFeature.hasTimestamp = false;
+        tempogramACTFeature.hasTimestamp = false;
+        featureSet[1].push_back(tempogramDFTFeature);
+        featureSet[2].push_back(tempogramACTFeature);
     }
     
     //Calculate cyclic tempogram
@@ -482,7 +507,7 @@
             float sum = 0;
             
             for (int j = 0; j < m_cyclicTempogramNumberOfOctaves; j++){
-                sum += tempogram[block][logBins[j][i]];
+                sum += tempogramDFT[block][logBins[j][i]];
             }
             cyclicTempogramFeature.values.push_back(sum/m_cyclicTempogramNumberOfOctaves);
             assert(!isnan(cyclicTempogramFeature.values.back()));
@@ -551,7 +576,7 @@
     m_tempogramMinBin = (max(floor(((m_tempogramMinBPM/60)/tempogramInputSampleRate)*m_tempogramFftLength), (float)0.0));
     m_tempogramMaxBin = (min(ceil(((m_tempogramMaxBPM/60)/tempogramInputSampleRate)*m_tempogramFftLength), (float)m_tempogramFftLength/2));
     
-    if (m_tempogramMinBPM > m_cyclicTempogramMinBPM) m_cyclicTempogramMinBPM = m_tempogramMinBPM;
+    if (m_tempogramMinBPM > m_cyclicTempogramMinBPM) m_cyclicTempogramMinBPM = m_tempogramMinBPM; //m_cyclicTempogram can't be less than default = 30
     float cyclicTempogramMaxBPM = 480;
     if (m_tempogramMaxBPM < cyclicTempogramMaxBPM) cyclicTempogramMaxBPM = m_tempogramMaxBPM;
     
diff -r 957b83524c06 -r fe23998968b4 TempogramPlugin.h
--- a/TempogramPlugin.h	Tue Aug 19 18:20:26 2014 +0100
+++ b/TempogramPlugin.h	Wed Aug 20 16:00:37 2014 +0100
@@ -24,12 +24,15 @@
 #include "WindowFunction.h"
 #include "NoveltyCurveProcessor.h"
 #include "SpectrogramProcessor.h"
+#include "AutocorrelationProcessor.h"
 #include <vamp-sdk/FFT.h>
 
 #include <cmath>
 #include <fstream>
 #include <cassert>
 #include <string>
+#include <sstream>
+#include <stdexcept>
 
 using std::string;
 using std::vector;