Sub Cam_rec()
    'Cameron's recursive algorithm for finding F(w) (first published in ESA Journal 1982)
    'Copyright D. Pelz 1/2011
    Dim a(20) As Single, b(20) As Single
    Dim pa(20) As Single, pb(20) As Single, pc(20) As Single, pd(20) As Single
    Dim qa(20) As Single, qb(20) As Single, qc(20) As Single, qd(20) As Single
    
    n = Val(Filpar1.n_res_.Text)
    f1 = Val(Filpar1.f_low_.Text): f2 = Val(Filpar1.f_high_.Text)
    fc = Sqr(f1 * f2)
    BW = f2 - f1
    om = freq * (1 - f1 * f2 / freq ^ 2) / BW
    RL = Val(Filpar1.pbrl_.Text)
    eps = 1 / Sqr(10 ^ (RL / 10) - 1)
    
    ' establish nftz and read FTZs
    nftz = 0
    If Filpar2.ftz1_.Text <> "" Then
        nftz = nftz + 1
        FTZ(nftz) = Val(Filpar2.ftz1_.Text)
    End If
    If Filpar2.ftz2_.Text <> "" Then
        nftz = nftz + 1
        FTZ(nftz) = Val(Filpar2.ftz2_.Text)
    End If
    If Filpar2.ftz3_.Text <> "" Then
        nftz = nftz + 1
        FTZ(nftz) = Val(Filpar2.ftz3_.Text)
    End If
        
    'set TZs at inf
    For i = nftz + 1 To n
        psi(i) = 10 ^ 7
    Next i
    
    'set normalised FTZs
    For i = 1 To nftz
        psi(i) = FTZ(i) * (1 - f1 * f2 / FTZ(i) ^ 2) / BW
    Next i
    'For i = 1 To n: Debug.Print n, nftz, i, psi(i): Next i
    
    
    
    
    For i = 1 To n: a(i) = 1 / psi(i): b(i) = Sqr(1 - 1 / psi(i) ^ 2): Next i
    
    
    p(0) = -a(1): p(1) = 1                          'set initial polynomial P1
    q(0) = b(1)                                     'dito for Q1
    Np = 1: Nq = 0                                  'initial degrees of p and q-polynomial
    
    For i = 2 To n                                  'major iteration loop
        For j = 1 To Np + 1
            pa(j) = p(j - 1)                        'multiply p by omega = shift coeffs. by one
        Next j
        pa(0) = 0                                   'destroy coeff.(0)
        
        For j = 0 To Np
            pb(j) = -p(j) * a(i)                    'multiply p by -a(i)
        Next j
        
        For j = 2 To Nq + 2
            pc(j) = q(j - 2) * b(i)                 ' multiply p by omega^2 = shift coeffs. by 2
        Next j                                      ' lowest coeff. of pc is now 2
        pc(0) = 0: pc(1) = 0                        ' destroy coeff.(0) and coeff.(1)
        
        For j = 0 To Nq
            pd(j) = -q(j) * b(i)                    ' multiply q by -b(i)
        Next j
        
        '-----------------------------------------------------------------------------------------
        
        For j = 1 To Nq + 1
            qa(j) = q(j - 1)                        'multiply q by omega = shift coeffs. by one
        Next j
        pa(0) = 0                                   'destroy coeff.(0)
        
        
        For j = 0 To Nq
            qb(j) = -q(j) * a(i)                    'multiply q by -a(i)
        Next j
        
        For j = 0 To Np
            qc(j) = p(j) * b(i)                     ' multiply p by omega^2 = shift coeffs. by 2
        Next j                                      ' lowest coeff. of pc is now 2
        
        
        Np = Np + 1                                 ' set new degree of p, q
        Nq = Nq + 1
        
        For j = 0 To Np                             ' build new polynomial p
            p(j) = pa(j) + pb(j) + pc(j) + pd(j)
        Next j
        
        For j = 0 To Nq                             ' build new polynomial q
            q(j) = qa(j) + qb(j) + qc(j)
        Next j
        
    Next i
    
    For i = 0 To n
        p(i) = pa(i) + pb(i) + pc(i) + pd(i)        ' build final p, q from its components
        q(i) = qa(i) + qb(i) + qc(i)
    Next i
    
    'Debug.Print p(0), p(1), p(2), p(3), p(4)
    ' ran 1st time without errors/bugs ! very sus.... :)
    poly = 1
End Sub