0
$\begingroup$

Using formula from Audio EQ Cookbook I implemented biquad IIR filters in Excel.

I can now plot every transfer function and total of 8 band EQ with shelvings and peakings.

My problem is, I can not implement any of the phase response formula in excel. I'm not very good at these kinds of functions and especially implementing complex number functions in excel.

Here is the transfer function.

enter image description here

I have $\omega_0, \alpha, A_x, \omega$ of the frequency, $\phi$ of the frequency $( = 4 * \sin{(\omega/2)}^2$, I don't know what this is), $b_0,b_1,b_2,a_1,a_2$ coefficients of every filter and dB values at each frequency ($H(s)$).

How can I write the phase response function in excel?

Improve this question
$\endgroup$
2

1 Answer 1

Reset to default
2
$\begingroup$

Excel is an awkward language for this type of thing so we can make easier if we do a little math upfront. Doing complex math in Excel is just about as much fun as slamming your fingers in the car door, so we shall try to avoid it.

We have

$$H(z) = \frac{b_0 + b_1z^{-1}+ b_2z^{-2}}{a_0 + a_1z^{-1}+ a_2z^{-2}}= \frac{B(z)}{A(z)}$$

If we want to evaluate this at normalized frequency $z = e^{j\omega}$ we get

$$B(z) = b_0 + b_1\cos(\omega)+b_2\cos(2\omega)-j\left[b_1\sin(\omega)+b_2\sin(2\omega)\right] = B_r+jB_i$$ with $$B_r = b_0 + b_1\cos(\omega)+b_2\cos(2\omega) \\ B_i = -b_1\sin(\omega)-b_2\sin(2\omega) $$

So the whole thing becomes

$$H(\omega) = \frac{B_r+jB_i}{A_r+jA_i}$$

With this we can calculate both magnitude and phase

$$|H(\omega)| = \sqrt{\frac{B_r^2+B_i^2}{A_r^2+A_i^2}} $$

$$\angle{H(\omega)} = \tan^{-1}\left(\frac{B_i}{B_r}\right)-\tan^{-1}\left(\frac{A_i}{A_r}\right)$$

We can implement this in excel by first calculating $B_r$, $B_i$, $A_r$, and $A_i$ and then using the last two formulas to calculate magnitude and phase. The phase can be done using the Excel function atan2() which in typical Excel fashion uses non-standard conventions and reverses the argument order.

You can find an example with a second order Butterworth lowpass at https://docs.google.com/spreadsheets/d/1ye7-WzmwNRpTglJsGk6Dtvb8kge9PICc/edit?usp=sharing&ouid=110367350423053204811&rtpof=true&sd=true

EDIT based on comments

It's been rightfully pointed out that the subtraction of the two phases can run into numerical problems and wrapping problems.

A better way to calculate the phase is the following one

$$H(\omega) = \frac{B_r+jB_i}{A_r+jA_i} = \frac{B_r+jB_i}{A_r+jA_i} \cdot \frac{A_r-jA_i}{A_r-jA_i} = \cdot \frac{B_rA_r+B_iA_i+j(B_iA_r-B_rA_i)}{X}, x \in \mathbb{R}$$

So we can get the phase of a single biquad as

$$\angle{H(\omega)} = \text{atan2}(B_iA_r-B_rA_i, B_rA_r+B_iA_i) $$ where $\text{atan2}()$ is the quadrant correct version of the inverse tangent. The spreadsheet has an added column for the new phase calculation.

Improve this answer
$\endgroup$
8
  • $\begingroup$ That's amazing, you're very helpful. It works. Thanks a lot. I was trying to do this with =COMPLEX(), IMEXP() and =IMARGUMENT() functions but i guess more things get wrong. Just one thing, there is a very little shift at high frequencies, maybe about 0,2°. It's not much important but i wonder why? Image link: [link]im.ge/i/1cTVFG $\endgroup$
    – BugraKezan
    Sep 30, 2022 at 14:12
  • $\begingroup$ Hard to tell. Could be numerical: Sigma Studio generally works for fixed point processors. Could also be coefficient quantization (number of digits). If you want to know which one is correct, you need to use a known-good reference (Matlab, Octave or Python perhaps) or use filter where you can calculate the result by hand. $\endgroup$
    – Hilmar
    Sep 30, 2022 at 14:45
  • $\begingroup$ @BugraKezan if this helped you please also check this off as the right answer to close the open question, thank you! $\endgroup$
    – Dan Boschen
    Sep 30, 2022 at 16:10
  • 1
    $\begingroup$ I know the two arctans are split up, but I am concerned that there still will be a disparity between $\angle{H(\omega)} = \tan^{-1}\left(\frac{B_i}{B_r}\right)-\tan^{-1}\left(\frac{A_i}{A_r}\right)$ and the general complex argument. Also, if the OP wants the phase unwrapped, we should be computing the delta of phase between adjacent points of the phase response. $\endgroup$
    – robert bristow-johnson
    Sep 30, 2022 at 19:53
  • 1
    $\begingroup$ @BugraKezan , please take a look at this answer to a related question. $\endgroup$
    – robert bristow-johnson
    Oct 9, 2022 at 21:40

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.