Update qm-dsp library (v1.7.1-20-g4d15479)
This commit is contained in:
parent
c05e6b2069
commit
b6768b4616
@ -1,7 +1,6 @@
|
||||
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
|
||||
|
||||
---
|
||||
|
||||
|
@ -81,7 +81,7 @@ public:
|
||||
}
|
||||
|
||||
const double *getWindow() const {
|
||||
return &m_window[0];
|
||||
return m_window.data();
|
||||
}
|
||||
|
||||
void cut(double *src) const {
|
||||
|
@ -37,7 +37,7 @@ public:
|
||||
}
|
||||
|
||||
const double *getWindow() const {
|
||||
return &m_window[0];
|
||||
return m_window.data();
|
||||
}
|
||||
|
||||
void cut(double *src) const {
|
||||
|
@ -34,7 +34,7 @@ int Chromagram::initialise( ChromaConfig Config )
|
||||
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 ];
|
||||
@ -112,7 +112,7 @@ void Chromagram::unityNormalise(double *src)
|
||||
|
||||
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;
|
||||
|
||||
@ -153,19 +153,17 @@ double* Chromagram::process( const double *real, const double *imag )
|
||||
}
|
||||
|
||||
// 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 ]);
|
||||
}
|
||||
|
@ -21,10 +21,10 @@
|
||||
#include "ConstantQ.h"
|
||||
|
||||
struct ChromaConfig{
|
||||
unsigned int FS;
|
||||
int FS;
|
||||
double min;
|
||||
double max;
|
||||
unsigned int BPO;
|
||||
int BPO;
|
||||
double CQThresh;
|
||||
MathUtilities::NormaliseType normalise;
|
||||
};
|
||||
@ -44,10 +44,10 @@ public:
|
||||
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();
|
||||
@ -58,13 +58,13 @@ private:
|
||||
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;
|
||||
|
@ -112,7 +112,6 @@ void DecimatorB::doProcess()
|
||||
{
|
||||
int filteridx = 0;
|
||||
int factorDone = 1;
|
||||
int factorRemaining = m_decFactor;
|
||||
|
||||
while (factorDone < m_decFactor) {
|
||||
|
||||
|
@ -59,13 +59,11 @@ void DFProcess::initialise( DFProcConfig Config )
|
||||
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;
|
||||
@ -193,7 +191,7 @@ void DFProcess::removeDCNormalize( double *src, double*dst )
|
||||
|
||||
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;
|
||||
}
|
||||
|
@ -81,8 +81,6 @@ private:
|
||||
double* m_filtScratchIn;
|
||||
double* m_filtScratchOut;
|
||||
|
||||
FilterConfig m_FilterConfigParams;
|
||||
|
||||
FiltFilt* m_FiltFilt;
|
||||
|
||||
bool m_isMedianPositive;
|
||||
|
@ -13,43 +13,22 @@
|
||||
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;
|
||||
@ -57,7 +36,6 @@ void FiltFilt::process(double *src, double *dst, unsigned int length)
|
||||
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];
|
||||
}
|
||||
@ -68,14 +46,13 @@ void FiltFilt::process(double *src, double *dst, unsigned int length)
|
||||
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
|
||||
@ -85,56 +62,56 @@ void FiltFilt::process(double *src, double *dst, unsigned int length)
|
||||
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;
|
||||
|
||||
}
|
||||
|
||||
|
@ -20,30 +20,21 @@
|
||||
|
||||
/**
|
||||
* 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
|
||||
|
@ -4,7 +4,6 @@
|
||||
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
|
||||
@ -15,73 +14,112 @@
|
||||
|
||||
#include "Filter.h"
|
||||
|
||||
//////////////////////////////////////////////////////////////////////
|
||||
// Construction/Destruction
|
||||
//////////////////////////////////////////////////////////////////////
|
||||
#include <stdexcept>
|
||||
|
||||
Filter::Filter( FilterConfig Config )
|
||||
using namespace std;
|
||||
|
||||
Filter::Filter(Parameters params)
|
||||
{
|
||||
m_ord = 0;
|
||||
m_outBuffer = NULL;
|
||||
m_inBuffer = NULL;
|
||||
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;
|
||||
|
||||
initialise( Config );
|
||||
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 ];
|
||||
if (!m_fir) {
|
||||
m_bufa.assign(m_bufa.size(), 0.0);
|
||||
}
|
||||
|
||||
reset();
|
||||
m_bufb.assign(m_bufb.size(), 0.0);
|
||||
}
|
||||
|
||||
void Filter::deInitialise()
|
||||
void
|
||||
Filter::process(const double *const __restrict__ in,
|
||||
double *const __restrict__ out,
|
||||
const int n)
|
||||
{
|
||||
delete[] m_inBuffer;
|
||||
delete[] m_outBuffer;
|
||||
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::reset()
|
||||
{
|
||||
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; }
|
||||
}
|
||||
|
||||
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 */
|
||||
}
|
||||
|
||||
|
||||
|
||||
|
@ -4,7 +4,6 @@
|
||||
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
|
||||
@ -16,46 +15,53 @@
|
||||
#ifndef FILTER_H
|
||||
#define FILTER_H
|
||||
|
||||
#ifndef NULL
|
||||
#define NULL 0
|
||||
#endif
|
||||
#include <vector>
|
||||
|
||||
/**
|
||||
* 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;
|
||||
};
|
||||
|
||||
/**
|
||||
* 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
|
||||
|
@ -70,10 +70,6 @@ void TempoTrack::initialise( TTParams Params )
|
||||
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;
|
||||
@ -120,9 +116,9 @@ void TempoTrack::deInitialise()
|
||||
|
||||
}
|
||||
|
||||
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 )
|
||||
{
|
||||
@ -147,15 +143,15 @@ double TempoTrack::tempoMM(double* ACF, double* weight, int tsig)
|
||||
|
||||
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++ )
|
||||
@ -476,7 +472,7 @@ void TempoTrack::constDetect( double* periodP, int currentIdx, int* flag )
|
||||
}
|
||||
}
|
||||
|
||||
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 );
|
||||
@ -491,7 +487,7 @@ int TempoTrack::findMeter(double *ACF, unsigned int len, double 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 )
|
||||
{
|
||||
@ -548,7 +544,7 @@ int TempoTrack::findMeter(double *ACF, unsigned int len, double period)
|
||||
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;
|
||||
@ -584,7 +580,7 @@ void TempoTrack::createPhaseExtractor(double *Filter, unsigned int winLength, do
|
||||
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++ )
|
||||
@ -604,7 +600,7 @@ void TempoTrack::createPhaseExtractor(double *Filter, unsigned int winLength, do
|
||||
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
|
||||
@ -624,7 +620,7 @@ void TempoTrack::createPhaseExtractor(double *Filter, unsigned int winLength, do
|
||||
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 );
|
||||
@ -667,7 +663,7 @@ int TempoTrack::phaseMM(double *DF, double *weighting, unsigned int winLength, d
|
||||
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;
|
||||
|
||||
@ -712,39 +708,39 @@ vector<int> TempoTrack::process( vector <double> DF,
|
||||
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();
|
||||
|
||||
@ -752,7 +748,7 @@ vector<int> TempoTrack::process( vector <double> DF,
|
||||
|
||||
int TTLoopIndex = 0;
|
||||
|
||||
for( unsigned int i = 0; i < TTFrames; i++ )
|
||||
for( int i = 0; i < TTFrames; i++ )
|
||||
{
|
||||
m_DFFramer.getFrame( m_rawDFFrame );
|
||||
|
||||
|
@ -30,16 +30,16 @@ using std::vector;
|
||||
|
||||
struct WinThresh
|
||||
{
|
||||
unsigned int pre;
|
||||
unsigned int post;
|
||||
int pre;
|
||||
int post;
|
||||
};
|
||||
|
||||
struct TTParams
|
||||
{
|
||||
unsigned int winLength; //Analysis window length
|
||||
unsigned int lagLength; //Lag & Stride size
|
||||
unsigned int alpha; //alpha-norm parameter
|
||||
unsigned int LPOrd; // low-pass Filter order
|
||||
int winLength; //Analysis window length
|
||||
int lagLength; //Lag & Stride size
|
||||
int alpha; //alpha-norm parameter
|
||||
int LPOrd; // low-pass Filter order
|
||||
double* LPACoeffs; //low pass Filter den coefficients
|
||||
double* LPBCoeffs; //low pass Filter num coefficients
|
||||
WinThresh WinT;//window size in frames for adaptive thresholding [pre post]:
|
||||
@ -59,22 +59,22 @@ private:
|
||||
void initialise( TTParams Params );
|
||||
void deInitialise();
|
||||
|
||||
int beatPredict( unsigned int FSP, double alignment, double period, unsigned int step);
|
||||
int phaseMM( double* DF, double* weighting, unsigned int winLength, double period );
|
||||
void createPhaseExtractor( double* Filter, unsigned int winLength, double period, unsigned int fsp, unsigned int lastBeat );
|
||||
int findMeter( double* ACF, unsigned int len, double period );
|
||||
int beatPredict( int FSP, double alignment, double period, int step);
|
||||
int phaseMM( double* DF, double* weighting, int winLength, double period );
|
||||
void createPhaseExtractor( double* Filter, int winLength, double period, int fsp, int lastBeat );
|
||||
int findMeter( double* ACF, int len, double period );
|
||||
void constDetect( double* periodP, int currentIdx, int* flag );
|
||||
void stepDetect( double* periodP, double* periodG, int currentIdx, int* flag );
|
||||
void createCombFilter( double* Filter, unsigned int winLength, unsigned int TSig, double beatLag );
|
||||
void createCombFilter( double* Filter, int winLength, int TSig, double beatLag );
|
||||
double tempoMM( double* ACF, double* weight, int sig );
|
||||
|
||||
unsigned int m_dataLength;
|
||||
unsigned int m_winLength;
|
||||
unsigned int m_lagLength;
|
||||
int m_dataLength;
|
||||
int m_winLength;
|
||||
int m_lagLength;
|
||||
|
||||
double m_rayparam;
|
||||
double m_sigma;
|
||||
double m_DFWVNnorm;
|
||||
double m_rayparam;
|
||||
double m_sigma;
|
||||
double m_DFWVNnorm;
|
||||
|
||||
vector<int> m_beats; // Vector of detected beats
|
||||
|
||||
|
@ -36,7 +36,7 @@ void TCSGram::getTCSVector(int iPosition, TCSVector& rTCSVector) const
|
||||
{
|
||||
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;
|
||||
|
@ -32,7 +32,7 @@ public:
|
||||
|
||||
void printDebug()
|
||||
{
|
||||
for (int i = 0; i < size(); i++)
|
||||
for (int i = 0; i < int(size()); i++)
|
||||
{
|
||||
std::cout << (*this)[i] << ";";
|
||||
}
|
||||
@ -68,7 +68,7 @@ public:
|
||||
|
||||
void printDebug()
|
||||
{
|
||||
for (int i = 0; i < size(); i++)
|
||||
for (int i = 0; i < int(size()); i++)
|
||||
{
|
||||
std::cout << (*this)[i] << ";";
|
||||
}
|
||||
|
@ -4,6 +4,12 @@
|
||||
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
|
||||
|
@ -1843,7 +1843,10 @@ Wavelet::createDecompositionFilters(Type wavelet,
|
||||
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()));
|
||||
}
|
||||
|
||||
|
@ -275,7 +275,8 @@ void kf_work(
|
||||
|
||||
if (m==1) {
|
||||
do{
|
||||
*Fout = *f;
|
||||
Fout->r = f->r;
|
||||
Fout->i = f->i;
|
||||
f += fstride*in_stride;
|
||||
}while(++Fout != Fout_end );
|
||||
}else{
|
||||
|
1
libs/qm-dsp/gitrev.txt
Normal file
1
libs/qm-dsp/gitrev.txt
Normal file
@ -0,0 +1 @@
|
||||
qm-vamp-plugins-v1.7.1-20-g4d15479
|
@ -34,7 +34,7 @@ double CosineDistance::distance(const vector<double> &v1,
|
||||
}
|
||||
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];
|
||||
|
@ -20,6 +20,7 @@
|
||||
#include <vector>
|
||||
#include <cmath>
|
||||
|
||||
using namespace std;
|
||||
|
||||
double MathUtilities::mod(double x, double y)
|
||||
{
|
||||
@ -38,9 +39,9 @@ double MathUtilities::princarg(double ang)
|
||||
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;
|
||||
|
||||
@ -56,17 +57,16 @@ void MathUtilities::getAlphaNorm(const double *data, unsigned int len, unsigned
|
||||
*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;
|
||||
@ -84,13 +84,13 @@ double MathUtilities::round(double x)
|
||||
}
|
||||
}
|
||||
|
||||
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) {
|
||||
@ -100,9 +100,9 @@ double MathUtilities::median(const double *src, unsigned int len)
|
||||
}
|
||||
}
|
||||
|
||||
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++)
|
||||
@ -113,7 +113,7 @@ double MathUtilities::sum(const double *src, unsigned int len)
|
||||
return retVal;
|
||||
}
|
||||
|
||||
double MathUtilities::mean(const double *src, unsigned int len)
|
||||
double MathUtilities::mean(const double *src, int len)
|
||||
{
|
||||
double retVal =0.0;
|
||||
|
||||
@ -126,9 +126,9 @@ double MathUtilities::mean(const double *src, unsigned int len)
|
||||
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.;
|
||||
|
||||
@ -142,9 +142,9 @@ double MathUtilities::mean(const std::vector<double> &src,
|
||||
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) {
|
||||
@ -171,10 +171,10 @@ void MathUtilities::getFrameMinMax(const double *data, unsigned int len, double
|
||||
}
|
||||
}
|
||||
|
||||
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];
|
||||
@ -197,15 +197,15 @@ int MathUtilities::getMax( double* pData, unsigned int Length, double* pMax )
|
||||
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 ];
|
||||
|
||||
@ -286,7 +286,7 @@ void MathUtilities::normalise(double *data, int length, NormaliseType type)
|
||||
}
|
||||
}
|
||||
|
||||
void MathUtilities::normalise(std::vector<double> &data, NormaliseType type)
|
||||
void MathUtilities::normalise(vector<double> &data, NormaliseType type)
|
||||
{
|
||||
switch (type) {
|
||||
|
||||
@ -317,20 +317,46 @@ void MathUtilities::normalise(std::vector<double> &data, NormaliseType 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);
|
||||
}
|
||||
|
@ -35,32 +35,32 @@ public:
|
||||
* 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
|
||||
@ -73,14 +73,21 @@ public:
|
||||
*/
|
||||
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 );
|
||||
/**
|
||||
* 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 void circShift( double* data, int length, int shift);
|
||||
|
||||
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,
|
||||
@ -94,12 +101,35 @@ public:
|
||||
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.
|
||||
*/
|
||||
|
@ -80,22 +80,36 @@ private:
|
||||
|
||||
// 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
|
||||
|
||||
@ -113,9 +127,9 @@ double TPolyFit::PolyFit2 (const vector<double> &x,
|
||||
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;
|
||||
@ -148,13 +162,13 @@ double TPolyFit::PolyFit2 (const vector<double> &x,
|
||||
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;
|
||||
@ -229,8 +243,8 @@ bool TPolyFit::GaussJordan (Matrix &b,
|
||||
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;
|
||||
@ -278,7 +292,7 @@ bool TPolyFit::GaussJordan2(Matrix &b,
|
||||
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)
|
||||
@ -355,53 +369,6 @@ bool TPolyFit::GaussJordan2(Matrix &b,
|
||||
} // { 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
|
||||
|
||||
|
@ -252,7 +252,7 @@ void tqli(double* d, double* e, int n, double** z)
|
||||
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
|
||||
|
||||
|
Loading…
Reference in New Issue
Block a user