source: trunk/external/atnf/PKSIO/PKSMS2writer.cc @ 1325

Last change on this file since 1325 was 1325, checked in by mar637, 17 years ago

Changes to use casacore instead of casa_asap/aips++\nAdded atnf PKSIO library snapshot to external and linking against this local copy
File size: 33.3 KB
Line 
1//#---------------------------------------------------------------------------
2//# PKSMS2writer.cc: Class to write Parkes multibeam data to a measurementset.
3//#---------------------------------------------------------------------------
4//# Copyright (C) 2000-2006
5//# Associated Universities, Inc. Washington DC, USA.
6//#
7//# This library is free software; you can redistribute it and/or modify it
8//# under the terms of the GNU Library General Public License as published by
9//# the Free Software Foundation; either version 2 of the License, or (at your
10//# option) any later version.
11//#
12//# This library is distributed in the hope that it will be useful, but
13//# WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
14//# or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Library General Public
15//# License for more details.
16//#
17//# You should have received a copy of the GNU Library General Public License
18//# along with this library; if not, write to the Free Software Foundation,
19//# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
20//#
21//# Correspondence concerning AIPS++ should be addressed as follows:
22//#        Internet email: aips2-request@nrao.edu.
23//#        Postal address: AIPS++ Project Office
24//#                        National Radio Astronomy Observatory
25//#                        520 Edgemont Road
26//#                        Charlottesville, VA 22903-2475 USA
27//#
28//# $Id: PKSMS2writer.cc,v 19.11 2006/11/06 22:25:22 mmarquar Exp $
29//#---------------------------------------------------------------------------
30
31#include <atnf/PKSIO/PKSMS2writer.h>
32
33#include <casa/Arrays/ArrayUtil.h>
34#include <casa/Arrays/ArrayMath.h>
35#include <casa/Arrays/ArrayLogical.h>
36#include <casa/Arrays/Cube.h>
37#include <casa/BasicSL/Complex.h>
38#include <casa/BasicSL/Constants.h>
39#include <casa/Quanta/QC.h>
40#include <measures/Measures/Stokes.h>
41#include <tables/Tables/ArrColDesc.h>
42#include <tables/Tables/IncrementalStMan.h>
43#include <tables/Tables/ScaColDesc.h>
44#include <tables/Tables/SetupNewTab.h>
45#include <tables/Tables/StandardStMan.h>
46#include <tables/Tables/Table.h>
47#include <tables/Tables/TableDesc.h>
48#include <tables/Tables/TiledShapeStMan.h>
49
50//------------------------------------------------- PKSMS2writer::PKSMS2writer
51
52// Default constructor.
53
54PKSMS2writer::PKSMS2writer()
55{
56}
57
58//------------------------------------------------ PKSMS2writer::~PKSMS2writer
59
60// Destructor.
61
62PKSMS2writer::~PKSMS2writer()
63{
64  close();
65}
66
67//------------------------------------------------------- PKSMS2writer::create
68
69// Create the output MS and and write static data.
70
71Int PKSMS2writer::create(
72        const String msName,
73        const String observer,
74        const String project,
75        const String antName,
76        const Vector<Double> antPosition,
77        const String obsMode,
78        const Float  equinox,
79        const String dopplerFrame,
80        const Vector<uInt> nChan,
81        const Vector<uInt> nPol,
82        const Vector<Bool> haveXPol,
83        const Bool   haveBase)
84{
85  // Open a MS table.
86  TableDesc pksDesc = MS::requiredTableDesc();
87
88  cNChan.assign(nChan);
89  cNPol.assign(nPol);
90  cHaveXPol.assign(haveXPol);
91
92  Int maxNPol = max(cNPol);
93
94
95  // Add the non-standard CALFCTR column.
96  pksDesc.addColumn(ArrayColumnDesc<Float>("CALFCTR", "Calibration factors",
97              IPosition(1,maxNPol), ColumnDesc::Direct));
98
99  // Add the optional FLOAT_DATA column.
100  MS::addColumnToDesc(pksDesc, MS::FLOAT_DATA, 2);
101  pksDesc.rwColumnDesc(MS::columnName(MS::FLOAT_DATA)).rwKeywordSet().
102                define("UNIT", String("Jy"));
103  pksDesc.rwColumnDesc(MS::columnName(MS::FLOAT_DATA)).rwKeywordSet().
104                define("MEASURE_TYPE", "");
105
106  if ((cHaveBase = haveBase)) {
107    // Add the non-standard BASELIN and BASESUB columns.
108    pksDesc.addColumn(ArrayColumnDesc<Float>("BASELIN", "Linear baseline fit",
109                IPosition(2,2,maxNPol), ColumnDesc::Direct));
110    pksDesc.addColumn(ArrayColumnDesc<Float>("BASESUB", "Baseline subtracted",
111                IPosition(2,9,maxNPol), ColumnDesc::Direct));
112  }
113
114  // Add the optional DATA column if cross-polarizations are to be recorded.
115  if (ntrue(cHaveXPol)) {
116    // Add the non-standard XCALFCTR column.
117    pksDesc.addColumn(ScalarColumnDesc<Complex>("XCALFCTR",
118                "Cross-polarization calibration factor"));
119
120    MS::addColumnToDesc(pksDesc, MS::DATA, 2);
121    pksDesc.rwColumnDesc(MS::columnName(MS::DATA)).rwKeywordSet().
122                define("UNIT", "Jy");
123    pksDesc.rwColumnDesc(MS::columnName(MS::DATA)).rwKeywordSet().
124                define("MEASURE_TYPE", "");
125  }
126
127  // Define hypercube for the float data (without coordinates).
128  pksDesc.defineHypercolumn("TiledData", 3,
129                stringToVector(MS::columnName(MS::FLOAT_DATA)));
130
131  SetupNewTable newtab(msName, pksDesc, Table::New);
132
133  // Set Incremental Storage Manager as the default.
134  IncrementalStMan incrStMan("ISMData");
135  newtab.bindAll(incrStMan, True);
136
137  // Use TiledShapeStMan for the FLOAT_DATA hypercube with tile size 16 kiB.
138  TiledShapeStMan tiledStMan("TiledData", IPosition(3,1,128,32));
139  newtab.bindColumn(MS::columnName(MS::FLOAT_DATA), tiledStMan);
140
141  // Use Standard Storage Manager to handle columns that change for each row.
142  StandardStMan stdStMan;
143  newtab.bindColumn(MS::columnName(MS::SCAN_NUMBER), stdStMan);
144  newtab.bindColumn(MS::columnName(MS::TIME), stdStMan);
145  newtab.bindColumn(MS::columnName(MS::SIGMA), stdStMan);
146  if (maxNPol > 2) {
147    newtab.bindColumn(MS::columnName(MS::DATA), stdStMan);
148  }
149
150  // Create the measurementset.
151  cPKSMS = new MeasurementSet(newtab, 0);
152
153  // Create subtables.
154  TableDesc antennaDesc = MSAntenna::requiredTableDesc();
155  SetupNewTable antennaSetup(cPKSMS->antennaTableName(), antennaDesc,
156                Table::New);
157  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::ANTENNA),
158                Table(antennaSetup));
159
160  TableDesc dataDescDesc = MSDataDescription::requiredTableDesc();
161  SetupNewTable dataDescSetup(cPKSMS->dataDescriptionTableName(), dataDescDesc,
162                Table::New);
163  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::DATA_DESCRIPTION),
164                Table(dataDescSetup));
165
166  TableDesc dopplerDesc = MSDoppler::requiredTableDesc();
167  SetupNewTable dopplerSetup(cPKSMS->dopplerTableName(), dopplerDesc,
168                Table::New);
169  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::DOPPLER),
170                Table(dopplerSetup));
171
172  TableDesc feedDesc = MSFeed::requiredTableDesc();
173  MSFeed::addColumnToDesc(feedDesc, MSFeedEnums::FOCUS_LENGTH);
174  SetupNewTable feedSetup(cPKSMS->feedTableName(), feedDesc, Table::New);
175  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::FEED),
176                Table(feedSetup));
177
178  TableDesc fieldDesc = MSField::requiredTableDesc();
179  SetupNewTable fieldSetup(cPKSMS->fieldTableName(), fieldDesc, Table::New);
180  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::FIELD),
181                Table(fieldSetup));
182
183  TableDesc flagCmdDesc = MSFlagCmd::requiredTableDesc();
184  SetupNewTable flagCmdSetup(cPKSMS->flagCmdTableName(), flagCmdDesc,
185                Table::New);
186  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::FLAG_CMD),
187                Table(flagCmdSetup));
188
189  TableDesc historyDesc = MSHistory::requiredTableDesc();
190  SetupNewTable historySetup(cPKSMS->historyTableName(), historyDesc,
191                Table::New);
192  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::HISTORY),
193                Table(historySetup));
194
195  TableDesc observationDesc = MSObservation::requiredTableDesc();
196  SetupNewTable observationSetup(cPKSMS->observationTableName(),
197                observationDesc, Table::New);
198  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::OBSERVATION),
199                Table(observationSetup));
200
201  TableDesc pointingDesc = MSPointing::requiredTableDesc();
202  SetupNewTable pointingSetup(cPKSMS->pointingTableName(), pointingDesc,
203                Table::New);
204  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::POINTING),
205                Table(pointingSetup));
206
207  TableDesc polarizationDesc = MSPolarization::requiredTableDesc();
208  SetupNewTable polarizationSetup(cPKSMS->polarizationTableName(),
209                polarizationDesc, Table::New);
210  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::POLARIZATION),
211                Table(polarizationSetup));
212
213  TableDesc processorDesc = MSProcessor::requiredTableDesc();
214  SetupNewTable processorSetup(cPKSMS->processorTableName(), processorDesc,
215                Table::New);
216  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::PROCESSOR),
217                Table(processorSetup));
218
219  TableDesc sourceDesc = MSSource::requiredTableDesc();
220  MSSource::addColumnToDesc(sourceDesc, MSSourceEnums::TRANSITION, 1);
221  MSSource::addColumnToDesc(sourceDesc, MSSourceEnums::REST_FREQUENCY,
222                1);
223  MSSource::addColumnToDesc(sourceDesc, MSSourceEnums::SYSVEL, 1);
224  SetupNewTable sourceSetup(cPKSMS->sourceTableName(), sourceDesc,
225                Table::New);
226  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::SOURCE),
227                Table(sourceSetup));
228
229  TableDesc spectralWindowDesc = MSSpectralWindow::requiredTableDesc();
230  MSSpectralWindow::addColumnToDesc(spectralWindowDesc,
231                MSSpectralWindowEnums::DOPPLER_ID);
232  SetupNewTable spectralWindowSetup(cPKSMS->spectralWindowTableName(),
233                spectralWindowDesc, Table::New);
234  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::SPECTRAL_WINDOW),
235                Table(spectralWindowSetup));
236
237  TableDesc stateDesc = MSState::requiredTableDesc();
238  SetupNewTable stateSetup(cPKSMS->stateTableName(), stateDesc, Table::New);
239  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::STATE),
240                Table(stateSetup));
241
242  TableDesc sysCalDesc = MSSysCal::requiredTableDesc();
243  MSSysCal::addColumnToDesc(sysCalDesc, MSSysCalEnums::TCAL, 1);
244  MSSysCal::addColumnToDesc(sysCalDesc, MSSysCalEnums::TSYS, 1);
245  SetupNewTable sysCalSetup(cPKSMS->sysCalTableName(), sysCalDesc,
246                Table::New);
247  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::SYSCAL),
248                Table(sysCalSetup));
249
250  TableDesc weatherDesc = MSWeather::requiredTableDesc();
251  MSWeather::addColumnToDesc(weatherDesc, MSWeatherEnums::PRESSURE);
252  MSWeather::addColumnToDesc(weatherDesc, MSWeatherEnums::REL_HUMIDITY);
253  MSWeather::addColumnToDesc(weatherDesc, MSWeatherEnums::TEMPERATURE);
254  SetupNewTable weatherSetup(cPKSMS->weatherTableName(), weatherDesc,
255                Table::New);
256  cPKSMS->rwKeywordSet().defineTable(MS::keywordName(MS::WEATHER),
257                Table(weatherSetup));
258
259  cPKSMS->initRefs();
260
261  // Measurementset subtables.
262  cAntenna         = cPKSMS->antenna();
263  cDataDescription = cPKSMS->dataDescription();
264  cDoppler         = cPKSMS->doppler();
265  cFeed            = cPKSMS->feed();
266  cField           = cPKSMS->field();
267  cFlagCmd         = cPKSMS->flagCmd();
268  cHistory         = cPKSMS->history();
269  cObservation     = cPKSMS->observation();
270  cPointing        = cPKSMS->pointing();
271  cPolarization    = cPKSMS->polarization();
272  cProcessor       = cPKSMS->processor();
273  cSource          = cPKSMS->source();
274  cSpectralWindow  = cPKSMS->spectralWindow();
275  cState           = cPKSMS->state();
276  cSysCal          = cPKSMS->sysCal();
277  cWeather         = cPKSMS->weather();
278
279  // Measurementset table columns;
280  cMSCols           = new MSColumns(*cPKSMS);
281  cAntennaCols      = new MSAntennaColumns(cAntenna);
282  cDataDescCols     = new MSDataDescColumns(cDataDescription);
283  cDopplerCols      = new MSDopplerColumns(cDoppler);
284  cFeedCols         = new MSFeedColumns(cFeed);
285  cFieldCols        = new MSFieldColumns(cField);
286  cFlagCmdCols      = new MSFlagCmdColumns(cFlagCmd);
287  cHistoryCols      = new MSHistoryColumns(cHistory);
288  cObservationCols  = new MSObservationColumns(cObservation);
289  cPointingCols     = new MSPointingColumns(cPointing);
290  cPolarizationCols = new MSPolarizationColumns(cPolarization);
291  cProcessorCols    = new MSProcessorColumns(cProcessor);
292  cSourceCols       = new MSSourceColumns(cSource);
293  cSpWindowCols     = new MSSpWindowColumns(cSpectralWindow);
294  cStateCols        = new MSStateColumns(cState);
295  cSysCalCols       = new MSSysCalColumns(cSysCal);
296  cWeatherCols      = new MSWeatherColumns(cWeather);
297
298  cCalFctrCol  = new ArrayColumn<Float>(*cPKSMS, "CALFCTR");
299  if (cHaveBase) {
300    cBaseLinCol = new ArrayColumn<Float>(*cPKSMS, "BASELIN");
301    cBaseSubCol = new ArrayColumn<Float>(*cPKSMS, "BASESUB");
302  }
303  if (ntrue(cHaveXPol)) {
304    cXCalFctrCol = new ScalarColumn<Complex>(*cPKSMS, "XCALFCTR");
305  }
306
307
308  // Define Measure references.
309  Vector<String> flagCat(1, "BAD");
310  cMSCols->sigma().rwKeywordSet().define("UNIT", "K");
311  cMSCols->flagCategory().rwKeywordSet().define("CATEGORY", flagCat);
312
313  String dirref;
314  if (equinox == 1950.0f) {
315    dirref = "B1950";
316  } else {
317    dirref = "J2000";
318  }
319
320  cFieldCols->delayDir().rwKeywordSet().asrwRecord("MEASINFO").
321                 define("Ref", dirref);
322  cFieldCols->phaseDir().rwKeywordSet().asrwRecord("MEASINFO").
323                 define("Ref", dirref);
324  cFieldCols->referenceDir().rwKeywordSet().asrwRecord("MEASINFO").
325                 define("Ref", dirref);
326
327  cPointingCols->direction().rwKeywordSet().asrwRecord("MEASINFO").
328                 define("Ref", dirref);
329  cPointingCols->target().rwKeywordSet().asrwRecord("MEASINFO").
330                 define("Ref", dirref);
331
332  cSourceCols->direction().rwKeywordSet().asrwRecord("MEASINFO").
333                 define("Ref", dirref);
334  cSourceCols->restFrequency().rwKeywordSet().asrwRecord("MEASINFO").
335                 define("Ref", "REST");
336
337  // Translate Doppler frame name.
338  if (dopplerFrame == "TOPOCENT") {
339    MFrequency::getType(cDopplerFrame, "TOPO");
340  } else if (dopplerFrame == "GEOCENTR") {
341    MFrequency::getType(cDopplerFrame, "GEO");
342  } else if (dopplerFrame == "BARYCENT") {
343    MFrequency::getType(cDopplerFrame, "BARY");
344  } else if (dopplerFrame == "GALACTOC") {
345    MFrequency::getType(cDopplerFrame, "GALACTO");
346  } else if (dopplerFrame == "LOCALGRP") {
347    MFrequency::getType(cDopplerFrame, "LGROUP");
348  } else if (dopplerFrame == "CMBDIPOL") {
349    MFrequency::getType(cDopplerFrame, "CMB");
350  } else if (dopplerFrame == "SOURCE") {
351    MFrequency::getType(cDopplerFrame, "REST");
352  }
353
354
355  // Store static data.
356  addAntennaEntry(antName, antPosition);
357  addDopplerEntry();
358  addFeedEntry();
359  addObservationEntry(observer, project);
360  addProcessorEntry();
361
362  return 0;
363}
364
365//-------------------------------------------------------- PKSMS2writer::write
366
367// Write the next data record.
368
369Int PKSMS2writer::write(
370        const Int             scanNo,
371        const Int             cycleNo,
372        const Double          mjd,
373        const Double          interval,
374        const String          fieldName,
375        const String          srcName,
376        const Vector<Double>  srcDir,
377        const Vector<Double>  srcPM,
378        const Double          srcVel,
379        const String          obsMode,
380        const Int             IFno,
381        const Double          refFreq,
382        const Double          bandwidth,
383        const Double          freqInc,
384        const Double          restFreq,
385        const Vector<Float>   tcal,
386        const String          tcalTime,
387        const Float           azimuth,
388        const Float           elevation,
389        const Float           parAngle,
390        const Float           focusAxi,
391        const Float           focusTan,
392        const Float           focusRot,
393        const Float           temperature,
394        const Float           pressure,
395        const Float           humidity,
396        const Float           windSpeed,
397        const Float           windAz,
398        const Int             refBeam,
399        const Int             beamNo,
400        const Vector<Double>  direction,
401        const Vector<Double>  scanRate,
402        const Vector<Float>   tsys,
403        const Vector<Float>   sigma,
404        const Vector<Float>   calFctr,
405        const Matrix<Float>   baseLin,
406        const Matrix<Float>   baseSub,
407        const Matrix<Float>   &spectra,
408        const Matrix<uChar>   &flagged,
409        const Complex         xCalFctr,
410        const Vector<Complex> &xPol)
411{
412  // Extend the time range in the OBSERVATION subtable.
413  Vector<Double> timerange(2);
414  cObservationCols->timeRange().get(0, timerange);
415  Double time = mjd*86400.0;
416  if (timerange(0) == 0.0) {
417    timerange(0) = time;
418  }
419  timerange(1) = time;
420  cObservationCols->timeRange().put(0, timerange);
421
422  Int iIF = IFno - 1;
423  Int nChan = cNChan(iIF);
424  Int nPol  = cNPol(iIF);
425
426  // IFno is the 1-relative row number in the DATA_DESCRIPTION,
427  // SPECTRAL_WINDOW, and POLARIZATION subtables.
428  if (Int(cDataDescription.nrow()) < IFno) {
429    // Add a new entry to each subtable.
430    addDataDescriptionEntry(IFno);
431    addSpectralWindowEntry(IFno, nChan, refFreq, bandwidth, freqInc);
432    addPolarizationEntry(IFno, nPol);
433  }
434
435  // Find or add the source to the SOURCE subtable.
436  Int srcId = addSourceEntry(srcName, srcDir, srcPM, restFreq, srcVel);
437
438  // Find or add the obsMode to the STATE subtable.
439  Int stateId = addStateEntry(obsMode);
440
441  // FIELD subtable.
442  Int fieldId = addFieldEntry(fieldName, time, direction, scanRate, srcId);
443
444  // POINTING subtable.
445  addPointingEntry(time, interval, fieldName, direction, scanRate);
446
447  // SYSCAL subtable.
448  addSysCalEntry(beamNo, iIF, time, interval, tcal, tsys);
449
450  // Handle weather information.
451  ROScalarColumn<Double> wTime(cWeather, "TIME");
452  Int nWeather = wTime.nrow();
453  if (nWeather == 0 || time > wTime(nWeather-1)) {
454    addWeatherEntry(time, interval, pressure, humidity, temperature);
455  }
456
457
458  // Extend the main table.
459  cPKSMS->addRow();
460  Int irow = cPKSMS->nrow() - 1;
461
462  // Keys.
463  cMSCols->time().put(irow, time);
464  cMSCols->antenna1().put(irow, 0);
465  cMSCols->antenna2().put(irow, 0);
466  cMSCols->feed1().put(irow, beamNo-1);
467  cMSCols->feed2().put(irow, beamNo-1);
468  cMSCols->dataDescId().put(irow, iIF);
469  cMSCols->processorId().put(irow, 0);
470  cMSCols->fieldId().put(irow, fieldId);
471
472  // Non-key attributes.
473  cMSCols->interval().put(irow, interval);
474  cMSCols->exposure().put(irow, interval);
475  cMSCols->timeCentroid().put(irow, time);
476  cMSCols->scanNumber().put(irow, scanNo);
477  cMSCols->arrayId().put(irow, 0);
478  cMSCols->observationId().put(irow, 0);
479  cMSCols->stateId().put(irow, stateId);
480
481  Vector<Double> uvw(3, 0.0);
482  cMSCols->uvw().put(irow, uvw);
483
484  // Baseline fit parameters.
485  if (cHaveBase) {
486    cBaseLinCol->put(irow, baseLin);
487
488    if (baseSub.nrow() == 9) {
489      cBaseSubCol->put(irow, baseSub);
490
491    } else {
492      Matrix<Float> tmp(9, 2, 0.0f);
493      for (Int ipol = 0; ipol < nPol; ipol++) {
494        for (uInt j = 0; j < baseSub.nrow(); j++) {
495          tmp(j,ipol) = baseSub(j,ipol);
496        }
497      }
498      cBaseSubCol->put(irow, tmp);
499    }
500  }
501
502  // Transpose spectra.
503  Matrix<Float> tmpData(nPol, nChan);
504  Matrix<Bool>  tmpFlag(nPol, nChan);
505  for (Int ipol = 0; ipol < nPol; ipol++) {
506    for (Int ichan = 0; ichan < nChan; ichan++) {
507      tmpData(ipol,ichan) = spectra(ichan,ipol);
508      tmpFlag(ipol,ichan) = flagged(ichan,ipol);
509    }
510  }
511
512  cCalFctrCol->put(irow, calFctr);
513  cMSCols->floatData().put(irow, tmpData);
514  cMSCols->flag().put(irow, tmpFlag);
515
516  // Cross-polarization spectra.
517  if (cHaveXPol(iIF)) {
518    cXCalFctrCol->put(irow, xCalFctr);
519    cMSCols->data().put(irow, xPol);
520  }
521
522  cMSCols->sigma().put(irow, sigma);
523
524  Vector<Float> weight(1, 1.0f);
525  cMSCols->weight().put(irow, weight);
526
527  // Flag information.
528  Cube<Bool> flags(nPol, nChan, 1, False);
529  cMSCols->flag().put(irow, flags.xyPlane(0));
530  cMSCols->flagCategory().put(irow, flags);
531  cMSCols->flagRow().put(irow, False);
532
533  return 0;
534}
535
536//-------------------------------------------------------- PKSMS2writer::close
537
538// Close the measurementset, flushing all associated tables.
539
540void PKSMS2writer::close()
541{
542  // Delete table column accessors.
543  delete cMSCols; cMSCols=0;
544  delete cAntennaCols; cAntennaCols=0;
545  delete cDataDescCols; cDataDescCols=0;
546  delete cDopplerCols; cDopplerCols=0;
547  delete cFeedCols; cFeedCols=0;
548  delete cFieldCols; cFieldCols=0;
549  delete cFlagCmdCols; cFlagCmdCols=0;
550  delete cHistoryCols; cHistoryCols=0;
551  delete cObservationCols; cObservationCols=0;
552  delete cPointingCols; cPointingCols=0;
553  delete cPolarizationCols; cPolarizationCols=0;
554  delete cProcessorCols; cProcessorCols=0;
555  delete cSourceCols; cSourceCols=0;
556  delete cSpWindowCols; cSpWindowCols=0;
557  delete cStateCols; cStateCols=0;
558  delete cSysCalCols; cSysCalCols=0;
559  delete cWeatherCols; cWeatherCols=0;
560
561  delete cCalFctrCol; cCalFctrCol=0;
562  if (cHaveBase) {
563    delete cBaseLinCol; cBaseLinCol=0;
564    delete cBaseSubCol; cBaseSubCol=0;
565  }
566  if (ntrue(cHaveXPol)) {
567    delete cXCalFctrCol; cXCalFctrCol=0;
568  }
569 
570  // Release all subtables.
571  cAntenna         = MSAntenna();
572  cDataDescription = MSDataDescription();
573  cDoppler         = MSDoppler();
574  cFeed            = MSFeed();
575  cField           = MSField();
576  cFlagCmd         = MSFlagCmd();
577  cHistory         = MSHistory();
578  cObservation     = MSObservation();
579  cPointing        = MSPointing();
580  cPolarization    = MSPolarization();
581  cProcessor       = MSProcessor();
582  cSource          = MSSource();
583  cSpectralWindow  = MSSpectralWindow();
584  cState           = MSState();
585  cSysCal          = MSSysCal();
586  cWeather         = MSWeather();
587  // Release the main table.
588  delete cPKSMS; cPKSMS=0;
589}
590
591//---------------------------------------------- PKSMS2writer::addAntennaEntry
592
593// Add an entry to the ANTENNA subtable.
594
595Int PKSMS2writer::addAntennaEntry(
596        const String antName,
597        const Vector<Double> &antPosition)
598{
599  // Extend the ANTENNA subtable.
600  cAntenna.addRow();
601  Int n = cAntenna.nrow() - 1;
602
603  // Data.
604  cAntennaCols->name().put(n, antName);
605  cAntennaCols->station().put(n, "ATNF_PARKES");
606  cAntennaCols->type().put(n, "GROUND-BASED");
607  cAntennaCols->mount().put(n, "ALT-AZ");
608  cAntennaCols->position().put(n, antPosition);
609  Vector<Double> antOffset(3, 0.0);
610  cAntennaCols->offset().put(n, antOffset);
611  cAntennaCols->dishDiameter().put(n, 64.0);
612
613  // Flags.
614  cAntennaCols->flagRow().put(n, False);
615
616  return n;
617}
618
619//-------------------------------------- PKSMS2writer::addDataDescriptionEntry
620
621// Add an entry to the DATA_DESCRIPTION subtable.
622
623Int PKSMS2writer::addDataDescriptionEntry(
624        const Int IFno)
625{
626  // Extend the DATA_DESCRIPTION subtable.
627  while (Int(cDataDescription.nrow()) < IFno) {
628    cDataDescription.addRow();
629  }
630  Int n = IFno - 1;
631
632  // Data.
633  cDataDescCols->spectralWindowId().put(n, n);
634  cDataDescCols->polarizationId().put(n, n);
635
636  // Flags.
637  cDataDescCols->flagRow().put(n, False);
638
639  return n;
640}
641
642//---------------------------------------------- PKSMS2writer::addDopplerEntry
643
644// Add an entry to the DOPPLER subtable.
645
646Int PKSMS2writer::addDopplerEntry()
647{
648  // Extend the DOPPLER subtable.
649  cDoppler.addRow();
650  Int n = cDoppler.nrow() - 1;
651
652  // Keys.
653  cDopplerCols->dopplerId().put(n, n);
654  cDopplerCols->sourceId().put(n, 0);
655
656  // Data.
657  cDopplerCols->transitionId().put(n, 0);
658
659  return n;
660}
661
662//------------------------------------------------- PKSMS2writer::addFeedEntry
663
664// Add an entry to the FEED subtable.
665
666Int PKSMS2writer::addFeedEntry()
667{
668  Int n = cFeed.nrow() - 1;
669  for (Int iBeam = 0; iBeam < 13; iBeam++) {
670    // Extend the FEED subtable.
671    cFeed.addRow();
672    n++;
673
674    // Keys.
675    cFeedCols->antennaId().put(n, cAntenna.nrow()-1);
676    cFeedCols->feedId().put(n, iBeam);
677    cFeedCols->spectralWindowId().put(n, -1);
678    cFeedCols->time().put(n, 0.0);
679    cFeedCols->interval().put(n, -1.0);
680
681    // Data description.
682    cFeedCols->numReceptors().put(n, 2);
683
684    // Data.
685    cFeedCols->beamId().put(n, -1);
686
687    Matrix<Double> beamOffset(2, 2, 0.0);
688    cFeedCols->beamOffset().put(n, beamOffset);
689
690    cFeedCols->focusLength().put(n, 26.0);
691
692    Vector<String> polarizationType(2);
693    polarizationType(0) = "X";
694    polarizationType(1) = "Y";
695    cFeedCols->polarizationType().put(n, polarizationType);
696
697    Matrix<Complex> polResponse(2, 2, Complex(0.0));
698    for (Int i = 0; i < 2; i++) {
699      polResponse(i,i) = Complex(1.0, 0.0);
700    }
701    cFeedCols->polResponse().put(n, polResponse);
702
703    Vector<Double> position(3, 0.0);
704    cFeedCols->position().put(n, position);
705
706    Vector<Double> receptorAngle(2, C::pi_4);
707    receptorAngle(1) += C::pi_2;
708    cFeedCols->receptorAngle().put(n, receptorAngle);
709  }
710
711  return n;
712}
713
714//------------------------------------------------ PKSMS2writer::addFieldEntry
715
716// Add an entry to the FIELD subtable.
717
718Int PKSMS2writer::addFieldEntry(
719        const String fieldName,
720        const Double time,
721        const Vector<Double> direction,
722        const Vector<Double> scanRate,
723        const Int srcId)
724{
725  // Extend the FIELD subtable.
726  cField.addRow();
727  Int n = cField.nrow() - 1;
728
729  // Data.
730  cFieldCols->name().put(n, fieldName);
731  cFieldCols->code().put(n, "DRIFT");
732  cFieldCols->time().put(n, time);
733
734  Matrix<Double> track(2, 2);
735  track.column(0) = direction;
736  track.column(1) = scanRate;
737  cFieldCols->numPoly().put(n, 1);
738  cFieldCols->delayDir().put(n, track);
739  cFieldCols->phaseDir().put(n, track);
740  cFieldCols->referenceDir().put(n, track);
741  cFieldCols->sourceId().put(n, srcId);
742
743  // Flags.
744  cFieldCols->flagRow().put(n, False);
745
746  return n;
747}
748
749//------------------------------------------ PKSMS2writer::addObservationEntry
750
751// Add an entry to the OBSERVATION subtable.
752
753Int PKSMS2writer::addObservationEntry(
754        const String observer,
755        const String project)
756{
757  // Extend the OBSERVATION subtable.
758  cObservation.addRow();
759  Int n = cObservation.nrow() - 1;
760
761  // Data.
762  cObservationCols->telescopeName().put(n, "Parkes");
763  Vector<Double> timerange(2, 0.0);
764  cObservationCols->timeRange().put(n, timerange);
765  cObservationCols->observer().put(n, observer);
766  Vector<String> log(1, "none");
767  cObservationCols->log().put(n, log);
768  cObservationCols->scheduleType().put(n, "ATNF");
769  Vector<String> schedule(1, "Not available");
770  cObservationCols->schedule().put(n, schedule);
771  cObservationCols->project().put(n, project);
772  cObservationCols->releaseDate().put(n, 0.0);
773
774  // Flags.
775  cObservationCols->flagRow().put(n, False);
776
777  return n;
778}
779
780//--------------------------------------------- PKSMS2writer::addPointingEntry
781
782// Add an entry to the POINTING subtable.  This compulsory subtable simply
783// duplicates information in the FIELD subtable.
784
785Int PKSMS2writer::addPointingEntry(
786        const Double time,
787        const Double interval,
788        const String fieldName,
789        const Vector<Double> direction,
790        const Vector<Double> scanRate)
791{
792  // Extend the POINTING subtable.
793  cPointing.addRow();
794  Int n = cPointing.nrow() - 1;
795
796  // Keys.
797  cPointingCols->antennaId().put(n, 0);
798  cPointingCols->time().put(n, time);
799  cPointingCols->interval().put(n, interval);
800
801  // Data.
802  cPointingCols->name().put(n, fieldName);
803  cPointingCols->numPoly().put(n, 1);
804  cPointingCols->timeOrigin().put(n, time);
805
806  Matrix<Double> track(2, 2);
807  track.column(0) = direction;
808  track.column(1) = scanRate;
809  cPointingCols->direction().put(n, track);
810  cPointingCols->target().put(n, track);
811  cPointingCols->tracking().put(n, True);
812
813  return n;
814}
815
816//----------------------------------------- PKSMS2writer::addPolarizationEntry
817
818// Add an entry to the POLARIZATION subtable.
819
820Int PKSMS2writer::addPolarizationEntry(
821        const Int IFno,
822        const Int nPol)
823{
824  // Extend the POLARIZATION subtable.
825  while (Int(cPolarization.nrow()) < IFno) {
826    cPolarization.addRow();
827  }
828  Int n = IFno - 1;
829
830  // Data description.
831  cPolarizationCols->numCorr().put(n, nPol);
832
833  // Data.
834  Vector<Int> corrType(2);
835  corrType(0) = Stokes::XX;
836  corrType(1) = Stokes::YY;
837  cPolarizationCols->corrType().put(n, corrType);
838
839  Matrix<Int> corrProduct(2,2,1);
840  if (nPol == 2) {
841    corrProduct(1,0) = 0;
842    corrProduct(0,1) = 0;
843  }
844  cPolarizationCols->corrProduct().put(n, corrProduct);
845
846  // Flags.
847  cPolarizationCols->flagRow().put(n, False);
848
849  return n;
850}
851
852
853//-------------------------------------------- PKSMS2writer::addProcessorEntry
854
855// Add an entry to the PROCESSOR subtable.
856
857Int PKSMS2writer::addProcessorEntry()
858{
859  // Extend the PROCESSOR subtable.
860  cProcessor.addRow();
861  Int n = cProcessor.nrow() - 1;
862
863  // Data.
864  cProcessorCols->type().put(n, "SPECTROMETER");
865  cProcessorCols->subType().put(n, "MULTIBEAM");
866  cProcessorCols->typeId().put(n, -1);
867  cProcessorCols->modeId().put(n, -1);
868
869  // Flags.
870  cProcessorCols->flagRow().put(n, False);
871
872  return n;
873}
874
875//----------------------------------------------- PKSMS2writer::addSourceEntry
876
877// Add an entry to the SOURCE subtable.
878
879Int PKSMS2writer::addSourceEntry(
880        const String name,
881        const Vector<Double> direction,
882        const Vector<Double> properMotion,
883        const Double restFreq,
884        const Double radialVelocity)
885{
886  // Look for an entry in the SOURCE subtable.
887  ROScalarColumn<String> sources(cSource, "NAME");
888  Int n;
889  Int nSrc = sources.nrow();
890  for (n = 0; n < nSrc; n++) {
891    if (sources(n) == name) {
892      break;
893    }
894  }
895
896  if (n == nSrc) {
897    // Not found, add a new entry to the SOURCE subtable.
898    cSource.addRow();
899
900    // Keys.
901    cSourceCols->sourceId().put(n, n);
902    cSourceCols->time().put(n, 0.0);
903    cSourceCols->interval().put(n, -1.0);
904    cSourceCols->spectralWindowId().put(n, -1);
905
906    // Data description.
907    cSourceCols->numLines().put(n, 1);
908
909    // Data.
910    cSourceCols->name().put(n, name);
911    cSourceCols->calibrationGroup().put(n, 0);
912    cSourceCols->code().put(n, "");
913    cSourceCols->direction().put(n, direction);
914//  Vector<Double> position(3, 0.0);
915//  cSourceCols->position().put(n, position);
916    cSourceCols->properMotion().put(n, properMotion);
917    Vector<Double> restFrequency(1, restFreq);
918    cSourceCols->restFrequency().put(n, restFrequency);
919    Vector<Double> sysvel(1, radialVelocity);
920    cSourceCols->sysvel().put(n, sysvel);
921  }
922
923  return n;
924}
925
926//--------------------------------------- PKSMS2writer::addSpectralWindowEntry
927
928// Add an entry to the SPECTRAL_WINDOW subtable.
929
930Int PKSMS2writer::addSpectralWindowEntry(
931        const Int IFno,
932        const Int nChan,
933        const Double refFreq,
934        const Double bandwidth,
935        const Double freqInc)
936{
937  // Extend the SPECTRAL_WINDOW subtable.
938  while (Int(cSpectralWindow.nrow()) < IFno) {
939    cSpectralWindow.addRow();
940  }
941  Int n = IFno - 1;
942
943  // Data description.
944  cSpWindowCols->numChan().put(n, nChan);
945
946  // Data.
947  cSpWindowCols->name().put(n, "L-band");
948  cSpWindowCols->refFrequency().put(n, refFreq);
949
950  // 0-relative reference channel number.
951  Double refChan = nChan / 2;
952  Vector<Double> freqs(nChan);
953  for (Int i = 0; i < nChan; i++) {
954    freqs(i) = refFreq + (i - refChan)*freqInc;
955  }
956  cSpWindowCols->chanFreq().put(n, freqs);
957
958  Vector<Double> chanWidths(nChan, freqInc);
959  cSpWindowCols->chanWidth().put(n, chanWidths);
960
961  cSpWindowCols->measFreqRef().put(n, cDopplerFrame);
962  cSpWindowCols->effectiveBW().put(n, chanWidths);
963
964  Vector<Double> resolution(nChan, fabs(freqInc));
965  cSpWindowCols->resolution().put(n, resolution);
966
967  cSpWindowCols->totalBandwidth().put(n, bandwidth);
968  cSpWindowCols->netSideband().put(n, 0);
969  cSpWindowCols->ifConvChain().put(n, -1);
970  cSpWindowCols->freqGroup().put(n, 0);
971  cSpWindowCols->freqGroupName().put(n, " ");
972  cSpWindowCols->dopplerId().put(n, 0);
973
974  // Flags.
975  cSpWindowCols->flagRow().put(n, False);
976
977  return n;
978}
979
980//------------------------------------------------ PKSMS2writer::addStateEntry
981
982// Add an entry to the STATE subtable.
983
984Int PKSMS2writer::addStateEntry(
985        const String obsMode)
986{
987  // Look for an entry in the STATE subtable.
988  for (uInt n = 0; n < cStateCols->nrow(); n++) {
989    if (cStateCols->obsMode()(n) == obsMode) {
990      return n;
991    }
992  }
993
994  // Not found, extend the STATE subtable.
995  cState.addRow();
996  uInt n = cStateCols->nrow() - 1;
997
998  // Data.
999  if (obsMode.contains("RF")) {
1000    cStateCols->sig().put(n, False);
1001    cStateCols->ref().put(n, True);
1002  } else if (!obsMode.contains("PA")) {
1003    // Signal and reference are both false for "paddle" data.
1004    cStateCols->sig().put(n, True);
1005    cStateCols->ref().put(n, False);
1006  }
1007
1008  cStateCols->load().put(n, 0.0);
1009  cStateCols->cal().put(n, 0.0);
1010  cStateCols->subScan().put(n, 0);
1011  cStateCols->obsMode().put(n, obsMode);
1012
1013  // Flags.
1014  cStateCols->flagRow().put(n, False);
1015
1016  return n;
1017}
1018
1019//----------------------------------------------- PKSMS2writer::addSysCalEntry
1020
1021// Add an entry to the SYSCAL subtable.
1022
1023Int PKSMS2writer::addSysCalEntry(
1024        const Int beamNo,
1025        const Int spWinId,
1026        const Double time,
1027        const Double interval,
1028        const Vector<Float> tcal,
1029        const Vector<Float> tsys)
1030{
1031  // Extend the SYSCAL subtable.
1032  cSysCal.addRow();
1033  Int n = cSysCal.nrow() - 1;
1034
1035  // Keys.
1036  cSysCalCols->antennaId().put(n, 0);
1037  cSysCalCols->feedId().put(n, beamNo-1);
1038  cSysCalCols->spectralWindowId().put(n, spWinId);
1039  cSysCalCols->time().put(n, time);
1040  cSysCalCols->interval().put(n, interval);
1041
1042  // Data.
1043  cSysCalCols->tcal().put(n, tcal);
1044  cSysCalCols->tsys().put(n, tsys);
1045
1046  return n;
1047}
1048
1049//---------------------------------------------- PKSMS2writer::addWeatherEntry
1050
1051// Add an entry to the WEATHER subtable.
1052
1053Int PKSMS2writer::addWeatherEntry(
1054        const Double time,
1055        const Double interval,
1056        const Double pressure,
1057        const Double relHumidity,
1058        const Double temperature)
1059{
1060  // Extend the WEATHER subtable.
1061  cWeather.addRow();
1062  Int n = cWeather.nrow() - 1;
1063
1064  // Keys.
1065  cWeatherCols->antennaId().put(n, 0);
1066  cWeatherCols->time().put(n, time);
1067  cWeatherCols->interval().put(n, interval);
1068
1069  // Data.
1070  cWeatherCols->pressure().put(n, pressure);
1071  cWeatherCols->relHumidity().put(n, relHumidity);
1072  cWeatherCols->temperature().put(n, temperature);
1073
1074  return n;
1075}
Note: See TracBrowser for help on using the repository browser.