Version 1.9 Build 1367

News	FAQ
Search	Home

Next: The inverse problem Up: Recommendations for the AIPS++ Imaging Model Previous: Introduction

A linear algebra model for Imaging

Many instruments are linear and can be represented by an operator A operating on an input sky I, to produce an data set D. If we use pixellated images then both I and D can be represented by vectors and the operator A can be represented by a matrix. We then have the matrix equation:

AI = D

(1)

The linearity means that images I^a and I^b which separately give rise to data D^a and D^b, together give rise to data D^a + D^b. The measurement matrix A is nearly always non-square and singular. Examples of A and D are as follows:

Convolution

A represents convolution by a shift-invariant Point Spread Function. Optical Imaging and the standard interferometry convolution equation differ in the statistics of the noise: in the former, the noise is independent from pixel to pixel, whereas in the latter the noise is correlated with the shape of the dirty beam.

A_{i, j} = B( $\displaystyle \bf x_{i}^{}$ - $\displaystyle \bf x_{j}^{}$ )

(2)

In both cases, A is a Toeplitz matrix (i.e. matrix element A_{i, j} is a function of $\bf x_{i}^{}$ - $\bf x_{j}^{}$ only) and D is a vector representing the image.

Interferometer

A is a matrix performing a Fourier transform, and D represents the visibility data. If the u, v plane vector is $\bf u$ and the image plane vector is $\bf x$ , then the i, j element of A is given by:

A_{i, j} = e^j2.

(3)

(actually one should separate into real and imaginary components but that is not important here).

Total Power

A represents convolution of the brightness by the primary beam B, and D is the vector of total power samples:

A_{i, j} = B( $\displaystyle \bf x_{i}^{}$ , $\displaystyle \bf x_{j}^{}$ )

(4)

Beam Switched Total Power

A is PSF(on)-PSF(off), and D is the vector of beam switched total power samples (Emerson et al., 1979):

A_{i, j} = B( $\displaystyle \bf x_{i}^{}$ - $\displaystyle \bf x_{j}^{}$ ) - B( $\displaystyle \bf x_{i}^{}$ - $\displaystyle \bf x_{j}^{}$ - $\displaystyle \Delta$ $\displaystyle \bf x$ )

(5)

where $\Delta$ $\bf x$ is the beam throw.

Mosaic

A is a matrix representing multiplication by a primary beam centered at a specific fixed pointing position $\bf x^{p}_{}$ , followed by Fourier transformation:

A_{i, j} = B( $\displaystyle \bf x_{j}^{}$ - $\displaystyle \bf x^{p}_{}$ )e^j2.

(6)

Mosaic with variable pointing

A is a matrix representing multiplication by a primary beam centered at a variable pointing position $\bf x^{p}_{i}$ , followed by Fourier transformation:

A_{i, j} = B( $\displaystyle \bf x_{j}^{}$ - $\displaystyle \bf x^{p}_{i}$ )e^j2.

(7)

Wide field transform

A is a matrix like that for the interferometer but including the w phase term.

Primary beam tapering

is a strange case which is throws up some interesting points. The idea is to represent only the primary beam tapering of an interferometric array. The synthesized beam is ignored and so D is actually a collection of dirty images. We return to this case below in the discussion of linear mosaic.

Overall, this is a clean way of writing down the operation of a Telescope, but it does look rather unfamiliar in some aspects. The key insight which connects it to most of the algorithms described in the literature is that in practice we hardly ever actually use linear algebra to calculate AI. Instead we use special symmetries of A to make shortcuts. For example, multiplication by a Toeplitz matrix is best performed using the circulant approximation relying upon zero-padding to twice the number of pixels on each axis, following by FFT convolution (see Andrews and Hunts, 1977 for a description of this trick). Similarly, we approximate AI for the interferometer by the familiar degridding step, and we use either the 3D or the polyhedron approach (Cornwell and Perley, 1992) for wide-field imaging. It is also fruitful to consider more complicated algorithms as using various approximations for A. For example, the Clark CLEAN algorithm (Clark 1980) can be regarded as the result of alternately using a sparse and a circulant approximation for the Toeplitz dirty beam matrix (more on this below).

This formulation becomes even more useful when applied to the inverse problem of determining the sky brightness from the data. I discuss this next.

Next: The inverse problem Up: Recommendations for the AIPS++ Imaging Model Previous: Introduction