diff options
Diffstat (limited to 'libs/qm-dsp/maths/MathUtilities.cpp')
-rw-r--r-- | libs/qm-dsp/maths/MathUtilities.cpp | 401 |
1 files changed, 401 insertions, 0 deletions
diff --git a/libs/qm-dsp/maths/MathUtilities.cpp b/libs/qm-dsp/maths/MathUtilities.cpp new file mode 100644 index 0000000000..809874121e --- /dev/null +++ b/libs/qm-dsp/maths/MathUtilities.cpp @@ -0,0 +1,401 @@ +/* -*- c-basic-offset: 4 indent-tabs-mode: nil -*- vi:set ts=8 sts=4 sw=4: */ + +/* + 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 + 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. +*/ + +#include "MathUtilities.h" + +#include <iostream> +#include <cmath> + + +double MathUtilities::mod(double x, double y) +{ + double a = floor( x / y ); + + double b = x - ( y * a ); + return b; +} + +double MathUtilities::princarg(double ang) +{ + double ValOut; + + ValOut = mod( ang + M_PI, -2 * M_PI ) + M_PI; + + return ValOut; +} + +void MathUtilities::getAlphaNorm(const double *data, unsigned int len, unsigned int alpha, double* ANorm) +{ + unsigned int i; + 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; + a = ::pow( a, ( 1.0 / (double) alpha ) ); + + *ANorm = a; +} + +double MathUtilities::getAlphaNorm( const std::vector <double> &data, unsigned int alpha ) +{ + unsigned int i; + unsigned 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; + a = ::pow( a, ( 1.0 / (double) alpha ) ); + + return a; +} + +double MathUtilities::round(double x) +{ + double val = (double)floor(x + 0.5); + + return val; +} + +double MathUtilities::median(const double *src, unsigned int len) +{ + unsigned int i, j; + double tmp = 0.0; + double tempMedian; + double medianVal; + + double* scratch = new double[ len ];//Vector < double > sortedX = Vector < double > ( size ); + + for ( i = 0; i < len; i++ ) + { + scratch[i] = src[i]; + } + + for ( i = 0; i < len - 1; i++ ) + { + for ( j = 0; j < len - 1 - i; j++ ) + { + if ( scratch[j + 1] < scratch[j] ) + { + // compare the two neighbors + tmp = scratch[j]; // swap a[j] and a[j+1] + scratch[j] = scratch[j + 1]; + scratch[j + 1] = tmp; + } + } + } + int middle; + if ( len % 2 == 0 ) + { + middle = len / 2; + tempMedian = ( scratch[middle] + scratch[middle - 1] ) / 2; + } + else + { + middle = ( int )floor( len / 2.0 ); + tempMedian = scratch[middle]; + } + + medianVal = tempMedian; + + delete [] scratch; + return medianVal; +} + +double MathUtilities::sum(const double *src, unsigned int len) +{ + unsigned int i ; + double retVal =0.0; + + for( i = 0; i < len; i++) + { + retVal += src[ i ]; + } + + return retVal; +} + +double MathUtilities::mean(const double *src, unsigned int len) +{ + double retVal =0.0; + + double s = sum( src, len ); + + retVal = s / (double)len; + + return retVal; +} + +double MathUtilities::mean(const std::vector<double> &src, + unsigned int start, + unsigned int count) +{ + double sum = 0.; + + for (int i = 0; i < count; ++i) + { + sum += src[start + i]; + } + + return sum / count; +} + +void MathUtilities::getFrameMinMax(const double *data, unsigned int len, double *min, double *max) +{ + unsigned int i; + double temp = 0.0; + double a=0.0; + + if (len == 0) { + *min = *max = 0; + return; + } + + *min = data[0]; + *max = data[0]; + + for( i = 0; i < len; i++) + { + temp = data[ i ]; + + if( temp < *min ) + { + *min = temp ; + } + if( temp > *max ) + { + *max = temp ; + } + + } +} + +int MathUtilities::getMax( double* pData, unsigned int Length, double* pMax ) +{ + unsigned int index = 0; + unsigned int i; + double temp = 0.0; + + double max = pData[0]; + + for( i = 0; i < Length; i++) + { + temp = pData[ i ]; + + if( temp > max ) + { + max = temp ; + index = i; + } + + } + + if (pMax) *pMax = max; + + + return index; +} + +int MathUtilities::getMax( const std::vector<double> & data, double* pMax ) +{ + unsigned int index = 0; + unsigned int i; + double temp = 0.0; + + double max = data[0]; + + for( i = 0; i < data.size(); i++) + { + temp = data[ i ]; + + if( temp > max ) + { + max = temp ; + index = i; + } + + } + + if (pMax) *pMax = max; + + + return index; +} + +void MathUtilities::circShift( double* pData, int length, int shift) +{ + shift = shift % length; + double temp; + int i,n; + + for( i = 0; i < shift; i++) + { + temp=*(pData + length - 1); + + for( n = length-2; n >= 0; n--) + { + *(pData+n+1)=*(pData+n); + } + + *pData = temp; + } +} + +int MathUtilities::compareInt (const void * a, const void * b) +{ + return ( *(int*)a - *(int*)b ); +} + +void MathUtilities::normalise(double *data, int length, NormaliseType type) +{ + switch (type) { + + case NormaliseNone: return; + + case NormaliseUnitSum: + { + double sum = 0.0; + for (int i = 0; i < length; ++i) { + sum += data[i]; + } + if (sum != 0.0) { + for (int i = 0; i < length; ++i) { + data[i] /= sum; + } + } + } + break; + + case NormaliseUnitMax: + { + double max = 0.0; + for (int i = 0; i < length; ++i) { + if (fabs(data[i]) > max) { + max = fabs(data[i]); + } + } + if (max != 0.0) { + for (int i = 0; i < length; ++i) { + data[i] /= max; + } + } + } + break; + + } +} + +void MathUtilities::normalise(std::vector<double> &data, NormaliseType type) +{ + switch (type) { + + case NormaliseNone: return; + + case NormaliseUnitSum: + { + double sum = 0.0; + for (int i = 0; i < data.size(); ++i) sum += data[i]; + if (sum != 0.0) { + for (int i = 0; i < data.size(); ++i) data[i] /= sum; + } + } + break; + + case NormaliseUnitMax: + { + double max = 0.0; + for (int i = 0; i < data.size(); ++i) { + if (fabs(data[i]) > max) max = fabs(data[i]); + } + if (max != 0.0) { + for (int i = 0; i < data.size(); ++i) data[i] /= max; + } + } + break; + + } +} + +void MathUtilities::adaptiveThreshold(std::vector<double> &data) +{ + int sz = int(data.size()); + if (sz == 0) return; + + std::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); + + smoothed[i] = mean(data, first, last - first + 1); + } + + for (int i = 0; i < sz; i++) { + data[i] -= smoothed[i]; + if (data[i] < 0.0) data[i] = 0.0; + } +} + +bool +MathUtilities::isPowerOfTwo(int x) +{ + if (x < 2) return false; + if (x & (x-1)) return false; + return true; +} + +int +MathUtilities::nextPowerOfTwo(int x) +{ + if (isPowerOfTwo(x)) return x; + int n = 1; + while (x) { x >>= 1; n <<= 1; } + return n; +} + +int +MathUtilities::previousPowerOfTwo(int x) +{ + if (isPowerOfTwo(x)) return x; + int n = 1; + x >>= 1; + while (x) { x >>= 1; n <<= 1; } + return n; +} + +int +MathUtilities::nearestPowerOfTwo(int x) +{ + if (isPowerOfTwo(x)) return x; + int n0 = previousPowerOfTwo(x), n1 = nearestPowerOfTwo(x); + if (x - n0 < n1 - x) return n0; + else return n1; +} + |