Changeset 2703 for trunk/src/STCalTsys.cpp

Timestamp:

12/19/12 19:19:53 (11 years ago)

Author:

Takeshi Nakazato

Message:

New Development: Yes

JIRA Issue: Yes CAS-4770, 4771, 4772

Ready for Test: Yes

Interface Changes: No

What Interface Changed: Please list interface changes

Test Programs: List test programs

Put in Release Notes: Yes/No?

Module(s): Module Names change impacts.

Description: Describe your changes here...

First version of calibration classes for sky calibration (STCalSkyPSAlma)
and Tsys calibration (STCalTsys), which are derived from base class
(STCalibration).

The caltable classes, STCalSky and STCalTsys, are renamed as STCalSkyTable
and STCalTsysTable, respectively.

File:

• trunk/src/STCalTsys.cpp (modified) (1 diff)

trunk/src/STCalTsys.cpp

@@ -5 +5 @@
 //
 //
-// Author: Takeshi Nakazato
+// Author: Takeshi Nakazato <takeshi.nakazato@nao.ac.jp> (C) 2012
 //
 // Copyright: See COPYING file that comes with this distribution
 //
 //
-#include <casa/Exceptions/Error.h>
-#include <tables/Tables/TableDesc.h>
-#include <tables/Tables/SetupNewTab.h>
-#include <tables/Tables/ArrColDesc.h>
-#include <tables/Tables/ScaColDesc.h>
-#include <tables/Tables/TableRecord.h>
-#include <measures/TableMeasures/TableMeasDesc.h>
-#include <measures/TableMeasures/TableMeasRefDesc.h>
-#include <measures/TableMeasures/TableMeasValueDesc.h>
 
-#include "Scantable.h"
+#include <vector>
+#include "STSelector.h"
 #include "STCalTsys.h"
+#include "RowAccumulator.h"
+#include "STIdxIter.h"
+#include "STDefs.h"
+#include <atnf/PKSIO/SrcType.h>
 
+#include <casa/Arrays/ArrayMath.h>
 
