# Bilinear Transform (Tustin's Method) applied to the Derivative

I hope that I have not misunderstood something terribly wrong, but the continuous derivative $$D=d/dt$$ can be considered a transfer function in Laplace space $$D(s) = s$$, right?

So when I try to discretize it using the bilinear transform (Tustin's method) I trivially get

$$D(z) = \frac{2}{T} \frac{1-z^{-1}}{1+z^{-1}}$$

When I apply this to a series containing one discrete impulse, the response oscillates at the Nyquist frequency. Even worse, the spectrum around $$\omega=0$$ is quadratic and not $$\sim i\omega$$ like it would be expected from the derivative. (Edit: the latter was just due to roundoff error, because the low-freq amplitude got swamped by the Nyquist-peak)

Although I know of course how well the bilinear transform works for discretizing all kinds of filters, it is hard for me to understand, why it is considered superior if it seems to fail so miserably for something as simple as the derivative, which can otherwise be easily represented by first order finite differences:

$$D(z)=\frac{1-z^{-1}}{T}$$

or even second order (symmetric) finite differences

$$D(z)=\frac{1-z^{-2}}{2Tz^{-1}}$$

I am sure it all has a very simple explanation, but I can't see it.

PS: What is all the more confusing: when I apply the bilinear D(z) to a step function, the result is (correctly) a single peak. Consequently the inverse of the bilinear D(z) applied to an impulse yields the step function, like it has to be. What is going on there?

It's a matter of perspective if the bilinear transform "fails so miserably" when trying to approximate a derivative. First of all, it's not true that the approximation is quadratic for small $$\omega$$, it's linear as it should be:

$$D(e^{j\omega})=\frac{2}{T}\frac{1-e^{-j\omega}}{1+e^{-j\omega}}=\frac{2j}{T}\tan\left(\frac{\omega}{2}\right)\approx \frac{j\omega}{T}\quad\textrm{for}\quad |\omega|\ll 1$$

So for small frequencies the discrete-time system behaves as a (discrete-time) differentiator. For larger frequencies we expect the effect of frequency warping, which is exactly what happens here: the infinite continuous-time frequency range is mapped to the interval $$\omega\in[0,\pi]$$.

In sum, the bilinear transform behaves as expected, but it is also clear that there are better ways to approximate a derivative in discrete-time than applying the bilinear transform to the continuous-time transfer function $$H(s)=s$$.

• Okay, I will have to check the part with the low-frequency behavior again in my numerics. Strange... – oliver May 20 at 11:09
• The apparent hyperbolic behavior was simply due to me choosing an odd FFT size, which excludes the Nyquist-frequency from the spectrum. But with even FFT I still can't reproduce the expected low-frequency behavior ~iw. Maybe roundoff error, because the Nyquist-peak is so big against the rest of the spectrum? – oliver May 20 at 13:03
• Yeah, near w->0 my spectrum was somehow corrupted by finite precision. If I prefilter the impulse with a butterworth at Nyquist/2, then the low-frequency amplitude in the spectrum is proportional to iw. Sorry for the false alarm. – oliver May 20 at 13:12
• One last thing I would like to understand: what specific feature of the derivative operator allows me to tell before trying, whether it is one that is not so ideal for applying the bilinear transform to? Or, how could I have predicted the oscillating behavior (is that already called an instability?) of the impulse response? – oliver May 20 at 13:20
• @oliver: The filter is unstable, because it has a pole at $z=-1$. Poles in the region $|z|\ge 1$ make a (causal) filter unstable. Note that the continuous-time filter is also unstable. The bilinear transform preserves (in-)stability. – Matt L. May 20 at 13:23