function lucy

Performs a lucy deconvolution of X using PSF as dirty beam and N iterations. THRESHOLD defaults to 0.000003.

LUCY DECONVOLUTION ALGORITHM. "Lucy has shown that his iterative scheme is related to the maximum likelihood solution of the deconvolution problem and would converge to that solution after an infinite number of iterations."

   iteration number : k
   dirty beam : dbeam
   dirty map : dmap
   estimated object : o(k)
   reconvolved estimated
   object : oc(k)
   correctin function : t(k)
   multiplication : *
   convolution : (*)

start: o(0) = dmap

   iteration:
   oc(k) = o(k) (*) dbeam
   t(k) = (dmap/oc(k)) (*) dbeam
   o(k+1)= o(k) * t(k)

Due to calculation of the quotient for t(k) the algorithm is sensitive to noise. Therefore a weighting function has been implemented which is unity if the dmap signal is larger than "thresh" and zero elsewhere. t(k) is only calculated in that area and hence the deconvolution is only effective there. Threshold equals input if positive. If input is negative thresh is determined via mean and rms in the area below abs(thresh). If thresh is zero thresh is determined in the area below -0.1 times the maximum flux in dmap. The redetermined threshold is mean value plus 3 times rms.

To check on convergence of algorithm the rms between input dmap and current lucy deconvolved map reconvolved with dbeam is calculated.

Syntax
result = lucy(X, PSF, N [, THRESHOLD])

Arguments

X:   A matrix.

PSF:   A matrix.

N:   An integer number.

THRESHOLD:   A real number.

See also
function clean
function wien

Notes
The memory required for lucy deconvolution (in bytes) calculates as 44*(Number of pixels in the image).

X:	A matrix.
PSF:	A matrix.
N:	An integer number.
THRESHOLD:	A real number.