+using namespace std;
 using namespace casa;
 
 namespace asap {
-
-const String STCalTsys::name_ = "APPLY_TSYS";
-
-STCalTsys::STCalTsys(const Scantable& parent)
-  : STApplyTable(parent, name_)
+  STCalTsys::STCalTsys(CountedPtr<Scantable> &s, vector<int> &iflist)
+    : STCalibration(s),
+      iflist_(iflist)
 {
-  setup();
+  applytable_ = new STCalTsysTable(*s);
 }
 
-STCalTsys::~STCalTsys()
+void STCalTsys::calibrate()
 {
+  STSelector selOrg = scantable_->getSelection();
+  STSelector sel;
+  sel.setIFs(iflist_);
+  scantable_->setSelection(sel);
+  
+  fillCalTable();
+
+  scantable_->setSelection(selOrg);
 }
 
-void STCalTsys::setup()
+void STCalTsys::fillCalTable()
 {
-  table_.addColumn(ArrayColumnDesc<Float>("TSYS"));
-  table_.addColumn(ScalarColumnDesc<Float>("ELEVATION"));
+  RowAccumulator acc(W_TINT);
+  
+  vector<string> cols(3);
+  cols[0] = "IFNO";
+  cols[1] = "POLNO";
+  cols[2] = "BEAMNO";
+  STIdxIterAcc iter(scantable_, cols);
 
-  table_.rwKeywordSet().define("ApplyType", "TSYS");
+  ROScalarColumn<Double> *tcol = new ROScalarColumn<Double>(scantable_->table(), "TIME");
+  Vector<Double> timeSec = tcol->getColumn() * 86400.0;
+  tcol->attach(scantable_->table(), "INTERVAL");
+  Vector<Double> intervalSec = tcol->getColumn();
+  delete tcol;
+  ROScalarColumn<Float> *ecol = new ROScalarColumn<Float>(scantable_->table(), "ELEVATION");
+  Vector<Float> elevation = ecol->getColumn();
+  delete ecol;
 
-  attachOptionalColumns();
+  ROArrayColumn<Float> specCol(scantable_->table(), "TSYS");
+  ROArrayColumn<uChar> flagCol(scantable_->table(), "FLAGTRA");
+
+  // dummy Tsys: the following process doesn't need Tsys but RowAccumulator 
+  //             requires to set it with spectral data
+  Vector<Float> tsys(1, 1.0);
+
+  Double timeCen = 0.0;
+  Float elCen = 0.0;
+  uInt count = 0;
+
+  while(!iter.pastEnd()) {
+    Vector<uInt> rows = iter.getRows(SHARE);
+    Vector<uInt> current = iter.current();
+    uInt len = rows.nelements();
+    if (len == 0) {
+      iter.next();
+      continue;
+    }
+    
+    uInt nchan = scantable_->nchan(scantable_->getIF(rows[0]));
+    Vector<uChar> flag(nchan);
+    Vector<Bool> bflag(nchan);
+    Vector<Float> spec(nchan);
+
+    Vector<Double> timeSep(len); 
+    for (uInt i = 0; i < len-1; i++) {
+      timeSep[i] = timeSec[rows[i+1]] - timeSec[rows[i]] ;
+    }
+    Double tMedian = median(timeSep(IPosition(1,0), IPosition(1,len-2)));
+    timeSep[len-1] = tMedian * 10000.0 ; // any large value
+
+    uInt irow ;
+    uInt jrow ;
+    for (uInt i = 0; i < len; i++) {
+      irow = rows[i];
+      jrow = (i < len-1) ? rows[i+1] : rows[i];
+      // accumulate data
+      flagCol.get(irow, flag);
+      convertArray(bflag, flag);
+      specCol.get(irow, spec);
+      if ( !allEQ(bflag,True) ) 
+        acc.add( spec, !bflag, tsys, intervalSec[irow], timeSec[irow] ) ;
+      timeCen += timeSec[irow];
+      elCen += elevation[irow];
+      count++;
+
+      // check time gap
+      double gap = 2.0 * timeSep[i] / (intervalSec[jrow] + intervalSec[irow]);
+      if ( gap > 5.0 ) {
+        if ( acc.state() ) {
+          acc.replaceNaN() ;
+//           const Vector<Bool> &msk = acc.getMask();
+//           convertArray(flag, !msk);
+//           for (uInt k = 0; k < nchan; ++k) {
+//             uChar userFlag = 1 << 7;
+//             if (msk[k]==True) userFlag = 0 << 7;
+//             flag(k) = userFlag;
+//           }
+          timeCen /= (Double)count * 86400.0; // sec->day
+          elCen /= (Float)count;
+          STCalTsysTable *p = dynamic_cast<STCalTsysTable *>(&(*applytable_));
+          p->appenddata(0, 0, current[2], current[0], current[1],
+                        timeCen, elCen, acc.getSpectrum());
+        }
+        acc.reset() ;
+        timeCen = 0.0;
+        elCen = 0.0;
+        count = 0;
+      }
+    }
+
+    iter.next() ;
+  }
 }
 
-void STCalTsys::attachOptionalColumns()
-{
-  tsysCol_.attach(table_, "TSYS");
-  elCol_.attach(table_,"ELEVATION");
-  
 }
-
-void STCalTsys::setdata(uInt irow, uInt scanno, uInt cycleno, 
-                        uInt beamno, uInt ifno, uInt polno, 
-                        Double time, Float elevation, Vector<Float> tsys)
-{
-  if (irow >= (uInt)nrow()) {
-    throw AipsError("row index out of range");
-  }
-
-  if (!sel_.empty()) {
-    os_.origin(LogOrigin("STCalTsys","setdata",WHERE));
-    os_ << LogIO::WARN << "Data selection is effective. Specified row index may be wrong." << LogIO::POST;
-  }  
-
-  setbasedata(irow, scanno, cycleno, beamno, ifno, polno, time);
-  elCol_.put(irow, elevation);
-  tsysCol_.put(irow, tsys);
-}
-
-void STCalTsys::appenddata(uInt scanno, uInt cycleno, 
-                           uInt beamno, uInt ifno, uInt polno, 
-                           Double time, Float elevation, Vector<Float> tsys)
-{
-  uInt irow = nrow();
-  table_.addRow(1, True);
-  setdata(irow, scanno, cycleno, beamno, ifno, polno, time, elevation, tsys);
-}
-}