P =8
    S  =int( log( P) /log( 2) +0.9999)

    Pi =3.14159265
    R1 =2^S

    R =R1 -1
    R2 =div( R1,  2)
    R4 =div( R1,  4)
    R3 =R4 +R2

    Dim Re( R1), Im( R1), Co( R3)

    for N =0 to P -1
        read dummy: Re( N) =dummy
        read dummy: Im( N) =dummy
    next N

    data    1, 0,      1, 0,      1, 0,      1, 0,      0, 0,     0, 0,      0, 0,       0, 0

    S2 =div( S, 2)
    S1 =S -S2
    P1 =2^S1
    P2 =2^S2

    dim V( P1 -1)
    V( 0) =0
    DV =1
    DP =P1

    for J =1 to S1
        HA =div( DP, 2)
        PT =P1 -HA
        for I =HA to PT step DP
            V( I) =V( I -HA) +DV
        next I
        DV =DV +DV
        DP =HA
    next J

    K =2 *Pi /R1

    for X =0 to R4
        COX =cos( K *X)
        Co( X) =COX
        Co( R2 -X) =0 -COX
        Co( R2 +X) =0 -COX
    next X

    print "FFT: bit reversal"

    for I =0 to P1 -1
        IP =I *P2
        for J =0 to P2 -1
            H =IP +J
            G =V( J) *P2 +V( I)
            if G >H then temp =Re( G): Re( G) =Re( H): Re( H) =temp
            if G >H then temp =Im( G): Im( G) =Im( H): Im( H) =temp
        next J
    next I

    T =1

    for stage =0 to S -1
        print "  Stage:- "; stage
        D =div( R2, T)
        for Z =0 to T -1
            L   =D *Z
            LS  =L +R4
            for I =0 to D -1
                A      =2 *I *T +Z
                B      =A +T
                F1     =Re( A)
                F2     =Im( A)
                P1     =Co( L)  *Re( B)
                P2     =Co( LS) *Im( B)
                P3     =Co( LS) *Re( B)
                P4     =Co( L)  *Im( B)
                Re( A) =F1 +P1 -P2
                Im( A) =F2 +P3 +P4
                Re( B) =F1 -P1 +P2
                Im( B) =F2 -P3 -P4
            next I
        next Z
        T =T +T
    next stage

    print "   M          Re( M)       Im( M)"

    for M =0 to R
        if abs( Re( M)) <10^-5 then Re( M) =0
        if abs( Im( M)) <10^-5 then Im( M) =0
        print "   "; M, Re( M), Im( M)
    next M

    end


    wait

    function div( a, b)
        div =int( a /b)
    end function

    end