casacore
LatticeCleaner.h
Go to the documentation of this file.
1//# Cleaner.h: this defines Cleaner a class for doing convolution
2//# Copyright (C) 1996,1997,1998,1999,2000,2001,2002,2003
3//# Associated Universities, Inc. Washington DC, USA.
4//#
5//# This library is free software; you can redistribute it and/or modify it
6//# under the terms of the GNU Library General Public License as published by
7//# the Free Software Foundation; either version 2 of the License, or (at your
8//# option) any later version.
9//#
10//# This library is distributed in the hope that it will be useful, but WITHOUT
11//# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12//# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
13//# License for more details.
14//#
15//# You should have received a copy of the GNU Library General Public License
16//# along with this library; if not, write to the Free Software Foundation,
17//# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
18//#
19//# Correspondence concerning AIPS++ should be addressed as follows:
20//# Internet email: aips2-request@nrao.edu.
21//# Postal address: AIPS++ Project Office
22//# National Radio Astronomy Observatory
23//# 520 Edgemont Road
24//# Charlottesville, VA 22903-2475 USA
25//#
26//#
27//# $Id$
28
29#ifndef LATTICES_LATTICECLEANER_H
30#define LATTICES_LATTICECLEANER_H
31
32//# Includes
33#include <casacore/casa/aips.h>
34#include <casacore/casa/Quanta/Quantum.h>
35#include <casacore/lattices/Lattices/TempLattice.h>
36#include <casacore/casa/Arrays/IPosition.h>
37#include <casacore/casa/Arrays/Vector.h>
38#include <casacore/casa/Containers/Block.h>
39
40namespace casacore { //# NAMESPACE CASACORE - BEGIN
41
42//# Forward Declarations
43class LatticeCleanProgress;
44template <class T> class TempLattice;
45
46// <summary>Lists the different types of Convolutions that can be done</summary>
47// <synopsis>This enumerator is brought out as a separate class because g++
48// currently cannot handle enumerators in a templated class. When it can this
49// class will go away and this enumerator moved into the Cleaner
50// class</synopsis>
51class CleanEnums {
52public:
53 enum CleanType {
54 // Hogbom
55 HOGBOM,
56 // Multi-scale
58 // Clark
59 CLARK
60 };
61};
62
63// <summary>A class for doing multi-dimensional cleaning</summary>
64
65// <use visibility=export>
66
67// <reviewed reviewer="" date="yyyy/mm/dd" tests="tLatticeCleaner">
68// </reviewed>
69
70// <prerequisite>
71// <li> The mathematical concept of deconvolution
72// </prerequisite>
73//
74// <etymology>
75
76// The LatticeCleaner class will deconvolve Lattices.
77
78// </etymology>
79//
80// <synopsis>
81// This class will perform various types of Clean deconvolution
82// on Lattices.
83//
84// </synopsis>
85//
86// <example>
87// <srcblock>
88// </srcblock>
89// </example>
90//
91// <motivation>
92// </motivation>
93//
94// <thrown>
95// <li> AipsError: if psf has more dimensions than the model.
96// </thrown>
97//
98// <todo asof="yyyy/mm/dd">
99// <li> Allow the psf to be specified with a
100// <linkto class=Function>Function</linkto>.
101// </todo>
102
103template<class T> class LatticeCleaner
104{
105public:
106
107 // Create a cleaner : default constructor
109
110 // Create a cleaner for a specific dirty image and PSF
111 LatticeCleaner(const Lattice<T> & psf, const Lattice<T> & dirty);
112
113 // The copy constructor uses reference semantics
115
116 // The assignment operator also uses reference semantics
118
119 // The destructor does nothing special.
121
122 // Update the dirty image only
123 void update(const Lattice<T> & dirty);
124
125 // Set a number of scale sizes. The units of the scale are pixels.
126 Bool setscales(const Int nscales, const Float scaleInc=1.0);
127
128 // Set a specific set of scales
129 Bool setscales(const Vector<Float> & scales);
130
131 // Set up control parameters
132 // cleanType - type of the cleaning algorithm to use (HOGBOM, MULTISCALE)
133 // niter - number of iterations
134 // gain - loop gain used in cleaning (a fraction of the maximum
135 // subtracted at every iteration)
136 // aThreshold - absolute threshold to stop iterations
137 // fThreshold - fractional threshold (i.e. given w.r.t. maximum residual)
138 // to stop iterations. This parameter is specified as
139 // Quantity so it can be given in per cents.
140 // choose - unused at the moment, specify False. Original meaning is
141 // to allow interactive decision on whether to continue iterations.
142 // This method always returns True.
143 Bool setcontrol(CleanEnums::CleanType cleanType, const Int niter,
144 const Float gain, const Quantity& aThreshold,
145 const Quantity& fThreshold,
146 const Bool choose=True);
147
148 // This version of the method disables stopping on fractional threshold
149 Bool setcontrol(CleanEnums::CleanType cleanType, const Int niter,
150 const Float gain, const Quantity& threshold,
151 const Bool choose=True);
152
153 // return how many iterations we did do
154 Int iteration() const { return itsIteration; }
156
157 // what iteration number to start on
158 void startingIteration(const Int starting = 0) {itsStartingIter = starting; }
159
160 // Clean an image.
161 //return value gives you a hint of what's happening
162 // 1 = converged
163 // 0 = not converged but behaving normally
164 // -1 = not converged and stopped on cleaning consecutive smallest scale
165 // -2 = not converged and either large scale hit negative or diverging
166 // -3 = clean is diverging rather than converging
168
169 // Set the mask
170 // mask - input mask lattice
171 // maskThreshold - if positive, the value is treated as a threshold value to determine
172 // whether a pixel is good (mask value is greater than the threshold) or has to be
173 // masked (mask value is below the threshold). Negative threshold switches mask clipping
174 // off. The mask value is used to weight the flux during cleaning. This mode is used
175 // to implement cleaning based on the signal-to-noise as opposed to the standard cleaning
176 // based on the flux. The default threshold value is 0.9, which ensures the behavior of the
177 // code is exactly the same as before this parameter has been introduced.
178 void setMask(Lattice<T> & mask, const T& maskThreshold = T(0.9));
179
180 // Tell the algorithm to NOT clean just the inner quarter
181 // (This is useful when multiscale clean is being used
182 // inside a major cycle for MF or WF algorithms)
183 // if True, the full image deconvolution will be attempted
185
186 // Consider the case of a point source:
187 // the flux on all scales is the same, and the first scale will be chosen.
188 // Now, consider the case of a point source with a *little* bit of extended structure:
189 // thats right, the largest scale will be chosen. In this case, we should provide some
190 // bias towards the small scales, or against the large scales. We do this in
191 // an ad hoc manner, multiplying the maxima found at each scale by
192 // 1.0 - itsSmallScaleBias * itsScaleSizes(scale)/itsScaleSizes(nScalesToClean-1);
193 // Typical bias values range from 0.2 to 1.0.
194 void setSmallScaleBias(const Float x=0.5) { itsSmallScaleBias = x; }
195
196 // During early iterations of a cycled MS Clean in mosaicing, it common
197 // to come across an ocsilatory pattern going between positive and
198 // negative in the large scale. If this is set, we stop at the first
199 // negative in the largest scale.
201
202 // Some algorithms require that the cycles be terminated when the image
203 // is dominated by point sources; if we get nStopPointMode of the
204 // smallest scale components in a row, we terminate the cycles
205 void stopPointMode(Int nStopPointMode) {itsStopPointMode = nStopPointMode; }
206
207 // After completion of cycle, querry this to find out if we stopped because
208 // of stopPointMode
210
211 // speedup() will speed the clean iteration by raising the
212 // threshold. This may be required if the threshold is
213 // accidentally set too low (ie, lower than can be achieved
214 // given errors in the approximate PSF).
215 //
216 // threshold(iteration) = threshold(0)
217 // * ( exp( (iteration - startingiteration)/Ndouble )/ 2.718 )
218 // If speedup() is NOT invoked, no effect on threshold
219 void speedup(const Float Ndouble);
220
221 // Look at what WE think the residuals look like
222 // Assumes the first scale is zero-sized
224
225 // Method to return threshold, including any speedup factors
227
228 // Method to return the strength optimum achieved at the last clean iteration
229 // The output of this method makes sense only if it is called after clean
231
232 // Helper function to optimize adding
233 static void addTo(Lattice<T>& to, const Lattice<T>& add);
234
235protected:
236 // Make sure that the peak of the Psf is within the image
238
239 // Make an lattice of the specified scale
240 void makeScale(Lattice<T>& scale, const Float& scaleSize);
241
242 // Make Spheroidal function for scale images
244
245 // Find the Peak of the Lattice
246 static Bool findMaxAbsLattice(const Lattice<T>& lattice,
247 T& maxAbs, IPosition& posMax);
248
249 // Find the Peak of the lattice, applying a mask
251 T& maxAbs, IPosition& posMax);
252
253 // Helper function to reduce the box sizes until the have the same
254 // size keeping the centers intact
255 static void makeBoxesSameSize(IPosition& blc1, IPosition& trc1,
256 IPosition &blc2, IPosition& trc2);
257
258
261 Int itsMaxNiter; // maximum possible number of iterations
265private:
266
267 //# The following functions are used in various places in the code and are
268 //# documented in the .cc file. Static functions are used when the functions
269 //# do not modify the object state. They ensure that implicit assumptions
270 //# about the current state and implicit side-effects are not possible
271 //# because all information must be supplied in the input arguments
272
273
276
279
285
287
288 Int itsIteration; // what iteration did we get to?
289 Int itsStartingIter; // what iteration did we get to?
291
294
295
298
299 // Memory to be allocated per TempLattice
301
302 // Let the user choose whether to stop
304
305 // Threshold speedup factors:
306 Bool itsDoSpeedup; // if false, threshold does not change with iteration
308
309 //# Stop now?
310 //#// Bool stopnow(); Removed on 8-Apr-2004 by GvD
311
312 // Calculate index into PsfConvScales
313 Int index(const Int scale, const Int otherscale);
314
317
318
326
327 // threshold for masks. If negative, mask values are used as weights and no pixels are
328 // discarded (although effectively they would be discarded if the mask value is 0.)
330
331};
332
333} //# NAMESPACE CASACORE - END
334
335#ifndef CASACORE_NO_AUTO_TEMPLATES
336#include <casacore/lattices/LatticeMath/LatticeCleaner.tcc>
337#endif //# CASACORE_NO_AUTO_TEMPLATES
338#endif
@ MULTISCALE
Multi-scale.
A class for doing multi-dimensional cleaning.
Bool findMaxAbsMaskLattice(const Lattice< T > &lattice, const Lattice< T > &mask, T &maxAbs, IPosition &posMax)
Find the Peak of the lattice, applying a mask.
PtrBlock< TempLattice< T > * > itsScales
T itsMaskThreshold
threshold for masks.
Vector< Float > itsScaleSizes
~LatticeCleaner()
The destructor does nothing special.
Bool itsChoose
Let the user choose whether to stop.
void stopPointMode(Int nStopPointMode)
Some algorithms require that the cycles be terminated when the image is dominated by point sources; i...
Int clean(Lattice< T > &model, LatticeCleanProgress *progress=0)
Clean an image.
Bool setscales(const Vector< Float > &scales)
Set a specific set of scales.
LatticeCleaner(const Lattice< T > &psf, const Lattice< T > &dirty)
Create a cleaner for a specific dirty image and PSF.
Double itsMemoryMB
Memory to be allocated per TempLattice.
Bool setcontrol(CleanEnums::CleanType cleanType, const Int niter, const Float gain, const Quantity &aThreshold, const Quantity &fThreshold, const Bool choose=True)
Set up control parameters cleanType - type of the cleaning algorithm to use (HOGBOM,...
static void addTo(Lattice< T > &to, const Lattice< T > &add)
Helper function to optimize adding.
void makeScale(Lattice< T > &scale, const Float &scaleSize)
Make an lattice of the specified scale.
PtrBlock< TempLattice< T > * > itsDirtyConvScales
void setMask(Lattice< T > &mask, const T &maskThreshold=T(0.9))
Set the mask mask - input mask lattice maskThreshold - if positive, the value is treated as a thresho...
Quantum< Double > itsFracThreshold
Bool setcontrol(CleanEnums::CleanType cleanType, const Int niter, const Float gain, const Quantity &threshold, const Bool choose=True)
This version of the method disables stopping on fractional threshold.
static Bool findMaxAbsLattice(const Lattice< T > &lattice, T &maxAbs, IPosition &posMax)
Find the Peak of the Lattice.
Bool setscales(const Int nscales, const Float scaleInc=1.0)
Set a number of scale sizes.
Bool validatePsf(const Lattice< T > &psf)
Make sure that the peak of the Psf is within the image.
void speedup(const Float Ndouble)
speedup() will speed the clean iteration by raising the threshold.
PtrBlock< TempLattice< Complex > * > itsScaleXfrs
void ignoreCenterBox(Bool huh)
Tell the algorithm to NOT clean just the inner quarter (This is useful when multiscale clean is being...
Float threshold() const
Method to return threshold, including any speedup factors.
LatticeCleaner< T > & operator=(const LatticeCleaner< T > &other)
The assignment operator also uses reference semantics.
void startingIteration(const Int starting=0)
what iteration number to start on
TempLattice< T > * itsMask
void setSmallScaleBias(const Float x=0.5)
Consider the case of a point source: the flux on all scales is the same, and the first scale will be ...
LatticeCleaner()
Create a cleaner : default constructor.
CleanEnums::CleanType itsCleanType
void update(const Lattice< T > &dirty)
Update the dirty image only.
Lattice< T > * residual()
Look at what WE think the residuals look like Assumes the first scale is zero-sized.
Bool itsDoSpeedup
Threshold speedup factors:
Vector< Float > itsTotalFluxScale
T strengthOptimum() const
Method to return the strength optimum achieved at the last clean iteration The output of this method ...
Int index(const Int scale, const Int otherscale)
Calculate index into PsfConvScales.
Bool queryStopPointMode() const
After completion of cycle, querry this to find out if we stopped because of stopPointMode.
static void makeBoxesSameSize(IPosition &blc1, IPosition &trc1, IPosition &blc2, IPosition &trc2)
Helper function to reduce the box sizes until the have the same size keeping the centers intact
TempLattice< Complex > * itsXfr
Float spheroidal(Float nu)
Make Spheroidal function for scale images.
Int iteration() const
return how many iterations we did do
Quantum< Double > itsThreshold
TempLattice< T > * itsDirty
PtrBlock< TempLattice< T > * > itsPsfConvScales
PtrBlock< TempLattice< T > * > itsScaleMasks
LatticeCleaner(const LatticeCleaner< T > &other)
The copy constructor uses reference semantics.
void stopAtLargeScaleNegative()
During early iterations of a cycled MS Clean in mosaicing, it common to come across an ocsilatory pat...
A drop-in replacement for Block<T*>.
Definition: Block.h:814
this file contains all the compiler specific defines
Definition: mainpage.dox:28
LatticeExprNode mask(const LatticeExprNode &expr)
This function returns the mask of the given expression.
float Float
Definition: aipstype.h:54
int Int
Definition: aipstype.h:50
bool Bool
Define the standard types used by Casacore.
Definition: aipstype.h:42
const Bool True
Definition: aipstype.h:43
double Double
Definition: aipstype.h:55