This is a stripped down version of the Queen Mary DSP library
-Version 1.7.1 from Sept 2015.
-https://code.soundsoftware.ac.uk/attachments/download/1582/qm-dsp-1.7.1.tar.gz
+https://github.com/c4dm/qm-dsp -- see gitrev.txt for version
---
}
const double *getWindow() const {
- return &m_window[0];
+ return m_window.data();
}
void cut(double *src) const {
}
const double *getWindow() const {
- return &m_window[0];
+ return m_window.data();
}
void cut(double *src) const {
m_normalise = Config.normalise; // if frame normalisation is required
// No. of constant Q bins
- m_uK = ( unsigned int ) ceil( m_BPO * log(m_FMax/m_FMin)/log(2.0));
+ m_uK = (int) ceil( m_BPO * log(m_FMax/m_FMin)/log(2.0));
// Create array for chroma result
m_chromadata = new double[ m_BPO ];
MathUtilities::getFrameMinMax( src, m_BPO, & min, &max );
- for( unsigned int i = 0; i < m_BPO; i++ )
+ for (int i = 0; i < m_BPO; i++)
{
val = src[ i ] / max;
}
// initialise chromadata to 0
- for (unsigned i = 0; i < m_BPO; i++) m_chromadata[i] = 0;
+ for (int i = 0; i < m_BPO; i++) m_chromadata[i] = 0;
- double cmax = 0.0;
- double cval = 0;
// Calculate ConstantQ frame
m_ConstantQ->process( real, imag, m_CQRe, m_CQIm );
// add each octave of cq data into Chromagram
- const unsigned octaves = (int)floor(double( m_uK/m_BPO))-1;
- for (unsigned octave = 0; octave <= octaves; octave++)
+ const int octaves = (int)floor(double( m_uK/m_BPO))-1;
+ for (int octave = 0; octave <= octaves; octave++)
{
- unsigned firstBin = octave*m_BPO;
- for (unsigned i = 0; i < m_BPO; i++)
+ int firstBin = octave*m_BPO;
+ for (int i = 0; i < m_BPO; i++)
{
m_chromadata[i] += kabs( m_CQRe[ firstBin + i ], m_CQIm[ firstBin + i ]);
}
#include "ConstantQ.h"
struct ChromaConfig{
- unsigned int FS;
+ int FS;
double min;
double max;
- unsigned int BPO;
+ int BPO;
double CQThresh;
MathUtilities::NormaliseType normalise;
};
double kabs( double real, double imag );
// Results
- unsigned int getK() { return m_uK;}
- unsigned int getFrameSize() { return m_frameSize; }
- unsigned int getHopSize() { return m_hopSize; }
-
+ int getK() { return m_uK;}
+ int getFrameSize() { return m_frameSize; }
+ int getHopSize() { return m_hopSize; }
+
private:
int initialise( ChromaConfig Config );
int deInitialise();
double* m_chromadata;
double m_FMin;
double m_FMax;
- unsigned int m_BPO;
- unsigned int m_uK;
+ int m_BPO;
+ int m_uK;
MathUtilities::NormaliseType m_normalise;
- unsigned int m_frameSize;
- unsigned int m_hopSize;
+ int m_frameSize;
+ int m_hopSize;
FFTReal* m_FFT;
ConstantQ* m_ConstantQ;
{
int filteridx = 0;
int factorDone = 1;
- int factorRemaining = m_decFactor;
while (factorDone < m_decFactor) {
filtSrc = new double[ m_length ];
filtDst = new double[ m_length ];
-
- //Low Pass Smoothing Filter Config
- m_FilterConfigParams.ord = Config.LPOrd;
- m_FilterConfigParams.ACoeffs = Config.LPACoeffs;
- m_FilterConfigParams.BCoeffs = Config.LPBCoeffs;
-
- m_FiltFilt = new FiltFilt( m_FilterConfigParams );
+ Filter::Parameters params;
+ params.a = std::vector<double>(Config.LPACoeffs, Config.LPACoeffs + Config.LPOrd + 1);
+ params.b = std::vector<double>(Config.LPBCoeffs, Config.LPBCoeffs + Config.LPOrd + 1);
+
+ m_FiltFilt = new FiltFilt(params);
//add delta threshold
m_delta = Config.delta;
MathUtilities::getAlphaNorm( src, m_length, m_alphaNormParam, &DFAlphaNorm );
- for( unsigned int i = 0; i< m_length; i++)
+ for (int i = 0; i < m_length; i++)
{
dst[ i ] = ( src[ i ] - DFMin ) / DFAlphaNorm;
}
double* m_filtScratchIn;
double* m_filtScratchOut;
- FilterConfig m_FilterConfigParams;
-
FiltFilt* m_FiltFilt;
bool m_isMedianPositive;
COPYING included with this distribution for more information.
*/
-#include <stdio.h>
#include "FiltFilt.h"
//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////
-FiltFilt::FiltFilt( FilterConfig Config )
+FiltFilt::FiltFilt(Filter::Parameters parameters) :
+ m_filter(parameters)
{
- m_filtScratchIn = NULL;
- m_filtScratchOut = NULL;
- m_ord = 0;
-
- initialise( Config );
+ m_ord = m_filter.getOrder();
}
FiltFilt::~FiltFilt()
{
- deInitialise();
-}
-
-void FiltFilt::initialise( FilterConfig Config )
-{
- m_ord = Config.ord;
- m_filterConfig.ord = Config.ord;
- m_filterConfig.ACoeffs = Config.ACoeffs;
- m_filterConfig.BCoeffs = Config.BCoeffs;
-
- m_filter = new Filter( m_filterConfig );
-}
-
-void FiltFilt::deInitialise()
-{
- delete m_filter;
}
-
void FiltFilt::process(double *src, double *dst, unsigned int length)
{
unsigned int i;
if (length == 0) return;
if (length < 2) {
- fprintf (stderr, "FiltFilt::process called for %d samples\n", length);
for( i = 0; i < length; i++ ) {
dst[i] = src [i];
}
unsigned int nFact = 3 * ( nFilt - 1);
unsigned int nExt = length + 2 * nFact;
- m_filtScratchIn = new double[ nExt ];
- m_filtScratchOut = new double[ nExt ];
-
+ double *filtScratchIn = new double[ nExt ];
+ double *filtScratchOut = new double[ nExt ];
for( i = 0; i< nExt; i++ )
{
- m_filtScratchIn[ i ] = 0.0;
- m_filtScratchOut[ i ] = 0.0;
+ filtScratchIn[ i ] = 0.0;
+ filtScratchOut[ i ] = 0.0;
}
// Edge transients reflection
unsigned int index = 0;
for( i = nFact; i > 0; i-- )
{
- m_filtScratchIn[ index++ ] = sample0 - src[ i ];
+ filtScratchIn[ index++ ] = sample0 - src[ i ];
}
index = 0;
for( i = 0; i < nFact && i + 2 < length; i++ )
{
- m_filtScratchIn[ (nExt - nFact) + index++ ] = sampleN - src[ (length - 2) - i ];
+ filtScratchIn[ (nExt - nFact) + index++ ] = sampleN - src[ (length - 2) - i ];
}
for(; i < nFact; i++ )
{
- m_filtScratchIn[ (nExt - nFact) + index++ ] = 0;
+ filtScratchIn[ (nExt - nFact) + index++ ] = 0;
}
index = 0;
for( i = 0; i < length; i++ )
{
- m_filtScratchIn[ i + nFact ] = src[ i ];
+ filtScratchIn[ i + nFact ] = src[ i ];
}
////////////////////////////////
// Do 0Ph filtering
- m_filter->process( m_filtScratchIn, m_filtScratchOut, nExt);
+ m_filter.process( filtScratchIn, filtScratchOut, nExt);
// reverse the series for FILTFILT
for ( i = 0; i < nExt; i++)
{
- m_filtScratchIn[ i ] = m_filtScratchOut[ nExt - i - 1];
+ filtScratchIn[ i ] = filtScratchOut[ nExt - i - 1];
}
// do FILTER again
- m_filter->process( m_filtScratchIn, m_filtScratchOut, nExt);
+ m_filter.process( filtScratchIn, filtScratchOut, nExt);
// reverse the series back
for ( i = 0; i < nExt; i++)
{
- m_filtScratchIn[ i ] = m_filtScratchOut[ nExt - i - 1 ];
+ filtScratchIn[ i ] = filtScratchOut[ nExt - i - 1 ];
}
for ( i = 0;i < nExt; i++)
{
- m_filtScratchOut[ i ] = m_filtScratchIn[ i ];
+ filtScratchOut[ i ] = filtScratchIn[ i ];
}
index = 0;
for( i = 0; i < length; i++ )
{
- dst[ index++ ] = m_filtScratchOut[ i + nFact ];
+ dst[ index++ ] = filtScratchOut[ i + nFact ];
}
- delete [] m_filtScratchIn;
- delete [] m_filtScratchOut;
+ delete [] filtScratchIn;
+ delete [] filtScratchOut;
}
/**
* Zero-phase digital filter, implemented by processing the data
- * through a filter specified by the given FilterConfig structure (see
+ * through a filter specified by the given filter parameters (see
* Filter) and then processing it again in reverse.
*/
class FiltFilt
{
public:
- FiltFilt( FilterConfig Config );
+ FiltFilt(Filter::Parameters);
virtual ~FiltFilt();
void reset();
void process( double* src, double* dst, unsigned int length );
private:
- void initialise( FilterConfig Config );
- void deInitialise();
-
- unsigned int m_ord;
-
- Filter* m_filter;
-
- double* m_filtScratchIn;
- double* m_filtScratchOut;
-
- FilterConfig m_filterConfig;
+ Filter m_filter;
+ int m_ord;
};
#endif
QM DSP Library
Centre for Digital Music, Queen Mary, University of London.
- This file 2005-2006 Christian Landone.
This program is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public License as
#include "Filter.h"
-//////////////////////////////////////////////////////////////////////
-// Construction/Destruction
-//////////////////////////////////////////////////////////////////////
+#include <stdexcept>
-Filter::Filter( FilterConfig Config )
-{
- m_ord = 0;
- m_outBuffer = NULL;
- m_inBuffer = NULL;
+using namespace std;
- initialise( Config );
+Filter::Filter(Parameters params)
+{
+ if (params.a.empty()) {
+ m_fir = true;
+ if (params.b.empty()) {
+ throw logic_error("Filter must have at least one pair of coefficients");
+ }
+ } else {
+ m_fir = false;
+ if (params.a.size() != params.b.size()) {
+ throw logic_error("Inconsistent numbers of filter coefficients");
+ }
+ }
+
+ m_sz = int(params.b.size());
+ m_order = m_sz - 1;
+
+ m_a = params.a;
+ m_b = params.b;
+
+ // We keep some empty space at the start of the buffer, and
+ // encroach gradually into it as we add individual sample
+ // calculations at the start. Then when we run out of space, we
+ // move the buffer back to the end and begin again. This is
+ // significantly faster than moving the whole buffer along in
+ // 1-sample steps every time.
+
+ m_offmax = 20;
+ m_offa = m_offmax;
+ m_offb = m_offmax;
+
+ if (!m_fir) {
+ m_bufa.resize(m_order + m_offmax);
+ }
+
+ m_bufb.resize(m_sz + m_offmax);
}
Filter::~Filter()
{
- deInitialise();
}
-void Filter::initialise( FilterConfig Config )
+void
+Filter::reset()
{
- m_ord = Config.ord;
- m_ACoeffs = Config.ACoeffs;
- m_BCoeffs = Config.BCoeffs;
+ m_offb = m_offmax;
+ m_offa = m_offmax;
- m_inBuffer = new double[ m_ord + 1 ];
- m_outBuffer = new double[ m_ord + 1 ];
-
- reset();
-}
+ if (!m_fir) {
+ m_bufa.assign(m_bufa.size(), 0.0);
+ }
-void Filter::deInitialise()
-{
- delete[] m_inBuffer;
- delete[] m_outBuffer;
+ m_bufb.assign(m_bufb.size(), 0.0);
}
-void Filter::reset()
+void
+Filter::process(const double *const __restrict__ in,
+ double *const __restrict__ out,
+ const int n)
{
- for( unsigned int i = 0; i < m_ord+1; i++ ){ m_inBuffer[ i ] = 0.0; }
- for(unsigned int i = 0; i < m_ord+1; i++ ){ m_outBuffer[ i ] = 0.0; }
+ for (int s = 0; s < n; ++s) {
+
+ if (m_offb > 0) --m_offb;
+ else {
+ for (int i = m_sz - 2; i >= 0; --i) {
+ m_bufb[i + m_offmax + 1] = m_bufb[i];
+ }
+ m_offb = m_offmax;
+ }
+ m_bufb[m_offb] = in[s];
+
+ double b_sum = 0.0;
+ for (int i = 0; i < m_sz; ++i) {
+ b_sum += m_b[i] * m_bufb[i + m_offb];
+ }
+
+ double outval;
+
+ if (m_fir) {
+
+ outval = b_sum;
+
+ } else {
+
+ double a_sum = 0.0;
+ for (int i = 0; i < m_order; ++i) {
+ a_sum += m_a[i + 1] * m_bufa[i + m_offa];
+ }
+
+ outval = b_sum - a_sum;
+
+ if (m_offa > 0) --m_offa;
+ else {
+ for (int i = m_order - 2; i >= 0; --i) {
+ m_bufa[i + m_offmax + 1] = m_bufa[i];
+ }
+ m_offa = m_offmax;
+ }
+ m_bufa[m_offa] = outval;
+ }
+
+ out[s] = outval;
+ }
}
-void Filter::process( double *src, double *dst, unsigned int length )
-{
- unsigned int SP,i,j;
-
- double xin,xout;
-
- for (SP=0;SP<length;SP++)
- {
- xin=src[SP];
- /* move buffer */
- for ( i = 0; i < m_ord; i++) {m_inBuffer[ m_ord - i ]=m_inBuffer[ m_ord - i - 1 ];}
- m_inBuffer[0]=xin;
-
- xout=0.0;
- for (j=0;j< m_ord + 1; j++)
- xout = xout + m_BCoeffs[ j ] * m_inBuffer[ j ];
- for (j = 0; j < m_ord; j++)
- xout= xout - m_ACoeffs[ j + 1 ] * m_outBuffer[ j ];
-
- dst[ SP ] = xout;
- for ( i = 0; i < m_ord - 1; i++ ) { m_outBuffer[ m_ord - i - 1 ] = m_outBuffer[ m_ord - i - 2 ];}
- m_outBuffer[0]=xout;
-
- } /* end of SP loop */
-}
-
-
-
QM DSP Library
Centre for Digital Music, Queen Mary, University of London.
- This file 2005-2006 Christian Landone.
This program is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public License as
#ifndef FILTER_H
#define FILTER_H
-#ifndef NULL
-#define NULL 0
-#endif
-
-/**
- * Filter specification. For a filter of order ord, the ACoeffs and
- * BCoeffs arrays must point to ord+1 values each. ACoeffs provides
- * the denominator and BCoeffs the numerator coefficients of the
- * filter.
- */
-struct FilterConfig{
- unsigned int ord;
- double* ACoeffs;
- double* BCoeffs;
-};
+#include <vector>
-/**
- * Digital filter specified through FilterConfig structure.
- */
-class Filter
+class Filter
{
public:
- Filter( FilterConfig Config );
- virtual ~Filter();
+ struct Parameters {
+ std::vector<double> a;
+ std::vector<double> b;
+ };
+
+ /**
+ * Construct an IIR filter with numerators b and denominators
+ * a. The filter will have order b.size()-1. To make an FIR
+ * filter, leave the vector a in the param struct empty.
+ * Otherwise, a and b must have the same number of values.
+ */
+ Filter(Parameters params);
+
+ ~Filter();
void reset();
- void process( double *src, double *dst, unsigned int length );
+ /**
+ * Filter the input sequence \arg in of length \arg n samples, and
+ * write the resulting \arg n samples into \arg out. There must be
+ * enough room in \arg out for \arg n samples to be written.
+ */
+ void process(const double *const __restrict__ in,
+ double *const __restrict__ out,
+ const int n);
+ int getOrder() const { return m_order; }
+
private:
- void initialise( FilterConfig Config );
- void deInitialise();
+ int m_order;
+ int m_sz;
+ std::vector<double> m_a;
+ std::vector<double> m_b;
+ std::vector<double> m_bufa;
+ std::vector<double> m_bufb;
+ int m_offa;
+ int m_offb;
+ int m_offmax;
+ bool m_fir;
- unsigned int m_ord;
-
- double* m_inBuffer;
- double* m_outBuffer;
-
- double* m_ACoeffs;
- double* m_BCoeffs;
+ Filter(const Filter &); // not supplied
+ Filter &operator=(const Filter &); // not supplied
};
-
+
#endif
m_tempoScratch = new double[ m_lagLength ];
m_smoothRCF = new double[ m_lagLength ];
-
- unsigned int winPre = Params.WinT.pre;
- unsigned int winPost = Params.WinT.post;
-
m_DFFramer.configure( m_winLength, m_lagLength );
m_DFPParams.length = m_winLength;
}
-void TempoTrack::createCombFilter(double* Filter, unsigned int winLength, unsigned int TSig, double beatLag)
+void TempoTrack::createCombFilter(double* Filter, int winLength, int /* TSig */, double beatLag)
{
- unsigned int i;
+ int i;
if( beatLag == 0 )
{
double period = 0;
double maxValRCF = 0.0;
- unsigned int maxIndexRCF = 0;
+ int maxIndexRCF = 0;
double* pdPeaks;
- unsigned int maxIndexTemp;
- double maxValTemp;
- unsigned int count;
+ int maxIndexTemp;
+ double maxValTemp;
+ int count;
- unsigned int numelem,i,j;
+ int numelem,i,j;
int a, b;
for( i = 0; i < m_lagLength; i++ )
}
}
-int TempoTrack::findMeter(double *ACF, unsigned int len, double period)
+int TempoTrack::findMeter(double *ACF, int len, double period)
{
int i;
int p = (int)MathUtilities::round( period );
double temp4B = 0.0;
double* dbf = new double[ len ]; int t = 0;
- for( unsigned int u = 0; u < len; u++ ){ dbf[ u ] = 0.0; }
+ for( int u = 0; u < len; u++ ){ dbf[ u ] = 0.0; }
if( (double)len < 6 * p + 2 )
{
return tsig;
}
-void TempoTrack::createPhaseExtractor(double *Filter, unsigned int winLength, double period, unsigned int fsp, unsigned int lastBeat)
+void TempoTrack::createPhaseExtractor(double *Filter, int /* winLength */, double period, int fsp, int lastBeat)
{
int p = (int)MathUtilities::round( period );
int predictedOffset = 0;
double sigma = (double)p/8;
double PhaseMin = 0.0;
double PhaseMax = 0.0;
- unsigned int scratchLength = p*2;
+ int scratchLength = p*2;
double temp = 0.0;
for( int i = 0; i < scratchLength; i++ )
std::cerr << "predictedOffset = " << predictedOffset << std::endl;
#endif
- unsigned int index = 0;
+ int index = 0;
for (int i = p - ( predictedOffset - 1); i < p + ( p - predictedOffset) + 1; i++)
{
#ifdef DEBUG_TEMPO_TRACK
delete [] phaseScratch;
}
-int TempoTrack::phaseMM(double *DF, double *weighting, unsigned int winLength, double period)
+int TempoTrack::phaseMM(double *DF, double *weighting, int winLength, double period)
{
int alignment = 0;
int p = (int)MathUtilities::round( period );
return alignment;
}
-int TempoTrack::beatPredict(unsigned int FSP0, double alignment, double period, unsigned int step )
+int TempoTrack::beatPredict(int FSP0, double alignment, double period, int step )
{
int beat = 0;
causalDF = DF;
//Prepare Causal Extension DFData
- unsigned int DFCLength = m_dataLength + m_winLength;
+// int DFCLength = m_dataLength + m_winLength;
- for( unsigned int j = 0; j < m_winLength; j++ )
+ for( int j = 0; j < m_winLength; j++ )
{
causalDF.push_back( 0 );
}
double* RW = new double[ m_lagLength ];
- for( unsigned int clear = 0; clear < m_lagLength; clear++){ RW[ clear ] = 0.0;}
+ for (int clear = 0; clear < m_lagLength; clear++){ RW[ clear ] = 0.0;}
double* GW = new double[ m_lagLength ];
- for(unsigned int clear = 0; clear < m_lagLength; clear++){ GW[ clear ] = 0.0;}
+ for (int clear = 0; clear < m_lagLength; clear++){ GW[ clear ] = 0.0;}
double* PW = new double[ m_lagLength ];
- for(unsigned clear = 0; clear < m_lagLength; clear++){ PW[ clear ] = 0.0;}
+ for(int clear = 0; clear < m_lagLength; clear++){ PW[ clear ] = 0.0;}
m_DFFramer.setSource( &causalDF[0], m_dataLength );
- unsigned int TTFrames = m_DFFramer.getMaxNoFrames();
+ int TTFrames = m_DFFramer.getMaxNoFrames();
#ifdef DEBUG_TEMPO_TRACK
std::cerr << "TTFrames = " << TTFrames << std::endl;
#endif
double* periodP = new double[ TTFrames ];
- for(unsigned clear = 0; clear < TTFrames; clear++){ periodP[ clear ] = 0.0;}
+ for(int clear = 0; clear < TTFrames; clear++){ periodP[ clear ] = 0.0;}
double* periodG = new double[ TTFrames ];
- for(unsigned clear = 0; clear < TTFrames; clear++){ periodG[ clear ] = 0.0;}
+ for(int clear = 0; clear < TTFrames; clear++){ periodG[ clear ] = 0.0;}
double* alignment = new double[ TTFrames ];
- for(unsigned clear = 0; clear < TTFrames; clear++){ alignment[ clear ] = 0.0;}
+ for(int clear = 0; clear < TTFrames; clear++){ alignment[ clear ] = 0.0;}
m_beats.clear();
int TTLoopIndex = 0;
- for( unsigned int i = 0; i < TTFrames; i++ )
+ for( int i = 0; i < TTFrames; i++ )
{
m_DFFramer.getFrame( m_rawDFFrame );
\r
struct WinThresh\r
{\r
- unsigned int pre;\r
- unsigned int post;\r
+ int pre;\r
+ int post;\r
};\r
\r
struct TTParams\r
{\r
- unsigned int winLength; //Analysis window length\r
- unsigned int lagLength; //Lag & Stride size\r
- unsigned int alpha; //alpha-norm parameter\r
- unsigned int LPOrd; // low-pass Filter order\r
+ int winLength; //Analysis window length\r
+ int lagLength; //Lag & Stride size\r
+ int alpha; //alpha-norm parameter\r
+ int LPOrd; // low-pass Filter order\r
double* LPACoeffs; //low pass Filter den coefficients\r
double* LPBCoeffs; //low pass Filter num coefficients\r
WinThresh WinT;//window size in frames for adaptive thresholding [pre post]:\r
void initialise( TTParams Params );\r
void deInitialise();\r
\r
- int beatPredict( unsigned int FSP, double alignment, double period, unsigned int step);\r
- int phaseMM( double* DF, double* weighting, unsigned int winLength, double period );\r
- void createPhaseExtractor( double* Filter, unsigned int winLength, double period, unsigned int fsp, unsigned int lastBeat );\r
- int findMeter( double* ACF, unsigned int len, double period );\r
+ int beatPredict( int FSP, double alignment, double period, int step);\r
+ int phaseMM( double* DF, double* weighting, int winLength, double period );\r
+ void createPhaseExtractor( double* Filter, int winLength, double period, int fsp, int lastBeat );\r
+ int findMeter( double* ACF, int len, double period );\r
void constDetect( double* periodP, int currentIdx, int* flag );\r
void stepDetect( double* periodP, double* periodG, int currentIdx, int* flag );\r
- void createCombFilter( double* Filter, unsigned int winLength, unsigned int TSig, double beatLag );\r
+ void createCombFilter( double* Filter, int winLength, int TSig, double beatLag );\r
double tempoMM( double* ACF, double* weight, int sig );\r
\r
- unsigned int m_dataLength;\r
- unsigned int m_winLength;\r
- unsigned int m_lagLength;\r
+ int m_dataLength;\r
+ int m_winLength;\r
+ int m_lagLength;\r
\r
- double m_rayparam;\r
- double m_sigma;\r
- double m_DFWVNnorm;\r
+ double m_rayparam;\r
+ double m_sigma;\r
+ double m_DFWVNnorm;\r
\r
vector<int> m_beats; // Vector of detected beats\r
\r
{
if (iPosition < 0)
rTCSVector = TCSVector();
- else if (iPosition >= m_VectorList.size())
+ else if (iPosition >= int(m_VectorList.size()))
rTCSVector = TCSVector();
else
rTCSVector = m_VectorList[iPosition].second;
void printDebug()
{
- for (int i = 0; i < size(); i++)
+ for (int i = 0; i < int(size()); i++)
{
std::cout << (*this)[i] << ";";
}
void printDebug()
{
- for (int i = 0; i < size(); i++)
+ for (int i = 0; i < int(size()); i++)
{
std::cout << (*this)[i] << ";";
}
QM DSP Library
Centre for Digital Music, Queen Mary, University of London.
+
+ 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. See the file
+ COPYING included with this distribution for more information.
*/
#ifndef FFT_H
break;
}
- assert(flength == lpd.size());
- assert(flength == hpd.size());
+ // avoid compiler warning for unused value if assert is not compiled in:
+ (void)flength;
+
+ assert(flength == int(lpd.size()));
+ assert(flength == int(hpd.size()));
}
if (m==1) {
do{
- *Fout = *f;
+ Fout->r = f->r;
+ Fout->i = f->i;
f += fstride*in_stride;
}while(++Fout != Fout_end );
}else{
--- /dev/null
+qm-vamp-plugins-v1.7.1-20-g4d15479
}
else
{
- for(int i=0; i<v1.size(); i++)
+ for(int i=0; i<int(v1.size()); i++)
{
dSum1 += v1[i]*v2[i];
dDen1 += v1[i]*v1[i];
#include <vector>
#include <cmath>
+using namespace std;
double MathUtilities::mod(double x, double y)
{
return ValOut;
}
-void MathUtilities::getAlphaNorm(const double *data, unsigned int len, unsigned int alpha, double* ANorm)
+void MathUtilities::getAlphaNorm(const double *data, int len, int alpha, double* ANorm)
{
- unsigned int i;
+ int i;
double temp = 0.0;
double a=0.0;
*ANorm = a;
}
-double MathUtilities::getAlphaNorm( const std::vector <double> &data, unsigned int alpha )
+double MathUtilities::getAlphaNorm( const vector <double> &data, int alpha )
{
- unsigned int i;
- unsigned int len = data.size();
+ int i;
+ int len = data.size();
double temp = 0.0;
double a=0.0;
for( i = 0; i < len; i++)
{
temp = data[ i ];
-
a += ::pow( fabs(temp), double(alpha) );
}
a /= ( double )len;
}
}
-double MathUtilities::median(const double *src, unsigned int len)
+double MathUtilities::median(const double *src, int len)
{
if (len == 0) return 0;
- std::vector<double> scratch;
+ vector<double> scratch;
for (int i = 0; i < len; ++i) scratch.push_back(src[i]);
- std::sort(scratch.begin(), scratch.end());
+ sort(scratch.begin(), scratch.end());
int middle = len/2;
if (len % 2 == 0) {
}
}
-double MathUtilities::sum(const double *src, unsigned int len)
+double MathUtilities::sum(const double *src, int len)
{
- unsigned int i ;
+ int i ;
double retVal =0.0;
for( i = 0; i < len; i++)
return retVal;
}
-double MathUtilities::mean(const double *src, unsigned int len)
+double MathUtilities::mean(const double *src, int len)
{
double retVal =0.0;
return retVal;
}
-double MathUtilities::mean(const std::vector<double> &src,
- unsigned int start,
- unsigned int count)
+double MathUtilities::mean(const vector<double> &src,
+ int start,
+ int count)
{
double sum = 0.;
return sum / count;
}
-void MathUtilities::getFrameMinMax(const double *data, unsigned int len, double *min, double *max)
+void MathUtilities::getFrameMinMax(const double *data, int len, double *min, double *max)
{
- unsigned int i;
+ int i;
double temp = 0.0;
if (len == 0) {
}
}
-int MathUtilities::getMax( double* pData, unsigned int Length, double* pMax )
+int MathUtilities::getMax( double* pData, int Length, double* pMax )
{
- unsigned int index = 0;
- unsigned int i;
+ int index = 0;
+ int i;
double temp = 0.0;
double max = pData[0];
return index;
}
-int MathUtilities::getMax( const std::vector<double> & data, double* pMax )
+int MathUtilities::getMax( const vector<double> & data, double* pMax )
{
- unsigned int index = 0;
- unsigned int i;
+ int index = 0;
+ int i;
double temp = 0.0;
double max = data[0];
- for( i = 0; i < data.size(); i++)
+ for( i = 0; i < int(data.size()); i++)
{
temp = data[ i ];
}
}
-void MathUtilities::normalise(std::vector<double> &data, NormaliseType type)
+void MathUtilities::normalise(vector<double> &data, NormaliseType type)
{
switch (type) {
}
}
-void MathUtilities::adaptiveThreshold(std::vector<double> &data)
+double MathUtilities::getLpNorm(const vector<double> &data, int p)
+{
+ double tot = 0.0;
+ for (int i = 0; i < int(data.size()); ++i) {
+ tot += abs(pow(data[i], p));
+ }
+ return pow(tot, 1.0 / p);
+}
+
+vector<double> MathUtilities::normaliseLp(const vector<double> &data,
+ int p,
+ double threshold)
+{
+ int n = int(data.size());
+ if (n == 0 || p == 0) return data;
+ double norm = getLpNorm(data, p);
+ if (norm < threshold) {
+ return vector<double>(n, 1.0 / pow(n, 1.0 / p)); // unit vector
+ }
+ vector<double> out(n);
+ for (int i = 0; i < n; ++i) {
+ out[i] = data[i] / norm;
+ }
+ return out;
+}
+
+void MathUtilities::adaptiveThreshold(vector<double> &data)
{
int sz = int(data.size());
if (sz == 0) return;
- std::vector<double> smoothed(sz);
+ vector<double> smoothed(sz);
int p_pre = 8;
int p_post = 7;
for (int i = 0; i < sz; ++i) {
- int first = std::max(0, i - p_pre);
- int last = std::min(sz - 1, i + p_post);
+ int first = max(0, i - p_pre);
+ int last = min(sz - 1, i + p_post);
smoothed[i] = mean(data, first, last - first + 1);
}
* Return through min and max pointers the highest and lowest
* values in the given array of the given length.
*/
- static void getFrameMinMax( const double* data, unsigned int len, double* min, double* max );
+ static void getFrameMinMax( const double* data, int len, double* min, double* max );
/**
* Return the mean of the given array of the given length.
*/
- static double mean( const double* src, unsigned int len );
+ static double mean( const double* src, int len );
/**
* Return the mean of the subset of the given vector identified by
* start and count.
*/
static double mean( const std::vector<double> &data,
- unsigned int start, unsigned int count );
+ int start, int count );
/**
* Return the sum of the values in the given array of the given
* length.
*/
- static double sum( const double* src, unsigned int len );
+ static double sum( const double* src, int len );
/**
* Return the median of the values in the given array of the given
* length. If the array is even in length, the returned value will
* be half-way between the two values adjacent to median.
*/
- static double median( const double* src, unsigned int len );
+ static double median( const double* src, int len );
/**
* The principle argument function. Map the phase angle ang into
*/
static double mod( double x, double y);
- static void getAlphaNorm(const double *data, unsigned int len, unsigned int alpha, double* ANorm);
- static double getAlphaNorm(const std::vector <double> &data, unsigned int alpha );
-
- static void circShift( double* data, int length, int shift);
+ /**
+ * The alpha norm is the alpha'th root of the mean alpha'th power
+ * magnitude. For example if alpha = 2 this corresponds to the RMS
+ * of the input data, and when alpha = 1 this is the mean
+ * magnitude.
+ */
+ static void getAlphaNorm(const double *data, int len, int alpha, double* ANorm);
- static int getMax( double* data, unsigned int length, double* max = 0 );
- static int getMax( const std::vector<double> &data, double* max = 0 );
- static int compareInt(const void * a, const void * b);
+ /**
+ * The alpha norm is the alpha'th root of the mean alpha'th power
+ * magnitude. For example if alpha = 2 this corresponds to the RMS
+ * of the input data, and when alpha = 1 this is the mean
+ * magnitude.
+ */
+ static double getAlphaNorm(const std::vector <double> &data, int alpha );
enum NormaliseType {
NormaliseNone,
static void normalise(std::vector<double> &data,
NormaliseType n = NormaliseUnitMax);
+ /**
+ * Calculate the L^p norm of a vector. Equivalent to MATLAB's
+ * norm(data, p).
+ */
+ static double getLpNorm(const std::vector<double> &data,
+ int p);
+
+ /**
+ * Normalise a vector by dividing through by its L^p norm. If the
+ * norm is below the given threshold, the unit vector for that
+ * norm is returned. p may be 0, in which case no normalisation
+ * happens and the data is returned unchanged.
+ */
+ static std::vector<double> normaliseLp(const std::vector<double> &data,
+ int p,
+ double threshold = 1e-6);
+
/**
* Threshold the input/output vector data against a moving-mean
* average filter.
*/
static void adaptiveThreshold(std::vector<double> &data);
+ static void circShift( double* data, int length, int shift);
+
+ static int getMax( double* data, int length, double* max = 0 );
+ static int getMax( const std::vector<double> &data, double* max = 0 );
+ static int compareInt(const void * a, const void * b);
+
/**
* Return true if x is 2^n for some integer n >= 0.
*/
// some utility functions
-namespace NSUtility
+struct NSUtility
{
- inline void swap(double &a, double &b) {double t = a; a = b; b = t;}
- void zeroise(vector<double> &array, int n);
- void zeroise(vector<int> &array, int n);
- void zeroise(vector<vector<double> > &matrix, int m, int n);
- void zeroise(vector<vector<int> > &matrix, int m, int n);
- inline double sqr(const double &x) {return x * x;}
+ static void swap(double &a, double &b) {double t = a; a = b; b = t;}
+ // fills a vector with zeros.
+ static void zeroise(vector<double> &array, int n) {
+ array.clear();
+ for(int j = 0; j < n; ++j) array.push_back(0);
+ }
+ // fills a vector with zeros.
+ static void zeroise(vector<int> &array, int n) {
+ array.clear();
+ for(int j = 0; j < n; ++j) array.push_back(0);
+ }
+ // fills a (m by n) matrix with zeros.
+ static void zeroise(vector<vector<double> > &matrix, int m, int n) {
+ vector<double> zero;
+ zeroise(zero, n);
+ matrix.clear();
+ for(int j = 0; j < m; ++j) matrix.push_back(zero);
+ }
+ // fills a (m by n) matrix with zeros.
+ static void zeroise(vector<vector<int> > &matrix, int m, int n) {
+ vector<int> zero;
+ zeroise(zero, n);
+ matrix.clear();
+ for(int j = 0; j < m; ++j) matrix.push_back(zero);
+ }
+ static double sqr(const double &x) {return x * x;}
};
-//---------------------------------------------------------------------------
-// Implementation
-//---------------------------------------------------------------------------
-using namespace NSUtility;
-//------------------------------------------------------------------------------------------
-
// main PolyFit routine
const int npoints(x.size());
const int nterms(coefs.size());
double correl_coef;
- zeroise(g, nterms);
- zeroise(a, nterms, nterms);
- zeroise(xmatr, npoints, nterms);
+ NSUtility::zeroise(g, nterms);
+ NSUtility::zeroise(a, nterms, nterms);
+ NSUtility::zeroise(xmatr, npoints, nterms);
if (nterms < 1) {
std::cerr << "ERROR: PolyFit called with less than one term" << std::endl;
return 0;
yc = 0.0;
for(j = 0; j < nterms; ++j)
yc += coefs [j] * xmatr [i][j];
- srs += sqr (yc - yi);
+ srs += NSUtility::sqr (yc - yi);
sum_y += yi;
sum_y2 += yi * yi;
}
// If all Y values are the same, avoid dividing by zero
- correl_coef = sum_y2 - sqr (sum_y) / npoints;
+ correl_coef = sum_y2 - NSUtility::sqr (sum_y) / npoints;
// Either return 0 or the correct value of correlation coefficient
if (correl_coef != 0)
correl_coef = srs / correl_coef;
vector<vector<int> >index;
Matrix w;
- zeroise(w, ncol, ncol);
- zeroise(index, ncol, 3);
+ NSUtility::zeroise(w, ncol, ncol);
+ NSUtility::zeroise(index, ncol, 3);
if(!GaussJordan2(b, y, w, index))
return false;
double big, t;
double pivot;
double determ;
- int irow, icol;
+ int irow = 0, icol = 0;
int ncol(b.size());
int nv = 1; // single constant vector
for(int i = 0; i < ncol; ++i)
} // { i-loop }
return true;
}
-//----------------------------------------------------------------------------------------------
-
-//------------------------------------------------------------------------------------
-
-// Utility functions
-//--------------------------------------------------------------------------
-
-// fills a vector with zeros.
-void NSUtility::zeroise(vector<double> &array, int n)
-{
- array.clear();
- for(int j = 0; j < n; ++j)
- array.push_back(0);
-}
-//--------------------------------------------------------------------------
-
-// fills a vector with zeros.
-void NSUtility::zeroise(vector<int> &array, int n)
-{
- array.clear();
- for(int j = 0; j < n; ++j)
- array.push_back(0);
-}
-//--------------------------------------------------------------------------
-
-// fills a (m by n) matrix with zeros.
-void NSUtility::zeroise(vector<vector<double> > &matrix, int m, int n)
-{
- vector<double> zero;
- zeroise(zero, n);
- matrix.clear();
- for(int j = 0; j < m; ++j)
- matrix.push_back(zero);
-}
-//--------------------------------------------------------------------------
-
-// fills a (m by n) matrix with zeros.
-void NSUtility::zeroise(vector<vector<int> > &matrix, int m, int n)
-{
- vector<int> zero;
- zeroise(zero, n);
- matrix.clear();
- for(int j = 0; j < m; ++j)
- matrix.push_back(zero);
-}
-//--------------------------------------------------------------------------
-
#endif
void pca_project(double** data, int n, int m, int ncomponents)
{
int i, j, k, k2;
- double **symmat, **symmat2, *evals, *interm;
+ double **symmat, /* **symmat2, */ *evals, *interm;
//TODO: assert ncomponents < m
projects / ardour.git / commitdiff