| OLD | NEW |
|---|
| (Empty) | |
| 1 <!doctype html> |
| 2 <html> |
| 3 <head> |
| 4 <title>Test Tail Time for IIRFilter</title> |
| 5 <script src="../../resources/testharness.js"></script> |
| 6 <script src="../../resources/testharnessreport.js"></script> |
| 7 <script src="../resources/audit-util.js"></script> |
| 8 <script src="../resources/audit.js"></script> |
| 9 </head> |
| 10 |
| 11 <body> |
| 12 <script> |
| 13 let audit = Audit.createTaskRunner(); |
| 14 let renderQuantumFrames = 128; |
| 15 |
| 16 // Must be a power of two to eliminate round-off differences between thsi |
| 17 // JS code and the WebAudio implementation. Otherwise, the sample rate is |
| 18 // arbitrary. |
| 19 let sampleRate = 16384; |
| 20 |
| 21 // Fairly arbitrary, but should be long enough so that the node propagates |
| 22 // silence before the end of the offline context. |
| 23 let renderDuration = 1; |
| 24 |
| 25 |
| 26 audit.define('1-pole tail', (task, should) => { |
| 27 let pole = 0.99; |
| 28 let IIROptions = {feedforward: [1], feedback: [1, -pole]}; |
| 29 // For the given filter, we can actually compute where the tail |
| 30 // begins. The impulse response for the 1-pole filter is h(n) = |
| 31 // a^n, where a = 0.9. The tail here starts when a^n < eps = |
| 32 // 1/32768. So n > log(eps)/log(a), or 98.7. Round that up to the |
| 33 // nearest render quantum frames. |
| 34 let tail = Math.ceil(Math.log(1 / 32768) / Math.log(pole)); |
| 35 |
| 36 runTest(should, IIROptions, tail, '1-pole').then(() => task.done()); |
| 37 }); |
| 38 |
| 39 audit.define('2 real pole test', (task, should) => { |
| 40 // Simple example of a 2-pole IIR filter where both poles are real. |
| 41 // We arbitrarily select a pole at 9.99 and one at -0.5. The IIRFilter |
| 42 // is then |
| 43 // 1 / ((z-0.99) * (z + 0.5)) |
| 44 // = 1/(z^2-0.49z-0.495) |
| 45 // = z^-2/(1-0.49/z-0.495/z^2) |
| 46 let IIROptions = {feedforward: [0, 0, 1], feedback: [1, -0.49, -0.495]}; |
| 47 |
| 48 // For this particular filter, we can analytically compute the impulse |
| 49 // response using partical fractios: |
| 50 // |
| 51 // 1 / ((z-0.99) * (z + 0.5)) |
| 52 // = 1/(-0.5-0.99)/(z + 0.5) - 1/(-0.5-0.99)/(z - 0.99) |
| 53 // = 1/1.49*(1/(z-0.99) - 1/(z+0.5)) |
| 54 // = 1/1.49*[1/z*sum(.99^n/z^n,n,0,inf) |
| 55 // - 1/z*sum((-0.5)^n/z^n,n,0,inf)] |
| 56 // = 1/1.49/z*sum((0.99^n-(-0.5)^n)/z^n) |
| 57 // |
| 58 // So the tail begins when 1/1.49*(0.99^n-(-0.5)^n) < 1/32768. This can |
| 59 // be solved numerically to give n = 995. |
| 60 let tail = 995; |
| 61 tail = renderQuantumFrames * Math.ceil(tail / renderQuantumFrames); |
| 62 |
| 63 runTest(should, IIROptions, tail, '2 real poles') |
| 64 .then(() => task.done()); |
| 65 }); |
| 66 |
| 67 audit.define('2 complex poles', (task, should) => { |
| 68 // Simple example of a 2-pole IIR filter where both poles are complex |
| 69 // conjugates. In this case, the poles will be r*exp(+/-i*theta) where |
| 70 // r = 0.99 and theta = 0.01. The filter is then |
| 71 // |
| 72 // 1/(z^2-2*r*cos(theta) + r^2) |
| 73 // = z^(-2)/(1-2*r*cos(theta)/z + r^2/z^2) |
| 74 let r = 0.99; |
| 75 let theta = 0.01; |
| 76 let IIROptions = { |
| 77 feedforward: [0, 0, 1], |
| 78 feedback: [1, -2 * r * Math.cos(theta), r * r] |
| 79 }; |
| 80 |
| 81 // Again, we can use partial fractions as for 2 real pole case to get an |
| 82 // analytically solution for the impulse response. For simplicity, let |
| 83 // p1 = r*exp(i*theta), p2 = r*exp(-i*theta). Then: |
| 84 // |
| 85 // 1/(z^2-2*r*cos(theta) + r^2) |
| 86 // = 1/(z-p1)/(z-p2) |
| 87 // = 1/(p2-p1)*[1/(z-p2) - 1/(z-p1)] |
| 88 // = 1/(p2-p1)*[1/z*sum(p2^n/z^n) - 1/z*sum(p1^n/z^n)] |
| 89 // = 1/(p2-p1)/z*sum((p2^n-p1^n)/z^n) |
| 90 // |
| 91 // So the tail begins when |
| 92 // 1/32768 > |1/(p2-p1)*(p2^n-p1^n)| |
| 93 // = 1/(r*sin(theta))*|r^n*(exp(-i*theta*n)-exp(i*theta*n))| |
| 94 // = 1/(2*r*sin(theta))*(2*r^n*|sin(theta*n)|); |
| 95 // = r^(n-1)*|sin(theta*n)|/sin(theta) |
| 96 // |
| 97 // This can be solved numerically to for n; |
| 98 let tail = 1474.256; |
| 99 tail = renderQuantumFrames * Math.ceil(tail / renderQuantumFrames); |
| 100 |
| 101 runTest(should, IIROptions, tail, '2 complex poles') |
| 102 .then(() => task.done()); |
| 103 }); |
| 104 |
| 105 audit.define('repeated poles', (task, should) => { |
| 106 // Two repeated roots. Let p be the repeated pole. Then the filter is |
| 107 // |
| 108 // 1/(z-p)^2 |
| 109 // = z^(-2)/(1-p/z)^2 |
| 110 // = z^(-2)/(1-2*p/z+p*p/z^2) |
| 111 |
| 112 let pole = 0.99; |
| 113 let IIROptions = { |
| 114 feedforward: [0,0,1], |
| 115 feedback: [1, -2*pole, pole*pole] |
| 116 }; |
| 117 |
| 118 // We can analytically compute the impulse response of this filter to be |
| 119 // |
| 120 // 1/z^2*sum(p^n*(n+1)/z^n, n, 0, inf) |
| 121 // = sum(p^n*(n+1)/z^(n+2), n, 0, inf) |
| 122 // = 1/p^2*sum((p^k*(k-1))/z^k,k,2,inf)) |
| 123 // |
| 124 // Therefore the tail starts when p^(k-2)*(k-1) < 1/32768. We can solve |
| 125 // this numerically to be 1781.213; |
| 126 |
| 127 let tail = 1781.213; |
| 128 runTest(should, IIROptions, tail, '2 repeated poles') |
| 129 .then(() => task.done()); |
| 130 |
| 131 }); |
| 132 |
| 133 audit.define('4-th order', (task, should) => { |
| 134 // Test consistency of tail times between a 4-th order direct IIR filter |
| 135 // and the equivalent cascade of second-order sections. The first |
| 136 // channel of the output is the cascaded biquad, and the second channel |
| 137 // is the 4-th order equivalent. |
| 138 let context = |
| 139 new OfflineAudioContext(2, renderDuration * sampleRate, sampleRate); |
| 140 |
| 141 let src = new AudioBufferSourceNode( |
| 142 context, {buffer: createImpulseBuffer(context, 1)}); |
| 143 |
| 144 // This is a 4-th order lowpass elliptic filter designed using |
| 145 // http://rtoy.github.io/webaudio-hacks/more/filter-design/filter-design
.html. |
| 146 // The sample rate is 16384 Hz with a passband at 3600 Hz with a 0.25 dB |
| 147 // attenuation, and a stopband at 4800 Hz, with a stopband attenuation |
| 148 // of 30 dB. (Nothing really special except that this gives a 4-th order |
| 149 // filter). |
| 150 |
| 151 let f0 = context.createIIRFilter( |
| 152 [0.6410686464424084, 0.2607836369670137, 0.6410686464424084], |
| 153 [1, -0.2287413068432929, 0.7716622366951231]); |
| 154 let f1 = context.createIIRFilter( |
| 155 [0.21283904239866536, 0.3184888523034876, 0.21283904239866536], |
| 156 [1, -0.4686913542990081, 0.21285829139982618]); |
| 157 |
| 158 // The poles for f0 are 0.1143706534216465 +/- 0.8709658950447078*i or |
| 159 // 0.8784430753868592*%e^(+/-1.440228658066206*%i), |
| 160 // |
| 161 // The poles for f1 are 0.2343456771495041 +/- 0.3974171548903829*i or |
| 162 // 0.4613656807780854*%e^(+/-1.038005727602151*%i. |
| 163 // |
| 164 // Thus, the tail time for f0 is approximately 80, but this is an |
| 165 // approximation since we didn't include the affect of the numerator. |
| 166 // Round this up to the next render to get an actual tail time of 128. |
| 167 // |
| 168 // Similarly, for f0, the tail time is 14.3. Thus, the actual tail time |
| 169 // is alos 128 for this filter. |
| 170 // |
| 171 // Since these biquads are cascaded, the total tail time for both is the |
| 172 // sum or 256 frames. However, the tail actually ends two render quanta |
| 173 // after this for a total of 512 frames. |
| 174 |
| 175 let biquadTailEnd = 512; |
| 176 |
| 177 // The equivalent 4-th order filter created multiplying the f0 and f1 |
| 178 // coefficients together appropriately. |
| 179 let f = context.createIIRFilter( |
| 180 [ |
| 181 0.136444436820611, 0.259678157018493, 0.355945554878375, |
| 182 0.259678157018493, 0.136444436820611 |
| 183 ], |
| 184 [ |
| 185 1.000000000000000, -0.697432661142301, 1.091729600983457, |
| 186 -0.410360902525266, 0.164254705240692 |
| 187 ]); |
| 188 |
| 189 let merger = context.createChannelMerger(2); |
| 190 merger.connect(context.destination); |
| 191 |
| 192 src.connect(f0).connect(f1).connect(merger, 0, 0); |
| 193 src.connect(f).connect(merger, 0, 1); |
| 194 |
| 195 src.start(); |
| 196 |
| 197 context.startRendering() |
| 198 .then(renderedBuffer => { |
| 199 // c0 = cascaded biquads |
| 200 // c1 = 4-th order filter |
| 201 let c0 = renderedBuffer.getChannelData(0); |
| 202 let c1 = renderedBuffer.getChannelData(1); |
| 203 |
| 204 // Sanity check: The two filters should have the same output |
| 205 // within some rounding error. |
| 206 should( |
| 207 c0.slice(0, biquadTailEnd), |
| 208 'Filter outputs[0:' + (biquadTailEnd - 1) + ']') |
| 209 .beCloseToArray( |
| 210 c1.slice(0, biquadTailEnd), |
| 211 {absoluteThreshold: 1.4902e-8}); |
| 212 should( |
| 213 c0.slice(biquadTailEnd), |
| 214 'Filter outputs[' + biquadTailEnd + ':]') |
| 215 .beEqualToArray(c1.slice(biquadTailEnd)); |
| 216 |
| 217 // Verify that after the tail time, the outputs are zero, and not |
| 218 // before for both the biquads and 4-th order filters. |
| 219 should( |
| 220 c0.slice(0, biquadTailEnd), |
| 221 'cascaded biquad output[0:' + (biquadTailEnd - 1) + ']') |
| 222 .notBeConstantValueOf(0); |
| 223 should( |
| 224 c0.slice(biquadTailEnd), |
| 225 'cascaded biquad output[' + biquadTailEnd + ':]') |
| 226 .beConstantValueOf(0); |
| 227 |
| 228 should( |
| 229 c1.slice(0, biquadTailEnd), |
| 230 '4-th order output[0:' + (biquadTailEnd - 1) + ']') |
| 231 .notBeConstantValueOf(0); |
| 232 should( |
| 233 c1.slice(biquadTailEnd), |
| 234 '4-th order output[' + biquadTailEnd + ':]') |
| 235 .beConstantValueOf(0); |
| 236 }) |
| 237 .then(() => task.done()); |
| 238 }); |
| 239 |
| 240 function runTest(should, IIROptions, tailFrames, prefix) { |
| 241 let context = |
| 242 new OfflineAudioContext(1, renderDuration * sampleRate, sampleRate); |
| 243 |
| 244 let src = new AudioBufferSourceNode( |
| 245 context, {buffer: createImpulseBuffer(context, 1)}); |
| 246 |
| 247 let iir = new IIRFilterNode(context, IIROptions); |
| 248 |
| 249 src.connect(iir).connect(context.destination); |
| 250 |
| 251 src.start(); |
| 252 |
| 253 return context.startRendering().then(renderedBuffer => { |
| 254 let audio = renderedBuffer.getChannelData(0); |
| 255 |
| 256 // Round up the tailFrames to the nearest render quantum. |
| 257 let tailQuantum = |
| 258 renderQuantumFrames * Math.ceil(tailFrames / renderQuantumFrames); |
| 259 let tailEndFrame = tailQuantum + 2 * renderQuantumFrames; |
| 260 |
| 261 should(tailEndFrame, prefix + ': tail end frame') |
| 262 .beLessThanOrEqualTo(context.length); |
| 263 |
| 264 // Clamp to the render duration so we don't go off the end. |
| 265 tailEndFrame = Math.min(tailEndFrame, context.length); |
| 266 |
| 267 for (let k = 0; k < tailEndFrame; k += renderQuantumFrames) { |
| 268 should( |
| 269 audio.slice(k, k + renderQuantumFrames), |
| 270 prefix + ': output[' + k + ':' + (k + renderQuantumFrames - 1) + |
| 271 ']') |
| 272 .notBeConstantValueOf(0); |
| 273 } |
| 274 |
| 275 if (tailEndFrame < context.length) { |
| 276 // All frames after should be zero because we're propagating |
| 277 // silence. |
| 278 should( |
| 279 audio.slice(tailEndFrame), |
| 280 'output[' + tailEndFrame + ':' + (context.length - 1) + ']') |
| 281 .beConstantValueOf(0); |
| 282 } |
| 283 }); |
| 284 } |
| 285 |
| 286 audit.run(); |
| 287 </script> |
| 288 </body> |
| 289 </html> |
| OLD | NEW |
|---|
Chromium Code Reviews
Issue 2851873003:
Compute the tail time for an IIRFilter from its coefficients (Closed)
