Mercurial > hg > qm-dsp

comparison dsp/signalconditioning/FiltFilt.cpp @ 483:fdaa63607c15

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Untabify, indent, tidy
author Chris Cannam <cannam@all-day-breakfast.com>
date Fri, 31 May 2019 11:54:32 +0100
parents fa851e147e3f
children b1f72e469ec8
comparison
equal deleted inserted replaced
482:cbe668c7d724 483:fdaa63607c15
28 FiltFilt::~FiltFilt() 28 FiltFilt::~FiltFilt()
29 { 29 {
30 } 30 }
31 31
32 void FiltFilt::process(double *src, double *dst, unsigned int length) 32 void FiltFilt::process(double *src, double *dst, unsigned int length)
33 { 33 {
34 unsigned int i; 34 unsigned int i;
35 35
36 if (length == 0) return; 36 if (length == 0) return;
37 37
38 unsigned int nFilt = m_ord + 1; 38 unsigned int nFilt = m_ord + 1;
39 unsigned int nFact = 3 * ( nFilt - 1); 39 unsigned int nFact = 3 * ( nFilt - 1);
40 unsigned int nExt = length + 2 * nFact; 40 unsigned int nExt = length + 2 * nFact;
41 41
42 double *filtScratchIn = new double[ nExt ]; 42 double *filtScratchIn = new double[ nExt ];
43 double *filtScratchOut = new double[ nExt ]; 43 double *filtScratchOut = new double[ nExt ];
44 44
45 for( i = 0; i< nExt; i++ ) 45 for( i = 0; i< nExt; i++ ) {
46 { 46 filtScratchIn[ i ] = 0.0;
47 filtScratchIn[ i ] = 0.0; 47 filtScratchOut[ i ] = 0.0;
48 filtScratchOut[ i ] = 0.0;
49 } 48 }
50 49
51 // Edge transients reflection 50 // Edge transients reflection
52 double sample0 = 2 * src[ 0 ]; 51 double sample0 = 2 * src[ 0 ];
53 double sampleN = 2 * src[ length - 1 ]; 52 double sampleN = 2 * src[ length - 1 ];
54 53
55 unsigned int index = 0; 54 unsigned int index = 0;
56 for( i = nFact; i > 0; i-- ) 55 for( i = nFact; i > 0; i-- ) {
57 { 56 filtScratchIn[ index++ ] = sample0 - src[ i ];
58 filtScratchIn[ index++ ] = sample0 - src[ i ];
59 } 57 }
60 index = 0; 58 index = 0;
61 for( i = 0; i < nFact; i++ ) 59 for( i = 0; i < nFact; i++ ) {
62 { 60 filtScratchIn[ (nExt - nFact) + index++ ] = sampleN - src[ (length - 2) - i ];
63 filtScratchIn[ (nExt - nFact) + index++ ] = sampleN - src[ (length - 2) - i ];
64 } 61 }
65 62
66 index = 0; 63 index = 0;
67 for( i = 0; i < length; i++ ) 64 for( i = 0; i < length; i++ ) {
68 { 65 filtScratchIn[ i + nFact ] = src[ i ];
69 filtScratchIn[ i + nFact ] = src[ i ];
70 } 66 }
71 67
72 //////////////////////////////// 68 ////////////////////////////////
73 // Do 0Ph filtering 69 // Do 0Ph filtering
74 m_filter.process( filtScratchIn, filtScratchOut, nExt); 70 m_filter.process( filtScratchIn, filtScratchOut, nExt);
75 71
76 // reverse the series for FILTFILT 72 // reverse the series for FILTFILT
77 for ( i = 0; i < nExt; i++) 73 for ( i = 0; i < nExt; i++) {
78 { 74 filtScratchIn[ i ] = filtScratchOut[ nExt - i - 1];
79 filtScratchIn[ i ] = filtScratchOut[ nExt - i - 1];
80 } 75 }
81 76
82 // do FILTER again 77 // do FILTER again
83 m_filter.process( filtScratchIn, filtScratchOut, nExt); 78 m_filter.process( filtScratchIn, filtScratchOut, nExt);
84 79
85 // reverse the series back 80 // reverse the series back
86 for ( i = 0; i < nExt; i++) 81 for ( i = 0; i < nExt; i++) {
87 { 82 filtScratchIn[ i ] = filtScratchOut[ nExt - i - 1 ];
88 filtScratchIn[ i ] = filtScratchOut[ nExt - i - 1 ];
89 } 83 }
90 for ( i = 0;i < nExt; i++) 84 for ( i = 0;i < nExt; i++) {
91 { 85 filtScratchOut[ i ] = filtScratchIn[ i ];
92 filtScratchOut[ i ] = filtScratchIn[ i ];
93 } 86 }
94 87
95 index = 0; 88 index = 0;
96 for( i = 0; i < length; i++ ) 89 for( i = 0; i < length; i++ ) {
97 { 90 dst[ index++ ] = filtScratchOut[ i + nFact ];
98 dst[ index++ ] = filtScratchOut[ i + nFact ]; 91 }
99 }
100 92
101 delete [] filtScratchIn; 93 delete [] filtScratchIn;
102 delete [] filtScratchOut; 94 delete [] filtScratchOut;
103
104 } 95 }
105 96
106 void FiltFilt::reset() 97 void FiltFilt::reset()
107 { 98 {
108 99