  The CLEAN algorithm, widely used in radio interferometry for the
deconvolution of radio images, performs well only if the raw radio image (dirty
image) is, to good approximation, a simple convolution between the instrumental
point-spread function (dirty beam) and the true distribution of emission across
the sky. An important case in which this approximation breaks down is during
frequency synthesis if the observing bandwidth is wide enough for variations in
the spectrum of the sky to become significant. The convolution assumption also
breaks down, in any situation but snapshot observations, if sources in the
field vary significantly in flux density over the duration of the observation.
Such time-variation can even be instrumental in nature, for example due to
jitter or rotation of the primary beam pattern on the sky during an
observation. An algorithm already exists for dealing with the spectral
variation encountered in wide-band frequency synthesis interferometry. This
algorithm is an extension of CLEAN in which, at each iteration, a set of N
`dirty beams' are fitted and subtracted in parallel, instead of just a single
dirty beam as in standard CLEAN. In the wide-band algorithm the beams are
obtained by expanding a nominal source spectrum in a Taylor series, each term
of the series generating one of the beams. In the present paper this algorithm
is extended to images which contain sources which vary over both frequency and
time. Different expansion schemes (or bases) on the time and frequency axes are
compared, and issues such as Gibbs ringing and non-orthogonality are discussed.
It is shown that practical considerations make it often desirable to
orthogonalize the set of beams before commencing the cleaning. This is easily
accomplished via a Gram-Schmidt technique.
