void eq_magplot(unsigned x_res, unsigned long from, unsigned long to,
                long *coefs, long *graph, unsigned intbits)
{
    unsigned long w = from;
    unsigned i;
    unsigned long delta = (to - from)/x_res;
    const long one = 1 << (32 - 2*intbits + 1);
    const long b0 = coefs[0], b1 = coefs[1], b2 = coefs[2];
    const long a1 = coefs[3], a2 = coefs[4];
    const long num1 = FRACMUL(b0, b0) + FRACMUL(b1, b1) + FRACMUL(b2, b2);
    const long num2 = 2*(FRACMUL(b0, b1) + FRACMUL(b1, b2));
    const long num3 = 2*FRACMUL(b0, b2);
    const long den1 = one + FRACMUL(a1, a1) + FRACMUL(a2, a2);
    const long den2 = 2*(-(a1 >> (intbits - 1)) + FRACMUL(a1, a2));
    const long den3 = -2*(a2 >> (intbits - 1));

    for (i = 0; i < x_res; i++, w += delta) {
        long cosw = fcos(w), cos2w = fcos(2*w);
        long num = num1 + FRACMUL(num2, cosw) + FRACMUL(num3, cos2w);
        long den = den1 + FRACMUL(den2, cosw) + FRACMUL(den3, cos2w);
        //printf("%f %f\n", MAKE_FLOAT(num, 25), MAKE_FLOAT(den, 25));
        long div = fsqrt(DIV64(num, den, 23/*4*8 - intbits*2 + 1*/), 23);
        printf("%d %d %f\n", num, den, MAKE_FLOAT(div, 23));
        graph[i] = div;
    }
}