Line data Source code
1 :
2 : // -----------------------------------------------------------------------------
3 :
4 : /*
5 :
6 : CalStats.h
7 :
8 : Description:
9 : ------------
10 : This header file contains definitions for the CalStats class.
11 :
12 : Classes:
13 : --------
14 : CalStats - This class calculates statistics of new CASA caltables.
15 :
16 : Modification history:
17 : ---------------------
18 : 2011 Nov 11 - Nick Elias, NRAO
19 : Initial version.
20 : 2012 Jan 25 - Nick Elias, NRAO
21 : Logging capability added. Error checking added.
22 :
23 : */
24 :
25 : // -----------------------------------------------------------------------------
26 : // Start of define macro to prevent multiple loading
27 : // -----------------------------------------------------------------------------
28 :
29 : #ifndef CAL_STATS_H
30 : #define CAL_STATS_H
31 :
32 : // -----------------------------------------------------------------------------
33 : // Includes
34 : // -----------------------------------------------------------------------------
35 :
36 : #include <cstring>
37 :
38 : #define _USE_MATH_DEFINES
39 : #include <cmath>
40 :
41 : #include <casacore/casa/BasicSL/String.h>
42 :
43 : #include <casacore/casa/aips.h>
44 :
45 : #include <casacore/casa/Exceptions/Error.h>
46 : #include <casacore/casa/Logging/LogIO.h>
47 :
48 : #include <casacore/casa/Arrays/IPosition.h>
49 : #include <casacore/casa/Arrays/Array.h>
50 : #include <casacore/casa/Arrays/Vector.h>
51 : #include <casacore/casa/Arrays/Cube.h>
52 : #include <casacore/casa/Arrays/ArrayIter.h>
53 : #include <casacore/casa/Arrays/ArrayMath.h>
54 : #include <casacore/casa/Arrays/ArrayLogical.h>
55 :
56 : #include <calanalysis/CalAnalysis/CalStatsFitter.h>
57 :
58 : // -----------------------------------------------------------------------------
59 : // Start of casa namespace
60 : // -----------------------------------------------------------------------------
61 :
62 : namespace casa {
63 :
64 : // -----------------------------------------------------------------------------
65 : // Start of CalStats class definition
66 : // -----------------------------------------------------------------------------
67 :
68 : /*
69 :
70 : CalStats
71 :
72 : Description:
73 : ------------
74 : This class calculates statistics of new CASA caltables.
75 :
76 : NB: At present this class gets data and calculates fit statistics, but other
77 : things such as histogram statistics will be added later.
78 :
79 : In a nutshell:
80 : --------------
81 : * This class can get data (no statistics calculated) and calculate fit
82 : statistics from these data.
83 : * Hooks are present in the code for calculating histogram statistics in the
84 : future.
85 : * This class employs internal iterators to move uniformly through the data
86 : cubes. This process is invisible to the user.
87 : * The input data are cubes whose axes are feed, frequency, and time. The other
88 : axes, such as antenna 1, antenna 2, etc. are handled in another class. There
89 : are two iteration axes and one non-iteration axis, which means that this class
90 : returns ONE-dimensional quantities (data, fits, or histograms) for each
91 : iteration. This class does not deal with multi-dimensional fits and
92 : histograms.
93 : * The feed axis is always an iteration axis (fits cannot be calculated along
94 : it). The user-defined iteration axis is either frequency or time, which means
95 : that the other axis is either time or frequency.
96 : * The input data are obtained from the NewCalTable and NewCalTabIter classes
97 : that iterate along antenna1, antenna2, spw, etc.
98 : * Once an instance is created, the stats<T> member function is the main user
99 : interface to calculate statistics. The choice of T determines which function
100 : is used (getting data, fit statistics, histogram statistics).
101 : * To minimize the amount of duplicated code, some of the member functions are
102 : templated and some of the templated member functions are specialized. The
103 : class is not templated.
104 :
105 : Nested classes:
106 : ---------------
107 : AXES - This nested class contains the axes for the CalStats class.
108 : DATA - This nested class contains the data for the CalStats class.
109 : ARG<T> - This nested template class contains the arguments for the
110 : CalStats::stats<T>() template member function.
111 : OUT<T> - This nested template class contains the outputs for the
112 : CalStats::stats<T>() template member function.
113 :
114 : Class public member functions:
115 : ------------------------------
116 : CalStats - This constructor saves input abscissae and data cubes to internal
117 : copies so that statistics can be calculated.
118 : CalStats - This copy constructor is unused by this class and unavailable when
119 : an instance is created.
120 : operator= - This operator= function is unused by this class and unavailable when
121 : an instance is created.
122 : ~CalStats - This destructor deallocates the internal memory of an instance.
123 :
124 : Class public state member functions:
125 : ------------------------------------
126 : axisIterID - This member function returns the iteration axis IDs.
127 : axisNonIterID - This member function returns the non-iteration axis ID.
128 : axisIterFeed - This member function returns the feed iteration axis values.
129 : axisIterUser - This member function returns the user-defined iteration axis
130 : values.
131 : axisNonIter - This member function returns the non-iteration axis values.
132 : statsShape - This member function returns the shape of the output statistics
133 : cube.
134 : value - This member function returns the input value.
135 : valueErr - This member function returns the input value errors.
136 : flag - This member function returns the input flags.
137 :
138 : Class static public member functions:
139 : -------------------------------------
140 : axisName - This function returns the string corresponding to the CalStats::AXIS
141 : enum.
142 :
143 : Class template public stats member functions:
144 : ---------------------------------------------
145 : stats<T> - This member function is the main user interface for calculating the
146 : the statistics for all iterations. Allowed T: CalStats::NONE only
147 : returns the input data, CalStatsFitter::FIT calculates fit
148 : statistics, and CalStatsHist::HIST calculates histogram statistics.
149 :
150 : Class specialized template public stats member functions:
151 : ---------------------------------------------------------
152 : statsWrap<T> - This member function wraps statistics functions and provides a
153 : uniform interface to stats<T>(). Allowed T: CalStats::NONE only
154 : returns the input data, CalStatsFitter::FIT calculates fit
155 : statistics, and CalStatsHist::HIST calculates histogram
156 : statistics.
157 :
158 : Class protected member functions:
159 : ---------------------------------
160 : CalStats - This default constructor is unused by this class and unavailable when
161 : an instance is created.
162 : next - This member function simultaneously iterates all of the internal
163 : copies of the input data cubes.
164 : reset - This member function simultaneously resets the iterators of all of
165 : the internal copies of the input data cubes.
166 :
167 : Modification history:
168 : ---------------------
169 : 2011 Nov 11 - Nick Elias, NRAO
170 : Initial version created with public member functions CalStats()
171 : (generic) and ~CalStats(); protected member functions CalStats()
172 : (default), CalStats() (copy), operator=(), next() and reset().
173 : 2011 Nov 15 - Nick Elias, NRAO
174 : Moved the CalTypes namespace and its members, originally defined
175 : in this file, to the CalStatsFitter class. Also, the
176 : CalTypes::AXIS typedef was replaced by CalStats::AXIS.
177 : 2011 Dec 08 - Nick Elias, NRAO
178 : Added the fit axis vector to the internal variables.
179 : 2011 Dec 11 - Nick Elias, NRAO
180 : The structure CalStats::FIT was added and replaces the
181 : CalStatsFitter::FIT structure output of the calcFit() public stats
182 : member function (the latter is now part of the former). Added
183 : init() and dealloc() public member functions.
184 : 2011 Dec 14 - Nick Elias, NRAO
185 : The structures CalStats::AXES, CalStats::DATA, and CalStats::NONE
186 : added. The nested class CalStats::OUT<T> added (C++ does not
187 : allow templated typedefs, so a nested class is used instead).
188 : Public member function getData() added. Public static member
189 : functions initAxes(), initGet(), initResultNone(),
190 : initResultFit(), and dealloc() added (overloaded for
191 : CalStats::DATA and CalStats::OUT<CalStatsFitter::FIT>). Removed
192 : public member functions init() and dealloc() (non-overloaded
193 : version).
194 : 2011 Dec 15 - Nick Elias, NRAO
195 : Private member functions next() and reset() now protected member
196 : functions. State public member functions axisIterID(),
197 : axisNonIterID(), axisIterFeed(), axisIterUser(), axisNonIter(),
198 : statsShape(), value(), valueErr(), and flag() added.
199 : 2011 Dec 16 - Nick Elias, NRAO
200 : Public member functions getData() and calcFit() replaced by
201 : stats<T>() template public member function. Specialized template
202 : public member function statsWrap<T>() added. Public static member
203 : functions initAxes(), initData(), initResultNone(), and
204 : initResultFit() replaced by template public static member function
205 : init<T>().
206 : 2011 Dec 21 - Nick Elias, NRAO
207 : casacore::Template public static member functions init<T>() and dealloc<T>
208 : removed because all of their duties are subsumed by the nested
209 : classes AXES, DATA, ARG, and OUT (they were previously
210 : structures).
211 : 2012 Jan 25 - Nick Elias, NRAO
212 : Created working versions of CalStats() (copy) and operator=() and
213 : turned them into public member functions.
214 : 2012 Mar 05 - Nick Elias, NRAO
215 : Static public member function axisName() added.
216 :
217 : */
218 :
219 : // -----------------------------------------------------------------------------
220 :
221 : class CalStats {
222 :
223 : public:
224 :
225 : // Axis enums. There are always two iteration axes. The FEED axis is
226 : // always the first interation axis. Either the FREQUENCY or TIME axis is
227 : // the other (user-defined) iteration axis. The remaining axis (TIME or
228 : // FREQUENCY) is therefore the non-iteration axis. NB: If additional enums
229 : // are added, additional names must be added to the axisName() static
230 : // private member function.
231 : typedef enum AXIS {
232 : INIT=-1, FEED=0, FREQUENCY, TIME
233 : } AXIS;
234 :
235 : // AXES nested class
236 : class AXES {
237 : public:
238 : AXIS eAxisIterFeedID; // FEED iteration axis ID
239 : AXIS eAxisIterUserID; // User-defined iteration axis ID
240 : AXIS eAxisNonIterID; // Non-iteration axis ID
241 : casacore::String sFeed; // FEED axis value
242 : casacore::Double dAxisIterUser; // User-defined iteration axis value
243 : AXES( void );
244 : AXES( const AXES& oAxes );
245 : ~AXES( void );
246 : AXES& operator=( const AXES& oAxes );
247 : };
248 :
249 : // DATA nested class
250 : class DATA {
251 : public:
252 : casacore::Vector<casacore::Double> oAbs; // The abscissae (non-iteration axis values)
253 : casacore::Vector<casacore::Double> oValue; // The values
254 : casacore::Vector<casacore::Double> oValueErr; // The value errors
255 : casacore::Vector<casacore::Bool> oFlag; // The flags
256 : DATA( void );
257 : DATA( const DATA& oDataIn );
258 : ~DATA( void );
259 : DATA& operator=( const DATA& oDataIn );
260 : };
261 :
262 : // Statistics ARG nested class (allowed T: CalStats::NONE,
263 : // CalStatsFitter::FIT, or CalStatsHist::HIST), used as an input to
264 : // stats<T>() and statsWrap<T>(). C++ also does not allow explicit template
265 : // specialization of nested classes within the parent class, so they are
266 : // defined immediately after this class.
267 : template <typename T> class ARG {};
268 :
269 : // NONE nested class
270 : class NONE {};
271 :
272 : // Statistics OUT nested class (allowed T: CalStats::NONE,
273 : // CalStatsFitter::FIT, or CalStatsHist::HIST), used to hold the output of
274 : // statsWrap<T>().
275 : template <typename T> class OUT {
276 : public:
277 : AXES oAxes;
278 : DATA oData;
279 : T oT;
280 : OUT( void );
281 : OUT( const OUT& oOut );
282 : ~OUT( void );
283 : OUT& operator=( const OUT& oOut );
284 : };
285 :
286 : // Generic constructor
287 : CalStats( const casacore::Cube<casacore::Double>& oValue, const casacore::Cube<casacore::Double>& oValueErr,
288 : const casacore::Cube<casacore::Bool>& oFlag, const casacore::Vector<casacore::String>& oFeed,
289 : const casacore::Vector<casacore::Double>& oFrequency, const casacore::Vector<casacore::Double>& oTime,
290 : const AXIS& eAxisIterUser );
291 :
292 : // Copy constructor and operator=() function
293 : CalStats( const CalStats& oCalStats );
294 : CalStats& operator=( const CalStats& oCalStats );
295 :
296 : // Destructor
297 : virtual ~CalStats( void );
298 :
299 : // Axis ID states
300 : casacore::IPosition axisIterID( void ) const;
301 : AXIS axisNonIterID( void ) const;
302 :
303 : // Axis value states
304 : casacore::Vector<casacore::String> axisIterFeed( void ) const;
305 : casacore::Vector<casacore::Double> axisIterUser( void ) const;
306 : casacore::Vector<casacore::Double> axisNonIter( void ) const;
307 :
308 : // Output statistics cube shape state
309 : casacore::IPosition statsShape( void ) const;
310 :
311 : // casacore::Input data states
312 : casacore::Cube<casacore::Double>& value( void ) const;
313 : casacore::Cube<casacore::Double>& valueErr( void ) const;
314 : casacore::Cube<casacore::Bool>& flag( void ) const;
315 :
316 : // The axis names
317 : static casacore::String axisName( const AXIS& eAxis );
318 :
319 : // Calculate statistics (allowed T: CalStats::NONE gets data without
320 : // calculating statistics, CalStatsFitter::FIT calculates fits, and
321 : // CalStatsHist::HIST calculates histogram statistics). Member function
322 : // stats() is the main user interface and statsWrap() is the supporting
323 : // wrapper.
324 : template <typename T> casacore::Matrix<OUT<T> > stats( const ARG<T>& oArg );
325 : template <typename T> T statsWrap( const casacore::Vector<casacore::Double>& oAbs,
326 : const casacore::Vector<casacore::Double>& oValue, const casacore::Vector<casacore::Double>& oValueErr,
327 : casacore::Vector<casacore::Bool>& oFlag, const ARG<T>& oArg );
328 :
329 : protected:
330 :
331 : // The axis IDs. The two iteration axes are FEED (always) and either TIME
332 : // or FREQUENCY (user defined). The non-iteration axis is either FREQUENCY
333 : // or TIME (the opposite of the user-defined iteration axis).
334 : casacore::IPosition oAxisIterID;
335 : AXIS eAxisNonIterID;
336 :
337 : // Internal copies of the iteration and non-iteration axis values
338 : casacore::Vector<casacore::String> oAxisIterFeed; // Feed axis iteration axis values
339 : casacore::Vector<casacore::Double> oAxisIterUser; // User-defined iteration axis values
340 : casacore::Vector<casacore::Double> oAxisNonIter; // Non-iteration axis values
341 :
342 : // Shape of the output statistics cubes
343 : casacore::IPosition oStatsShape;
344 :
345 : // Internal copies of input parameter cubes
346 : casacore::Cube<casacore::Double>* poValue;
347 : casacore::Cube<casacore::Double>* poValueErr;
348 : casacore::Cube<casacore::Bool>* poFlag;
349 :
350 : // casacore::Input parameter cube iterators
351 : casacore::ArrayIterator<casacore::Double>* poValueIter;
352 : casacore::ArrayIterator<casacore::Double>* poValueErrIter;
353 : casacore::ArrayIterator<casacore::Bool>* poFlagIter;
354 :
355 : // Unused constructor
356 : CalStats( void );
357 :
358 : // Simultaneously increment and reset all input parameter cube iterators
359 : void next( void );
360 : void reset( void );
361 :
362 : };
363 :
364 : // -----------------------------------------------------------------------------
365 : // End of CalStats class definition
366 : // -----------------------------------------------------------------------------
367 :
368 : // -----------------------------------------------------------------------------
369 : // Start of ARG<T> specialized class templates
370 : // -----------------------------------------------------------------------------
371 :
372 : // Specialization for the CalStats::NONE() class. It tells the
373 : // CalStats::stats<CalStats::NONE>() method just to return the data with no
374 : // processing.
375 : template <> class CalStats::ARG<CalStats::NONE> {};
376 :
377 : // Specialization for the CalStatsFitter::FIT() class. It tells the
378 : // CalStats::stats<CalStatsFitter::FIT>() method to perform fits and return the
379 : // data and fit parameters.
380 : template <> class CalStats::ARG<CalStatsFitter::FIT> {
381 : public:
382 : CalStatsFitter::ORDER eOrder;
383 : CalStatsFitter::TYPE eType;
384 : CalStatsFitter::WEIGHT eWeight;
385 : ARG( void ) {
386 : eOrder = CalStatsFitter::ORDER_INIT;
387 : eType = CalStatsFitter::TYPE_INIT;
388 : eWeight = CalStatsFitter::WEIGHT_INIT;
389 : return;
390 : }
391 : };
392 :
393 : // -----------------------------------------------------------------------------
394 : // End of ARG<T> specialized class templates
395 : // -----------------------------------------------------------------------------
396 :
397 : // -----------------------------------------------------------------------------
398 : // Start of OUT<T> specialized class template public member functions
399 : // -----------------------------------------------------------------------------
400 :
401 : // Default constructor
402 : template <typename T>
403 400 : CalStats::OUT<T>::OUT( void ) {
404 400 : oAxes = CalStats::AXES();
405 400 : oData = CalStats::DATA();
406 400 : oT = T();
407 400 : }
408 :
409 : // Copy constructor
410 : template <typename T>
411 200 : CalStats::OUT<T>::OUT( const CalStats::OUT<T>& oOut ) {
412 200 : oAxes = CalStats::AXES( oOut.oAxes );
413 200 : oData = CalStats::DATA( oOut.oData );
414 200 : oT = T( oOut.oT );
415 200 : }
416 :
417 : // Destructor
418 : template <typename T>
419 600 : CalStats::OUT<T>::~OUT( void ) {}
420 :
421 : // operator=
422 : template <typename T>
423 200 : CalStats::OUT<T>& CalStats::OUT<T>::operator=( const CalStats::OUT<T>& oOut ) {
424 200 : if ( this != &oOut ) {
425 200 : oAxes = CalStats::AXES( oOut.oAxes );
426 200 : oData = CalStats::DATA( oOut.oData );
427 200 : oT = T( oOut.oT );
428 : }
429 200 : return( *this );
430 : }
431 :
432 : // -----------------------------------------------------------------------------
433 : // End of OUT<T> specialized class template public member functions
434 : // -----------------------------------------------------------------------------
435 :
436 : // -----------------------------------------------------------------------------
437 : // Start of CalStats::stats<T> template public statistics member function
438 : // -----------------------------------------------------------------------------
439 :
440 : /*
441 :
442 : CalStats::stats<T>
443 :
444 : Description:
445 : ------------
446 : This member function calculates the desired statistics. The allowed templates
447 : are CalStats::NONE (no statistics, just the input data), CalStatsFitter::FIT
448 : (fit statistics), and CalStatsHist::HIST (histogram statistics).
449 :
450 : Inputs:
451 : -------
452 : oArg - This reference to a CalStats::ARG<T> instance contains the extra input
453 : parameters.
454 :
455 : Outputs:
456 : --------
457 : The reference to the casacore::Matrix<CalStats::OUT<T> > instance containing the
458 : statistics, returned via the function value.
459 :
460 : Modification history:
461 : ---------------------
462 : 2011 Nov 11 - Nick Elias, NRAO
463 : Initial version. This template member function replaces the
464 : getData() and calcFit() member functions.
465 : 2012 Jan 25 - Nick Elias, NRAO
466 : Logging capability added.
467 :
468 : */
469 :
470 : // -----------------------------------------------------------------------------
471 :
472 : template <typename T>
473 1 : casacore::Matrix<CalStats::OUT<T>> CalStats::stats( const CalStats::ARG<T>& oArg ) {
474 :
475 : // Initialize the CalStats::OUT<T> array and its iterator
476 :
477 1 : casacore::Array<CalStats::OUT<T> > out( oStatsShape );
478 :
479 1 : casacore::ArrayIterator<CalStats::OUT<T> > oOutIter( out, oAxisIterID, false );
480 :
481 :
482 : // For each iteration, convert the resulting arrays to vectors and feed them
483 : // to the CalStatsFitter::fit() function
484 :
485 401 : while ( !poValueIter->pastEnd() ) {
486 :
487 200 : casacore::IPosition oPos( poValueIter->pos() );
488 :
489 200 : casacore::uInt uiLength = poValueIter->array().nelements();
490 200 : casacore::IPosition oShape( 1, uiLength );
491 :
492 200 : casacore::Vector<casacore::Double> oAbs( oAxisNonIter.copy() );
493 400 : casacore::Vector<casacore::Double> oValue( poValueIter->array().copy().reform(oShape) );
494 400 : casacore::Vector<casacore::Double> oValueErr( poValueErrIter->array().copy().reform(oShape) );
495 400 : casacore::Vector<casacore::Bool> oFlag( poFlagIter->array().copy().reform(oShape) );
496 :
497 200 : CalStats::OUT<T> oOut;
498 :
499 200 : oOut.oAxes.eAxisIterFeedID = (CalStats::AXIS) oAxisIterID[0];
500 200 : oOut.oAxes.eAxisIterUserID = (CalStats::AXIS) oAxisIterID[1];
501 200 : oOut.oAxes.eAxisNonIterID = eAxisNonIterID;
502 200 : oOut.oAxes.sFeed = casacore::String( oAxisIterFeed[oPos[0]] );
503 200 : oOut.oAxes.dAxisIterUser = oAxisIterUser[oPos[oAxisIterID[1]]];
504 :
505 200 : oOut.oData.oAbs = casacore::Vector<casacore::Double>( oAbs );
506 200 : oOut.oData.oValue = casacore::Vector<casacore::Double>( oValue );
507 200 : oOut.oData.oValueErr = casacore::Vector<casacore::Double>( oValueErr );
508 :
509 : try {
510 200 : oOut.oT = statsWrap<T>( oAbs, oValue, oValueErr, oFlag, oArg );
511 : }
512 :
513 0 : catch ( casacore::AipsError oAE ) {
514 0 : casacore::LogIO log( casacore::LogOrigin( "CalStats", "stats<T>()", WHERE ) );
515 : log << casacore::LogIO::NORMAL << oAE.getMesg() << ", iteration: "
516 0 : << oPos.asVector() << ", continuing ..." << casacore::LogIO::POST;
517 0 : oOut.oT = T();
518 0 : }
519 :
520 : // The flag output vector is set here because robust fitting can change them
521 200 : oOut.oData.oFlag = casacore::Vector<casacore::Bool>( oFlag );
522 :
523 200 : oOutIter.array() = casacore::Vector<CalStats::OUT<T> >( 1, oOut );
524 :
525 200 : next(); oOutIter.next();
526 :
527 : }
528 :
529 :
530 : // Reset the input parameter iterators
531 :
532 1 : reset();
533 :
534 :
535 : // Return the reference to the casacore::Matrix<CalStats::OUT<T> > instance
536 :
537 1 : out.removeDegenerate();
538 :
539 1 : casacore::Matrix<CalStats::OUT<T> > outMatrix = out;
540 2 : return outMatrix;
541 :
542 1 : }
543 :
544 : // -----------------------------------------------------------------------------
545 : // End of CalStats::stats<T> template public statistics member function
546 : // -----------------------------------------------------------------------------
547 :
548 : };
549 :
550 : // -----------------------------------------------------------------------------
551 : // End of casa namespace
552 : // -----------------------------------------------------------------------------
553 :
554 : #endif
555 :
556 : // -----------------------------------------------------------------------------
557 : // End of define macro to prevent multiple loading
558 : // -----------------------------------------------------------------------------
|
