Line data Source code
1 : //# MomentsBase.h: base class for moment generator
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 : //# $Id: MomentsBase.h 20299 2008-04-03 05:56:44Z gervandiepen $
27 :
28 : #ifndef IMAGES_MOMENTSBASE_H
29 : #define IMAGES_MOMENTSBASE_H
30 :
31 : #include <casacore/casa/aips.h>
32 : #include <casacore/coordinates/Coordinates/CoordinateSystem.h>
33 : #include <casacore/casa/Quanta/Quantum.h>
34 : #include <casacore/measures/Measures/MDoppler.h>
35 : #include <casacore/casa/Logging/LogIO.h>
36 : #include <casacore/casa/Arrays/Vector.h>
37 : #include <iosfwd>
38 :
39 : namespace casacore{
40 :
41 : class IPosition;
42 : class String;
43 : class Unit;
44 : }
45 :
46 : namespace casa {
47 :
48 : template <class T> class MomentCalcBase;
49 :
50 : // <summary>
51 : // This class is a base class for generating moments from an image or a spectral data.
52 : // </summary>
53 :
54 : // <use visibility=export>
55 :
56 : // <reviewed reviewer="" date="yyyy/mm/dd" tests="" demos="">
57 : // </reviewed>
58 :
59 : // <prerequisite>
60 : // <li> <linkto class="ImageMoments">ImageMoments</linkto>
61 : // <li> <linkto class="MSMoments">MSMoments</linkto>
62 : // <li> <linkto class="casacore::LatticeApply">casacore::LatticeApply</linkto>
63 : // <li> <linkto class="MomentCalcBase">MomentCalcBase</linkto>
64 : // </prerequisite>
65 :
66 : // <etymology>
67 : // This class is an abstract class to compute moments from images or spectral data.
68 : // </etymology>
69 :
70 : // <synopsis>
71 : // The primary goal of MSMoments, ImageMoments, and MSMoments are to help
72 : // spectral-line astronomers analyze their multi-dimensional images or
73 : // spectral data (in the form of casacore::MeasurementSet) by generating moments of
74 : // a specified axis. ImageMoments is a specialized class to generate moments
75 : // from images, while MSMoments is designed properly for casacore::MeasurementSet input.
76 : // MSMoments class is an abstract class that is inherited by the above two
77 : // concrete classes.
78 : // MomentsBase provides interface layer to the MomentCalculators so that
79 : // functionalities in MomentCalculators can be shared with ImageMoments and
80 : // MSMoments.
81 : // The word "moment" is used loosely here. It refers to collapsing an axis
82 : // to one pixel and putting the value of that pixel (for all of the other
83 : // non-collapsed axes) to something computed from the data values along
84 : // the moment axis. For example, take an RA-DEC-Velocity cube, collapse
85 : // the velocity axis by computing the mean intensity at each RA-DEC
86 : // pixel. This class and its inheritances offer many different moments and a variety of
87 : // interactive and automatic ways to compute them.
88 : //
89 :
90 : // <motivation>
91 : // MSMoments is defined so that moments can be generated from both images
92 : // and spectral data (in the form of casacore::MeasurementSet).
93 : // </motivation>
94 :
95 :
96 : template <class T> class MomentsBase {
97 : public:
98 :
99 : // Note that if I don't put MomentCalcBase as a forward declaration
100 : // and use instead "friend class MomentCalcBase<T>" The gnu compiler
101 : // fails with a typedef problem. But I can't solve it with say
102 : // <src>typedef MomentCalcBase<T> gpp_type;</src> because you can only do a
103 : // typedef with an actual type, not a <tt>T</tt> !
104 : friend class MomentCalcBase<T>;
105 :
106 : // The <src>enum MethodTypes</src> is provided for use with the
107 : // <src>setWinFitMethod</src> function. It gives the allowed moment
108 : // methods which are available with this function. The use of these
109 : // methods is described further with the description of this function
110 : // as well as in the general documentation earlier.
111 : enum MethodTypes {
112 : // Invokes the spectral windowing method
113 : WINDOW,
114 : // Invokes Gaussian fitting
115 : FIT,
116 : NMETHODS
117 : };
118 :
119 : // This <src>enum MomentTypes</src> is provided for use with the
120 : // <src>setMoments</src> function. It gives the allowed moment
121 : // types that you can ask for.
122 :
123 : enum MomentTypes {
124 : // The average intensity
125 : AVERAGE,
126 : // The integrated intensity
127 : INTEGRATED,
128 : // The intensity weighted mean coordinate (usually velocity)
129 : WEIGHTED_MEAN_COORDINATE,
130 : // The intensity weighted coordinate (usually velocity) dispersion
131 : WEIGHTED_DISPERSION_COORDINATE,
132 : // The median intensity
133 : MEDIAN,
134 : // The median coordinate (usually velocity). Treat the spectrum as
135 : // a probability distribution, generate the cumulative distribution,
136 : // and find the coordinate corresponding to the 50% value.
137 : MEDIAN_COORDINATE,
138 : // The standard deviation about the mean of the intensity
139 : STANDARD_DEVIATION,
140 : // The rms of the intensity
141 : RMS,
142 : // The absolute mean deviation of the intensity
143 : ABS_MEAN_DEVIATION,
144 : // The maximum value of the intensity
145 : MAXIMUM,
146 : // The coordinate (usually velocity) of the maximum value of the intensity
147 : MAXIMUM_COORDINATE,
148 : // The minimum value of the intensity
149 : MINIMUM,
150 : // The coordinate (usually velocity) of the minimum value of the intensity
151 : MINIMUM_COORDINATE,
152 : // Total number
153 : NMOMENTS,
154 : // Default value is the integrated intensity
155 : DEFAULT = INTEGRATED
156 : };
157 :
158 : // Destructor
159 : virtual ~MomentsBase();
160 :
161 : // Set the desired moments via an <src>casacore::Int</src> array. Each <src>casacore::Int</src>
162 : // specifies a different moment; the allowed values and their meanings
163 : // are given by the <src>enum MomentTypes</src>. A return value
164 : // of <src>false</src> indicates you asked for an out of range
165 : // moment. If you don't call this function, the default state of the
166 : // class is to request the integrated intensity.
167 : casacore::Bool setMoments (const casacore::Vector<casacore::Int>& moments);
168 :
169 : // Set the moment axis (0 relative). A return value of <src>false</src>
170 : // indicates that the axis was not contained in the image. If you don't
171 : // call this function, the default state of the class is to set the
172 : // moment axis to the spectral axis if it can find one. Otherwise
173 : // an error will result.
174 : virtual casacore::Bool setMomentAxis (casacore::Int) = 0;
175 :
176 : // The method by which you compute the moments is specified by calling
177 : // (or not calling) the <src>setWinFitMethod</src> and
178 : // <src>setSmoothMethod</src> functions. The default state of the class
179 : // is to compute directly on all (or some according to <src>setInExClude</src>)
180 : // of the pixels in the spectrum. Calling these functions modifies the
181 : // computational state to something more complicated.
182 : //
183 : // The <src>setWinMethod</src> function requires an <src>casacore::Int</src> array
184 : // as its argument. Each <src>casacore::Int</src> specifies a different method
185 : // that you can invoke (either separately or in combination). The
186 : // allowed values and their meanings are given by the
187 : // <src>enum MethodTypes</src>.
188 : //
189 : // Both the windowing and fitting methods have interactive modes. The
190 : // windowing method also has a fitting flavour, so if you set both
191 : // MomentsBase::WINDOW and MomentsBase::FIT, you would be invoking the
192 : // windowing method but determining the window by fitting Gaussians
193 : // automatically (as MomentsBase::INTERACTIVE) was not given.
194 : //
195 : // If you don't call this function, then neither the windowing nor fitting
196 : // methods are invoked. A return value of <src>false</src> indicates
197 : // that you asked for an illegal method.
198 : casacore::Bool setWinFitMethod(const casacore::Vector<casacore::Int>& method);
199 :
200 : // This function invokes smoothing of the input image. Give <src>casacore::Int</src>
201 : // arrays for the axes (0 relative) to be smoothed and the smoothing kernel
202 : // types (use the <src>enum KernelTypes</src>) for each axis. Give a
203 : // <src>casacore::Double</src> array for the widths (full width for BOXCAR and full
204 : // width at half maximum for GAUSSIAN) in pixels of the smoothing kernels for
205 : // each axis. For HANNING smoothing, you always get the quarter-half-quarter
206 : // kernel (no matter what you might ask for). A return value of <src>false</src>
207 : // indicates that you have given an inconsistent or invalid set of smoothing
208 : // parameters. If you don't call this function the default state of the
209 : // class is to do no smoothing. The kernel types are specified with
210 : // the casacore::VectorKernel::KernelTypes enum
211 : // <group>
212 : virtual casacore::Bool setSmoothMethod(
213 : const casacore::Vector<casacore::Int>&, const casacore::Vector<casacore::Int>&,
214 : const casacore::Vector<casacore::Quantum<casacore::Double> >&
215 : ) = 0;
216 :
217 : casacore::Bool setSmoothMethod(
218 : const casacore::Vector<casacore::Int>& smoothAxes,
219 : const casacore::Vector<casacore::Int>& kernelTypes,
220 : const casacore::Vector<casacore::Double>& kernelWidths
221 : );
222 : // </group>
223 :
224 : // You may specify a pixel intensity range as either one for which
225 : // all pixels in that range are included in the moment calculation,
226 : // or one for which all pixels in that range are excluded from the moment
227 : // calculations. One or the other of <src>include</src> and <src>exclude</src>
228 : // must therefore be a zero length vector if you call this function.
229 : // A return value of <src>false</src> indicates that you have given both
230 : // an <src>include</src> and an <src>exclude</src> range. If you don't call
231 : // this function, the default state of the class is to include all pixels.
232 : void setInExCludeRange(
233 : const casacore::Vector<T>& include, const casacore::Vector<T>& exclude
234 : );
235 :
236 : // This function is used to help assess whether a spectrum along the moment
237 : // axis is all noise or not. If it is all noise, there is not a lot of point
238 : // to trying to computing the moment.
239 : // <src>peakSNR</src> is the signal-to-noise ratio of the peak value in the
240 : // spectrum below which the spectrum is considered pure noise.
241 : // <src>stdDeviation</src> is the standard deviation of the noise for the
242 : // input image.
243 : //
244 : // Default values for one or the other parameter are indicated by giving zero.
245 : //
246 : // The default state of the class then is to set <src>peakSNR=3</src>
247 : // and/or to work out the noise level from a Gaussian fit to a histogram
248 : // (above 25%) of the entire image (it is very hard to get an accurate
249 : // estimate of the noise a single spectrum).
250 : void setSnr(const T& peakSNR, const T& stdDeviation);
251 :
252 : // This is the output file name for the smoothed image. It can be useful
253 : // to have access this to this image when trying to get the pixel
254 : // <src>include</src> or <src>exclude</src> range correct for the smooth-clip
255 : // method. The default state of the class is to not output the smoothed image.
256 : casacore::Bool setSmoothOutName(const casacore::String& smOut);
257 :
258 : // Set Velocity type. This is used for moments for which the moment axis is
259 : // a spectral axis for which the coordinate is traditionally presented in
260 : // km/s You can select the velocity definition. The default state of the
261 : // class is to use the radio definition.
262 : void setVelocityType (casacore::MDoppler::Types type);
263 :
264 : // Reset argument error condition. If you specify invalid arguments to
265 : // one of the above functions, an internal flag will be set which will
266 : // prevent the <src>createMoments</src> function, which is defined in its inheritances,
267 : // from doing anything
268 : // (should you have chosen to igmore the Boolean return values of the functions).
269 : // This function allows you to reset that internal state to good.
270 : void resetError () {goodParameterStatus_p = true; error_p = "";};
271 :
272 : // Recover last error message
273 0 : casacore::String errorMessage() const {return error_p;};
274 :
275 : // Get CoordinateSystem
276 : virtual const casacore::CoordinateSystem& coordinates() = 0;
277 :
278 : // Helper function to convert a string containing a list of desired methods to
279 : // the correct <src>casacore::Vector<casacore::Int></src> required for the <src>setWinFitMethod</src> function.
280 : // This may be usful if your user interface involves strings rather than integers.
281 : // A new value is added to the output vector (which is resized appropriately) if any of the
282 : // substrings "window", "fit" or "interactive" (actually "win", "box" and
283 : // "inter" will do) is present.
284 : static casacore::Vector<casacore::Int> toMethodTypes (const casacore::String& methods);
285 :
286 : virtual casacore::IPosition getShape() const = 0;
287 :
288 0 : casacore::Bool shouldConvertToVelocity() const { return convertToVelocity_p; }
289 :
290 : protected:
291 :
292 : // Constructor takes an image and a <src>casacore::LogIO</src> object for logging purposes.
293 : // You specify whether output images are automatically overwritten if pre-existing,
294 : // or whether an intercative choice dialog widget appears (overWriteOutput=F)
295 : // You may also determine whether a progress meter is displayed or not.
296 : MomentsBase(
297 : casacore::LogIO &os, casacore::Bool overWriteOutput=false,
298 : casacore::Bool showProgress=true
299 : );
300 :
301 : casacore::LogIO os_p;
302 : casacore::Bool showProgress_p;
303 : casacore::Int momentAxisDefault_p {-10};
304 : T peakSNR_p {T(3)};
305 : T stdDeviation_p {T(0)};
306 : T yMin_p {T(0)};
307 : T yMax_p {T(0)};
308 : casacore::String out_p;
309 : casacore::String smoothOut_p {};
310 : casacore::Bool goodParameterStatus_p {true};
311 : casacore::Bool doWindow_p {false};
312 : casacore::Bool doFit_p {false};
313 : casacore::Bool doSmooth_p {false};
314 : casacore::Bool noInclude_p {true};
315 : casacore::Bool noExclude_p {true};
316 : casacore::Bool fixedYLimits_p {false};
317 :
318 : casacore::Int momentAxis_p {momentAxisDefault_p};
319 : casacore::Int worldMomentAxis_p;
320 : casacore::Vector<casacore::Int> kernelTypes_p {};
321 : casacore::Vector<casacore::Quantity> kernelWidths_p {};
322 : casacore::Vector<casacore::Int> moments_p {1, INTEGRATED};
323 : casacore::Vector<T> selectRange_p {};
324 : casacore::Vector<casacore::Int> smoothAxes_p {};
325 : casacore::Bool overWriteOutput_p;
326 : casacore::String error_p {};
327 : casacore::Bool convertToVelocity_p {false};
328 : casacore::MDoppler::Types velocityType_p {casacore::MDoppler::RADIO};
329 :
330 : // Check that the combination of methods that the user has requested is valid
331 : // casacore::List a handy dandy table if not.
332 : void _checkMethod();
333 :
334 : // Take the user's data inclusion and exclusion data ranges and
335 : // generate the range and Booleans to say what sort it is
336 : void _setIncludeExclude (
337 : casacore::Vector<T>& range, casacore::Bool& noInclude,
338 : casacore::Bool& noExclude, const casacore::Vector<T>& include,
339 : const casacore::Vector<T>& exclude
340 : );
341 :
342 : // Set the output image suffixes and units
343 : casacore::Bool _setOutThings(
344 : casacore::String& suffix, casacore::Unit& momentUnits,
345 : const casacore::Unit& imageUnits, const casacore::String& momentAxisUnits,
346 : const casacore::Int moment, casacore::Bool convertToVelocity
347 : );
348 :
349 : // Make output Coordinates
350 : casacore::CoordinateSystem _makeOutputCoordinates(
351 : casacore::IPosition& outShape, const casacore::CoordinateSystem& cSysIn,
352 : const casacore::IPosition& inShape, casacore::Int momentAxis,
353 : casacore::Bool removeAxis
354 : );
355 : };
356 :
357 : }
358 :
359 : #ifndef AIPS_NO_TEMPLATE_SRC
360 : #include <imageanalysis/ImageAnalysis/MomentsBase.tcc>
361 : #endif
362 :
363 : #endif
|