qm-dsp 1.8
DetectionFunction.cpp
Go to the documentation of this file.
1/* -*- c-basic-offset: 4 indent-tabs-mode: nil -*- vi:set ts=8 sts=4 sw=4: */
2
3/*
4 QM DSP Library
5
6 Centre for Digital Music, Queen Mary, University of London.
7 This file 2005-2006 Christian Landone.
8
9 This program is free software; you can redistribute it and/or
10 modify it under the terms of the GNU General Public License as
11 published by the Free Software Foundation; either version 2 of the
12 License, or (at your option) any later version. See the file
13 COPYING included with this distribution for more information.
14*/
15
16#include "DetectionFunction.h"
17#include <cstring>
18
20// Construction/Destruction
22
24 m_window(0)
25{
30
31 initialise( Config );
32}
33
38
39
41{
44
45 m_DFType = Config.DFType;
46 m_stepSize = Config.stepSize;
47 m_dbRise = Config.dbRise;
48
52 if (m_whitenRelaxCoeff < 0) m_whitenRelaxCoeff = 0.9997;
53 if (m_whitenFloor < 0) m_whitenFloor = 0.01;
54
55 m_magHistory = new double[ m_halfLength ];
56 memset(m_magHistory,0, m_halfLength*sizeof(double));
57
58 m_phaseHistory = new double[ m_halfLength ];
59 memset(m_phaseHistory,0, m_halfLength*sizeof(double));
60
61 m_phaseHistoryOld = new double[ m_halfLength ];
62 memset(m_phaseHistoryOld,0, m_halfLength*sizeof(double));
63
64 m_magPeaks = new double[ m_halfLength ];
65 memset(m_magPeaks,0, m_halfLength*sizeof(double));
66
68
69 m_magnitude = new double[ m_halfLength ];
70 m_thetaAngle = new double[ m_halfLength ];
71 m_unwrapped = new double[ m_halfLength ];
72
74 m_windowed = new double[ m_dataLength ];
75}
76
78{
79 delete [] m_magHistory ;
80 delete [] m_phaseHistory ;
81 delete [] m_phaseHistoryOld ;
82 delete [] m_magPeaks ;
83
84 delete m_phaseVoc;
85
86 delete [] m_magnitude;
87 delete [] m_thetaAngle;
88 delete [] m_windowed;
89 delete [] m_unwrapped;
90
91 delete m_window;
92}
93
94double DetectionFunction::processTimeDomain(const double *samples)
95{
96 m_window->cut(samples, m_windowed);
97
100
101 if (m_whiten) whiten();
102
103 return runDF();
104}
105
107 const double *imags)
108{
111
112 if (m_whiten) whiten();
113
114 return runDF();
115}
116
118{
119 for (unsigned int i = 0; i < m_halfLength; ++i) {
120 double m = m_magnitude[i];
121 if (m < m_magPeaks[i]) {
122 m = m + (m_magPeaks[i] - m) * m_whitenRelaxCoeff;
123 }
124 if (m < m_whitenFloor) m = m_whitenFloor;
125 m_magPeaks[i] = m;
126 m_magnitude[i] /= m;
127 }
128}
129
131{
132 double retVal = 0;
133
134 switch( m_DFType )
135 {
136 case DF_HFC:
137 retVal = HFC( m_halfLength, m_magnitude);
138 break;
139
140 case DF_SPECDIFF:
142 break;
143
144 case DF_PHASEDEV:
145 // Using the instantaneous phases here actually provides the
146 // same results (for these calculations) as if we had used
147 // unwrapped phases, but without the possible accumulation of
148 // phase error over time
150 break;
151
152 case DF_COMPLEXSD:
154 break;
155
156 case DF_BROADBAND:
158 break;
159 }
160
161 return retVal;
162}
163
164double DetectionFunction::HFC(unsigned int length, double *src)
165{
166 unsigned int i;
167 double val = 0;
168
169 for( i = 0; i < length; i++)
170 {
171 val += src[ i ] * ( i + 1);
172 }
173 return val;
174}
175
176double DetectionFunction::specDiff(unsigned int length, double *src)
177{
178 unsigned int i;
179 double val = 0.0;
180 double temp = 0.0;
181 double diff = 0.0;
182
183 for( i = 0; i < length; i++)
184 {
185 temp = fabs( (src[ i ] * src[ i ]) - (m_magHistory[ i ] * m_magHistory[ i ]) );
186
187 diff= sqrt(temp);
188
189 // (See note in phaseDev below.)
190
191 val += diff;
192
193 m_magHistory[ i ] = src[ i ];
194 }
195
196 return val;
197}
198
199
200double DetectionFunction::phaseDev(unsigned int length, double *srcPhase)
201{
202 unsigned int i;
203 double tmpPhase = 0;
204 double tmpVal = 0;
205 double val = 0;
206
207 double dev = 0;
208
209 for( i = 0; i < length; i++)
210 {
211 tmpPhase = (srcPhase[ i ]- 2*m_phaseHistory[ i ]+m_phaseHistoryOld[ i ]);
212 dev = MathUtilities::princarg( tmpPhase );
213
214 // A previous version of this code only counted the value here
215 // if the magnitude exceeded 0.1. My impression is that
216 // doesn't greatly improve the results for "loud" music (so
217 // long as the peak picker is reasonably sophisticated), but
218 // does significantly damage its ability to work with quieter
219 // music, so I'm removing it and counting the result always.
220 // Same goes for the spectral difference measure above.
221
222 tmpVal = fabs(dev);
223 val += tmpVal ;
224
226 m_phaseHistory[ i ] = srcPhase[ i ];
227 }
228
229 return val;
230}
231
232
233double DetectionFunction::complexSD(unsigned int length, double *srcMagnitude, double *srcPhase)
234{
235 unsigned int i;
236 double val = 0;
237 double tmpPhase = 0;
238 double tmpReal = 0;
239 double tmpImag = 0;
240
241 double dev = 0;
242 ComplexData meas = ComplexData( 0, 0 );
243 ComplexData j = ComplexData( 0, 1 );
244
245 for( i = 0; i < length; i++)
246 {
247 tmpPhase = (srcPhase[ i ]- 2*m_phaseHistory[ i ]+m_phaseHistoryOld[ i ]);
248 dev= MathUtilities::princarg( tmpPhase );
249
250 meas = m_magHistory[i] - ( srcMagnitude[ i ] * exp( j * dev) );
251
252 tmpReal = real( meas );
253 tmpImag = imag( meas );
254
255 val += sqrt( (tmpReal * tmpReal) + (tmpImag * tmpImag) );
256
258 m_phaseHistory[ i ] = srcPhase[ i ];
259 m_magHistory[ i ] = srcMagnitude[ i ];
260 }
261
262 return val;
263}
264
265double DetectionFunction::broadband(unsigned int length, double *src)
266{
267 double val = 0;
268 for (unsigned int i = 0; i < length; ++i) {
269 double sqrmag = src[i] * src[i];
270 if (m_magHistory[i] > 0.0) {
271 double diff = 10.0 * log10(sqrmag / m_magHistory[i]);
272 if (diff > m_dbRise) val = val + 1;
273 }
274 m_magHistory[i] = sqrmag;
275 }
276 return val;
277}
278
283
#define DF_PHASEDEV
#define DF_HFC
#define DF_BROADBAND
#define DF_SPECDIFF
#define DF_COMPLEXSD
#define NULL
Definition Filter.h:20
complex< double > ComplexData
Definition MathAliases.h:23
@ HanningWindow
Definition Window.h:27
unsigned int m_stepSize
PhaseVocoder * m_phaseVoc
unsigned int m_dataLength
double broadband(unsigned int length, double *srcMagnitude)
Window< double > * m_window
double complexSD(unsigned int length, double *srcMagnitude, double *srcPhase)
unsigned int m_halfLength
double specDiff(unsigned int length, double *src)
DetectionFunction(DFConfig Config)
double phaseDev(unsigned int length, double *srcPhase)
void initialise(DFConfig Config)
double processFrequencyDomain(const double *reals, const double *imags)
Process a single frequency-domain frame, provided as frameLength/2+1 real and imaginary component val...
double HFC(unsigned int length, double *src)
double processTimeDomain(const double *samples)
Process a single time-domain frame of audio, provided as frameLength samples.
static double princarg(double ang)
The principle argument function.
void processFrequencyDomain(const double *reals, const double *imags, double *mag, double *phase, double *unwrapped)
Given one frame of frequency-domain samples, return the magnitudes, instantaneous phases,...
void processTimeDomain(const double *src, double *mag, double *phase, double *unwrapped)
Given one frame of time-domain samples, FFT and return the magnitudes, instantaneous phases,...
Various shaped windows for sample frame conditioning, including cosine windows (Hann etc) and triangu...
Definition Window.h:41
void cut(T *src) const
Definition Window.h:61
double whiteningRelaxCoeff
unsigned int stepSize
bool adaptiveWhitening
double whiteningFloor
unsigned int frameLength