source: trunk/src/MSWriter.cpp@ 2023

Last change on this file since 2023 was 2022, checked in by Takeshi Nakazato, 14 years ago

New Development: No

JIRA Issue: Yes CAS-2718

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...

Copy ASDM subtable that are optionally created by importasdm task.
For filler, their root names are stored as String in table keyword,
while for writer, existing tables are copied and are stored as Table
in table keyword.


File size: 61.9 KB
RevLine 
[1974]1//
2// C++ Interface: MSWriter
3//
4// Description:
5//
6// This class is specific writer for MS format
7//
8// Takeshi Nakazato <takeshi.nakazato@nao.ac.jp>, (C) 2010
9//
10// Copyright: See COPYING file that comes with this distribution
11//
12//
13
14#include <casa/OS/File.h>
15#include <casa/OS/RegularFile.h>
16#include <casa/OS/Directory.h>
17#include <casa/OS/SymLink.h>
18#include <casa/BasicSL/String.h>
[2016]19#include <casa/Containers/RecordField.h>
20#include <casa/Arrays/Cube.h>
[1974]21
[1975]22#include <tables/Tables/ExprNode.h>
[1974]23#include <tables/Tables/TableDesc.h>
24#include <tables/Tables/SetupNewTab.h>
25#include <tables/Tables/TableIter.h>
26#include <tables/Tables/RefRows.h>
[2016]27#include <tables/Tables/TableRow.h>
[1974]28
29#include <ms/MeasurementSets/MeasurementSet.h>
30#include <ms/MeasurementSets/MSColumns.h>
31#include <ms/MeasurementSets/MSPolIndex.h>
32#include <ms/MeasurementSets/MSDataDescIndex.h>
[1975]33#include <ms/MeasurementSets/MSSourceIndex.h>
[1974]34
35#include "MSWriter.h"
36#include "STHeader.h"
37#include "STFrequencies.h"
[1975]38#include "STMolecules.h"
[1977]39#include "STTcal.h"
[1974]40
[2016]41#include <ctime>
42#include <sys/time.h>
43
44// Boost
45#include <boost/pool/object_pool.hpp>
46
[1974]47using namespace casa ;
[2016]48using namespace std ;
[1974]49
[2016]50namespace asap {
51double MSWriter::gettimeofday_sec()
[1974]52{
[2016]53 struct timeval tv ;
54 gettimeofday( &tv, NULL ) ;
55 return tv.tv_sec + (double)tv.tv_usec*1.0e-6 ;
56}
[1974]57
58MSWriter::MSWriter(CountedPtr<Scantable> stable)
[2016]59 : table_(stable),
[2019]60 isTcal_(False),
61 isWeather_(False),
62 tcalSpec_(False),
63 tsysSpec_(False),
[2016]64 ptTabName_("")
[1974]65{
66 os_ = LogIO() ;
67 os_.origin( LogOrigin( "MSWriter", "MSWriter()", WHERE ) ) ;
68 os_ << "MSWriter::MSWriter()" << LogIO::POST ;
69
70 // initialize writer
71 init() ;
72}
73
74MSWriter::~MSWriter()
75{
76 os_.origin( LogOrigin( "MSWriter", "~MSWriter()", WHERE ) ) ;
77 os_ << "MSWriter::~MSWriter()" << LogIO::POST ;
[2016]78
79 if ( mstable_ != 0 )
80 delete mstable_ ;
[1974]81}
82
83bool MSWriter::write(const string& filename, const Record& rec)
84{
85 os_.origin( LogOrigin( "MSWriter", "write()", WHERE ) ) ;
[2016]86 double startSec = gettimeofday_sec() ;
87 os_ << "start MSWriter::write() startSec=" << startSec << LogIO::POST ;
[1974]88
89 filename_ = filename ;
90
91 // parsing MS options
92 Bool overwrite = False ;
93 if ( rec.isDefined( "ms" ) ) {
94 Record msrec = rec.asRecord( "ms" ) ;
95 if ( msrec.isDefined( "overwrite" ) ) {
96 overwrite = msrec.asBool( "overwrite" ) ;
97 }
98 }
99
100 os_ << "Parsing MS options" << endl ;
101 os_ << " overwrite = " << overwrite << LogIO::POST ;
102
103 File file( filename_ ) ;
104 if ( file.exists() ) {
105 if ( overwrite ) {
106 os_ << filename_ << " exists. Overwrite existing data... " << LogIO::POST ;
107 if ( file.isRegular() ) RegularFile(file).remove() ;
108 else if ( file.isDirectory() ) Directory(file).removeRecursive() ;
109 else SymLink(file).remove() ;
110 }
111 else {
112 os_ << LogIO::SEVERE << "ERROR: " << filename_ << " exists..." << LogIO::POST ;
113 return False ;
114 }
115 }
116
[2016]117 // TEST
118 // memory allocation by boost::object_pool
119// boost::object_pool<ROTableColumn> *tpoolr = new boost::object_pool<ROTableColumn> ;
120// ROTableColumn *tcol = 0 ;
121 //
122
[1974]123 // set up MS
124 setupMS() ;
125
126 // subtables
127 // OBSERVATION
128 fillObservation() ;
129
130 // ANTENNA
131 fillAntenna() ;
132
[1975]133 // PROCESSOR
134 fillProcessor() ;
135
136 // SOURCE
137 fillSource() ;
138
[1977]139 // WEATHER
[2019]140 if ( isWeather_ )
141 fillWeather() ;
[1977]142
[1974]143 // MAIN
144 // Iterate over several ids
145 Vector<uInt> processedFreqId( 0 ) ;
146 Int defaultFieldId = 0 ;
[2016]147
148 // row based
149 TableRow row( *mstable_ ) ;
150 TableRecord &trec = row.record() ;
151 NoticeTarget *dataRF = 0 ;
152 if ( useFloatData_ )
153 dataRF = new RecordFieldPtr< Array<Float> >( trec, "FLOAT_DATA" ) ;
154 else if ( useData_ )
155 dataRF = new RecordFieldPtr< Array<Complex> >( trec, "DATA" ) ;
156 RecordFieldPtr< Array<Bool> > flagRF( trec, "FLAG" ) ;
157 RecordFieldPtr<Bool> flagrowRF( trec, "FLAG_ROW" ) ;
158 RecordFieldPtr<Double> timeRF( trec, "TIME" ) ;
159 RecordFieldPtr<Double> timecRF( trec, "TIME_CENTROID" ) ;
160 RecordFieldPtr<Double> intervalRF( trec, "INTERVAL" ) ;
161 RecordFieldPtr<Double> exposureRF( trec, "EXPOSURE" ) ;
162 RecordFieldPtr< Array<Float> > weightRF( trec, "WEIGHT" ) ;
163 RecordFieldPtr< Array<Float> > sigmaRF( trec, "SIGMA" ) ;
164 RecordFieldPtr<Int> ddidRF( trec, "DATA_DESC_ID" ) ;
165 RecordFieldPtr<Int> stateidRF( trec, "STATE_ID" ) ;
[2020]166 RecordFieldPtr< Array<Bool> > flagcatRF( trec, "FLAG_CATEGORY" ) ;
[2016]167
168 // OBSERVATION_ID is always 0
169 RecordFieldPtr<Int> intRF( trec, "OBSERVATION_ID" ) ;
170 *intRF = 0 ;
171
172 // ANTENNA1 and ANTENNA2 are always 0
173 intRF.attachToRecord( trec, "ANTENNA1" ) ;
174 *intRF = 0 ;
175 intRF.attachToRecord( trec, "ANTENNA2" ) ;
176 *intRF = 0 ;
177
178 // ARRAY_ID is tentatively set to 0
179 intRF.attachToRecord( trec, "ARRAY_ID" ) ;
180 *intRF = 0 ;
181
182 // PROCESSOR_ID is tentatively set to 0
183 intRF.attachToRecord( trec, "PROCESSOR_ID" ) ;
184 *intRF = 0 ;
185
186 // UVW is always [0,0,0]
187 RecordFieldPtr< Array<Double> > uvwRF( trec, "UVW" ) ;
188 *uvwRF = Vector<Double>( 3, 0.0 ) ;
189
[1974]190 //
191 // ITERATION: FIELDNAME
192 //
193 TableIterator iter0( table_->table(), "FIELDNAME" ) ;
194 while( !iter0.pastEnd() ) {
[2016]195 //Table t0( iter0.table() ) ;
196 Table t0 = iter0.table() ;
197 ROTableColumn sharedCol( t0, "FIELDNAME" ) ;
198 String fieldName = sharedCol.asString( 0 ) ;
199// tcol = tpoolr->construct( t0, "FIELDNAME" ) ;
200// String fieldName = tcol->asString( 0 ) ;
201// tpoolr->destroy( tcol ) ;
202 sharedCol.attach( t0, "SRCNAME" ) ;
203 String srcName = sharedCol.asString( 0 ) ;
204// tcol = tpoolr->construct( t0, "SRCNAME" ) ;
205// String srcName = tcol->asString( 0 ) ;
206// tpoolr->destroy( tcol ) ;
207 sharedCol.attach( t0, "TIME" ) ;
208 Double minTime = (Double)sharedCol.asdouble( 0 ) * 86400.0 ; // day->sec
209// tcol = tpoolr->construct( t0, "TIME" ) ;
210// Double minTime = (Double)(tcol->asdouble( 0 )) * 86400.0 ; // day -> sec
211// tpoolr->destroy( tcol ) ;
[1975]212 ROArrayColumn<Double> scanRateCol( t0, "SCANRATE" ) ;
[2016]213 //Vector<Double> scanRate = scanRateCol.getColumn()[0] ;
214 Vector<Double> scanRate = scanRateCol( 0 ) ;
[1974]215 String::size_type pos = fieldName.find( "__" ) ;
216 Int fieldId = -1 ;
217 if ( pos != String::npos ) {
218 os_ << "fieldName.substr( pos+2 )=" << fieldName.substr( pos+2 ) << LogIO::POST ;
219 fieldId = String::toInt( fieldName.substr( pos+2 ) ) ;
220 fieldName = fieldName.substr( 0, pos ) ;
221 }
222 else {
223 os_ << "use default field id" << LogIO::POST ;
224 fieldId = defaultFieldId ;
225 defaultFieldId++ ;
226 }
227 os_ << "fieldId" << fieldId << ": " << fieldName << LogIO::POST ;
[2016]228
229 // FIELD_ID
230 intRF.attachToRecord( trec, "FIELD_ID" ) ;
231 *intRF = fieldId ;
232
[1974]233 //
234 // ITERATION: BEAMNO
235 //
236 TableIterator iter1( t0, "BEAMNO" ) ;
237 while( !iter1.pastEnd() ) {
[2016]238 Table t1 = iter1.table() ;
239// tcol = tpoolr->construct( t1, "BEAMNO" ) ;
240// uInt beamNo = tcol->asuInt( 0 ) ;
241// tpoolr->destroy( tcol ) ;
242 sharedCol.attach( t1, "BEAMNO" ) ;
243 uInt beamNo = sharedCol.asuInt( 0 ) ;
[1974]244 os_ << "beamNo = " << beamNo << LogIO::POST ;
[2016]245
246 // FEED1 and FEED2
247 intRF.attachToRecord( trec, "FEED1" ) ;
248 *intRF = beamNo ;
249 intRF.attachToRecord( trec, "FEED2" ) ;
250 *intRF = beamNo ;
251
[1974]252 //
253 // ITERATION: SCANNO
254 //
255 TableIterator iter2( t1, "SCANNO" ) ;
256 while( !iter2.pastEnd() ) {
[2016]257 Table t2 = iter2.table() ;
258// tcol = tpoolr->construct( t2, "SCANNO" ) ;
259// uInt scanNo = tcol->asuInt( 0 ) ;
260// tpoolr->destroy( tcol ) ;
261 sharedCol.attach( t2, "SCANNO" ) ;
262 uInt scanNo = sharedCol.asuInt( 0 ) ;
[1974]263 os_ << "scanNo = " << scanNo << LogIO::POST ;
[2016]264
265 // SCAN_NUMBER
266 // MS: 1-based
267 // Scantable: 0-based
268 intRF.attachToRecord( trec, "SCAN_NUMBER" ) ;
269 *intRF = scanNo + 1 ;
270
[1974]271 //
[1977]272 // ITERATION: IFNO
[1974]273 //
[1977]274 TableIterator iter3( t2, "IFNO" ) ;
[1974]275 while( !iter3.pastEnd() ) {
[2016]276 Table t3 = iter3.table() ;
277// tcol = tpoolr->construct( t3, "IFNO" ) ;
278// uInt ifNo = tcol->asuInt( 0 ) ;
279// tpoolr->destroy( tcol ) ;
280 sharedCol.attach( t3, "IFNO" ) ;
281 uInt ifNo = sharedCol.asuInt( 0 ) ;
[1977]282 os_ << "ifNo = " << ifNo << LogIO::POST ;
[2016]283// tcol = tpoolr->construct( t3, "FREQ_ID" ) ;
284// uInt freqId = tcol->asuInt( 0 ) ;
285// tpoolr->destroy( tcol ) ;
286 sharedCol.attach( t3, "FREQ_ID" ) ;
287 uInt freqId = sharedCol.asuInt( 0 ) ;
[1977]288 os_ << "freqId = " << freqId << LogIO::POST ;
[2016]289 Int subscan = 1 ; // 1-base
[1974]290 //
[2016]291 // ITERATION: SRCTYPE
[1974]292 //
[2016]293 TableIterator iter4( t3, "SRCTYPE" ) ;
[1974]294 while( !iter4.pastEnd() ) {
[2016]295 Table t4 = iter4.table() ;
296// tcol = tpoolr->construct( t4, "SRCTYPE" ) ;
297// uInt srcType = tcol->asInt( 0 ) ;
298// tpoolr->destroy( tcol ) ;
299 sharedCol.attach( t4, "SRCTYPE" ) ;
300 Int srcType = sharedCol.asInt( 0 ) ;
301 Int stateId = addState( srcType, subscan ) ;
302 *stateidRF = stateId ;
[1974]303 //
[2016]304 // ITERATION: CYCLENO and TIME
[1974]305 //
[2016]306 Block<String> cols( 2 ) ;
307 cols[0] = "CYCLENO" ;
308 cols[1] = "TIME" ;
309 TableIterator iter5( t4, cols ) ;
[1974]310 while( !iter5.pastEnd() ) {
[2019]311 Table t5 = iter5.table().sort("POLNO") ;
[2016]312 //sharedCol.attach( t5, "CYCLENO" ) ;
313 //uInt cycleNo = sharedCol.asuInt( 0 ) ;
[1974]314 Int nrow = t5.nrow() ;
315 os_ << "nrow = " << nrow << LogIO::POST ;
316
317 Vector<Int> polnos( nrow ) ;
318 indgen( polnos, 0 ) ;
319 Int polid = addPolarization( polnos ) ;
320 os_ << "polid = " << polid << LogIO::POST ;
[2016]321
322 // DATA/FLOAT_DATA
[1977]323 ROArrayColumn<Float> specCol( t5, "SPECTRA" ) ;
324 ROArrayColumn<uChar> flagCol( t5, "FLAGTRA" ) ;
325 uInt nchan = specCol( 0 ).size() ;
326 IPosition cellshape( 2, nrow, nchan ) ;
327 if ( useFloatData_ ) {
328 // FLOAT_DATA
[2016]329 Matrix<Float> dataArr( cellshape ) ;
330 Matrix<Bool> flagArr( cellshape ) ;
331 Vector<Bool> tmpB ;
[1977]332 for ( Int ipol = 0 ; ipol < nrow ; ipol++ ) {
[2016]333 dataArr.row( ipol ) = specCol( ipol ) ;
[2019]334 tmpB.reference( flagArr.row( ipol ) ) ;
335 convertArray( tmpB, flagCol( ipol ) ) ;
[1977]336 }
[2019]337 //*(*((RecordFieldPtr< Array<Float> > *)dataRF)) = dataArr ;
338 ((RecordFieldPtr< Array<Float> > *)dataRF)->define( dataArr ) ;
[1977]339
340 // FLAG
[2019]341 //*flagRF = flagArr ;
342 flagRF.define( flagArr ) ;
[1974]343 }
[1977]344 else if ( useData_ ) {
345 // DATA
346 // assume nrow = 4
[2016]347 Matrix<Complex> dataArr( cellshape ) ;
[1977]348 Vector<Float> zeroIm( nchan, 0 ) ;
[2016]349 Matrix<Float> dummy( IPosition( 2, 2, nchan ) ) ;
350 dummy.row( 0 ) = specCol( 0 ) ;
351 dummy.row( 1 ) = zeroIm ;
352 dataArr.row( 0 ) = RealToComplex( dummy ) ;
353 dummy.row( 0 ) = specCol( 1 ) ;
354 dataArr.row( 3 ) = RealToComplex( dummy ) ;
355 dummy.row( 0 ) = specCol( 2 ) ;
356 dummy.row( 1 ) = specCol( 3 ) ;
357 dataArr.row( 1 ) = RealToComplex( dummy ) ;
358 dataArr.row( 2 ) = conj( dataArr.row( 1 ) ) ;
[2019]359 //*(*((RecordFieldPtr< Array<Complex> > *)dataRF)) = dataArr ;
360 ((RecordFieldPtr< Array<Complex> > *)dataRF)->define( dataArr ) ;
[2016]361
362
[1977]363 // FLAG
[2016]364 Matrix<Bool> flagArr( cellshape ) ;
[2019]365 Vector<Bool> tmpB ;
366 tmpB.reference( flagArr.row( 0 ) ) ;
367 convertArray( tmpB, flagCol( 0 ) ) ;
368 tmpB.reference( flagArr.row( 3 ) ) ;
369 convertArray( tmpB, flagCol( 3 ) ) ;
370 tmpB.reference( flagArr.row( 1 ) ) ;
371 convertArray( tmpB, ( flagCol( 2 ) | flagCol( 3 ) ) ) ;
372 flagArr.row( 2 ) = flagArr.row( 1 ) ;
373 //*flagRF = flagArr ;
374 flagRF.define( flagArr ) ;
[1977]375 }
[2019]376
[1977]377 // FLAG_ROW
[2016]378// tcol = tpoolr->construct( t5, "FLAGROW" ) ;
379 sharedCol.attach( t5, "FLAGROW" ) ;
380 Vector<uInt> flagRowArr( nrow ) ;
381 for ( Int irow = 0 ; irow < nrow ; irow++ )
382 flagRowArr[irow] = sharedCol.asuInt( irow ) ;
383// flagRowArr[irow] = tcol->asuInt( irow ) ;
384// tpoolr->destroy( tcol ) ;
385 *flagrowRF = anyNE( flagRowArr, (uInt)0 ) ;
[2019]386
[1974]387 // TIME and TIME_CENTROID
[2016]388// tcol = tpoolr->construct( t5, "TIME" ) ;
389// Double mTimeV = (Double)(tcol->asdouble(0)) * 86400.0 ; // day->sec
390// tpoolr->destroy( tcol ) ;
391 sharedCol.attach( t5, "TIME" ) ;
392 Double mTimeV = (Double)sharedCol.asdouble( 0 ) * 86400.0 ; // day -> sec
393 *timeRF = mTimeV ;
394 *timecRF = mTimeV ;
[2019]395
[1974]396 // INTERVAL and EXPOSURE
[2016]397// tcol = tpoolr->construct( t5, "INTERVAL" ) ;
398// Double interval = (Double)(tcol->asdouble( 0 )) ;
399// tpoolr->destroy( tcol ) ;
400 sharedCol.attach( t5, "INTERVAL" ) ;
401 Double interval = (Double)sharedCol.asdouble( 0 ) ;
402 *intervalRF = interval ;
403 *exposureRF = interval ;
404
[1977]405 // WEIGHT and SIGMA
406 // always 1 at the moment
407 Vector<Float> wArr( nrow, 1.0 ) ;
[2021]408// *weightRF = wArr ;
409// *sigmaRF = wArr ;
410 weightRF.define( wArr ) ;
411 sigmaRF.define( wArr ) ;
[1974]412
413 // add DATA_DESCRIPTION row
414 Int ddid = addDataDescription( polid, ifNo ) ;
415 os_ << "ddid = " << ddid << LogIO::POST ;
[2016]416 *ddidRF = ddid ;
[1974]417
[1977]418 // for SYSCAL table
[2016]419// tcol = tpoolr->construct( t5, "TCAL_ID" ) ;
420 sharedCol.attach( t5, "TCAL_ID" ) ;
421 Vector<uInt> tcalIdArr( nrow ) ;
422 for ( Int irow = 0 ; irow < nrow ; irow++ )
423 tcalIdArr[irow] = sharedCol.asuInt( irow ) ;
424// tcalIdArr[irow] = tcol->asuInt( irow ) ;
425// tpoolr->destroy( tcol ) ;
[1977]426 os_ << "tcalIdArr = " << tcalIdArr << LogIO::POST ;
427 String key = String::toString( tcalIdArr[0] ) ;
[2018]428 if ( !tcalIdRec_.isDefined( key ) ) {
[1977]429 tcalIdRec_.define( key, tcalIdArr ) ;
[2018]430 tcalRowRec_.define( key, t5.rowNumbers() ) ;
431 }
432 else {
433 Vector<uInt> pastrows = tcalRowRec_.asArrayuInt( key ) ;
434 tcalRowRec_.define( key, concatenateArray( pastrows, t5.rowNumbers() ) ) ;
435 }
[1977]436
437 // fill STATE_ID
[2016]438 //ROScalarColumn<Int> srcTypeCol( t5, "SRCTYPE" ) ;
439 //Int srcType = srcTypeCol( 0 ) ;
440 //Int stateId = addState( srcType, subscan ) ;
441 //*stateidRF = stateId ;
442
[1977]443 // for POINTING table
[2016]444 if ( ptTabName_ == "" ) {
445 ROArrayColumn<Double> dirCol( t5, "DIRECTION" ) ;
446 Vector<Double> dir = dirCol( 0 ) ;
447 dirCol.attach( t5, "SCANRATE" ) ;
448 Vector<Double> rate = dirCol( 0 ) ;
449 Matrix<Double> msDir( 2, 1 ) ;
450 msDir.column( 0 ) = dir ;
451 if ( anyNE( rate, 0.0 ) ) {
452 msDir.resize( 2, 2 ) ;
453 msDir.column( 1 ) = rate ;
454 }
455 addPointing( fieldName, mTimeV, interval, msDir ) ;
456 }
457
458 // FLAG_CATEGORY is tentatively set
[2019]459 //*flagcatRF = Cube<Bool>( nrow, nchan, 1, False ) ;
460 flagcatRF.define( Cube<Bool>( nrow, nchan, 1, False ) ) ;
[2016]461
462 // add row
463 mstable_->addRow( 1, True ) ;
464 row.put( mstable_->nrow()-1 ) ;
465
[1974]466 iter5.next() ;
467 }
468
469 iter4.next() ;
470 }
[1977]471
472 // add SPECTRAL_WINDOW row
473 if ( allNE( processedFreqId, freqId ) ) {
474 uInt vsize = processedFreqId.size() ;
475 processedFreqId.resize( vsize+1, True ) ;
476 processedFreqId[vsize] = freqId ;
477 addSpectralWindow( ifNo, freqId ) ;
478 }
479
[1974]480 iter3.next() ;
481 }
482
483 iter2.next() ;
484 }
485
486 // add FEED row
487 addFeed( beamNo ) ;
488
489 iter1.next() ;
490 }
491
[1975]492 // add FIELD row
493 addField( fieldId, fieldName, srcName, minTime, scanRate ) ;
494
[1974]495 iter0.next() ;
496 }
497
[2016]498// delete tpoolr ;
499 delete dataRF ;
[1974]500
[2016]501 // SYSCAL
[2019]502 if ( isTcal_ )
503 fillSysCal() ;
[1974]504
[2022]505 // ASDM tables
506 const TableRecord &stKeys = table_->table().keywordSet() ;
507 TableRecord &msKeys = mstable_->rwKeywordSet() ;
508 uInt nfields = stKeys.nfields() ;
509 for ( uInt ifield = 0 ; ifield < nfields ; ifield++ ) {
510 String kname = stKeys.name( ifield ) ;
511 if ( kname.find( "ASDM" ) != String::npos ) {
512 String asdmpath = stKeys.asString( ifield ) ;
513 os_ << "found ASDM table: " << asdmpath << LogIO::POST ;
514 if ( Table::isReadable( asdmpath ) ) {
515 Table newAsdmTab( asdmpath, Table::Old ) ;
516 newAsdmTab.copy( filename_+"/"+kname, Table::New ) ;
517 os_ << "add subtable: " << kname << LogIO::POST ;
518 msKeys.defineTable( kname, Table( filename_+"/"+kname, Table::Old ) ) ;
519 }
520 }
521 }
522
[2016]523 // replace POINTING table with original one if exists
524 if ( ptTabName_ != "" ) {
525 delete mstable_ ;
526 mstable_ = 0 ;
527 Table newPtTab( ptTabName_, Table::Old ) ;
528 newPtTab.copy( filename_+"/POINTING", Table::New ) ;
529 }
[1974]530
[2016]531 double endSec = gettimeofday_sec() ;
532 os_ << "end MSWriter::write() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1975]533
[1974]534 return True ;
535}
536
537void MSWriter::init()
538{
539// os_.origin( LogOrigin( "MSWriter", "init()", WHERE ) ) ;
[2016]540 double startSec = gettimeofday_sec() ;
541 os_ << "start MSWriter::init() startSec=" << startSec << LogIO::POST ;
[1974]542
543 // access to scantable
544 header_ = table_->getHeader() ;
545
546 // FLOAT_DATA? or DATA?
547 if ( header_.npol > 2 ) {
[1977]548 useFloatData_ = False ;
549 useData_ = True ;
[1974]550 }
551 else {
[1977]552 useFloatData_ = True ;
553 useData_ = False ;
[1974]554 }
555
556 // polarization type
557 polType_ = header_.poltype ;
558 if ( polType_ == "" )
559 polType_ = "stokes" ;
560 else if ( polType_.find( "linear" ) != String::npos )
561 polType_ = "linear" ;
562 else if ( polType_.find( "circular" ) != String::npos )
563 polType_ = "circular" ;
564 else if ( polType_.find( "stokes" ) != String::npos )
565 polType_ = "stokes" ;
566 else if ( polType_.find( "linpol" ) != String::npos )
567 polType_ = "linpol" ;
568 else
569 polType_ = "notype" ;
570
[2019]571 // Check if some subtables are exists
572 if ( table_->tcal().table().nrow() != 0 ) {
573 ROTableColumn col( table_->tcal().table(), "TCAL" ) ;
574 if ( col.isDefined( 0 ) ) {
575 os_ << "nrow=" << table_->tcal().table().nrow() << ": TCAL table exists" << LogIO::POST ;
576 isTcal_ = True ;
577 }
578 else {
579 os_ << "no TCAL rows" << LogIO::POST ;
580 }
581 }
582 else {
583 os_ << "no TCAL rows" << LogIO::POST ;
584 }
585 if ( table_->weather().table().nrow() != 0 ) {
586 ROTableColumn col( table_->weather().table(), "TEMPERATURE" ) ;
587 if ( col.isDefined( 0 ) ) {
588 os_ << "nrow=" << table_->weather().table().nrow() << ": WEATHER table exists" << LogIO::POST ;
589 isWeather_ =True ;
590 }
591 else {
592 os_ << "no WEATHER rows" << LogIO::POST ;
593 }
594 }
595 else {
596 os_ << "no WEATHER rows" << LogIO::POST ;
597 }
598
[1977]599 // Are TCAL_SPECTRUM and TSYS_SPECTRUM necessary?
[2019]600 if ( isTcal_ && header_.nchan != 1 ) {
[1977]601 // examine TCAL subtable
602 Table tcaltab = table_->tcal().table() ;
603 ROArrayColumn<Float> tcalCol( tcaltab, "TCAL" ) ;
604 for ( uInt irow = 0 ; irow < tcaltab.nrow() ; irow++ ) {
605 if ( tcalCol( irow ).size() != 1 )
606 tcalSpec_ = True ;
607 }
608 // examine spectral data
609 TableIterator iter0( table_->table(), "IFNO" ) ;
610 while( !iter0.pastEnd() ) {
611 Table t0( iter0.table() ) ;
612 ROArrayColumn<Float> sharedFloatArrCol( t0, "SPECTRA" ) ;
613 uInt len = sharedFloatArrCol( 0 ).size() ;
614 if ( len != 1 ) {
615 sharedFloatArrCol.attach( t0, "TSYS" ) ;
616 if ( sharedFloatArrCol( 0 ).size() != 1 )
617 tsysSpec_ = True ;
618 }
619 iter0.next() ;
620 }
621 }
[2016]622
623 // check if reference for POINTING table exists
624 const TableRecord &rec = table_->table().keywordSet() ;
625 if ( rec.isDefined( "POINTING" ) ) {
626 ptTabName_ = rec.asString( "POINTING" ) ;
627 if ( !Table::isReadable( ptTabName_ ) ) {
628 ptTabName_ = "" ;
629 }
630 }
631
632 double endSec = gettimeofday_sec() ;
633 os_ << "end MSWriter::init() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1974]634}
635
636void MSWriter::setupMS()
637{
638// os_.origin( LogOrigin( "MSWriter", "setupMS()", WHERE ) ) ;
[2016]639 double startSec = gettimeofday_sec() ;
640 os_ << "start MSWriter::setupMS() startSec=" << startSec << LogIO::POST ;
[1974]641
642 TableDesc msDesc = MeasurementSet::requiredTableDesc() ;
[1977]643 if ( useFloatData_ )
644 MeasurementSet::addColumnToDesc( msDesc, MSMainEnums::FLOAT_DATA, 2 ) ;
645 else if ( useData_ )
646 MeasurementSet::addColumnToDesc( msDesc, MSMainEnums::DATA, 2 ) ;
[1974]647
648 SetupNewTable newtab( filename_, msDesc, Table::New ) ;
649
650 mstable_ = new MeasurementSet( newtab ) ;
651
652 // create subtables
653 TableDesc antennaDesc = MSAntenna::requiredTableDesc() ;
654 SetupNewTable antennaTab( mstable_->antennaTableName(), antennaDesc, Table::New ) ;
655 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::ANTENNA ), Table( antennaTab ) ) ;
656
657 TableDesc dataDescDesc = MSDataDescription::requiredTableDesc() ;
658 SetupNewTable dataDescTab( mstable_->dataDescriptionTableName(), dataDescDesc, Table::New ) ;
659 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::DATA_DESCRIPTION ), Table( dataDescTab ) ) ;
660
661 TableDesc dopplerDesc = MSDoppler::requiredTableDesc() ;
662 SetupNewTable dopplerTab( mstable_->dopplerTableName(), dopplerDesc, Table::New ) ;
663 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::DOPPLER ), Table( dopplerTab ) ) ;
664
665 TableDesc feedDesc = MSFeed::requiredTableDesc() ;
666 SetupNewTable feedTab( mstable_->feedTableName(), feedDesc, Table::New ) ;
667 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::FEED ), Table( feedTab ) ) ;
668
669 TableDesc fieldDesc = MSField::requiredTableDesc() ;
670 SetupNewTable fieldTab( mstable_->fieldTableName(), fieldDesc, Table::New ) ;
671 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::FIELD ), Table( fieldTab ) ) ;
672
673 TableDesc flagCmdDesc = MSFlagCmd::requiredTableDesc() ;
674 SetupNewTable flagCmdTab( mstable_->flagCmdTableName(), flagCmdDesc, Table::New ) ;
675 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::FLAG_CMD ), Table( flagCmdTab ) ) ;
676
677 TableDesc freqOffsetDesc = MSFreqOffset::requiredTableDesc() ;
678 SetupNewTable freqOffsetTab( mstable_->freqOffsetTableName(), freqOffsetDesc, Table::New ) ;
679 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::FREQ_OFFSET ), Table( freqOffsetTab ) ) ;
680
681 TableDesc historyDesc = MSHistory::requiredTableDesc() ;
682 SetupNewTable historyTab( mstable_->historyTableName(), historyDesc, Table::New ) ;
683 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::HISTORY ), Table( historyTab ) ) ;
684
685 TableDesc observationDesc = MSObservation::requiredTableDesc() ;
686 SetupNewTable observationTab( mstable_->observationTableName(), observationDesc, Table::New ) ;
687 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::OBSERVATION ), Table( observationTab ) ) ;
688
689 TableDesc pointingDesc = MSPointing::requiredTableDesc() ;
690 SetupNewTable pointingTab( mstable_->pointingTableName(), pointingDesc, Table::New ) ;
691 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::POINTING ), Table( pointingTab ) ) ;
692
693 TableDesc polarizationDesc = MSPolarization::requiredTableDesc() ;
694 SetupNewTable polarizationTab( mstable_->polarizationTableName(), polarizationDesc, Table::New ) ;
695 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::POLARIZATION ), Table( polarizationTab ) ) ;
696
697 TableDesc processorDesc = MSProcessor::requiredTableDesc() ;
698 SetupNewTable processorTab( mstable_->processorTableName(), processorDesc, Table::New ) ;
699 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::PROCESSOR ), Table( processorTab ) ) ;
700
701 TableDesc sourceDesc = MSSource::requiredTableDesc() ;
[1975]702 MSSource::addColumnToDesc( sourceDesc, MSSourceEnums::TRANSITION, 1 ) ;
703 MSSource::addColumnToDesc( sourceDesc, MSSourceEnums::REST_FREQUENCY, 1 ) ;
704 MSSource::addColumnToDesc( sourceDesc, MSSourceEnums::SYSVEL, 1 ) ;
[1974]705 SetupNewTable sourceTab( mstable_->sourceTableName(), sourceDesc, Table::New ) ;
706 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::SOURCE ), Table( sourceTab ) ) ;
707
708 TableDesc spwDesc = MSSpectralWindow::requiredTableDesc() ;
709 SetupNewTable spwTab( mstable_->spectralWindowTableName(), spwDesc, Table::New ) ;
710 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::SPECTRAL_WINDOW ), Table( spwTab ) ) ;
711
712 TableDesc stateDesc = MSState::requiredTableDesc() ;
713 SetupNewTable stateTab( mstable_->stateTableName(), stateDesc, Table::New ) ;
714 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::STATE ), Table( stateTab ) ) ;
715
716 // TODO: add TCAL_SPECTRUM and TSYS_SPECTRUM if necessary
717 TableDesc sysCalDesc = MSSysCal::requiredTableDesc() ;
[1977]718 MSSysCal::addColumnToDesc( sysCalDesc, MSSysCalEnums::TCAL, 2 ) ;
719 MSSysCal::addColumnToDesc( sysCalDesc, MSSysCalEnums::TSYS, 2 ) ;
720 if ( tcalSpec_ )
721 MSSysCal::addColumnToDesc( sysCalDesc, MSSysCalEnums::TCAL_SPECTRUM, 2 ) ;
722 if ( tsysSpec_ )
723 MSSysCal::addColumnToDesc( sysCalDesc, MSSysCalEnums::TSYS_SPECTRUM, 2 ) ;
[1974]724 SetupNewTable sysCalTab( mstable_->sysCalTableName(), sysCalDesc, Table::New ) ;
725 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::SYSCAL ), Table( sysCalTab ) ) ;
726
727 TableDesc weatherDesc = MSWeather::requiredTableDesc() ;
[1977]728 MSWeather::addColumnToDesc( weatherDesc, MSWeatherEnums::TEMPERATURE ) ;
729 MSWeather::addColumnToDesc( weatherDesc, MSWeatherEnums::PRESSURE ) ;
730 MSWeather::addColumnToDesc( weatherDesc, MSWeatherEnums::REL_HUMIDITY ) ;
731 MSWeather::addColumnToDesc( weatherDesc, MSWeatherEnums::WIND_SPEED ) ;
732 MSWeather::addColumnToDesc( weatherDesc, MSWeatherEnums::WIND_DIRECTION ) ;
[1974]733 SetupNewTable weatherTab( mstable_->weatherTableName(), weatherDesc, Table::New ) ;
734 mstable_->rwKeywordSet().defineTable( MeasurementSet::keywordName( MeasurementSet::WEATHER ), Table( weatherTab ) ) ;
735
736 mstable_->initRefs() ;
737
[2016]738 double endSec = gettimeofday_sec() ;
739 os_ << "end MSWriter::setupMS() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1974]740}
741
742void MSWriter::fillObservation()
743{
[2016]744 double startSec = gettimeofday_sec() ;
745 os_ << "start MSWriter::fillObservation() startSec=" << startSec << LogIO::POST ;
746
[1974]747 // only 1 row
748 mstable_->observation().addRow( 1, True ) ;
749 MSObservationColumns msObsCols( mstable_->observation() ) ;
750 msObsCols.observer().put( 0, header_.observer ) ;
751 // tentatively put antennaname (from ANTENNA subtable)
752 String hAntennaName = header_.antennaname ;
753 String::size_type pos = hAntennaName.find( "//" ) ;
754 String telescopeName ;
755 if ( pos != String::npos ) {
756 telescopeName = hAntennaName.substr( 0, pos ) ;
757 }
758 else {
759 pos = hAntennaName.find( "@" ) ;
760 telescopeName = hAntennaName.substr( 0, pos ) ;
761 }
762 os_ << "telescopeName = " << telescopeName << LogIO::POST ;
763 msObsCols.telescopeName().put( 0, telescopeName ) ;
764 msObsCols.project().put( 0, header_.project ) ;
765 //ScalarMeasColumn<MEpoch> timeCol( table_->table().sort("TIME"), "TIME" ) ;
766 Table sortedtable = table_->table().sort("TIME") ;
767 ScalarMeasColumn<MEpoch> timeCol( sortedtable, "TIME" ) ;
768 Vector<MEpoch> trange( 2 ) ;
769 trange[0] = timeCol( 0 ) ;
770 trange[1] = timeCol( table_->nrow()-1 ) ;
771 msObsCols.timeRangeMeas().put( 0, trange ) ;
[2016]772
773 double endSec = gettimeofday_sec() ;
774 os_ << "end MSWriter::fillObservation() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1974]775}
776
777void MSWriter::fillAntenna()
778{
[2016]779 double startSec = gettimeofday_sec() ;
780 os_ << "start MSWriter::fillAntenna() startSec=" << startSec << LogIO::POST ;
781
[1974]782 // only 1 row
783 mstable_->antenna().addRow( 1, True ) ;
784 MSAntennaColumns msAntCols( mstable_->antenna() ) ;
785
786 String hAntennaName = header_.antennaname ;
787 String::size_type pos = hAntennaName.find( "//" ) ;
788 String antennaName ;
789 String stationName ;
790 if ( pos != String::npos ) {
791 hAntennaName = hAntennaName.substr( pos+2 ) ;
792 }
793 pos = hAntennaName.find( "@" ) ;
794 if ( pos != String::npos ) {
795 antennaName = hAntennaName.substr( 0, pos ) ;
796 stationName = hAntennaName.substr( pos+1 ) ;
797 }
798 else {
799 antennaName = hAntennaName ;
800 stationName = hAntennaName ;
801 }
802 os_ << "antennaName = " << antennaName << LogIO::POST ;
803 os_ << "stationName = " << stationName << LogIO::POST ;
804
805 msAntCols.name().put( 0, antennaName ) ;
806 msAntCols.station().put( 0, stationName ) ;
807
808 os_ << "antennaPosition = " << header_.antennaposition << LogIO::POST ;
809
810 msAntCols.position().put( 0, header_.antennaposition ) ;
[2016]811
812 double endSec = gettimeofday_sec() ;
813 os_ << "end MSWriter::fillAntenna() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1974]814}
815
[1975]816void MSWriter::fillProcessor()
817{
[2016]818 double startSec = gettimeofday_sec() ;
819 os_ << "start MSWriter::fillProcessor() startSec=" << startSec << LogIO::POST ;
[1975]820
821 // only add empty 1 row
822 MSProcessor msProc = mstable_->processor() ;
823 msProc.addRow( 1, True ) ;
[2016]824
825 double endSec = gettimeofday_sec() ;
826 os_ << "end MSWriter::fillProcessor() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1975]827}
828
829void MSWriter::fillSource()
830{
[2016]831 double startSec = gettimeofday_sec() ;
832 os_ << "start MSWriter::fillSource() startSec=" << startSec << LogIO::POST ;
833
[1975]834 // access to MS SOURCE subtable
835 MSSource msSrc = mstable_->source() ;
836
837 // access to MOLECULE subtable
838 STMolecules stm = table_->molecules() ;
839
840 Int srcId = 0 ;
841 Vector<Double> restFreq ;
842 Vector<String> molName ;
843 Vector<String> fMolName ;
[2016]844
845 // row based
846 TableRow row( msSrc ) ;
847 TableRecord &rec = row.record() ;
848 RecordFieldPtr<Int> srcidRF( rec, "SOURCE_ID" ) ;
849 RecordFieldPtr<String> nameRF( rec, "NAME" ) ;
850 RecordFieldPtr< Array<Double> > srcpmRF( rec, "PROPER_MOTION" ) ;
851 RecordFieldPtr< Array<Double> > srcdirRF( rec, "DIRECTION" ) ;
852 RecordFieldPtr<Int> numlineRF( rec, "NUM_LINES" ) ;
853 RecordFieldPtr< Array<Double> > restfreqRF( rec, "REST_FREQUENCY" ) ;
854 RecordFieldPtr< Array<Double> > sysvelRF( rec, "SYSVEL" ) ;
855 RecordFieldPtr< Array<String> > transitionRF( rec, "TRANSITION" ) ;
856 RecordFieldPtr<Double> timeRF( rec, "TIME" ) ;
857 RecordFieldPtr<Double> intervalRF( rec, "INTERVAL" ) ;
858 RecordFieldPtr<Int> spwidRF( rec, "SPECTRAL_WINDOW_ID" ) ;
859
[1975]860 //
861 // ITERATION: SRCNAME
862 //
863 TableIterator iter0( table_->table(), "SRCNAME" ) ;
864 while( !iter0.pastEnd() ) {
[2016]865 //Table t0( iter0.table() ) ;
866 Table t0 = iter0.table() ;
[1975]867
868 // get necessary information
869 ROScalarColumn<String> srcNameCol( t0, "SRCNAME" ) ;
870 String srcName = srcNameCol( 0 ) ;
871 ROArrayColumn<Double> sharedDArrRCol( t0, "SRCPROPERMOTION" ) ;
872 Vector<Double> srcPM = sharedDArrRCol( 0 ) ;
873 sharedDArrRCol.attach( t0, "SRCDIRECTION" ) ;
874 Vector<Double> srcDir = sharedDArrRCol( 0 ) ;
875 ROScalarColumn<Double> srcVelCol( t0, "SRCVELOCITY" ) ;
876 Double srcVel = srcVelCol( 0 ) ;
877
[2016]878 // NAME
879 *nameRF = srcName ;
880
881 // SOURCE_ID
882 *srcidRF = srcId ;
883
884 // PROPER_MOTION
885 *srcpmRF = srcPM ;
886
887 // DIRECTION
888 *srcdirRF = srcDir ;
889
[1975]890 //
891 // ITERATION: MOLECULE_ID
892 //
893 TableIterator iter1( t0, "MOLECULE_ID" ) ;
894 while( !iter1.pastEnd() ) {
[2016]895 //Table t1( iter1.table() ) ;
896 Table t1 = iter1.table() ;
[1975]897
898 // get necessary information
899 ROScalarColumn<uInt> molIdCol( t1, "MOLECULE_ID" ) ;
900 uInt molId = molIdCol( 0 ) ;
901 stm.getEntry( restFreq, molName, fMolName, molId ) ;
902
[2016]903 uInt numFreq = restFreq.size() ;
904
905 // NUM_LINES
906 *numlineRF = numFreq ;
907
908 // REST_FREQUENCY
909 *restfreqRF = restFreq ;
910
911 // TRANSITION
[2019]912 //*transitionRF = fMolName ;
913 Vector<String> transition ;
914 if ( fMolName.size() != 0 ) {
915 transition = fMolName ;
916 }
917 else if ( molName.size() != 0 ) {
918 transition = molName ;
919 }
920 else {
921 transition.resize( numFreq ) ;
922 transition = "" ;
923 }
924 *transitionRF = transition ;
[2016]925
926 // SYSVEL
927 Vector<Double> sysvelArr( numFreq, srcVel ) ;
928 *sysvelRF = sysvelArr ;
929
[1975]930 //
931 // ITERATION: IFNO
932 //
933 TableIterator iter2( t1, "IFNO" ) ;
934 while( !iter2.pastEnd() ) {
[2016]935 //Table t2( iter2.table() ) ;
936 Table t2 = iter2.table() ;
937 uInt nrow = msSrc.nrow() ;
[1975]938
939 // get necessary information
940 ROScalarColumn<uInt> ifNoCol( t2, "IFNO" ) ;
941 uInt ifno = ifNoCol( 0 ) ; // IFNO = SPECTRAL_WINDOW_ID
[2016]942 Double midTime ;
[1975]943 Double interval ;
[2016]944 getValidTimeRange( midTime, interval, t2 ) ;
[1975]945
946 // fill SPECTRAL_WINDOW_ID
[2016]947 *spwidRF = ifno ;
[1975]948
949 // fill TIME, INTERVAL
[2016]950 *timeRF = midTime ;
951 *intervalRF = interval ;
[1975]952
[2016]953 // add row
954 msSrc.addRow( 1, True ) ;
955 row.put( nrow ) ;
956
[1975]957 iter2.next() ;
958 }
959
960 iter1.next() ;
961 }
962
963 // increment srcId if SRCNAME changed
964 srcId++ ;
965
966 iter0.next() ;
967 }
[2016]968
969 double endSec = gettimeofday_sec() ;
970 os_ << "end MSWriter::fillSource() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1975]971}
972
[1977]973void MSWriter::fillWeather()
974{
[2016]975 double startSec = gettimeofday_sec() ;
976 os_ << "start MSWriter::fillWeather() startSec=" << startSec << LogIO::POST ;
[1977]977
978 // access to MS WEATHER subtable
979 MSWeather msw = mstable_->weather() ;
980
981 // access to WEATHER subtable
982 Table stw = table_->weather().table() ;
983 uInt nrow = stw.nrow() ;
984
985 if ( nrow == 0 )
986 return ;
987
988 msw.addRow( nrow, True ) ;
989 MSWeatherColumns mswCols( msw ) ;
990
991 // ANTENNA_ID is always 0
992 Vector<Int> antIdArr( nrow, 0 ) ;
993 mswCols.antennaId().putColumn( antIdArr ) ;
994
995 // fill weather status
996 ROScalarColumn<Float> sharedFloatCol( stw, "TEMPERATURE" ) ;
997 mswCols.temperature().putColumn( sharedFloatCol ) ;
998 sharedFloatCol.attach( stw, "PRESSURE" ) ;
999 mswCols.pressure().putColumn( sharedFloatCol ) ;
1000 sharedFloatCol.attach( stw, "HUMIDITY" ) ;
1001 mswCols.relHumidity().putColumn( sharedFloatCol ) ;
1002 sharedFloatCol.attach( stw, "WINDSPEED" ) ;
1003 mswCols.windSpeed().putColumn( sharedFloatCol ) ;
1004 sharedFloatCol.attach( stw, "WINDAZ" ) ;
1005 mswCols.windDirection().putColumn( sharedFloatCol ) ;
1006
1007 // fill TIME and INTERVAL
[2016]1008 Double midTime ;
[1977]1009 Double interval ;
1010 Vector<Double> intervalArr( nrow, 0.0 ) ;
1011 TableIterator iter( table_->table(), "WEATHER_ID" ) ;
1012 while( !iter.pastEnd() ) {
[2016]1013 //Table tab( iter.table() ) ;
1014 Table tab = iter.table() ;
[1977]1015
1016 ROScalarColumn<uInt> widCol( tab, "WEATHER_ID" ) ;
1017 uInt wid = widCol( 0 ) ;
1018
[2016]1019 getValidTimeRange( midTime, interval, tab ) ;
1020 mswCols.time().put( wid, midTime ) ;
[1977]1021 intervalArr[wid] = interval ;
1022
1023 iter.next() ;
1024 }
1025 mswCols.interval().putColumn( intervalArr ) ;
[2016]1026
1027 double endSec = gettimeofday_sec() ;
1028 os_ << "end MSWriter::fillWeather() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1977]1029}
1030
1031void MSWriter::fillSysCal()
1032{
[2016]1033 double startSec = gettimeofday_sec() ;
1034 os_ << "start MSWriter::fillSysCal() startSec=" << startSec << LogIO::POST ;
[1977]1035
[2016]1036 //tcalIdRec_.print( cout ) ;
[1977]1037
1038 // access to MS SYSCAL subtable
1039 MSSysCal mssc = mstable_->sysCal() ;
1040
1041 // access to TCAL subtable
[2016]1042 Table stt = table_->tcal().table() ;
1043 uInt nrow = stt.nrow() ;
[1977]1044
[2016]1045 // access to MAIN table
[2018]1046 Block<String> cols( 6 ) ;
[2016]1047 cols[0] = "TIME" ;
1048 cols[1] = "TCAL_ID" ;
1049 cols[2] = "TSYS" ;
1050 cols[3] = "BEAMNO" ;
1051 cols[4] = "IFNO" ;
[2018]1052 cols[5] = "INTERVAL" ;
[2016]1053 Table tab = table_->table().project( cols ) ;
1054
[1977]1055 if ( nrow == 0 )
1056 return ;
1057
1058 nrow = tcalIdRec_.nfields() ;
1059
[2016]1060 Double midTime ;
[1977]1061 Double interval ;
1062 String timeStr ;
1063
[2016]1064 // row base
1065 TableRow row( mssc ) ;
1066 TableRecord &trec = row.record() ;
1067 RecordFieldPtr<Int> antennaRF( trec, "ANTENNA_ID" ) ;
1068 RecordFieldPtr<Int> feedRF( trec, "FEED_ID" ) ;
1069 RecordFieldPtr<Int> spwRF( trec, "SPECTRAL_WINDOW_ID" ) ;
1070 RecordFieldPtr<Double> timeRF( trec, "TIME" ) ;
1071 RecordFieldPtr<Double> intervalRF( trec, "INTERVAL" ) ;
1072 RecordFieldPtr< Array<Float> > tsysRF( trec, "TSYS" ) ;
1073 RecordFieldPtr< Array<Float> > tcalRF( trec, "TCAL" ) ;
1074 RecordFieldPtr< Array<Float> > tsysspRF ;
1075 RecordFieldPtr< Array<Float> > tcalspRF ;
1076 if ( tsysSpec_ )
1077 tsysspRF.attachToRecord( trec, "TSYS_SPECTRUM" ) ;
1078 if ( tcalSpec_ )
1079 tcalspRF.attachToRecord( trec, "TCAL_SPECTRUM" ) ;
[1977]1080
[2016]1081 // ANTENNA_ID is always 0
1082 *antennaRF = 0 ;
[1977]1083
[2016]1084 Table sortedstt = stt.sort( "ID" ) ;
1085 ROArrayColumn<Float> tcalCol( sortedstt, "TCAL" ) ;
1086 ROTableColumn idCol( sortedstt, "ID" ) ;
[2018]1087 ROArrayColumn<Float> tsysCol( tab, "TSYS" ) ;
1088 ROTableColumn tcalidCol( tab, "TCAL_ID" ) ;
1089 ROTableColumn timeCol( tab, "TIME" ) ;
1090 ROTableColumn intervalCol( tab, "INTERVAL" ) ;
1091 ROTableColumn beamnoCol( tab, "BEAMNO" ) ;
1092 ROTableColumn ifnoCol( tab, "IFNO" ) ;
[2016]1093 for ( uInt irow = 0 ; irow < nrow ; irow++ ) {
1094 double t1 = gettimeofday_sec() ;
1095 Vector<uInt> ids = tcalIdRec_.asArrayuInt( irow ) ;
1096 os_ << "ids = " << ids << LogIO::POST ;
1097 uInt npol = ids.size() ;
[2018]1098 Vector<uInt> rows = tcalRowRec_.asArrayuInt( irow ) ;
1099 os_ << "rows = " << rows << LogIO::POST ;
1100 Vector<Double> atime( rows.nelements() ) ;
1101 Vector<Double> ainterval( rows.nelements() ) ;
1102 Vector<uInt> atcalid( rows.nelements() ) ;
1103 for( uInt jrow = 0 ; jrow < rows.nelements() ; jrow++ ) {
1104 atime[jrow] = (Double)timeCol.asdouble( rows[jrow] ) ;
1105 ainterval[jrow] = (Double)intervalCol.asdouble( rows[jrow] ) ;
1106 atcalid[jrow] = tcalidCol.asuInt( rows[jrow] ) ;
1107 }
1108 Vector<Float> dummy = tsysCol( rows[0] ) ;
1109 Matrix<Float> tsys( npol,dummy.nelements() ) ;
1110 tsys.row( 0 ) = dummy ;
1111 for ( uInt jrow = 1 ; jrow < npol ; jrow++ )
1112 tsys.row( jrow ) = tsysCol( rows[jrow] ) ;
[1977]1113
[2016]1114 // FEED_ID
[2018]1115 *feedRF = beamnoCol.asuInt( rows[0] ) ;
[1977]1116
[2016]1117 // SPECTRAL_WINDOW_ID
[2018]1118 *spwRF = ifnoCol.asuInt( rows[0] ) ;
[1977]1119
[2016]1120 // TIME and INTERVAL
[2018]1121 getValidTimeRange( midTime, interval, atime, ainterval ) ;
[2016]1122 *timeRF = midTime ;
1123 *intervalRF = interval ;
1124
1125 // TCAL and TSYS
1126 Matrix<Float> tcal ;
1127 Table t ;
1128 if ( idCol.asuInt( ids[0] ) == ids[0] ) {
1129 os_ << "sorted at irow=" << irow << " ids[0]=" << ids[0] << LogIO::POST ;
1130 Vector<Float> dummyC = tcalCol( ids[0] ) ;
1131 tcal.resize( npol, dummyC.size() ) ;
1132 tcal.row( 0 ) = dummyC ;
[1977]1133 }
[2016]1134 else {
1135 os_ << "NOT sorted at irow=" << irow << " ids[0]=" << ids[0] << LogIO::POST ;
1136 t = stt( stt.col("ID") == ids[0] ) ;
1137 Vector<Float> dummyC = tcalCol( 0 ) ;
1138 tcal.resize( npol, dummyC.size() ) ;
1139 tcal.row( 0 ) = dummyC ;
1140 }
1141 if ( npol == 2 ) {
1142 if ( idCol.asuInt( ids[1] ) == ids[1] ) {
1143 os_ << "sorted at irow=" << irow << " ids[1]=" << ids[1] << LogIO::POST ;
1144 tcal.row( 1 ) = tcalCol( ids[1] ) ;
1145 }
1146 else {
1147 os_ << "NOT sorted at irow=" << irow << " ids[1]=" << ids[1] << LogIO::POST ;
1148 t = stt( stt.col("ID") == ids[1] ) ;
1149 tcalCol.attach( t, "TCAL" ) ;
1150 tcal.row( 1 ) = tcalCol( 1 ) ;
1151 }
1152 }
1153 else if ( npol == 3 ) {
1154 if ( idCol.asuInt( ids[2] ) == ids[2] )
1155 tcal.row( 1 ) = tcalCol( ids[2] ) ;
1156 else {
1157 t = stt( stt.col("ID") == ids[2] ) ;
1158 tcalCol.attach( t, "TCAL" ) ;
1159 tcal.row( 1 ) = tcalCol( 0 ) ;
1160 }
1161 if ( idCol.asuInt( ids[1] ) == ids[1] )
1162 tcal.row( 2 ) = tcalCol( ids[1] ) ;
1163 else {
1164 t = stt( stt.col("ID") == ids[1] ) ;
1165 tcalCol.attach( t, "TCAL" ) ;
1166 tcal.row( 2 ) = tcalCol( 0 ) ;
1167 }
1168 }
1169 else if ( npol == 4 ) {
1170 if ( idCol.asuInt( ids[2] ) == ids[2] )
1171 tcal.row( 1 ) = tcalCol( ids[2] ) ;
1172 else {
1173 t = stt( stt.col("ID") == ids[2] ) ;
1174 tcalCol.attach( t, "TCAL" ) ;
1175 tcal.row( 1 ) = tcalCol( 0 ) ;
1176 }
1177 if ( idCol.asuInt( ids[3] ) == ids[3] )
1178 tcal.row( 2 ) = tcalCol( ids[3] ) ;
1179 else {
1180 t = stt( stt.col("ID") == ids[3] ) ;
1181 tcalCol.attach( t, "TCAL" ) ;
1182 tcal.row( 2 ) = tcalCol( 0 ) ;
1183 }
1184 if ( idCol.asuInt( ids[1] ) == ids[1] )
1185 tcal.row( 2 ) = tcalCol( ids[1] ) ;
1186 else {
1187 t = stt( stt.col("ID") == ids[1] ) ;
1188 tcalCol.attach( t, "TCAL" ) ;
1189 tcal.row( 3 ) = tcalCol( 0 ) ;
1190 }
1191 }
1192 if ( tcalSpec_ ) {
1193 // put TCAL_SPECTRUM
[2019]1194 //*tcalspRF = tcal ;
1195 tcalspRF.define( tcal ) ;
[2016]1196 // set TCAL (mean of TCAL_SPECTRUM)
1197 Matrix<Float> tcalMean( npol, 1 ) ;
1198 for ( uInt iid = 0 ; iid < npol ; iid++ ) {
1199 tcalMean( iid, 0 ) = mean( tcal.row(iid) ) ;
1200 }
1201 // put TCAL
1202 *tcalRF = tcalMean ;
1203 }
1204 else {
1205 // put TCAL
1206 *tcalRF = tcal ;
1207 }
[1977]1208
[2016]1209 if ( tsysSpec_ ) {
1210 // put TSYS_SPECTRUM
[2019]1211 //*tsysspRF = tsys ;
1212 tsysspRF.define( tsys ) ;
[2016]1213 // set TSYS (mean of TSYS_SPECTRUM)
1214 Matrix<Float> tsysMean( npol, 1 ) ;
1215 for ( uInt iid = 0 ; iid < npol ; iid++ ) {
1216 tsysMean( iid, 0 ) = mean( tsys.row(iid) ) ;
1217 }
1218 // put TSYS
1219 *tsysRF = tsysMean ;
1220 }
1221 else {
1222 // put TSYS
1223 *tsysRF = tsys ;
1224 }
[1977]1225
[2016]1226 // add row
1227 mssc.addRow( 1, True ) ;
1228 row.put( mssc.nrow()-1 ) ;
1229
1230 double t2 = gettimeofday_sec() ;
1231 os_ << irow << "th loop elapsed time = " << t2-t1 << "sec" << LogIO::POST ;
[1977]1232 }
1233
[2016]1234 double endSec = gettimeofday_sec() ;
1235 os_ << "end MSWriter::fillSysCal() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1977]1236}
1237
[1974]1238void MSWriter::addFeed( Int id )
1239{
[2016]1240 double startSec = gettimeofday_sec() ;
1241 os_ << "start MSWriter::addFeed() startSec=" << startSec << LogIO::POST ;
[1974]1242
1243 // add row
1244 MSFeed msFeed = mstable_->feed() ;
1245 msFeed.addRow( 1, True ) ;
1246 Int nrow = msFeed.nrow() ;
[2020]1247 Int numReceptors = 2 ;
1248 Vector<String> polType( numReceptors ) ;
1249 Matrix<Double> beamOffset( 2, numReceptors ) ;
1250 beamOffset = 0.0 ;
1251 Vector<Double> receptorAngle( numReceptors, 0.0 ) ;
1252 if ( polType_ == "linear" ) {
1253 polType[0] = "X" ;
1254 polType[1] = "Y" ;
1255 }
1256 else if ( polType_ == "circular" ) {
1257 polType[0] = "R" ;
1258 polType[1] = "L" ;
1259 }
1260 else {
1261 polType[0] = "X" ;
1262 polType[1] = "Y" ;
1263 }
1264 Matrix<Complex> polResponse( numReceptors, numReceptors, 0.0 ) ;
1265 for ( Int i = 0 ; i < numReceptors ; i++ )
1266 polResponse( i, i ) = 0.0 ;
[1974]1267
1268 MSFeedColumns msFeedCols( mstable_->feed() ) ;
1269
1270 msFeedCols.feedId().put( nrow-1, id ) ;
1271 msFeedCols.antennaId().put( nrow-1, 0 ) ;
[2020]1272 msFeedCols.numReceptors().put( nrow-1, numReceptors ) ;
1273 msFeedCols.polarizationType().put( nrow-1, polType ) ;
1274 msFeedCols.beamOffset().put( nrow-1, beamOffset ) ;
1275 msFeedCols.receptorAngle().put( nrow-1, receptorAngle ) ;
1276 msFeedCols.polResponse().put( nrow-1, polResponse ) ;
[2016]1277
1278 double endSec = gettimeofday_sec() ;
1279 os_ << "end MSWriter::addFeed() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1974]1280}
1281
1282void MSWriter::addSpectralWindow( Int spwid, Int freqid )
1283{
[2016]1284 double startSec = gettimeofday_sec() ;
1285 os_ << "start MSWriter::addSpectralWindow() startSec=" << startSec << LogIO::POST ;
[1974]1286
1287 // add row
1288 MSSpectralWindow msSpw = mstable_->spectralWindow() ;
1289 while( (Int)msSpw.nrow() <= spwid ) {
1290 msSpw.addRow( 1, True ) ;
1291 }
1292
1293 MSSpWindowColumns msSpwCols( msSpw ) ;
1294
1295 STFrequencies stf = table_->frequencies() ;
1296
1297 // MEAS_FREQ_REF
1298 msSpwCols.measFreqRef().put( spwid, stf.getFrame( True ) ) ;
1299
1300 Double refpix ;
1301 Double refval ;
1302 Double inc ;
1303 stf.getEntry( refpix, refval, inc, (uInt)freqid ) ;
1304
1305 // NUM_CHAN
[2016]1306 Int nchan = (Int)(refpix * 2) ;
[1974]1307 msSpwCols.numChan().put( spwid, nchan ) ;
1308
1309 // TOTAL_BANDWIDTH
1310 Double bw = nchan * inc ;
1311 msSpwCols.totalBandwidth().put( spwid, bw ) ;
1312
1313 // REF_FREQUENCY
1314 Double refFreq = refval - refpix * inc ;
1315 msSpwCols.refFrequency().put( spwid, refFreq ) ;
1316
1317 // NET_SIDEBAND
1318 // tentative: USB->0, LSB->1
1319 Int netSideband = 0 ;
1320 if ( inc < 0 )
1321 netSideband = 1 ;
1322 msSpwCols.netSideband().put( spwid, netSideband ) ;
1323
1324 // RESOLUTION, CHAN_WIDTH, EFFECTIVE_BW
1325 Vector<Double> sharedDoubleArr( nchan, inc ) ;
1326 msSpwCols.resolution().put( spwid, sharedDoubleArr ) ;
1327 msSpwCols.chanWidth().put( spwid, sharedDoubleArr ) ;
1328 msSpwCols.effectiveBW().put( spwid, sharedDoubleArr ) ;
1329
1330 // CHAN_FREQ
1331 indgen( sharedDoubleArr, refFreq, inc ) ;
1332 msSpwCols.chanFreq().put( spwid, sharedDoubleArr ) ;
[2016]1333
1334 double endSec = gettimeofday_sec() ;
1335 os_ << "end MSWriter::addSpectralWindow() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1974]1336}
1337
[1975]1338void MSWriter::addField( Int fid, String fieldname, String srcname, Double t, Vector<Double> rate )
1339{
[2016]1340 double startSec = gettimeofday_sec() ;
1341 os_ << "start MSWriter::addField() startSec=" << startSec << LogIO::POST ;
[1975]1342
1343 MSField msField = mstable_->field() ;
1344 while( (Int)msField.nrow() <= fid ) {
1345 msField.addRow( 1, True ) ;
1346 }
1347 MSFieldColumns msFieldCols( msField ) ;
1348
1349 // Access to SOURCE table
1350 MSSource msSrc = mstable_->source() ;
1351
1352 // fill target row
1353 msFieldCols.name().put( fid, fieldname ) ;
1354 msFieldCols.time().put( fid, t ) ;
1355 Int numPoly = 0 ;
1356 if ( anyNE( rate, 0.0 ) )
1357 numPoly = 1 ;
1358 msFieldCols.numPoly().put( fid, numPoly ) ;
1359 MSSourceIndex msSrcIdx( msSrc ) ;
1360 Int srcId = -1 ;
1361 Vector<Int> srcIdArr = msSrcIdx.matchSourceName( srcname ) ;
1362 if ( srcIdArr.size() != 0 ) {
1363 srcId = srcIdArr[0] ;
1364 MSSource msSrcSel = msSrc( msSrc.col("SOURCE_ID") == srcId ) ;
1365 ROMSSourceColumns msSrcCols( msSrcSel ) ;
1366 Vector<Double> srcDir = msSrcCols.direction()( 0 ) ;
[2016]1367 Matrix<Double> srcDirA( IPosition( 2, 2, 1+numPoly ) ) ;
[1975]1368 os_ << "srcDirA = " << srcDirA << LogIO::POST ;
[2016]1369 os_ << "sliced srcDirA = " << srcDirA.column( 0 ) << LogIO::POST ;
1370 srcDirA.column( 0 ) = srcDir ;
[1975]1371 os_ << "srcDirA = " << srcDirA << LogIO::POST ;
1372 if ( numPoly != 0 )
[2016]1373 srcDirA.column( 1 ) = rate ;
[1975]1374 msFieldCols.phaseDir().put( fid, srcDirA ) ;
1375 msFieldCols.referenceDir().put( fid, srcDirA ) ;
1376 msFieldCols.delayDir().put( fid, srcDirA ) ;
1377 }
1378 msFieldCols.sourceId().put( fid, srcId ) ;
[2016]1379
1380 double endSec = gettimeofday_sec() ;
1381 os_ << "end MSWriter::addField() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1975]1382}
1383
[2016]1384void MSWriter::addPointing( String &name, Double &me, Double &interval, Matrix<Double> &dir )
[1977]1385{
[2016]1386 double startSec = gettimeofday_sec() ;
1387 os_ << "start MSWriter::addPointing() startSec=" << startSec << LogIO::POST ;
[1977]1388
1389 // access to POINTING subtable
1390 MSPointing msp = mstable_->pointing() ;
1391 uInt nrow = msp.nrow() ;
1392
1393 // add row
1394 msp.addRow( 1, True ) ;
1395
1396 // fill row
[2016]1397 TableRow row( msp ) ;
1398 TableRecord &rec = row.record() ;
1399 RecordFieldPtr<Int> antennaRF( rec, "ANTENNA_ID" ) ;
1400 *antennaRF = 0 ;
1401 RecordFieldPtr<Int> numpolyRF( rec, "NUM_POLY" ) ;
1402 *numpolyRF = dir.ncolumn() ;
1403 RecordFieldPtr<Double> timeRF( rec, "TIME" ) ;
1404 *timeRF = me ;
1405 RecordFieldPtr<Double> toriginRF( rec, "TIME_ORIGIN" ) ;
1406 *toriginRF = me ;
1407 RecordFieldPtr<Double> intervalRF( rec, "INTERVAL" ) ;
1408 *intervalRF = interval ;
1409 RecordFieldPtr<String> nameRF( rec, "NAME" ) ;
1410 *nameRF = name ;
1411 RecordFieldPtr<Bool> trackRF( rec, "TRACKING" ) ;
1412 *trackRF = True ;
1413 RecordFieldPtr< Array<Double> > dirRF( rec, "DIRECTION" ) ;
1414 *dirRF = dir ;
1415 RecordFieldPtr< Array<Double> > targetRF( rec, "TARGET" ) ;
1416 *dirRF = dir ;
1417 row.put( nrow ) ;
1418
1419 double endSec = gettimeofday_sec() ;
1420 os_ << "end MSWriter::addPointing() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1977]1421}
1422
[1974]1423Int MSWriter::addPolarization( Vector<Int> polnos )
1424{
[2016]1425 double startSec = gettimeofday_sec() ;
1426 os_ << "start MSWriter::addPolarization() startSec=" << startSec << LogIO::POST ;
[1974]1427
1428 os_ << "polnos = " << polnos << LogIO::POST ;
1429 MSPolarization msPol = mstable_->polarization() ;
1430 uInt nrow = msPol.nrow() ;
[2016]1431
1432 // only 1 POLARIZATION row for 1 scantable
1433 if ( nrow > 0 )
[2019]1434 return 0 ;
[1974]1435
1436 Vector<Int> corrType = toCorrType( polnos ) ;
1437
[2016]1438 ROArrayColumn<Int> corrtCol( msPol, "CORR_TYPE" ) ;
1439 Matrix<Int> corrTypeArr = corrtCol.getColumn() ;
[1974]1440 Int polid = -1 ;
[2016]1441 for ( uInt irow = 0 ; irow < nrow ; irow++ ) {
1442 if ( allEQ( corrType, corrTypeArr.column( irow ) ) ) {
1443 polid = irow ;
1444 break ;
1445 }
1446 }
[1974]1447
[2016]1448 if ( polid == -1 ) {
1449 MSPolarizationColumns msPolCols( msPol ) ;
1450
[1974]1451 // add row
1452 msPol.addRow( 1, True ) ;
1453 polid = (Int)nrow ;
1454
1455 // CORR_TYPE
1456 msPolCols.corrType().put( nrow, corrType ) ;
1457
1458 // NUM_CORR
1459 uInt npol = corrType.size() ;
1460 msPolCols.numCorr().put( nrow, npol ) ;
1461
1462 // CORR_PRODUCT
1463 Matrix<Int> corrProd( 2, npol, -1 ) ;
1464 if ( npol == 1 ) {
1465 corrProd = 0 ;
1466 }
1467 else if ( npol == 2 ) {
[1975]1468 corrProd.column( 0 ) = 0 ;
1469 corrProd.column( 1 ) = 1 ;
[1974]1470 }
1471 else {
[1975]1472 corrProd.column( 0 ) = 0 ;
1473 corrProd.column( 3 ) = 1 ;
1474 corrProd( 0,1 ) = 0 ;
1475 corrProd( 1,1 ) = 1 ;
1476 corrProd( 0,2 ) = 1 ;
1477 corrProd( 1,2 ) = 0 ;
[1974]1478 }
[2016]1479 msPolCols.corrProduct().put( nrow, corrProd ) ;
[1974]1480 }
1481
[2016]1482 double endSec = gettimeofday_sec() ;
1483 os_ << "end MSWriter::addPolarization() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
1484
[1974]1485 return polid ;
1486}
1487
[1975]1488Int MSWriter::addDataDescription( Int polid, Int spwid )
1489{
[2016]1490 double startSec = gettimeofday_sec() ;
1491 os_ << "start MSWriter::addDataDescription() startSec=" << startSec << LogIO::POST ;
[1975]1492
1493 MSDataDescription msDataDesc = mstable_->dataDescription() ;
1494 uInt nrow = msDataDesc.nrow() ;
1495
[2016]1496 // only 1 POLARIZATION_ID for 1 scantable
1497 Int ddid = -1 ;
1498 ROScalarColumn<Int> spwCol( msDataDesc, "SPECTRAL_WINDOW_ID" ) ;
1499 Vector<Int> spwIds = spwCol.getColumn() ;
1500 //ROScalarColumn<Int> polCol( msDataDesc, "POLARIZATION_ID" ) ;
1501 //Vector<Int> polIds = polCol.getColumn() ;
1502 for ( uInt irow = 0 ; irow < nrow ; irow++ ) {
1503 //if ( spwid == spwIds[irow] && polid == polIds[irow] ) {
1504 if ( spwid == spwIds[irow] ) {
1505 ddid = irow ;
1506 break ;
1507 }
1508 }
1509 os_ << "ddid = " << ddid << LogIO::POST ;
1510
[1975]1511
[2016]1512 if ( ddid == -1 ) {
[1975]1513 msDataDesc.addRow( 1, True ) ;
1514 MSDataDescColumns msDataDescCols( msDataDesc ) ;
1515 msDataDescCols.polarizationId().put( nrow, polid ) ;
1516 msDataDescCols.spectralWindowId().put( nrow, spwid ) ;
1517 ddid = (Int)nrow ;
1518 }
1519
[2016]1520 double endSec = gettimeofday_sec() ;
1521 os_ << "end MSWriter::addDataDescription() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
1522
[1975]1523 return ddid ;
1524}
1525
[1977]1526Int MSWriter::addState( Int st, Int &subscan )
1527{
[2016]1528 double startSec = gettimeofday_sec() ;
1529 os_ << "start MSWriter::addState() startSec=" << startSec << LogIO::POST ;
[1977]1530
1531 // access to STATE subtable
1532 MSState msState = mstable_->state() ;
1533 uInt nrow = msState.nrow() ;
1534
1535 String obsMode ;
1536 Bool isSignal ;
[2016]1537 Double tnoise ;
1538 Double tload ;
1539 queryType( st, obsMode, isSignal, tnoise, tload ) ;
[1977]1540 os_ << "obsMode = " << obsMode << " isSignal = " << isSignal << LogIO::POST ;
1541
1542 Int idx = -1 ;
[2016]1543 ROScalarColumn<String> obsModeCol( msState, "OBS_MODE" ) ;
1544 //ROScalarColumn<Bool> sigCol( msState, "SIG" ) ;
1545 //ROScalarColumn<Bool> refCol( msState, "REF" ) ;
1546 ROScalarColumn<Int> subscanCol( msState, "SUB_SCAN" ) ;
1547// Vector<String> obsModeArr = obsModeCol.getColumn() ;
1548// Vector<Bool> sigArr = sigCol.getColumn() ;
1549// Vector<Bool> refArr = refCol.getColumn() ;
1550// Vector<Int> subscanArr = subscanCol.getColumn() ;
1551 for ( uInt irow = 0 ; irow < nrow ; irow++ ) {
1552 if ( obsModeCol(irow) == obsMode
1553 //&& sigCol(irow) == isSignal
1554 //&& refCol(irow) != isSignal
1555 && subscanCol(irow) == subscan ) {
1556// if ( obsModeArr[irow] == obsMode
1557// && sigArr[irow] == isSignal
1558// && refArr[irow] != isSignal
1559// && subscanArr[irow] == subscan ) {
1560 idx = irow ;
1561 break ;
1562 }
1563 }
1564 if ( idx == -1 ) {
[1977]1565 msState.addRow( 1, True ) ;
[2016]1566 TableRow row( msState ) ;
1567 TableRecord &rec = row.record() ;
1568 RecordFieldPtr<String> obsmodeRF( rec, "OBS_MODE" ) ;
1569 *obsmodeRF = obsMode ;
1570 RecordFieldPtr<Bool> sigRF( rec, "SIG" ) ;
1571 *sigRF = isSignal ;
1572 RecordFieldPtr<Bool> refRF( rec, "REF" ) ;
1573 *refRF = !isSignal ;
1574 RecordFieldPtr<Int> subscanRF( rec, "SUB_SCAN" ) ;
1575 *subscanRF = subscan ;
1576 RecordFieldPtr<Double> noiseRF( rec, "CAL" ) ;
1577 *noiseRF = tnoise ;
1578 RecordFieldPtr<Double> loadRF( rec, "LOAD" ) ;
1579 *loadRF = tload ;
1580 row.put( nrow ) ;
1581// ScalarColumn<String> obsModeCol( msState, "OBS_MODE" ) ;
1582// obsModeCol.put( nrow, obsMode ) ;
1583// ScalarColumn<Bool> sharedBCol( msState, "SIG" ) ;
1584// sharedBCol.put( nrow, isSignal ) ;
1585// sharedBCol.attach( msState, "REF" ) ;
1586// sharedBCol.put( nrow, !isSignal ) ;
1587// ScalarColumn<Int> subscanCol( msState, "SUB_SCAN" ) ;
1588// subscanCol.put( nrow, subscan ) ;
[1977]1589 idx = nrow ;
1590 }
1591 subscan++ ;
[2016]1592
1593 double endSec = gettimeofday_sec() ;
1594 os_ << "end MSWriter::addState() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
1595
[1977]1596 return idx ;
1597}
1598
[1974]1599Vector<Int> MSWriter::toCorrType( Vector<Int> polnos )
1600{
[2016]1601 double startSec = gettimeofday_sec() ;
1602 os_ << "start MSWriter::toCorrType() startSec=" << startSec << LogIO::POST ;
1603
[1974]1604 uInt npol = polnos.size() ;
1605 Vector<Int> corrType( npol, Stokes::Undefined ) ;
1606
1607 if ( npol == 4 ) {
1608 if ( polType_ == "linear" ) {
1609 for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) {
1610 if ( polnos[ipol] == 0 )
1611 corrType[ipol] = Stokes::XX ;
1612 else if ( polnos[ipol] == 1 )
1613 corrType[ipol] = Stokes::XY ;
1614 else if ( polnos[ipol] == 2 )
1615 corrType[ipol] = Stokes::YX ;
1616 else if ( polnos[ipol] == 3 )
1617 corrType[ipol] = Stokes::YY ;
1618 }
1619 }
1620 else if ( polType_ == "circular" ) {
1621 for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) {
1622 if ( polnos[ipol] == 0 )
1623 corrType[ipol] = Stokes::RR ;
1624 else if ( polnos[ipol] == 1 )
1625 corrType[ipol] = Stokes::RL ;
1626 else if ( polnos[ipol] == 2 )
1627 corrType[ipol] = Stokes::LR ;
1628 else if ( polnos[ipol] == 3 )
1629 corrType[ipol] = Stokes::LL ;
1630 }
1631 }
1632 else if ( polType_ == "stokes" ) {
1633 for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) {
1634 if ( polnos[ipol] == 0 )
1635 corrType[ipol] = Stokes::I ;
1636 else if ( polnos[ipol] == 1 )
1637 corrType[ipol] = Stokes::Q ;
1638 else if ( polnos[ipol] == 2 )
1639 corrType[ipol] = Stokes::U ;
1640 else if ( polnos[ipol] == 3 )
1641 corrType[ipol] = Stokes::V ;
1642 }
1643 }
1644 }
1645 else if ( npol == 2 ) {
1646 if ( polType_ == "linear" ) {
1647 for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) {
1648 if ( polnos[ipol] == 0 )
1649 corrType[ipol] = Stokes::XX ;
1650 else if ( polnos[ipol] == 1 )
1651 corrType[ipol] = Stokes::YY ;
1652 }
1653 }
1654 else if ( polType_ == "circular" ) {
1655 for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) {
1656 if ( polnos[ipol] == 0 )
1657 corrType[ipol] = Stokes::RR ;
1658 else if ( polnos[ipol] == 1 )
1659 corrType[ipol] = Stokes::LL ;
1660 }
1661 }
1662 else if ( polType_ == "stokes" ) {
1663 for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) {
1664 if ( polnos[ipol] == 0 )
1665 corrType[ipol] = Stokes::I ;
1666 else if ( polnos[ipol] == 1 )
1667 corrType[ipol] = Stokes::V ;
1668 }
1669 }
1670 else if ( polType_ == "linpol" ) {
1671 for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) {
1672 if ( polnos[ipol] == 1 )
1673 corrType[ipol] = Stokes::Plinear ;
1674 else if ( polnos[ipol] == 2 )
1675 corrType[ipol] = Stokes::Pangle ;
1676 }
1677 }
1678 }
1679 else if ( npol == 1 ) {
1680 if ( polType_ == "linear" )
1681 corrType[0] = Stokes::XX ;
1682 else if ( polType_ == "circular" )
1683 corrType[0] = Stokes::RR ;
1684 else if ( polType_ == "stokes" )
1685 corrType[0] = Stokes::I ;
1686 }
1687
[2016]1688 double endSec = gettimeofday_sec() ;
1689 os_ << "end MSWriter::toCorrType() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
1690
[1974]1691 return corrType ;
1692}
1693
[2016]1694void MSWriter::getValidTimeRange( Double &me, Double &interval, Table &tab )
[1974]1695{
[2016]1696 double startSec = gettimeofday_sec() ;
1697 os_ << "start MSWriter::getVaridTimeRange() startSec=" << startSec << LogIO::POST ;
1698
1699 // sort table
1700 //Table stab = tab.sort( "TIME" ) ;
1701
[1975]1702 ROScalarColumn<Double> timeCol( tab, "TIME" ) ;
1703 Vector<Double> timeArr = timeCol.getColumn() ;
1704 Double minTime ;
1705 Double maxTime ;
1706 minMax( minTime, maxTime, timeArr ) ;
[2016]1707 Double midTime = 0.5 * ( minTime + maxTime ) * 86400.0 ;
1708 // unit for TIME
1709 // Scantable: "d"
1710 // MS: "s"
1711 me = midTime ;
1712 interval = ( maxTime - minTime ) * 86400.0 ;
1713
1714 double endSec = gettimeofday_sec() ;
1715 os_ << "end MSWriter::getValidTimeRange() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1975]1716}
[1974]1717
[2018]1718void MSWriter::getValidTimeRange( Double &me, Double &interval, Vector<Double> &atime, Vector<Double> &ainterval )
1719{
1720 double startSec = gettimeofday_sec() ;
1721 os_ << "start MSWriter::getVaridTimeRange() startSec=" << startSec << LogIO::POST ;
1722
1723 // sort table
1724 //Table stab = tab.sort( "TIME" ) ;
1725
1726 Double minTime ;
1727 Double maxTime ;
1728 minMax( minTime, maxTime, atime ) ;
1729 Double midTime = 0.5 * ( minTime + maxTime ) * 86400.0 ;
1730 // unit for TIME
1731 // Scantable: "d"
1732 // MS: "s"
1733 me = midTime ;
[2020]1734 interval = ( maxTime - minTime ) * 86400.0 + mean( ainterval ) ;
[2018]1735
1736 double endSec = gettimeofday_sec() ;
1737 os_ << "end MSWriter::getValidTimeRange() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
1738}
1739
[2016]1740//void MSWriter::queryType( Int type, String &stype, Bool &b )
1741void MSWriter::queryType( Int type, String &stype, Bool &b, Double &t, Double &l )
[1977]1742{
[2016]1743 double startSec = gettimeofday_sec() ;
1744 os_ << "start MSWriter::queryType() startSec=" << startSec << LogIO::POST ;
1745
[1977]1746 switch ( type ) {
1747 case SrcType::PSON:
1748 stype = "POSITION_SWITCH.OBSERVE_TARGET.ON_SOURCE" ;
1749 b = True ;
[2016]1750 t = 0.0 ;
1751 l = 0.0 ;
[1977]1752 break ;
1753 case SrcType::PSOFF:
1754 stype = "POSITION_SWITCH.OBSERVE_TARGET.OFF_SOURCE" ;
1755 b = False ;
[2016]1756 t = 0.0 ;
1757 l = 0.0 ;
[1977]1758 break ;
1759 case SrcType::NOD:
1760 stype = "NOD.OBSERVE_TARGET.ON_SOURCE" ;
1761 b = True ;
[2016]1762 t = 0.0 ;
1763 l = 0.0 ;
[1977]1764 break ;
1765 case SrcType::FSON:
1766 stype = "FREQUENCY_SWITCH.OBSERVE_TARGET.ON_SOURCE" ;
1767 b = True ;
[2016]1768 t = 0.0 ;
1769 l = 0.0 ;
[1977]1770 break ;
1771 case SrcType::FSOFF:
1772 stype = "FREQUENCY_SWITCH.OBSERVE_TARGET.ON_SOURCE" ;
1773 b = False ;
[2016]1774 t = 0.0 ;
1775 l = 0.0 ;
[1977]1776 break ;
1777 case SrcType::SKY:
1778 stype = "UNSPECIFIED.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1779 b = False ;
[2016]1780 t = 0.0 ;
1781 l = 1.0 ;
[1977]1782 break ;
1783 case SrcType::HOT:
1784 stype = "UNSPECIFIED.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1785 b = False ;
[2016]1786 t = 0.0 ;
1787 l = 1.0 ;
[1977]1788 break ;
1789 case SrcType::WARM:
1790 stype = "UNSPECIFIED.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
[2016]1791 t = 0.0 ;
[1977]1792 b = False ;
[2016]1793 l = 1.0 ;
1794 break ;
[1977]1795 case SrcType::COLD:
1796 stype = "UNSPECIFIED.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1797 b = False ;
[2016]1798 t = 0.0 ;
1799 l = 1.0 ;
[1977]1800 break ;
1801 case SrcType::PONCAL:
1802 stype = "POSITION_SWITCH.CALIBRATE_TEMPERATURE.ON_SOURCE" ;
1803 b = True ;
[2016]1804 t = 1.0 ;
1805 l = 0.0 ;
[1977]1806 break ;
1807 case SrcType::POFFCAL:
1808 stype = "POSITION_SWITCH.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1809 b = False ;
[2016]1810 t = 1.0 ;
1811 l = 0.0 ;
[1977]1812 break ;
1813 case SrcType::NODCAL:
1814 stype = "NOD.CALIBRATE_TEMPERATURE.ON_SOURCE" ;
1815 b = True ;
[2016]1816 t = 1.0 ;
1817 l = 0.0 ;
[1977]1818 break ;
1819 case SrcType::FONCAL:
1820 stype = "FREQUENCY_SWITCH.CALIBRATE_TEMPERATURE.ON_SOURCE" ;
1821 b = True ;
[2016]1822 t = 1.0 ;
1823 l = 0.0 ;
[1977]1824 break ;
1825 case SrcType::FOFFCAL:
1826 stype = "FREQUENCY_SWITCH.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1827 b = False ;
[2016]1828 t = 1.0 ;
1829 l = 0.0 ;
[1977]1830 break ;
1831 case SrcType::FSLO:
1832 stype = "FREQUENCY_SWITCH.OBSERVE_TARGET.ON_SOURCE" ;
1833 b = True ;
[2016]1834 t = 0.0 ;
1835 l = 0.0 ;
[1977]1836 break ;
1837 case SrcType::FLOOFF:
1838 stype = "FREQUENCY_SWITCH.OBSERVE_TARGET.OFF_SOURCE" ;
1839 b = False ;
[2016]1840 t = 0.0 ;
1841 l = 0.0 ;
[1977]1842 break ;
1843 case SrcType::FLOSKY:
1844 stype = "FREQUENCY_SWITCH.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1845 b = False ;
[2016]1846 t = 0.0 ;
1847 l = 1.0 ;
[1977]1848 break ;
1849 case SrcType::FLOHOT:
1850 stype = "FREQUENCY_SWITCH.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1851 b = False ;
[2016]1852 t = 0.0 ;
1853 l = 1.0 ;
[1977]1854 break ;
1855 case SrcType::FLOWARM:
1856 stype = "FREQUENCY_SWITCH.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1857 b = False ;
[2016]1858 t = 0.0 ;
1859 l = 1.0 ;
[1977]1860 break ;
1861 case SrcType::FLOCOLD:
1862 stype = "FREQUENCY_SWITCH.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1863 b = False ;
[2016]1864 t = 0.0 ;
1865 l = 1.0 ;
[1977]1866 break ;
1867 case SrcType::FSHI:
1868 stype = "FREQUENCY_SWITCH.OBSERVE_TARGET.ON_SOURCE" ;
1869 b = True ;
[2016]1870 t = 0.0 ;
1871 l = 0.0 ;
[1977]1872 break ;
1873 case SrcType::FHIOFF:
1874 stype = "FREQUENCY_SWITCH.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1875 b = False ;
[2016]1876 t = 0.0 ;
1877 l = 0.0 ;
[1977]1878 break ;
1879 case SrcType::FHISKY:
1880 stype = "FREQUENCY_SWITCH.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1881 b = False ;
[2016]1882 t = 0.0 ;
1883 l = 1.0 ;
[1977]1884 break ;
1885 case SrcType::FHIHOT:
1886 stype = "FREQUENCY_SWITCH.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1887 b = False ;
[2016]1888 t = 0.0 ;
1889 l = 1.0 ;
[1977]1890 break ;
1891 case SrcType::FHIWARM:
1892 stype = "FREQUENCY_SWITCH.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1893 b = False ;
[2016]1894 t = 0.0 ;
1895 l = 1.0 ;
[1977]1896 break ;
1897 case SrcType::FHICOLD:
1898 stype = "FREQUENCY_SWITCH.CALIBRATE_TEMPERATURE.OFF_SOURCE" ;
1899 b = False ;
[2016]1900 t = 0.0 ;
1901 l = 1.0 ;
[1977]1902 break ;
1903 case SrcType::SIG:
1904 stype = "UNSPECIFIED.OBSERVE_TARGET.ON_SOURCE" ;
1905 b = True ;
[2016]1906 t = 0.0 ;
1907 l = 0.0 ;
[1977]1908 break ;
1909 case SrcType::REF:
1910 stype = "UNSPECIFIED.OBSERVE_TARGET.ON_SOURCE" ;
1911 b = False ;
[2016]1912 t = 0.0 ;
1913 l = 0.0 ;
[1977]1914 break ;
1915 default:
1916 stype = "UNSPECIFIED" ;
1917 b = True ;
[2016]1918 t = 0.0 ;
1919 l = 0.0 ;
[1977]1920 break ;
1921 }
[2016]1922
1923 double endSec = gettimeofday_sec() ;
1924 os_ << "end MSWriter::queryType() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ;
[1974]1925}
[1977]1926}
Note: See TracBrowser for help on using the repository browser.