source: trunk/src/Scantable.cpp @ 889

Last change on this file since 889 was 888, checked in by mar637, 18 years ago

Fixed all nXYZ() functions, to get row 0 NOT the number 0

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 24.8 KB
RevLine 
[805]1//
2// C++ Implementation: Scantable
3//
4// Description:
5//
6//
7// Author: Malte Marquarding <asap@atnf.csiro.au>, (C) 2005
8//
9// Copyright: See COPYING file that comes with this distribution
10//
11//
[206]12#include <map>
13
[125]14#include <casa/aips.h>
[80]15#include <casa/iostream.h>
16#include <casa/iomanip.h>
[805]17#include <casa/OS/Path.h>
18#include <casa/OS/File.h>
[80]19#include <casa/Arrays/Array.h>
20#include <casa/Arrays/ArrayMath.h>
21#include <casa/Arrays/MaskArrMath.h>
22#include <casa/Arrays/ArrayLogical.h>
23#include <casa/Arrays/ArrayAccessor.h>
[455]24#include <casa/Arrays/VectorSTLIterator.h>
[206]25#include <casa/Arrays/Vector.h>
[418]26#include <casa/BasicMath/Math.h>
[504]27#include <casa/BasicSL/Constants.h>
[286]28#include <casa/Quanta/MVAngle.h>
[805]29#include <casa/Containers/RecordField.h>
[2]30
[80]31#include <tables/Tables/TableParse.h>
32#include <tables/Tables/TableDesc.h>
[488]33#include <tables/Tables/TableCopy.h>
[80]34#include <tables/Tables/SetupNewTab.h>
35#include <tables/Tables/ScaColDesc.h>
36#include <tables/Tables/ArrColDesc.h>
[805]37#include <tables/Tables/TableRow.h>
38#include <tables/Tables/TableVector.h>
39#include <tables/Tables/TableIter.h>
[2]40
[80]41#include <tables/Tables/ExprNode.h>
42#include <tables/Tables/TableRecord.h>
43#include <measures/Measures/MFrequency.h>
[805]44#include <measures/Measures/MEpoch.h>
[80]45#include <measures/Measures/MeasTable.h>
[805]46#include <measures/Measures/MeasRef.h>
47#include <measures/TableMeasures/TableMeasRefDesc.h>
48#include <measures/TableMeasures/TableMeasValueDesc.h>
49#include <measures/TableMeasures/TableMeasDesc.h>
50#include <measures/TableMeasures/ScalarMeasColumn.h>
[105]51#include <coordinates/Coordinates/CoordinateUtil.h>
[80]52#include <casa/Quanta/MVTime.h>
[281]53#include <casa/Quanta/MVAngle.h>
[2]54
[805]55#include "Scantable.h"
[878]56#include "STAttr.h"
[2]57
[125]58using namespace casa;
[2]59
[805]60namespace asap {
61
62Scantable::Scantable(Table::TableType ttype) :
[852]63  type_(ttype)
[206]64{
[805]65  setupMainTable();
[852]66  freqTable_ = STFrequencies(*this);
[805]67  table_.rwKeywordSet().defineTable("FREQUENCIES", freqTable_.table());
[852]68  weatherTable_ = STWeather(*this);
[805]69  table_.rwKeywordSet().defineTable("WEATHER", weatherTable_.table());
[852]70  focusTable_ = STFocus(*this);
[805]71  table_.rwKeywordSet().defineTable("FOCUS", focusTable_.table());
[852]72  tcalTable_ = STTcal(*this);
[805]73  table_.rwKeywordSet().defineTable("TCAL", tcalTable_.table());
[852]74  moleculeTable_ = STMolecules(*this);
[805]75  table_.rwKeywordSet().defineTable("MOLECULES", moleculeTable_.table());
[860]76  historyTable_ = STHistory(*this);
77  table_.rwKeywordSet().defineTable("HISTORY", historyTable_.table());
[805]78  setupFitTable();
[852]79  fitTable_ = table_.keywordSet().asTable("FITS");
[805]80  originalTable_ = table_;
[322]81  attach();
[18]82}
[206]83
[805]84Scantable::Scantable(const std::string& name, Table::TableType ttype) :
[852]85  type_(ttype)
[206]86{
[865]87  Table tab(name, Table::Update);
[499]88  Int version;
[483]89  tab.keywordSet().get("VERSION", version);
90  if (version != version_) {
91    throw(AipsError("Unsupported version of ASAP file."));
92  }
[805]93  if ( type_ == Table::Memory )
[852]94    table_ = tab.copyToMemoryTable(generateName());
[805]95  else
96    table_ = tab;
[859]97  attachSubtables();
[805]98  originalTable_ = table_;
[329]99  attach();
[2]100}
101
[805]102Scantable::Scantable( const Scantable& other, bool clear )
[206]103{
[805]104  // with or without data
[859]105  String newname = String(generateName());
[865]106  type_ = other.table_.tableType();
[859]107  if ( other.table_.tableType() == Table::Memory ) {
108      if ( clear ) {
109        table_ = TableCopy::makeEmptyMemoryTable(newname,
110                                                 other.table_, True);
111      } else
112        table_ = other.table_.copyToMemoryTable(newname);
[16]113  } else {
[865]114      other.table_.deepCopy(newname, Table::New, False, Table::AipsrcEndian,
115                            Bool(clear));
116      table_ = Table(newname, Table::Update);
117      copySubtables(other);
118      table_.markForDelete();
119  }
[859]120
121  attachSubtables();
[805]122  originalTable_ = table_;
[322]123  attach();
[2]124}
125
[865]126void Scantable::copySubtables(const Scantable& other) {
127  Table t = table_.rwKeywordSet().asTable("FREQUENCIES");
128  TableCopy::copyRows(t, other.freqTable_.table());
129  t = table_.rwKeywordSet().asTable("FOCUS");
130  TableCopy::copyRows(t, other.focusTable_.table());
131  t = table_.rwKeywordSet().asTable("WEATHER");
132  TableCopy::copyRows(t, other.weatherTable_.table());
133  t = table_.rwKeywordSet().asTable("TCAL");
134  TableCopy::copyRows(t, other.tcalTable_.table());
135  t = table_.rwKeywordSet().asTable("MOLECULES");
136  TableCopy::copyRows(t, other.moleculeTable_.table());
137  t = table_.rwKeywordSet().asTable("HISTORY");
138  TableCopy::copyRows(t, other.historyTable_.table());
139}
140
[859]141void Scantable::attachSubtables()
142{
143  freqTable_ = STFrequencies(table_);
144  focusTable_ = STFocus(table_);
145  weatherTable_ = STWeather(table_);
146  tcalTable_ = STTcal(table_);
147  moleculeTable_ = STMolecules(table_);
[860]148  historyTable_ = STHistory(table_);
[859]149}
150
[805]151Scantable::~Scantable()
[206]152{
[805]153  cout << "~Scantable() " << this << endl;
[2]154}
155
[805]156void Scantable::setupMainTable()
[206]157{
[805]158  TableDesc td("", "1", TableDesc::Scratch);
159  td.comment() = "An ASAP Scantable";
160  td.rwKeywordSet().define("VERSION", Int(version_));
[2]161
[805]162  // n Cycles
163  td.addColumn(ScalarColumnDesc<uInt>("SCANNO"));
164  // new index every nBeam x nIF x nPol
165  td.addColumn(ScalarColumnDesc<uInt>("CYCLENO"));
[2]166
[805]167  td.addColumn(ScalarColumnDesc<uInt>("BEAMNO"));
168  td.addColumn(ScalarColumnDesc<uInt>("IFNO"));
169  td.rwKeywordSet().define("POLTYPE", String("linear"));
170  td.addColumn(ScalarColumnDesc<uInt>("POLNO"));
[138]171
[805]172  td.addColumn(ScalarColumnDesc<uInt>("FREQ_ID"));
173  td.addColumn(ScalarColumnDesc<uInt>("MOLECULE_ID"));
174  // linear, circular, stokes [I Q U V], stokes1 [I Plinear Pangle V]
175  td.addColumn(ScalarColumnDesc<Int>("REFBEAMNO"));
[80]176
[805]177  td.addColumn(ScalarColumnDesc<Double>("TIME"));
178  TableMeasRefDesc measRef(MEpoch::UTC); // UTC as default
179  TableMeasValueDesc measVal(td, "TIME");
180  TableMeasDesc<MEpoch> mepochCol(measVal, measRef);
181  mepochCol.write(td);
[483]182
[805]183  td.addColumn(ScalarColumnDesc<Double>("INTERVAL"));
184
[2]185  td.addColumn(ScalarColumnDesc<String>("SRCNAME"));
[805]186  // Type of source (on=0, off=1, other=-1)
187  td.addColumn(ScalarColumnDesc<Int>("SRCTYPE", Int(-1)));
188  td.addColumn(ScalarColumnDesc<String>("FIELDNAME"));
189
190  //The actual Data Vectors
[2]191  td.addColumn(ArrayColumnDesc<Float>("SPECTRA"));
192  td.addColumn(ArrayColumnDesc<uChar>("FLAGTRA"));
[89]193  td.addColumn(ArrayColumnDesc<Float>("TSYS"));
[805]194
195  td.addColumn(ArrayColumnDesc<Double>("DIRECTION",
196                                       IPosition(1,2),
197                                       ColumnDesc::Direct));
198  TableMeasRefDesc mdirRef(MDirection::J2000); // default
199  TableMeasValueDesc tmvdMDir(td, "DIRECTION");
200  // the TableMeasDesc gives the column a type
201  TableMeasDesc<MDirection> mdirCol(tmvdMDir, mdirRef);
202  // writing create the measure column
203  mdirCol.write(td);
204  td.addColumn(ScalarColumnDesc<Double>("AZIMUTH"));
205  td.addColumn(ScalarColumnDesc<Double>("ELEVATION"));
[105]206  td.addColumn(ScalarColumnDesc<Float>("PARANGLE"));
207
[805]208  td.addColumn(ScalarColumnDesc<uInt>("TCAL_ID"));
209  td.addColumn(ScalarColumnDesc<uInt>("FIT_ID"));
210
211  td.addColumn(ScalarColumnDesc<uInt>("FOCUS_ID"));
212  td.addColumn(ScalarColumnDesc<uInt>("WEATHER_ID"));
213
214  td.rwKeywordSet().define("OBSMODE", String(""));
215
[418]216  // Now create Table SetUp from the description.
[859]217  SetupNewTable aNewTab(generateName(), td, Table::Scratch);
[852]218  table_ = Table(aNewTab, type_, 0);
[805]219  originalTable_ = table_;
[745]220
[805]221}
[418]222
[805]223void Scantable::setupFitTable()
[322]224{
[39]225  TableDesc td("", "1", TableDesc::Scratch);
[805]226  td.addColumn(ScalarColumnDesc<uInt>("FIT_ID"));
[455]227  td.addColumn(ArrayColumnDesc<String>("FUNCTIONS"));
228  td.addColumn(ArrayColumnDesc<Int>("COMPONENTS"));
229  td.addColumn(ArrayColumnDesc<Double>("PARAMETERS"));
230  td.addColumn(ArrayColumnDesc<Bool>("PARMASK"));
[465]231  td.addColumn(ArrayColumnDesc<String>("FRAMEINFO"));
[859]232  SetupNewTable aNewTab("fits", td, Table::Scratch);
[455]233  Table aTable(aNewTab, Table::Memory);
234  table_.rwKeywordSet().defineTable("FITS", aTable);
235}
236
[805]237void Scantable::attach()
[455]238{
[805]239  timeCol_.attach(table_, "TIME");
240  srcnCol_.attach(table_, "SRCNAME");
241  specCol_.attach(table_, "SPECTRA");
242  flagsCol_.attach(table_, "FLAGTRA");
[865]243  tsysCol_.attach(table_, "TSYS");
[805]244  cycleCol_.attach(table_,"CYCLENO");
245  scanCol_.attach(table_, "SCANNO");
246  beamCol_.attach(table_, "BEAMNO");
[847]247  ifCol_.attach(table_, "IFNO");
248  polCol_.attach(table_, "POLNO");
[805]249  integrCol_.attach(table_, "INTERVAL");
250  azCol_.attach(table_, "AZIMUTH");
251  elCol_.attach(table_, "ELEVATION");
252  dirCol_.attach(table_, "DIRECTION");
253  paraCol_.attach(table_, "PARANGLE");
254  fldnCol_.attach(table_, "FIELDNAME");
255  rbeamCol_.attach(table_, "REFBEAMNO");
[455]256
[805]257  mfitidCol_.attach(table_,"FIT_ID");
[859]258  //fitidCol_.attach(fitTable_,"FIT_ID");
[465]259
[805]260  mfreqidCol_.attach(table_, "FREQ_ID");
[465]261
[805]262  mtcalidCol_.attach(table_, "TCAL_ID");
[465]263
[805]264  mfocusidCol_.attach(table_, "FOCUS_ID");
[455]265
[805]266  mmolidCol_.attach(table_, "MOLECULE_ID");
[455]267}
268
[845]269void Scantable::setHeader(const SDHeader& sdh)
[206]270{
[18]271  table_.rwKeywordSet().define("nIF", sdh.nif);
272  table_.rwKeywordSet().define("nBeam", sdh.nbeam);
273  table_.rwKeywordSet().define("nPol", sdh.npol);
274  table_.rwKeywordSet().define("nChan", sdh.nchan);
275  table_.rwKeywordSet().define("Observer", sdh.observer);
276  table_.rwKeywordSet().define("Project", sdh.project);
277  table_.rwKeywordSet().define("Obstype", sdh.obstype);
278  table_.rwKeywordSet().define("AntennaName", sdh.antennaname);
279  table_.rwKeywordSet().define("AntennaPosition", sdh.antennaposition);
280  table_.rwKeywordSet().define("Equinox", sdh.equinox);
281  table_.rwKeywordSet().define("FreqRefFrame", sdh.freqref);
282  table_.rwKeywordSet().define("FreqRefVal", sdh.reffreq);
283  table_.rwKeywordSet().define("Bandwidth", sdh.bandwidth);
284  table_.rwKeywordSet().define("UTC", sdh.utc);
[206]285  table_.rwKeywordSet().define("FluxUnit", sdh.fluxunit);
286  table_.rwKeywordSet().define("Epoch", sdh.epoch);
[50]287}
[21]288
[845]289SDHeader Scantable::getHeader() const
[206]290{
[21]291  SDHeader sdh;
292  table_.keywordSet().get("nBeam",sdh.nbeam);
293  table_.keywordSet().get("nIF",sdh.nif);
294  table_.keywordSet().get("nPol",sdh.npol);
295  table_.keywordSet().get("nChan",sdh.nchan);
296  table_.keywordSet().get("Observer", sdh.observer);
297  table_.keywordSet().get("Project", sdh.project);
298  table_.keywordSet().get("Obstype", sdh.obstype);
299  table_.keywordSet().get("AntennaName", sdh.antennaname);
300  table_.keywordSet().get("AntennaPosition", sdh.antennaposition);
301  table_.keywordSet().get("Equinox", sdh.equinox);
302  table_.keywordSet().get("FreqRefFrame", sdh.freqref);
303  table_.keywordSet().get("FreqRefVal", sdh.reffreq);
304  table_.keywordSet().get("Bandwidth", sdh.bandwidth);
305  table_.keywordSet().get("UTC", sdh.utc);
[206]306  table_.keywordSet().get("FluxUnit", sdh.fluxunit);
307  table_.keywordSet().get("Epoch", sdh.epoch);
[21]308  return sdh;
[18]309}
[805]310
[845]311bool Scantable::conformant( const Scantable& other )
312{
313  return this->getHeader().conformant(other.getHeader());
314}
315
316
[837]317int Scantable::nscan() const {
[805]318  int n = 0;
319  int previous = -1; int current = 0;
[322]320  for (uInt i=0; i< scanCol_.nrow();i++) {
321    scanCol_.getScalar(i,current);
[50]322    if (previous != current) {
[89]323      previous = current;
[50]324      n++;
325    }
326  }
327  return n;
328}
329
[805]330std::string Scantable::formatSec(Double x) const
[206]331{
[105]332  Double xcop = x;
333  MVTime mvt(xcop/24./3600.);  // make days
[365]334
[105]335  if (x < 59.95)
[281]336    return  String("      ") + mvt.string(MVTime::TIME_CLEAN_NO_HM, 7)+"s";
[745]337  else if (x < 3599.95)
[281]338    return String("   ") + mvt.string(MVTime::TIME_CLEAN_NO_H,7)+" ";
339  else {
340    ostringstream oss;
341    oss << setw(2) << std::right << setprecision(1) << mvt.hour();
342    oss << ":" << mvt.string(MVTime::TIME_CLEAN_NO_H,7) << " ";
343    return String(oss);
[745]344  }
[105]345};
346
[805]347std::string Scantable::formatDirection(const MDirection& md) const
[281]348{
349  Vector<Double> t = md.getAngle(Unit(String("rad"))).getValue();
350  Int prec = 7;
351
352  MVAngle mvLon(t[0]);
353  String sLon = mvLon.string(MVAngle::TIME,prec);
354  MVAngle mvLat(t[1]);
355  String sLat = mvLat.string(MVAngle::ANGLE+MVAngle::DIG2,prec);
[380]356  return sLon + String(" ") + sLat;
[281]357}
358
359
[805]360std::string Scantable::getFluxUnit() const
[206]361{
[847]362  return table_.keywordSet().asString("FluxUnit");
[206]363}
364
[805]365void Scantable::setFluxUnit(const std::string& unit)
[218]366{
367  String tmp(unit);
368  Unit tU(tmp);
369  if (tU==Unit("K") || tU==Unit("Jy")) {
370     table_.rwKeywordSet().define(String("FluxUnit"), tmp);
371  } else {
372     throw AipsError("Illegal unit - must be compatible with Jy or K");
373  }
374}
375
[805]376void Scantable::setInstrument(const std::string& name)
[236]377{
[805]378  bool throwIt = true;
[878]379  Instrument ins = STAttr::convertInstrument(name, throwIt);
[236]380  String nameU(name);
381  nameU.upcase();
382  table_.rwKeywordSet().define(String("AntennaName"), nameU);
383}
384
[805]385MPosition Scantable::getAntennaPosition () const
386{
387  Vector<Double> antpos;
388  table_.keywordSet().get("AntennaPosition", antpos);
389  MVPosition mvpos(antpos(0),antpos(1),antpos(2));
390  return MPosition(mvpos);
391}
[281]392
[805]393void Scantable::makePersistent(const std::string& filename)
394{
395  String inname(filename);
396  Path path(inname);
397  inname = path.expandedName();
398  table_.deepCopy(inname, Table::New);
399}
400
[837]401int Scantable::nbeam( int scanno ) const
[805]402{
403  if ( scanno < 0 ) {
404    Int n;
405    table_.keywordSet().get("nBeam",n);
406    return int(n);
[105]407  } else {
[805]408    // take the first POLNO,IFNO,CYCLENO as nbeam shouldn't vary with these
[888]409    Table t = table_(table_.col("SCANNO") == scanno);
410    ROTableRow row(t);
411    const TableRecord& rec = row.get(0);
412    Table subt = t( t.col("IFNO") == Int(rec.asuInt("IFNO"))
413                    && t.col("POLNO") == Int(rec.asuInt("POLNO"))
414                    && t.col("CYCLENO") == Int(rec.asuInt("CYCLENO")) );
415    ROTableVector<uInt> v(subt, "BEAMNO");
[805]416    return int(v.nelements());
[105]417  }
[805]418  return 0;
419}
[455]420
[837]421int Scantable::nif( int scanno ) const
[805]422{
423  if ( scanno < 0 ) {
424    Int n;
425    table_.keywordSet().get("nIF",n);
426    return int(n);
427  } else {
428    // take the first POLNO,BEAMNO,CYCLENO as nbeam shouldn't vary with these
[888]429    Table t = table_(table_.col("SCANNO") == scanno);
430    ROTableRow row(t);
431    const TableRecord& rec = row.get(0);
432    Table subt = t( t.col("BEAMNO") == Int(rec.asuInt("BEAMNO"))
433                    && t.col("POLNO") == Int(rec.asuInt("POLNO"))
434                    && t.col("CYCLENO") == Int(rec.asuInt("CYCLENO")) );
435    if ( subt.nrow() == 0 ) return 0;
436    ROTableVector<uInt> v(subt, "IFNO");
[805]437    return int(v.nelements());
[2]438  }
[805]439  return 0;
440}
[321]441
[837]442int Scantable::npol( int scanno ) const
[805]443{
444  if ( scanno < 0 ) {
445    Int n;
446    table_.keywordSet().get("nPol",n);
447    return n;
448  } else {
449    // take the first POLNO,IFNO,CYCLENO as nbeam shouldn't vary with these
[888]450    Table t = table_(table_.col("SCANNO") == scanno);
451    ROTableRow row(t);
452    const TableRecord& rec = row.get(0);
453    Table subt = t( t.col("BEAMNO") == Int(rec.asuInt("BEAMNO"))
454                    && t.col("IFNO") == Int(rec.asuInt("IFNO"))
455                    && t.col("CYCLENO") == Int(rec.asuInt("CYCLENO")) );
456    if ( subt.nrow() == 0 ) return 0;
457    ROTableVector<uInt> v(subt, "POLNO");
[805]458    return int(v.nelements());
[321]459  }
[805]460  return 0;
[2]461}
[805]462
[845]463int Scantable::ncycle( int scanno ) const
[206]464{
[805]465  if ( scanno < 0 ) {
[837]466    Block<String> cols(2);
467    cols[0] = "SCANNO";
468    cols[1] = "CYCLENO";
469    TableIterator it(table_, cols);
470    int n = 0;
471    while ( !it.pastEnd() ) {
472      ++n;
473    }
474    return n;
[805]475  } else {
[888]476    Table t = table_(table_.col("SCANNO") == scanno);
477    ROTableRow row(t);
478    const TableRecord& rec = row.get(0);
479    Table subt = t( t.col("BEAMNO") == Int(rec.asuInt("BEAMNO"))
480                    && t.col("POLNO") == Int(rec.asuInt("POLNO"))
481                    && t.col("IFNO") == Int(rec.asuInt("IFNO")) );
482    if ( subt.nrow() == 0 ) return 0;
483    return int(subt.nrow());
[805]484  }
485  return 0;
[18]486}
[455]487
488
[845]489int Scantable::nrow( int scanno ) const
[805]490{
[845]491  return int(table_.nrow());
492}
493
494int Scantable::nchan( int ifno ) const
495{
496  if ( ifno < 0 ) {
[805]497    Int n;
498    table_.keywordSet().get("nChan",n);
499    return int(n);
500  } else {
[845]501    // take the first SCANNO,POLNO,BEAMNO,CYCLENO as nbeam shouldn't vary with these
[888]502    Table t = table_(table_.col("IFNO") == ifno);
503    if ( t.nrow() == 0 ) return 0;
504    ROArrayColumn<Float> v(t, "SPECTRA");
[837]505    return v(0).nelements();
[805]506  }
507  return 0;
[18]508}
[455]509
[847]510
511int Scantable::getBeam(int whichrow) const
512{
513  return beamCol_(whichrow);
514}
515
516int Scantable::getIF(int whichrow) const
517{
518  return ifCol_(whichrow);
519}
520
521int Scantable::getPol(int whichrow) const
522{
523  return polCol_(whichrow);
524}
525
[805]526std::string Scantable::formatTime(const MEpoch& me, bool showdate) const
527{
528  MVTime mvt(me.getValue());
529  if (showdate)
530    mvt.setFormat(MVTime::YMD);
531  else
532    mvt.setFormat(MVTime::TIME);
533  ostringstream oss;
534  oss << mvt;
535  return String(oss);
536}
[488]537
[805]538void Scantable::calculateAZEL()
539{
540  MPosition mp = getAntennaPosition();
541  MEpoch::ROScalarColumn timeCol(table_, "TIME");
542  ostringstream oss;
543  oss << "Computed azimuth/elevation using " << endl
544      << mp << endl;
545  for (uInt i=0; i<nrow(); ++i) {
546    MEpoch me = timeCol(i);
547    MDirection md = dirCol_(i);
548    dirCol_.get(i,md);
549    oss  << " Time: " << formatTime(me,False) << " Direction: " << formatDirection(md)
550         << endl << "     => ";
551    MeasFrame frame(mp, me);
552    Vector<Double> azel =
553        MDirection::Convert(md, MDirection::Ref(MDirection::AZEL,
554                                                frame)
555                            )().getAngle("rad").getValue();
556    azCol_.put(i,azel[0]);
557    elCol_.put(i,azel[1]);
558    oss << "azel: " << azel[0]/C::pi*180.0 << " "
559        << azel[1]/C::pi*180.0 << " (deg)" << endl;
[16]560  }
[805]561  pushLog(String(oss));
562}
[89]563
[865]564void Scantable::flag()
565{
566  if ( selector_.empty() )
567    throw(AipsError("Trying to flag whole scantable. Aborted."));
568}
569
[805]570std::vector<bool> Scantable::getMask(int whichrow) const
571{
572  Vector<uChar> flags;
573  flagsCol_.get(uInt(whichrow), flags);
574  Vector<Bool> bflag(flags.shape());
575  convertArray(bflag, flags);
576  bflag = !bflag;
577  std::vector<bool> mask;
578  bflag.tovector(mask);
579  return mask;
580}
[89]581
[805]582std::vector<float> Scantable::getSpectrum(int whichrow) const
583{
584  Vector<Float> arr;
585  specCol_.get(whichrow, arr);
586  std::vector<float> out;
587  arr.tovector(out);
588  return out;
[89]589}
[212]590
[884]591
592void asap::Scantable::setSpectrum( const std::vector<float>& spec,
593                                   int whichrow )
594{
595  Vector<Float> spectrum(spec);
596  Vector<Float> arr;
597  specCol_.get(whichrow, arr);
598  if ( spectrum.nelements() != arr.nelements() )
599    throw AipsError("The specturm has incorrect number of channels.");
600  specCol_.put(whichrow, spectrum);
601}
602
603
[805]604String Scantable::generateName()
[745]605{
[805]606  return (File::newUniqueName("./","temp")).baseName();
[212]607}
608
[805]609const casa::Table& Scantable::table( ) const
[212]610{
[805]611  return table_;
[212]612}
613
[805]614casa::Table& Scantable::table( )
[386]615{
[805]616  return table_;
[386]617}
618
[805]619std::string Scantable::getPolarizationLabel(bool linear, bool stokes,
620                                            bool linpol, int polidx) const
[401]621{
[805]622  uInt idx = 0;
623  if (polidx >=0) idx = polidx;
624  return "";
625  //return SDPolUtil::polarizationLabel(idx, linear, stokes, linpol);
[401]626}
[386]627
[805]628void Scantable::unsetSelection()
[380]629{
[805]630  table_ = originalTable_;
[847]631  attach();
[805]632  selector_.reset();
[380]633}
[386]634
[805]635void Scantable::setSelection( const STSelector& selection )
[430]636{
[805]637  Table tab = const_cast<STSelector&>(selection).apply(originalTable_);
638  if ( tab.nrow() == 0 ) {
639    throw(AipsError("Selection contains no data. Not applying it."));
640  }
641  table_ = tab;
[847]642  attach();
[805]643  selector_ = selection;
[430]644}
645
[837]646std::string Scantable::summary( bool verbose )
[447]647{
[805]648  // Format header info
649  ostringstream oss;
650  oss << endl;
651  oss << asap::SEPERATOR << endl;
652  oss << " Scan Table Summary" << endl;
653  oss << asap::SEPERATOR << endl;
654  oss.flags(std::ios_base::left);
655  oss << setw(15) << "Beams:" << setw(4) << nbeam() << endl
656      << setw(15) << "IFs:" << setw(4) << nif() << endl
657      << setw(15) << "Polarisations:" << setw(4) << npol() << endl
658      << setw(15) << "Channels:"  << setw(4) << nchan() << endl;
659  oss << endl;
660  String tmp;
[860]661  oss << setw(15) << "Observer:"
662      << table_.keywordSet().asString("Observer") << endl;
[805]663  oss << setw(15) << "Obs Date:" << getTime(-1,true) << endl;
664  table_.keywordSet().get("Project", tmp);
665  oss << setw(15) << "Project:" << tmp << endl;
666  table_.keywordSet().get("Obstype", tmp);
667  oss << setw(15) << "Obs. Type:" << tmp << endl;
668  table_.keywordSet().get("AntennaName", tmp);
669  oss << setw(15) << "Antenna Name:" << tmp << endl;
670  table_.keywordSet().get("FluxUnit", tmp);
671  oss << setw(15) << "Flux Unit:" << tmp << endl;
672  Vector<Float> vec;
673  oss << setw(15) << "Rest Freqs:";
674  if (vec.nelements() > 0) {
675      oss << setprecision(10) << vec << " [Hz]" << endl;
676  } else {
677      oss << "none" << endl;
678  }
679  oss << setw(15) << "Abcissa:" << "channel" << endl;
680  oss << selector_.print() << endl;
681  oss << endl;
682  // main table
683  String dirtype = "Position ("
684                  + MDirection::showType(dirCol_.getMeasRef().getType())
685                  + ")";
686  oss << setw(5) << "Scan"
687      << setw(15) << "Source"
688//      << setw(24) << dirtype
689      << setw(10) << "Time"
690      << setw(18) << "Integration" << endl
691      << setw(5) << "" << setw(10) << "Beam" << dirtype << endl
692      << setw(15) << "" << setw(5) << "IF"
693      << setw(8) << "Frame" << setw(16)
694      << "RefVal" << setw(10) << "RefPix" << setw(12) << "Increment" <<endl;
695  oss << asap::SEPERATOR << endl;
696  TableIterator iter(table_, "SCANNO");
697  while (!iter.pastEnd()) {
698    Table subt = iter.table();
699    ROTableRow row(subt);
700    MEpoch::ROScalarColumn timeCol(subt,"TIME");
701    const TableRecord& rec = row.get(0);
702    oss << setw(4) << std::right << rec.asuInt("SCANNO")
703        << std::left << setw(1) << ""
704        << setw(15) << rec.asString("SRCNAME")
705        << setw(10) << formatTime(timeCol(0), false);
706    // count the cycles in the scan
707    TableIterator cyciter(subt, "CYCLENO");
708    int nint = 0;
709    while (!cyciter.pastEnd()) {
710      ++nint;
711      ++cyciter;
712    }
713    oss << setw(3) << std::right << nint  << setw(3) << " x " << std::left
714        << setw(6) <<  formatSec(rec.asFloat("INTERVAL")) << endl;
[447]715
[805]716    TableIterator biter(subt, "BEAMNO");
717    while (!biter.pastEnd()) {
718      Table bsubt = biter.table();
719      ROTableRow brow(bsubt);
720      MDirection::ROScalarColumn bdirCol(bsubt,"DIRECTION");
721      const TableRecord& brec = brow.get(0);
722      oss << setw(6) << "" <<  setw(10) << brec.asuInt("BEAMNO");
723      oss  << setw(24) << formatDirection(bdirCol(0)) << endl;
724      TableIterator iiter(bsubt, "IFNO");
725      while (!iiter.pastEnd()) {
726        Table isubt = iiter.table();
727        ROTableRow irow(isubt);
728        const TableRecord& irec = irow.get(0);
729        oss << std::right <<setw(8) << "" << std::left << irec.asuInt("IFNO");
730        oss << frequencies().print(irec.asuInt("FREQ_ID"));
[447]731
[805]732        ++iiter;
733      }
734      ++biter;
735    }
736    ++iter;
[447]737  }
[805]738  /// @todo implement verbose mode
739  return String(oss);
[447]740}
741
[805]742std::string Scantable::getTime(int whichrow, bool showdate) const
[777]743{
[805]744  MEpoch::ROScalarColumn timeCol(table_, "TIME");
745  MEpoch me;
746  if (whichrow > -1) {
747    me = timeCol(uInt(whichrow));
748  } else {
749    Double tm;
750    table_.keywordSet().get("UTC",tm);
751    me = MEpoch(MVEpoch(tm));
[777]752  }
[805]753  return formatTime(me, showdate);
[777]754}
[805]755
[865]756std::vector< double > asap::Scantable::getAbcissa( int whichrow ) const
757{
758  if ( whichrow > table_.nrow() ) throw(AipsError("Illegal ro number"));
759  std::vector<double> stlout;
760  int nchan = specCol_(whichrow).nelements();
761  String us = freqTable_.getUnitString();
762  if ( us == "" || us == "pixel" || us == "channel" ) {
763    for (int i=0; i<nchan; ++i) {
764      stlout.push_back(double(i));
765    }
766    return stlout;
767  }
768
769  const MPosition& mp = getAntennaPosition();
770  const MDirection& md = dirCol_(whichrow);
771  const MEpoch& me = timeCol_(whichrow);
772  Double rf = moleculeTable_.getRestFrequency(mmolidCol_(whichrow));
773  SpectralCoordinate spc =
774    freqTable_.getSpectralCoordinate(md, mp, me, rf, mfreqidCol_(whichrow));
775  Vector<Double> pixel(nchan);
776  Vector<Double> world;
777  indgen(pixel);
778  if ( Unit(us) == Unit("Hz") ) {
779    for ( int i=0; i < nchan; ++i) {
780      Double world;
781      spc.toWorld(world, pixel[i]);
782      stlout.push_back(double(world));
783    }
784  } else if ( Unit(us) == Unit("km/s") ) {
785    Vector<Double> world;
786    spc.pixelToVelocity(world, pixel);
787    world.tovector(stlout);
788  }
789  return stlout;
790}
791
[847]792std::string Scantable::getAbcissaLabel( int whichrow ) const
793{
794  if ( whichrow > table_.nrow() ) throw(AipsError("Illegal ro number"));
795  const MPosition& mp = getAntennaPosition();
796  const MDirection& md = dirCol_(whichrow);
797  const MEpoch& me = timeCol_(whichrow);
798  const Double& rf = mmolidCol_(whichrow);
799  SpectralCoordinate spc =
800    freqTable_.getSpectralCoordinate(md, mp, me, rf, mfreqidCol_(whichrow));
801
802  String s = "Channel";
803  Unit u = Unit(freqTable_.getUnitString());
804  if (u == Unit("km/s")) {
805    s = CoordinateUtil::axisLabel(spc,0,True,True,True);
806  } else if (u == Unit("Hz")) {
807    Vector<String> wau(1);wau = u.getName();
808    spc.setWorldAxisUnits(wau);
809
810    s = CoordinateUtil::axisLabel(spc,0,True,True,False);
811  }
812  return s;
813
814}
815
816void asap::Scantable::setRestFrequencies( double rf, const std::string& unit )
817{
818  ///@todo lookup in line table
819  Unit u(unit);
820  Quantum<Double> urf(rf, u);
821  uInt id = moleculeTable_.addEntry(urf.getValue("Hz"), "", "");
822  TableVector<uInt> tabvec(table_, "MOLECULE_ID");
823  tabvec = id;
824}
825
826void asap::Scantable::setRestFrequencies( const std::string& name )
827{
828  throw(AipsError("setRestFrequencies( const std::string& name ) NYI"));
829  ///@todo implement
830}
831
[852]832std::vector< unsigned int > asap::Scantable::rownumbers( ) const
833{
834  std::vector<unsigned int> stlout;
835  Vector<uInt> vec = table_.rowNumbers();
836  vec.tovector(stlout);
837  return stlout;
838}
839
[865]840
841
[805]842}//namespace asap
Note: See TracBrowser for help on using the repository browser.