00001 00002 // ----------------------------------------------------------------------------- 00003 00004 /* 00005 00006 CalStats.h 00007 00008 Description: 00009 ------------ 00010 This header file contains definitions for the CalStats class. 00011 00012 Classes: 00013 -------- 00014 CalStats - This class calculates statistics of new CASA caltables. 00015 00016 Modification history: 00017 --------------------- 00018 2011 Nov 11 - Nick Elias, NRAO 00019 Initial version. 00020 2012 Jan 25 - Nick Elias, NRAO 00021 Logging capability added. Error checking added. 00022 00023 */ 00024 00025 // ----------------------------------------------------------------------------- 00026 // Start of define macro to prevent multiple loading 00027 // ----------------------------------------------------------------------------- 00028 00029 #ifndef CAL_STATS_H 00030 #define CAL_STATS_H 00031 00032 // ----------------------------------------------------------------------------- 00033 // Includes 00034 // ----------------------------------------------------------------------------- 00035 00036 #include <cstring> 00037 00038 #define _USE_MATH_DEFINES 00039 #include <cmath> 00040 00041 #include <casa/BasicSL/String.h> 00042 00043 #include <casa/aips.h> 00044 00045 #include <casa/Exceptions/Error.h> 00046 #include <casa/Logging/LogIO.h> 00047 00048 #include <casa/Arrays/IPosition.h> 00049 #include <casa/Arrays/Array.h> 00050 #include <casa/Arrays/Vector.h> 00051 #include <casa/Arrays/Cube.h> 00052 #include <casa/Arrays/ArrayIter.h> 00053 #include <casa/Arrays/ArrayMath.h> 00054 #include <casa/Arrays/ArrayLogical.h> 00055 00056 #include <calanalysis/CalAnalysis/CalStatsFitter.h> 00057 00058 // ----------------------------------------------------------------------------- 00059 // Start of casa namespace 00060 // ----------------------------------------------------------------------------- 00061 00062 namespace casa { 00063 00064 // ----------------------------------------------------------------------------- 00065 // Start of CalStats class definition 00066 // ----------------------------------------------------------------------------- 00067 00068 /* 00069 00070 CalStats 00071 00072 Description: 00073 ------------ 00074 This class calculates statistics of new CASA caltables. 00075 00076 NB: At present this class gets data and calculates fit statistics, but other 00077 things such as histogram statistics will be added later. 00078 00079 In a nutshell: 00080 -------------- 00081 * This class can get data (no statistics calculated) and calculate fit 00082 statistics from these data. 00083 * Hooks are present in the code for calculating histogram statistics in the 00084 future. 00085 * This class employs internal iterators to move uniformly through the data 00086 cubes. This process is invisible to the user. 00087 * The input data are cubes whose axes are feed, frequency, and time. The other 00088 axes, such as antenna 1, antenna 2, etc. are handled in another class. There 00089 are two iteration axes and one non-iteration axis, which means that this class 00090 returns ONE-dimensional quantities (data, fits, or histograms) for each 00091 iteration. This class does not deal with multi-dimensional fits and 00092 histograms. 00093 * The feed axis is always an iteration axis (fits cannot be calculated along 00094 it). The user-defined iteration axis is either frequency or time, which means 00095 that the other axis is either time or frequency. 00096 * The input data are obtained from the NewCalTable and NewCalTabIter classes 00097 that iterate along antenna1, antenna2, spw, etc. 00098 * Once an instance is created, the stats<T> member function is the main user 00099 interface to calculate statistics. The choice of T determines which function 00100 is used (getting data, fit statistics, histogram statistics). 00101 * To minimize the amount of duplicated code, some of the member functions are 00102 templated and some of the templated member functions are specialized. The 00103 class is not templated. 00104 00105 Nested classes: 00106 --------------- 00107 AXES - This nested class contains the axes for the CalStats class. 00108 DATA - This nested class contains the data for the CalStats class. 00109 ARG<T> - This nested template class contains the arguments for the 00110 CalStats::stats<T>() template member function. 00111 OUT<T> - This nested template class contains the outputs for the 00112 CalStats::stats<T>() template member function. 00113 00114 Class public member functions: 00115 ------------------------------ 00116 CalStats - This constructor saves input abscissae and data cubes to internal 00117 copies so that statistics can be calculated. 00118 CalStats - This copy constructor is unused by this class and unavailable when 00119 an instance is created. 00120 operator= - This operator= function is unused by this class and unavailable when 00121 an instance is created. 00122 ~CalStats - This destructor deallocates the internal memory of an instance. 00123 00124 Class public state member functions: 00125 ------------------------------------ 00126 axisIterID - This member function returns the iteration axis IDs. 00127 axisNonIterID - This member function returns the non-iteration axis ID. 00128 axisIterFeed - This member function returns the feed iteration axis values. 00129 axisIterUser - This member function returns the user-defined iteration axis 00130 values. 00131 axisNonIter - This member function returns the non-iteration axis values. 00132 statsShape - This member function returns the shape of the output statistics 00133 cube. 00134 value - This member function returns the input value. 00135 valueErr - This member function returns the input value errors. 00136 flag - This member function returns the input flags. 00137 00138 Class static public member functions: 00139 ------------------------------------- 00140 axisName - This function returns the string corresponding to the CalStats::AXIS 00141 enum. 00142 00143 Class template public stats member functions: 00144 --------------------------------------------- 00145 stats<T> - This member function is the main user interface for calculating the 00146 the statistics for all iterations. Allowed T: CalStats::NONE only 00147 returns the input data, CalStatsFitter::FIT calculates fit 00148 statistics, and CalStatsHist::HIST calculates histogram statistics. 00149 00150 Class specialized template public stats member functions: 00151 --------------------------------------------------------- 00152 statsWrap<T> - This member function wraps statistics functions and provides a 00153 uniform interface to stats<T>(). Allowed T: CalStats::NONE only 00154 returns the input data, CalStatsFitter::FIT calculates fit 00155 statistics, and CalStatsHist::HIST calculates histogram 00156 statistics. 00157 00158 Class protected member functions: 00159 --------------------------------- 00160 CalStats - This default constructor is unused by this class and unavailable when 00161 an instance is created. 00162 next - This member function simultaneously iterates all of the internal 00163 copies of the input data cubes. 00164 reset - This member function simultaneously resets the iterators of all of 00165 the internal copies of the input data cubes. 00166 00167 Modification history: 00168 --------------------- 00169 2011 Nov 11 - Nick Elias, NRAO 00170 Initial version created with public member functions CalStats() 00171 (generic) and ~CalStats(); protected member functions CalStats() 00172 (default), CalStats() (copy), operator=(), next() and reset(). 00173 2011 Nov 15 - Nick Elias, NRAO 00174 Moved the CalTypes namespace and its members, originally defined 00175 in this file, to the CalStatsFitter class. Also, the 00176 CalTypes::AXIS typedef was replaced by CalStats::AXIS. 00177 2011 Dec 08 - Nick Elias, NRAO 00178 Added the fit axis vector to the internal variables. 00179 2011 Dec 11 - Nick Elias, NRAO 00180 The structure CalStats::FIT was added and replaces the 00181 CalStatsFitter::FIT structure output of the calcFit() public stats 00182 member function (the latter is now part of the former). Added 00183 init() and dealloc() public member functions. 00184 2011 Dec 14 - Nick Elias, NRAO 00185 The structures CalStats::AXES, CalStats::DATA, and CalStats::NONE 00186 added. The nested class CalStats::OUT<T> added (C++ does not 00187 allow templated typedefs, so a nested class is used instead). 00188 Public member function getData() added. Public static member 00189 functions initAxes(), initGet(), initResultNone(), 00190 initResultFit(), and dealloc() added (overloaded for 00191 CalStats::DATA and CalStats::OUT<CalStatsFitter::FIT>). Removed 00192 public member functions init() and dealloc() (non-overloaded 00193 version). 00194 2011 Dec 15 - Nick Elias, NRAO 00195 Private member functions next() and reset() now protected member 00196 functions. State public member functions axisIterID(), 00197 axisNonIterID(), axisIterFeed(), axisIterUser(), axisNonIter(), 00198 statsShape(), value(), valueErr(), and flag() added. 00199 2011 Dec 16 - Nick Elias, NRAO 00200 Public member functions getData() and calcFit() replaced by 00201 stats<T>() template public member function. Specialized template 00202 public member function statsWrap<T>() added. Public static member 00203 functions initAxes(), initData(), initResultNone(), and 00204 initResultFit() replaced by template public static member function 00205 init<T>(). 00206 2011 Dec 21 - Nick Elias, NRAO 00207 Template public static member functions init<T>() and dealloc<T> 00208 removed because all of their duties are subsumed by the nested 00209 classes AXES, DATA, ARG, and OUT (they were previously 00210 structures). 00211 2012 Jan 25 - Nick Elias, NRAO 00212 Created working versions of CalStats() (copy) and operator=() and 00213 turned them into public member functions. 00214 2012 Mar 05 - Nick Elias, NRAO 00215 Static public member function axisName() added. 00216 00217 */ 00218 00219 // ----------------------------------------------------------------------------- 00220 00221 class CalStats { 00222 00223 public: 00224 00225 // Axis enums. There are always two iteration axes. The FEED axis is 00226 // always the first interation axis. Either the FREQUENCY or TIME axis is 00227 // the other (user-defined) iteration axis. The remaining axis (TIME or 00228 // FREQUENCY) is therefore the non-iteration axis. NB: If additional enums 00229 // are added, additional names must be added to the axisName() static 00230 // private member function. 00231 typedef enum AXIS { 00232 INIT=-1, FEED=0, FREQUENCY, TIME 00233 } AXIS; 00234 00235 // AXES nested class 00236 class AXES { 00237 public: 00238 AXIS eAxisIterFeedID; // FEED iteration axis ID 00239 AXIS eAxisIterUserID; // User-defined iteration axis ID 00240 AXIS eAxisNonIterID; // Non-iteration axis ID 00241 String sFeed; // FEED axis value 00242 Double dAxisIterUser; // User-defined iteration axis value 00243 AXES( void ); 00244 AXES( const AXES& oAxes ); 00245 ~AXES( void ); 00246 AXES& operator=( const AXES& oAxes ); 00247 }; 00248 00249 // DATA nested class 00250 class DATA { 00251 public: 00252 Vector<Double> oAbs; // The abscissae (non-iteration axis values) 00253 Vector<Double> oValue; // The values 00254 Vector<Double> oValueErr; // The value errors 00255 Vector<Bool> oFlag; // The flags 00256 DATA( void ); 00257 DATA( const DATA& oDataIn ); 00258 ~DATA( void ); 00259 DATA& operator=( const DATA& oDataIn ); 00260 }; 00261 00262 // Statistics ARG nested class (allowed T: CalStats::NONE, 00263 // CalStatsFitter::FIT, or CalStatsHist::HIST), used as an input to 00264 // stats<T>() and statsWrap<T>(). C++ also does not allow explicit template 00265 // specialization of nested classes within the parent class, so they are 00266 // defined immediately after this class. 00267 template <typename T> class ARG {}; 00268 00269 // NONE nested class 00270 class NONE {}; 00271 00272 // Statistics OUT nested class (allowed T: CalStats::NONE, 00273 // CalStatsFitter::FIT, or CalStatsHist::HIST), used to hold the output of 00274 // statsWrap<T>(). 00275 template <typename T> class OUT { 00276 public: 00277 AXES oAxes; 00278 DATA oData; 00279 T oT; 00280 OUT( void ); 00281 OUT( const OUT& oOut ); 00282 ~OUT( void ); 00283 OUT& operator=( const OUT& oOut ); 00284 }; 00285 00286 // Generic constructor 00287 CalStats( const Cube<Double>& oValue, const Cube<Double>& oValueErr, 00288 const Cube<Bool>& oFlag, const Vector<String>& oFeed, 00289 const Vector<Double>& oFrequency, const Vector<Double>& oTime, 00290 const AXIS& eAxisIterUser ); 00291 00292 // Copy constructor and operator=() function 00293 CalStats( const CalStats& oCalStats ); 00294 CalStats& operator=( const CalStats& oCalStats ); 00295 00296 // Destructor 00297 virtual ~CalStats( void ); 00298 00299 // Axis ID states 00300 IPosition& axisIterID( void ) const; 00301 AXIS& axisNonIterID( void ) const; 00302 00303 // Axis value states 00304 Vector<String>& axisIterFeed( void ) const; 00305 Vector<Double>& axisIterUser( void ) const; 00306 Vector<Double>& axisNonIter( void ) const; 00307 00308 // Output statistics cube shape state 00309 IPosition& statsShape( void ) const; 00310 00311 // Input data states 00312 Cube<Double>& value( void ) const; 00313 Cube<Double>& valueErr( void ) const; 00314 Cube<Bool>& flag( void ) const; 00315 00316 // The axis names 00317 static String& axisName( const AXIS& eAxis ); 00318 00319 // Calculate statistics (allowed T: CalStats::NONE gets data without 00320 // calculating statistics, CalStatsFitter::FIT calculates fits, and 00321 // CalStatsHist::HIST calculates histogram statistics). Member function 00322 // stats() is the main user interface and statsWrap() is the supporting 00323 // wrapper. 00324 template <typename T> Matrix<OUT<T> >& stats( const ARG<T>& oArg ); 00325 template <typename T> T& statsWrap( const Vector<Double>& oAbs, 00326 const Vector<Double>& oValue, const Vector<Double>& oValueErr, 00327 Vector<Bool>& oFlag, const ARG<T>& oArg ); 00328 00329 protected: 00330 00331 // The axis IDs. The two iteration axes are FEED (always) and either TIME 00332 // or FREQUENCY (user defined). The non-iteration axis is either FREQUENCY 00333 // or TIME (the opposite of the user-defined iteration axis). 00334 IPosition oAxisIterID; 00335 AXIS eAxisNonIterID; 00336 00337 // Internal copies of the iteration and non-iteration axis values 00338 Vector<String> oAxisIterFeed; // Feed axis iteration axis values 00339 Vector<Double> oAxisIterUser; // User-defined iteration axis values 00340 Vector<Double> oAxisNonIter; // Non-iteration axis values 00341 00342 // Shape of the output statistics cubes 00343 IPosition oStatsShape; 00344 00345 // Internal copies of input parameter cubes 00346 Cube<Double>* poValue; 00347 Cube<Double>* poValueErr; 00348 Cube<Bool>* poFlag; 00349 00350 // Input parameter cube iterators 00351 ArrayIterator<Double>* poValueIter; 00352 ArrayIterator<Double>* poValueErrIter; 00353 ArrayIterator<Bool>* poFlagIter; 00354 00355 // Unused constructor 00356 CalStats( void ); 00357 00358 // Simultaneously increment and reset all input parameter cube iterators 00359 void next( void ); 00360 void reset( void ); 00361 00362 }; 00363 00364 // ----------------------------------------------------------------------------- 00365 // End of CalStats class definition 00366 // ----------------------------------------------------------------------------- 00367 00368 // ----------------------------------------------------------------------------- 00369 // Start of ARG<T> specialized class templates 00370 // ----------------------------------------------------------------------------- 00371 00372 // Specialization for the CalStats::NONE() class. It tells the 00373 // CalStats::stats<CalStats::NONE>() method just to return the data with no 00374 // processing. 00375 template <> class CalStats::ARG<CalStats::NONE> {}; 00376 00377 // Specialization for the CalStatsFitter::FIT() class. It tells the 00378 // CalStats::stats<CalStatsFitter::FIT>() method to perform fits and return the 00379 // data and fit parameters. 00380 template <> class CalStats::ARG<CalStatsFitter::FIT> { 00381 public: 00382 CalStatsFitter::ORDER eOrder; 00383 CalStatsFitter::TYPE eType; 00384 CalStatsFitter::WEIGHT eWeight; 00385 ARG( void ) { 00386 eOrder = CalStatsFitter::ORDER_INIT; 00387 eType = CalStatsFitter::TYPE_INIT; 00388 eWeight = CalStatsFitter::WEIGHT_INIT; 00389 return; 00390 } 00391 }; 00392 00393 // ----------------------------------------------------------------------------- 00394 // End of ARG<T> specialized class templates 00395 // ----------------------------------------------------------------------------- 00396 00397 // ----------------------------------------------------------------------------- 00398 // Start of OUT<T> specialized class template public member functions 00399 // ----------------------------------------------------------------------------- 00400 00401 // Default constructor 00402 template <typename T> 00403 CalStats::OUT<T>::OUT( void ) { 00404 oAxes = CalStats::AXES(); 00405 oData = CalStats::DATA(); 00406 oT = T(); 00407 return; 00408 } 00409 00410 // Copy constructor 00411 template <typename T> 00412 CalStats::OUT<T>::OUT( const CalStats::OUT<T>& oOut ) { 00413 oAxes = CalStats::AXES( oOut.oAxes ); 00414 oData = CalStats::DATA( oOut.oData ); 00415 oT = T( oOut.oT ); 00416 return; 00417 } 00418 00419 // Destructor 00420 template <typename T> 00421 CalStats::OUT<T>::~OUT( void ) {} 00422 00423 // operator= 00424 template <typename T> 00425 CalStats::OUT<T>& CalStats::OUT<T>::operator=( const CalStats::OUT<T>& oOut ) { 00426 if ( this != &oOut ) { 00427 oAxes = CalStats::AXES( oOut.oAxes ); 00428 oData = CalStats::DATA( oOut.oData ); 00429 oT = T( oOut.oT ); 00430 } 00431 return( *this ); 00432 } 00433 00434 // ----------------------------------------------------------------------------- 00435 // End of OUT<T> specialized class template public member functions 00436 // ----------------------------------------------------------------------------- 00437 00438 // ----------------------------------------------------------------------------- 00439 // Start of CalStats::stats<T> template public statistics member function 00440 // ----------------------------------------------------------------------------- 00441 00442 /* 00443 00444 CalStats::stats<T> 00445 00446 Description: 00447 ------------ 00448 This member function calculates the desired statistics. The allowed templates 00449 are CalStats::NONE (no statistics, just the input data), CalStatsFitter::FIT 00450 (fit statistics), and CalStatsHist::HIST (histogram statistics). 00451 00452 Inputs: 00453 ------- 00454 oArg - This reference to a CalStats::ARG<T> instance contains the extra input 00455 parameters. 00456 00457 Outputs: 00458 -------- 00459 The reference to the Matrix<CalStats::OUT<T> > instance containing the 00460 statistics, returned via the function value. 00461 00462 Modification history: 00463 --------------------- 00464 2011 Nov 11 - Nick Elias, NRAO 00465 Initial version. This template member function replaces the 00466 getData() and calcFit() member functions. 00467 2012 Jan 25 - Nick Elias, NRAO 00468 Logging capability added. 00469 00470 */ 00471 00472 // ----------------------------------------------------------------------------- 00473 00474 template <typename T> 00475 Matrix<CalStats::OUT<T> >& CalStats::stats( const CalStats::ARG<T>& oArg ) { 00476 00477 // Initialize the CalStats::OUT<T> array and its iterator 00478 00479 Array<CalStats::OUT<T> >* poOut = new Array<OUT<T> >( oStatsShape ); 00480 00481 ArrayIterator<CalStats::OUT<T> > oOutIter( *poOut, oAxisIterID, False ); 00482 00483 00484 // For each iteration, convert the resulting arrays to vectors and feed them 00485 // to the CalStatsFitter::fit() function 00486 00487 while ( !poValueIter->pastEnd() ) { 00488 00489 IPosition oPos( poValueIter->pos() ); 00490 00491 uInt uiLength = poValueIter->array().nelements(); 00492 IPosition oShape( 1, uiLength ); 00493 00494 Vector<Double> oAbs( oAxisNonIter.copy() ); 00495 Vector<Double> oValue( poValueIter->array().copy().reform(oShape) ); 00496 Vector<Double> oValueErr( poValueErrIter->array().copy().reform(oShape) ); 00497 Vector<Bool> oFlag( poFlagIter->array().copy().reform(oShape) ); 00498 00499 CalStats::OUT<T> oOut; 00500 00501 oOut.oAxes.eAxisIterFeedID = (CalStats::AXIS) oAxisIterID[0]; 00502 oOut.oAxes.eAxisIterUserID = (CalStats::AXIS) oAxisIterID[1]; 00503 oOut.oAxes.eAxisNonIterID = eAxisNonIterID; 00504 oOut.oAxes.sFeed = String( oAxisIterFeed[oPos[0]] ); 00505 oOut.oAxes.dAxisIterUser = oAxisIterUser[oPos[oAxisIterID[1]]]; 00506 00507 oOut.oData.oAbs = Vector<Double>( oAbs ); 00508 oOut.oData.oValue = Vector<Double>( oValue ); 00509 oOut.oData.oValueErr = Vector<Double>( oValueErr ); 00510 00511 try { 00512 oOut.oT = statsWrap<T>( oAbs, oValue, oValueErr, oFlag, oArg ); 00513 } 00514 00515 catch ( AipsError oAE ) { 00516 LogIO log( LogOrigin( "CalStats", "stats<T>()", WHERE ) ); 00517 log << LogIO::WARN << oAE.getMesg() << ", iteration: " 00518 << oPos.asVector() << ", continuing ..." << LogIO::POST; 00519 oOut.oT = T(); 00520 } 00521 00522 // The flag output vector is set here because robust fitting can change them 00523 oOut.oData.oFlag = Vector<Bool>( oFlag ); 00524 00525 oOutIter.array() = Vector<CalStats::OUT<T> >( 1, oOut ); 00526 00527 next(); oOutIter.next(); 00528 00529 } 00530 00531 00532 // Reset the input parameter iterators 00533 00534 reset(); 00535 00536 00537 // Return the reference to the Matrix<CalStats::OUT<T> > instance 00538 00539 poOut->removeDegenerate(); 00540 Matrix<CalStats::OUT<T> >* poMatrix = (Matrix<CalStats::OUT<T> >*) poOut; 00541 00542 return( *poMatrix ); 00543 00544 } 00545 00546 // ----------------------------------------------------------------------------- 00547 // End of CalStats::stats<T> template public statistics member function 00548 // ----------------------------------------------------------------------------- 00549 00550 }; 00551 00552 // ----------------------------------------------------------------------------- 00553 // End of casa namespace 00554 // ----------------------------------------------------------------------------- 00555 00556 #endif 00557 00558 // ----------------------------------------------------------------------------- 00559 // End of define macro to prevent multiple loading 00560 // -----------------------------------------------------------------------------