source: trunk/src/MSWriter.cpp@ 2016

Last change on this file since 2016 was 2016, 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...

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