Forum Rule: Always post complete source code & details to reproduce any issue!

Code:
```void setCoefficients(uint32_t stage, const int *coefficients);
void setCoefficients(uint32_t stage, const double *coefficients) {
int coef[5];
// 1073741824 bytes = 1 GB
// Google! Prob has something to do with memory allocation
coef[0] = coefficients[0] * 1073741824.0;
coef[1] = coefficients[1] * 1073741824.0;
coef[2] = coefficients[2] * 1073741824.0;
coef[3] = coefficients[3] * 1073741824.0;
coef[4] = coefficients[4] * 1073741824.0;
setCoefficients(stage, coef);
}

// Compute common filter functions
// http://www.musicdsp.org/files/Audio-EQ-Cookbook.txt
void setLowpass(uint32_t stage, float frequency, float q = 0.7071) {
int coef[5];
double w0 = frequency * (2 * 3.141592654 / AUDIO_SAMPLE_RATE_EXACT);
double sinW0 = sin(w0);
double alpha = sinW0 / ((double)q * 2.0);
double cosW0 = cos(w0);
double scale = 1073741824.0 / (1.0 + alpha);
/* b0 */ coef[0] = ((1.0 - cosW0) / 2.0) * scale;
/* b1 */ coef[1] = (1.0 - cosW0) * scale;
/* b2 */ coef[2] = coef[0];
/* a1 */ coef[3] = (-2.0 * cosW0) * scale;
/* a2 */ coef[4] = (1.0 - alpha) * scale;
setCoefficients(stage, coef);
}```
I understand that we are using bilinear transformation to obtain the coefficients. But I am not quite sure about the scaling factor here which is 1073741824.0. If someone could explain the reason that we include such number in the code, I would very much appreciate it. Thank you!!

2. Originally Posted by yihanhu
Code:
```void setCoefficients(uint32_t stage, const int *coefficients);
void setCoefficients(uint32_t stage, const double *coefficients) {
int coef[5];
// 1073741824 bytes = 1 GB
// Google! Prob has something to do with memory allocation
coef[0] = coefficients[0] * 1073741824.0;
coef[1] = coefficients[1] * 1073741824.0;
coef[2] = coefficients[2] * 1073741824.0;
coef[3] = coefficients[3] * 1073741824.0;
coef[4] = coefficients[4] * 1073741824.0;
setCoefficients(stage, coef);
}

// Compute common filter functions
// http://www.musicdsp.org/files/Audio-EQ-Cookbook.txt
void setLowpass(uint32_t stage, float frequency, float q = 0.7071) {
int coef[5];
double w0 = frequency * (2 * 3.141592654 / AUDIO_SAMPLE_RATE_EXACT);
double sinW0 = sin(w0);
double alpha = sinW0 / ((double)q * 2.0);
double cosW0 = cos(w0);
double scale = 1073741824.0 / (1.0 + alpha);
/* b0 */ coef[0] = ((1.0 - cosW0) / 2.0) * scale;
/* b1 */ coef[1] = (1.0 - cosW0) * scale;
/* b2 */ coef[2] = coef[0];
/* a1 */ coef[3] = (-2.0 * cosW0) * scale;
/* a2 */ coef[4] = (1.0 - alpha) * scale;
setCoefficients(stage, coef);
}```
I understand that we are using bilinear transformation to obtain the coefficients. But I am not quite sure about the scaling factor here which is 1073741824.0. If someone could explain the reason that we include such number in the code, I would very much appreciate it. Thank you!!
FWIW, is 1073741824 is 0x40000000 (or 0x1UL << 30). I don't know why it is used.

3. To get good performance on older Teensy's, most of the audio code uses fixed-point instead of floating-point math.

In this case, biquad coefs are converted to fixed-point with sign bit, one integer bit, and 30 fraction bits. This effectively represents [-2.0, +2.0] with enough fraction bits to handle low-cutoff high-Q filters where coefs are very sensitive to quantization. DSP people usually call this Q30.

You can see how it works by following the scaling of the biquad loop:
Code:
```            in2 = *data;
sum = signed_multiply_accumulate_32x16b(sum, b0, in2);
sum = signed_multiply_accumulate_32x16t(sum, b1, bprev);
sum = signed_multiply_accumulate_32x16b(sum, b2, bprev);
sum = signed_multiply_accumulate_32x16t(sum, a1, aprev);
sum = signed_multiply_accumulate_32x16b(sum, a2, aprev);
out2 = signed_saturate_rshift(sum, 16, 14);
sum &= 0x3FFF;```
You can consider the input samples to be Q0 (no fraction bits).
Each multiply-add does sum += (a * b) >> 16, so the scaling is (Q30 * Q0) >> 16 = Q14 and sum keeps 14 fraction bits.
The final shift discards the 14 fraction bits and saturates the output to 16 integer bits.

A problem with just throwing away the fraction bits, is output is reused by the feedback path (aprev in the code above). A low-freq high-Q filter recirculates for a long time, and the accumulated quantization noise ruins the SNR. One solution is to keep fraction bits in the feedback path, but that requires a 32x32 multiply instead of 32x16 SIMD multiply which requires fewer loads.

You get almost as good result by using first-order quantization-error feedback. The above code implements this via "fraction saving." It is important to note how a signed fixed-point shift (sum >> 14) implements "truncation towards -infinity". That is, 1.xxx truncates to 1 and -1.xxx truncates to -2 so quantization error is guaranteed to be a positive fraction. The (sum &= 0x3FFF) extracts this fraction and saves it for the next iteration. Then by replacing the first multiply with a multiply-add, it feeds back into the next sum essentially for free.

Hope this helps.

#### Posting Permissions

• You may not post new threads
• You may not post replies
• You may not post attachments
• You may not edit your posts
•