The Z Transform

We will start by defining a transform function that at first looks so trivial that it appears that it cannot possibly be useful. Let us assume that data is digitized at a uniform sampling rate. Say our data consisted of the sequence of points:

$$x_t = ( 10, 5, 0, 4, -3, -6, 0, ... )$$ (1)

We can represent this data as a polynomial:

$$X(Z) = 10 + 5 Z + 0 Z^2 + 4 Z^3 - 3 Z^4 - 6 Z^5 + 0 Z^6 + ...$$ (2)

This polynomial is called the Z transform of the time series. The trick with Z transforms is how we interpret what the $Z$ actually represents. The $Z$ is called the unit delay operator because of the following, if we take (2) and multiply it by $Z$ we get,

$$Y(Z) = Z X(Z) = 10 Z + 5 Z^2 + 0 Z^3 + 4 Z^4 - 3 Z^5 - 6 Z^6 + 0 Z^7 + ...$$ (3)

which becomes the untransformed series:

$$y_t = ( 0, 10, 5, 0, 4, -3, -6, 0, ... )$$

which is just the original time series delayed by one unit of time (Be careful if you try to look the Z transform up in a book. In a variation of the pi throwing contest that I mentioned a few columns back, different authors use different sign conventions on the Z transform polynomial.).

This innocent looking transform and delay operator turn out to be very powerful. Consider the following, suppose that $x_t$ represents the response of the ground to an explosion (which we will call the impulse response). If another explosion occurred 10 time units later, then the signal that we would see is

$$y_t = ( 10, 5, 0, 4, -3, -6, 0, 0, 0, 0, 10, 5, 0, 4, -3, -6, 0, ... )$$

which can be represented by

$$Y(Z) = X(Z) + Z^{10} X(Z)$$

If the second explosion was half the strength of the first, then we'd write,

$$Y(Z) = X(Z) + \frac{1}{2} Z^{10} X(Z)$$

If the second was a half strength implosion, then the sign is just reversed,

$$Y(Z) = X(Z) - \frac{1}{2} Z^{10} X(Z)$$

If the second signal came in overlapping the first (say at time 4) then the numbers in the series is a bit messy looking,

$$y_t = ( 10, 5, 0, 4, -8, -8.5, 0, -2, 1.5, 3, 0, ... )$$

but rest assured this series can be written out as having the transform,

$$Y(Z) = X(Z) - \frac{1}{2} Z^4 X(Z)$$ (4)

$$= (1 - \frac{1}{2} Z^4) X(Z)$$ (5)

A Forth function to do this kind of factorization is shown in listing 1. This code has been contributed to the Forth Scientific Library and is awaiting review for its acceptance. (Volunteers for doing code reviews are eagerly accepted!).

In general we can do this for an arbitrarily complicated sequence of explosions and implosions as long as the response of the earth remains linear (an example of a nonlinear response would be the state of the ground within a mile or so of a nuclear explosion). So suppose we have a sequence consisting of an explosion at $t=0$, a half strength implosion at time $t=2$ and a quarter strength explosion at time $t=3$,

$$B(Z) = 1 - \frac{1}{2} Z^2 + \frac{1}{4} Z^3$$ (6)

then,

$$Y(Z) = B(Z) X(Z)$$ (7)

Now without doing any integrals, and just manipulating polynomials we have done the convolution theorem,

$$y_t = \int x(t - \tau) b(\tau) d\tau$$ (8)

Both equations (7) and (8) state that the signal $y_t$ is the convolution of $x_t$ and $b_t$, but the Z transform version just involves simple manipulation of polynomials. Listing 2 provides an implementation of convolution algorithms, this code is also an unreviewed contribution to the Forth Scientific Library.