b=939491 avoid producing infinite stream of subnormals in BiquadFilterNode tail r=padenot

Using double precision for the processing means the main loop is unlikely to
encounter subnormals and provides suitable extra precision to minimize
accumulation of error.

If the tail output values are denormalized when converted to single precision
then flush them to zero to reduce downstream computation cost.

--HG--
extra : transplant_source : %08%09u%C48%E3i%AB%23%1B%D8_kz%E5%A8_%3D%C4%93
This commit is contained in:
Karl Tomlinson 2013-12-03 12:07:16 +13:00
parent 33afa7bd37
commit 6bef079a18
2 changed files with 31 additions and 13 deletions

View File

@ -26,9 +26,15 @@
* THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*/
#include "DenormalDisabler.h"
// For M_PI from cmath
#ifdef _MSC_VER
# define _USE_MATH_DEFINES
#endif
#include "Biquad.h"
#include <cmath>
#include <float.h>
#include <algorithm>
namespace WebCore {
@ -47,8 +53,6 @@ Biquad::~Biquad()
void Biquad::process(const float* sourceP, float* destP, size_t framesToProcess)
{
int n = framesToProcess;
// Create local copies of member variables
double x1 = m_x1;
double x2 = m_x2;
@ -61,12 +65,12 @@ void Biquad::process(const float* sourceP, float* destP, size_t framesToProcess)
double a1 = m_a1;
double a2 = m_a2;
while (n--) {
for (size_t i = 0; i < framesToProcess; ++i) {
// FIXME: this can be optimized by pipelining the multiply adds...
float x = *sourceP++;
float y = b0*x + b1*x1 + b2*x2 - a1*y1 - a2*y2;
double x = sourceP[i];
double y = b0*x + b1*x1 + b2*x2 - a1*y1 - a2*y2;
*destP++ = y;
destP[i] = y;
// Update state variables
x2 = x1;
@ -75,12 +79,22 @@ void Biquad::process(const float* sourceP, float* destP, size_t framesToProcess)
y1 = y;
}
// Local variables back to member. Flush denormals here so we
// don't slow down the inner loop above.
m_x1 = DenormalDisabler::flushDenormalFloatToZero(x1);
m_x2 = DenormalDisabler::flushDenormalFloatToZero(x2);
m_y1 = DenormalDisabler::flushDenormalFloatToZero(y1);
m_y2 = DenormalDisabler::flushDenormalFloatToZero(y2);
// Avoid introducing a stream of subnormals when input is silent and the
// tail approaches zero.
if (x1 == 0.0 && x2 == 0.0 && (y1 != 0.0 || y2 != 0.0) &&
fabs(y1) < FLT_MIN && fabs(y2) < FLT_MIN) {
// Flush future values to zero (until there is new input).
y1 = y2 = 0.0;
// Flush calculated values.
for (int i = framesToProcess; i-- && fabs(destP[i]) < FLT_MIN; ) {
destP[i] = 0.0f;
}
}
// Local variables back to member.
m_x1 = x1;
m_x2 = x2;
m_y1 = y1;
m_y2 = y2;
}
void Biquad::reset()

View File

@ -93,6 +93,10 @@ private:
double m_a2;
// Filter memory
//
// Double precision for the output values is valuable because errors can
// accumulate. Input values are also stored as double so they need not be
// converted again for computation.
double m_x1; // input delayed by 1 sample
double m_x2; // input delayed by 2 samples
double m_y1; // output delayed by 1 sample