def bilinear(b, a, fs=<NUM_LIT:1.0>):
fs = float(fs)<EOL>a, b = map(atleast_1d, (a, b))<EOL>D = len(a) - <NUM_LIT:1><EOL>N = len(b) - <NUM_LIT:1><EOL>artype = float<EOL>M = max([N, D])<EOL>Np = M<EOL>Dp = M<EOL>bprime = numpy.zeros(Np + <NUM_LIT:1>, artype)<EOL>aprime = numpy.zeros(Dp + <NUM_LIT:1>, artype)<EOL>for j in range(Np + <NUM_LIT:1>):<EOL><INDENT>val = <NUM_LIT:0.0><EOL>for i in range(N + <NUM_LIT:1>):<EOL><INDENT>for k in range(i + <NUM_LIT:1>):<EOL><INDENT>for l in range(M - i + <NUM_LIT:1>):<EOL><INDENT>if k + l == j:<EOL><INDENT>val += (comb(i, k) * comb(M - i, l) * b[N - i] *<EOL>pow(<NUM_LIT:2> * fs, i) * (-<NUM_LIT:1>) ** k)<EOL><DEDENT><DEDENT><DEDENT><DEDENT>bprime[j] = real(val)<EOL><DEDENT>for j in range(Dp + <NUM_LIT:1>):<EOL><INDENT>val = <NUM_LIT:0.0><EOL>for i in range(D + <NUM_LIT:1>):<EOL><INDENT>for k in range(i + <NUM_LIT:1>):<EOL><INDENT>for l in range(M - i + <NUM_LIT:1>):<EOL><INDENT>if k + l == j:<EOL><INDENT>val += (comb(i, k) * comb(M - i, l) * a[D - i] *<EOL>pow(<NUM_LIT:2> * fs, i) * (-<NUM_LIT:1>) ** k)<EOL><DEDENT><DEDENT><DEDENT><DEDENT>aprime[j] = real(val)<EOL><DEDENT>return normalize(bprime, aprime)<EOL>
Return a digital filter from an analog one using a bilinear transform.

    The bilinear transform substitutes ``(z-1) / (z+1)`` for ``s``.