Why does this manual bilinear transformation yield different results from Matlab's? I've got a first-order Butterworth filter with the cutoff frequency $\omega_c$. Its transfer function is then

$$H(s) = \frac{\omega_c}{s+\omega_c}$$

Using the bilinear transform to find an $H(z)$ (what is that function called?), I get

$$H(z)=\frac{\omega_c}{\frac{2}{T}\frac{z-1}{z+1} + \omega_c}
= \frac{\omega_c z + \omega_c}{\left(\frac{2}{T}+\omega_c\right)z + \omega_c-\frac{2}{T}}$$

However, I can't reconcile this result with what Matlab is doing. It seems wrong, no matter what value of $T$. I assume that B and A below are the coefficients of $H(z)$.

&gt;&gt; [B,A] = butter(1,0.5)
B = 0.5000    0.5000
A = 1.0000   -0.0000
&gt;&gt; [B,A] = butter(1,0.6)
B = 0.5792    0.5792
A = 1.0000    0.1584
&gt;&gt; [B,A] = butter(1,0.7)
B = 0.6625    0.6625
A = 1.0000    0.3249
&gt;&gt; [B,A] = butter(1,0.8)
B = 0.7548    0.7548
A = 1.0000    0.5095


What am I misunderstanding?
A couple things:

Before making the substitution $s = \frac{2}{T} \frac{z-1}{z+1}$, you need to prewarp the cutoff frequency by making the substitution:

$$
\omega_{c,w} = \frac{2}{T} \tan (\omega_c\frac{T}{2})
$$

where $\omega_{c,w}$ is the warped cutoff frequency. This is necessary because the bilinear transform maps the left-half plane in the Laplace domain (used in analog filter design) to the unit circle in the $z$-domain in a nonlinear fashion. Therefore, as you approach the Nyquist rate (digital frequencies of $\pm \pi$), the approximation to the analog filter prototype becomes inaccurate.

Also, the second parameter that you are passing to the butter function is the normalized cutoff frequency, not the sampling interval $T$. The normalized frequency used by that function is in the interval $(0,1)$ and is equal to the ratio of the desired cutoff frequency to the Nyquist rate:

$$
\omega_n = \frac{\omega_c}{2 \pi \frac{f_s}{2}}
$$

$$
\omega_n = \frac{\omega_c}{\pi f_s}
$$

$$
\omega_n = \frac{\omega_c T}{\pi}
$$