I have omitted a bit of theory here that is found in PVI.

For an FT matrix of arbitrary size *N*,
the desired shift is *N*/2, so values at alternate points
in the time axis are multiplied by -1.
A subroutine for that purpose is `fth()`.

# FT a vector in a matrix, with first omega = - pi # subroutine fth( adj,sign, m1, n12, cx) integer i, adj, m1, n12 real sign complex cx(m1,n12) temporary complex temp(n12) do i= 1, n12 temp(i) = cx(1,i) if( adj == 0) { do i= 2, n12, 2 temp(i) = -temp(i) call ftu( sign, n12, temp) } else { call ftu( -sign, n12, temp) do i= 2, n12, 2 temp(i) = -temp(i) } do i= 1, n12 cx(1,i) = temp(i) return; end

To Fourier transform a 1024-point complex vector `cx(1024)`
and then inverse transform it, you would

call fth( 0, 1., 1, 1024, cx) call fth( 1, 1., 1, 1024, cx)

You might wonder about the apparent redundancy of using both
the argument `conj` and the argument `sign`.
Having two arguments instead of one allows
us to define the *forward* transform for a *time* axis
with the opposite sign as
the forward transform for a *space* axis.

The subroutine `fth()` is somewhat cluttered by
the inclusion of a
frequently needed practical feature--namely,
the facility to extract vectors from a matrix,
transform the vectors, and then restore them into the matrix.

10/31/1997