Problem Statement

Digital filters are used to apply a mathematical operation to a sampled signal. One of the common formulations is the "direct form II transposed" which can represent both infinite impulse response (IIR) and finite impulse response (FIR) filters, as well as being more numerically stable than other forms. [1] Filter a signal using an order 3 low-pass Butterworth filter. The coefficients for the filter are a=[1.00000000, -2.77555756e-16, 3.33333333e-01, -1.85037171e-17] and b = [0.16666667, 0.5, 0.5, 0.16666667] The signal that needs filtering is the following vector: [-0.917843918645, 0.141984778794, 1.20536903482, 0.190286794412, -0.662370894973, -1.00700480494, -0.404707073677 ,0.800482325044, 0.743500089861, 1.01090520172, 0.741527555207, 0.277841675195, 0.400833448236, -0.2085993586, -0.172842103641, -0.134316096293, 0.0259303398477, 0.490105989562, 0.549391221511, 0.9047198589] [Wikipedia on Butterworth filters]

Let's start with the solution:

include ..\Utilitys.pmt

( 1.00000000  -2.77555756e-16  3.33333333e-01  -1.85037171e-17 ) var a
( 0.16666667   0.5             0.5              0.16666667 ) var b
( -0.917843918645   0.141984778794   1.20536903482    0.190286794412 -0.662370894973 
  -1.00700480494   -0.404707073677   0.800482325044   0.743500089861  1.01090520172 
   0.741527555207   0.277841675195   0.400833448236  -0.2085993586   -0.172842103641 
  -0.134316096293   0.0259303398477  0.490105989562   0.549391221511  0.9047198589 ) 
  
len dup 0 swap repeat >ps

for var i
    0 >ps
    b len i min for var j
        j get rot i j - 1 + get rot * ps> + >ps swap
    endfor
    drop
    a len i min for var j
        j get ps> tps i j - 1 + get nip rot * - >ps
    endfor
    drop
    ps> a 1 get nip / ps> swap i set >ps
endfor
drop ps> ?