# length 27 (-13,...,13) linear phase FIR filter param n := 14; param pi := 4*atan(1); param eps := 1.0e-4; param sinc_w {k in -2*n+2..2*n-2} := if k != 0 then (sin(2*pi*k*0.9)/(pi*k)) else (0.8); param sinc_m {k in -2*n+2..2*n-2} := if k != 0 then ( sin(2*pi*k*0.1)/(pi*k)+ sin(2*pi*k*0.8)/(pi*k)) else (0.8); param sinc_t {k in -2*n+2..2*n-2} := if k != 0 then (sin(2*pi*k*0.2)/(pi*k)) else (0.4); var rho >= 0; var hw {0..n-1}; var hm {0..n-1}; var ht {0..n-1}; minimize power_bnd: rho; subject to passband: ((hw[0]+hm[0]+ht[0]-1)^2 + 2*sum {k in 1..n-1} (hw[k]+hm[k]+ht[k])^2) <= (eps); subject to wooferband: ( sum {k in -n+1..n-1, kk in -n+1..n-1} hw[abs(k)] * hw[abs(kk)] * sinc_w[k+kk] ) /rho <= (0.6)*rho; subject to midrangeband: ( sum {k in -n+1..n-1, kk in -n+1..n-1} hm[abs(k)]*hm[abs(kk)] * sinc_m[k+kk] ) /rho <= (0.4)*rho; subject to tweeterband: ( sum {k in -n+1..n-1, kk in -n+1..n-1} ht[abs(k)]*ht[abs(kk)] * sinc_t[k+kk] ) /rho <= (0.6)*rho; let hw[0] := 2; let hm[0] := 2; let ht[0] := 2; let rho := 1; solve; printf "%3d \n", card({-n+1..n-1}) > hw; printf "%3d \n", card({-n+1..n-1}) > hm; printf "%3d \n", card({-n+1..n-1}) > ht; printf "%3d \n", card({-n+1..n-1}) > hwmt; printf {k in -n+1..n-1}: "%10.6f \n", hw[abs(k)] > hw; printf {k in -n+1..n-1}: "%10.6f \n", hm[abs(k)] > hm; printf {k in -n+1..n-1}: "%10.6f \n", ht[abs(k)] > ht; printf {k in -n+1..n-1}: "%10.6f \n", hw[abs(k)]+hm[abs(k)]+ht[abs(k)] > hwmt; printf {nu in 0..1 by 1/1000}: "%7.4f %10.3e \n", nu, 20*log10(abs( hw[0] + 2* sum {k in 1..n-1} (hw[k]*cos(-2*pi*k*nu)) )) > 3aw.out; printf {nu in 0..1 by 1/1000}: "%7.4f %10.3e \n", nu, 20*log10(abs( hm[0] + 2* sum {k in 1..n-1} (hm[k]*cos(-2*pi*k*nu)) )) > 3am.out; printf {nu in 0..1 by 1/1000}: "%7.4f %10.3e \n", nu, 20*log10(abs( ht[0] + 2* sum {k in 1..n-1} (ht[k]*cos(-2*pi*k*nu)) )) > 3at.out;