#include #include #include #include #include "Scantable.h" #include "STTcal.h" #include "STIdxIter.h" #include "STSelector.h" using namespace casa; using namespace asap; namespace { // Iteration Helper // template // class IterationHelper // { // public: // static void Iterate(T &calibrator) // { // vector cols( 3 ) ; // cols[0] = "BEAMNO" ; // cols[1] = "POLNO" ; // cols[2] = "IFNO" ; // STIdxIter2 calibrator.GetIterator(cols);//iter( out, cols ) ; // STSelector sel ; // while ( !iter.pastEnd() ) { // Record ids = iter.currentValue() ; // stringstream ss ; // ss << "SELECT FROM $1 WHERE " // << "BEAMNO==" << ids.asuInt(cols[0]) << "&&" // << "POLNO==" << ids.asuInt(cols[1]) << "&&" // << "IFNO==" << ids.asuInt(cols[2]) ; // //cout << "TaQL string: " << ss.str() << endl ; // sel.setTaQL( ss.str() ) ; // Vector rows = iter.getRows( SHARE ) ; // calibrator.Calibrate(sel, rows); // aoff->setSelection( sel ) ; // // out should be an exact copy of s except that SPECTRA column is empty // calibrateALMA( out, s, aoff, rows ) ; // aoff->unsetSelection() ; // sel.reset() ; // iter.next() ; // } // } // }; // Interpolation Helper class TcalData { public: TcalData(CountedPtr s) : table_(s) {} ~TcalData() {} const String method_name() const {return "getTcalFromTime";} uInt nrow() const {return table_->nrow();} Vector GetEntry(int idx) const { String time; uInt tcalid = table_->getTcalId(idx); Vector return_value; table_->tcal().getEntry(time, return_value, tcalid); return return_value; } private: CountedPtr table_; }; class TsysData { public: TsysData(CountedPtr s) : tsyscolumn_(s->table(), "TSYS") {} ~TsysData() {} const String method_name() const {return "getTsysFromTime";} uInt nrow() const {return tsyscolumn_.nrow();} Vector GetEntry(int idx) const { return tsyscolumn_(idx); } private: ROArrayColumn tsyscolumn_; }; class SpectralData { public: SpectralData(Matrix s) : data_(s) {} ~SpectralData() {} const String method_name() const {return "getSpectraFromTime";} uInt nrow() const {return data_.ncolumn();} Vector GetEntry(int idx) const { return data_.column(idx); } private: Matrix data_; }; vector getRowIdFromTime(double reftime, const Vector &t) { // double reft = reftime ; double dtmin = 1.0e100 ; double dtmax = -1.0e100 ; // vector dt ; int just_before = -1 ; int just_after = -1 ; Vector dt = t - reftime ; for ( unsigned int i = 0 ; i < dt.size() ; i++ ) { if ( dt[i] > 0.0 ) { // after reftime if ( dt[i] < dtmin ) { just_after = i ; dtmin = dt[i] ; } } else if ( dt[i] < 0.0 ) { // before reftime if ( dt[i] > dtmax ) { just_before = i ; dtmax = dt[i] ; } } else { // just a reftime just_before = i ; just_after = i ; dtmax = 0 ; dtmin = 0 ; break ; } } vector v(2) ; v[0] = just_before ; v[1] = just_after ; return v ; } template class SimpleInterpolationHelper { public: static Vector GetFromTime(double reftime, const Vector &timeVec, const vector &idx, const T &data, const string mode) { Vector return_value; LogIO os_; LogIO os( LogOrigin( "STMath", data.method_name(), WHERE ) ) ; if ( data.nrow() == 0 ) { os << LogIO::SEVERE << "No row in the input scantable. Return empty tcal." << LogIO::POST ; } else if ( data.nrow() == 1 ) { return_value = data.GetEntry(0); } else { if ( mode == "before" ) { int id = -1 ; if ( idx[0] != -1 ) { id = idx[0] ; } else if ( idx[1] != -1 ) { os << LogIO::WARN << "Failed to find a scan before reftime. return a spectrum just after the reftime." << LogIO::POST ; id = idx[1] ; } return_value = data.GetEntry(id); } else if ( mode == "after" ) { int id = -1 ; if ( idx[1] != -1 ) { id = idx[1] ; } else if ( idx[0] != -1 ) { os << LogIO::WARN << "Failed to find a scan after reftime. return a spectrum just before the reftime." << LogIO::POST ; id = idx[1] ; } return_value = data.GetEntry(id); } else if ( mode == "nearest" ) { int id = -1 ; if ( idx[0] == -1 ) { id = idx[1] ; } else if ( idx[1] == -1 ) { id = idx[0] ; } else if ( idx[0] == idx[1] ) { id = idx[0] ; } else { double t0 = timeVec[idx[0]] ; double t1 = timeVec[idx[1]] ; if ( abs( t0 - reftime ) > abs( t1 - reftime ) ) { id = idx[1] ; } else { id = idx[0] ; } } return_value = data.GetEntry(id); } else if ( mode == "linear" ) { if ( idx[0] == -1 ) { // use after os << LogIO::WARN << "Failed to interpolate. return a spectrum just after the reftime." << LogIO::POST ; int id = idx[1] ; return_value = data.GetEntry(id); } else if ( idx[1] == -1 ) { // use before os << LogIO::WARN << "Failed to interpolate. return a spectrum just before the reftime." << LogIO::POST ; int id = idx[0] ; return_value = data.GetEntry(id); } else if ( idx[0] == idx[1] ) { // use before //os << "No need to interporate." << LogIO::POST ; int id = idx[0] ; return_value = data.GetEntry(id); } else { // do interpolation double t0 = timeVec[idx[0]] ; double t1 = timeVec[idx[1]] ; Vector value0 = data.GetEntry(idx[0]); Vector value1 = data.GetEntry(idx[1]); double tfactor = (reftime - t0) / (t1 - t0) ; for ( unsigned int i = 0 ; i < value0.size() ; i++ ) { value1[i] = ( value1[i] - value0[i] ) * tfactor + value0[i] ; } return_value = value1; } } else { os << LogIO::SEVERE << "Unknown mode" << LogIO::POST ; } } return return_value ; } }; } // anonymous namespace