- Timestamp:
- 02/18/11 19:26:01 (14 years ago)
- Location:
- trunk/src
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/src/MSFiller.cpp
r1993 r2002 23 23 #include <tables/Tables/TableParse.h> 24 24 #include <tables/Tables/RefRows.h> 25 26 #include <casa/Containers/Block.h> 25 #include <tables/Tables/TableRow.h> 26 27 #include <casa/Containers/RecordField.h> 27 28 #include <casa/Logging/LogIO.h> 28 29 #include <casa/Arrays/Slicer.h> … … 197 198 colTsys_ = "TSYS" ; 198 199 } 200 //colTcal_ = "TCAL" ; 201 //colTsys_ = "TSYS" ; 199 202 MSPointing pointtab = mstable_.pointing() ; 200 203 if ( mstable_.weather().nrow() == 0 ) … … 206 209 // memory allocation by boost::object_pool 207 210 boost::object_pool<ROTableColumn> *tpoolr = new boost::object_pool<ROTableColumn> ; 208 boost::object_pool<TableColumn> *tpoolw = new boost::object_pool<TableColumn> ;209 211 // 210 212 … … 222 224 // SUBTABLES: TCAL 223 225 if ( isSysCal_ ) 224 fillTcal( tpoolr , tpoolw) ;226 fillTcal( tpoolr ) ; 225 227 226 228 // SUBTABLES: FIT … … 232 234 // shared pointers 233 235 ROTableColumn *tcolr ; 234 TableColumn *tcolw ; 235 236 // Scantable columns 237 Table stab = table_->table() ; 238 TableColumn *scannoCol = tpoolw->construct( stab, "SCANNO" ) ; 239 TableColumn *cyclenoCol = tpoolw->construct( stab, "CYCLENO" ) ; 240 TableColumn *beamnoCol = tpoolw->construct( stab, "BEAMNO" ) ; 241 TableColumn *ifnoCol = tpoolw->construct( stab, "IFNO" ) ; 242 TableColumn *polnoCol = tpoolw->construct( stab, "POLNO" ) ; 243 TableColumn *freqidCol = tpoolw->construct( stab, "FREQ_ID" ) ; 244 TableColumn *molidCol = tpoolw->construct( stab, "MOLECULE_ID" ) ; 245 TableColumn *flagrowCol = tpoolw->construct( stab, "FLAGROW" ) ; 246 ScalarMeasColumn<MEpoch> *timeCol = new ScalarMeasColumn<MEpoch>( stab, "TIME" ) ; 247 TableColumn *intervalCol = tpoolw->construct( stab, "INTERVAL" ) ; 248 TableColumn *srcnameCol = tpoolw->construct( stab, "SRCNAME" ) ; 249 TableColumn *srctypeCol = tpoolw->construct( stab, "SRCTYPE" ) ; 250 TableColumn *fieldnameCol = tpoolw->construct( stab, "SRCNAME" ) ; 251 ArrayColumn<Float> *spCol = new ArrayColumn<Float>( stab, "SPECTRA" ) ; 252 ArrayColumn<uChar> *flCol = new ArrayColumn<uChar>( stab, "FLAGTRA" ) ; 253 ArrayColumn<Float> *tsysCol = new ArrayColumn<Float>( stab, "TSYS" ) ; 254 ArrayColumn<Double> *dirCol = new ArrayColumn<Double>( stab, "DIRECTION" ) ; 255 TableColumn *azCol = tpoolw->construct( stab, "AZIMUTH" ) ; 256 TableColumn *elCol = tpoolw->construct( stab, "ELEVATION" ) ; 257 TableColumn *tcalidCol = tpoolw->construct( stab, "TCAL_ID" ) ; 258 TableColumn *focusidCol = tpoolw->construct( stab, "FOCUS_ID" ) ; 259 TableColumn *weatheridCol = tpoolw->construct( stab, "WEATHER_ID" ) ; 260 ArrayColumn<Double> *srcpmCol = new ArrayColumn<Double>( stab, "SRCPROPERMOTION" ) ; 261 ArrayColumn<Double> *srcdirCol = new ArrayColumn<Double>( stab, "SRCDIRECTION" ) ; 262 TableColumn *srcvelCol = tpoolw->construct( stab, "SRCVELOCITY" ) ; 263 ArrayColumn<Double> *scanrateCol = new ArrayColumn<Double>( stab, "SCANRATE" ) ; 264 236 265 237 // MAIN 266 238 // Iterate over several ids 267 Int oldnr = table_->nrow() ;268 239 map<Int, uInt> ifmap ; // (IFNO, FREQ_ID) pair 269 240 ROArrayQuantColumn<Double> *sharedQDArrCol = new ROArrayQuantColumn<Double>( anttab, "POSITION" ) ; … … 289 260 os_ << "end fill init: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 290 261 262 // row based 263 Table &stab = table_->table() ; 264 TableRow row( stab ) ; 265 TableRecord &trec = row.record() ; 266 RecordFieldPtr< Array<Float> > spRF( trec, "SPECTRA" ) ; 267 RecordFieldPtr< Array<uChar> > ucarrRF( trec, "FLAGTRA" ) ; 268 RecordFieldPtr<Double> timeRF( trec, "TIME" ) ; 269 RecordFieldPtr< Array<Float> > tsysRF( trec, "TSYS" ) ; 270 RecordFieldPtr<Double> intervalRF( trec, "INTERVAL" ) ; 271 RecordFieldPtr< Array<Double> > dirRF( trec, "DIRECTION" ) ; 272 RecordFieldPtr<Float> azRF( trec, "AZIMUTH" ) ; 273 RecordFieldPtr<Float> elRF( trec, "ELEVATION" ) ; 274 RecordFieldPtr< Array<Double> > scrRF( trec, "SCANRATE" ) ; 275 RecordFieldPtr<uInt> cycleRF( trec, "CYCLENO" ) ; 276 RecordFieldPtr<uInt> flrRF( trec, "FLAGROW" ) ; 277 RecordFieldPtr<uInt> tcalidRF( trec, "TCAL_ID" ) ; 278 RecordFieldPtr<uInt> widRF( trec, "WEATHER_ID" ) ; 279 RecordFieldPtr<uInt> polnoRF( trec, "POLNO" ) ; 280 281 282 // REFBEAMNO 283 RecordFieldPtr<Int> intRF( trec, "REFBEAMNO" ) ; 284 *intRF = 0 ; 285 286 // FIT_ID 287 intRF.attachToRecord( trec, "FIT_ID" ) ; 288 *intRF = -1 ; 289 290 // OPACITY 291 RecordFieldPtr<Float> floatRF( trec, "OPACITY" ) ; 292 *floatRF = 0.0 ; 293 291 294 // 292 295 // ITERATION: OBSERVATION_ID 293 296 // 294 Int added0 = 0 ;295 Int current0 = table_->nrow() ;296 297 TableIterator iter0( mstable_, "OBSERVATION_ID" ) ; 297 298 while( !iter0.pastEnd() ) { … … 329 330 // ITERATION: FEED1 330 331 // 331 Int added1 = 0 ;332 Int current1 = table_->nrow() ;333 332 TableIterator iter1( t0, "FEED1" ) ; 334 333 while( !iter1.pastEnd() ) { … … 341 340 tpoolr->destroy( tcolr ) ; 342 341 nbeam++ ; 342 343 // BEAMNO 344 RecordFieldPtr<uInt> uintRF( trec, "BEAMNO" ) ; 345 *uintRF = feedId ; 346 347 // FOCUS_ID 348 uintRF.attachToRecord( trec, "FOCUS_ID" ) ; 349 *uintRF = 0 ; 350 343 351 time1 = gettimeofday_sec() ; 344 352 os_ << "end 1st iteration init: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; … … 346 354 // ITERATION: FIELD_ID 347 355 // 348 Int added2 = 0 ;349 Int current2 = table_->nrow() ;350 356 TableIterator iter2( t1, "FIELD_ID" ) ; 351 357 while( !iter2.pastEnd() ) { … … 362 368 String fieldName = tcolr->asString( fieldId ) + "__" + String::toString(fieldId) ; 363 369 tpoolr->destroy( tcolr ) ; 370 371 // FIELDNAME 372 RecordFieldPtr<String> strRF( trec, "FIELDNAME" ) ; 373 *strRF = fieldName ; 374 375 364 376 time1 = gettimeofday_sec() ; 365 377 os_ << "end 2nd iteration init: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; … … 367 379 // ITERATION: DATA_DESC_ID 368 380 // 369 Int added3 = 0 ;370 Int current3 = table_->nrow() ;371 381 TableIterator iter3( t2, "DATA_DESC_ID" ) ; 372 382 while( !iter3.pastEnd() ) { … … 383 393 Int spwId = tcolr->asInt( ddId ) ; 384 394 tpoolr->destroy( tcolr ) ; 395 396 // IFNO 397 uintRF.attachToRecord( trec, "IFNO" ) ; 398 *uintRF = (uInt)spwId ; 399 385 400 // polarization information 386 401 tcolr = tpoolr->construct( poltab, "NUM_CORR" ) ; … … 396 411 // source information 397 412 MSSource srctabSel( srctab( srctab.col("SOURCE_ID") == srcId && srctab.col("SPECTRAL_WINDOW_ID") == spwId ) ) ; 398 //os_ << "srcId = " << srcId << " spwId = " << spwId << " nrow = " << srctabSel.nrow() << LogIO::POST ;413 os_ << "srcId = " << srcId << " spwId = " << spwId << " nrow = " << srctabSel.nrow() << LogIO::POST ; 399 414 tcolr = tpoolr->construct( srctabSel, "NAME" ) ; 400 415 String srcName = tcolr->asString( 0 ) ; 401 416 tpoolr->destroy( tcolr ) ; 417 418 // SRCNAME 419 strRF.attachToRecord( trec, "SRCNAME" ) ; 420 *strRF = srcName ; 421 402 422 //os_ << "srcName = " << srcName << LogIO::POST ; 403 423 ROArrayColumn<Double> *roArrDCol = new ROArrayColumn<Double>( srctabSel, "PROPER_MOTION" ) ; 404 424 Array<Double> srcPM = (*roArrDCol)( 0 ) ; 405 425 delete roArrDCol ; 426 427 // SRCPROPERMOTION 428 RecordFieldPtr< Array<Double> > darrRF( trec, "SRCPROPERMOTION" ) ; 429 *darrRF = srcPM ; 430 406 431 //os_ << "srcPM = " << srcPM << LogIO::POST ; 407 432 roArrDCol = new ROArrayColumn<Double>( srctabSel, "DIRECTION" ) ; 408 433 Array<Double> srcDir = (*roArrDCol)( 0 ) ; 409 434 delete roArrDCol ; 435 436 // SRCDIRECTION 437 darrRF.attachToRecord( trec, "SRCDIRECTION" ) ; 438 *darrRF = srcDir ; 439 410 440 //os_ << "srcDir = " << srcDir << LogIO::POST ; 411 441 Array<Double> sysVels ; … … 421 451 sysVel = sysVels( IPosition(1,0) ) ; 422 452 } 453 454 // SRCVELOCITY 455 RecordFieldPtr<Double> doubleRF( trec, "SRCVELOCITY" ) ; 456 *doubleRF = sysVel ; 457 423 458 //delete tmpArrCol ; 424 459 //os_ << "sysVel = " << sysVel << LogIO::POST ; … … 451 486 } 452 487 uInt molId = table_->molecules().addEntry( restFreqs, transitionName, transitionName ) ; 488 489 // MOLECULE_ID 490 uintRF.attachToRecord( trec, "MOLECULE_ID" ) ; 491 *uintRF = molId ; 492 453 493 // spectral setup 454 494 MeasFrame mf( me, mp, md ) ; … … 521 561 //os_ << "added to ifmap: (" << spwId << "," << freqId << ")" << LogIO::POST ; 522 562 } 563 564 // FREQ_ID 565 uintRF.attachToRecord( trec, "FREQ_ID" ) ; 566 *uintRF = freqId ; 567 523 568 // for TSYS and TCAL 524 569 MSSysCal caltabsel( caltab( caltab.col("ANTENNA_ID") == antenna_ && caltab.col("FEED_ID") == feedId && caltab.col("SPECTRAL_WINDOW_ID") == spwId ).sort("TIME") ) ; 570 ROScalarMeasColumn<MEpoch> scTimeCol( caltabsel, "TIME" ) ; 571 Block<MEpoch> scTime( caltabsel.nrow() ) ; 572 for ( uInt irow = 0 ; irow < caltabsel.nrow() ; irow++ ) 573 scTime[irow] = scTimeCol( irow ) ; 574 ROScalarColumn<Double> *scIntervalCol = new ROScalarColumn<Double>( caltabsel, "INTERVAL" ) ; 575 Vector<Double> scInterval = scIntervalCol->getColumn() ; 576 delete scIntervalCol ; 577 ROArrayColumn<Float> scTsysCol( caltabsel, "TSYS" ) ; 525 578 time1 = gettimeofday_sec() ; 526 579 os_ << "end 3rd iteration init: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; … … 528 581 // ITERATION: SCAN_NUMBER 529 582 // 530 Int added4 = 0 ;531 Int current4 = table_->nrow() ;532 583 TableIterator iter4( t3, "SCAN_NUMBER" ) ; 533 584 while( !iter4.pastEnd() ) { … … 538 589 Int scanNum = tcolr->asInt( 0 ) ; 539 590 tpoolr->destroy( tcolr ) ; 591 592 // SCANNO 593 uintRF.attachToRecord( trec, "SCANNO" ) ; 594 *uintRF = scanNum - 1 ; 595 540 596 uInt cycle = 0 ; 541 597 time1 = gettimeofday_sec() ; … … 544 600 // ITERATION: STATE_ID 545 601 // 546 Int added5 = 0 ;547 Int current5 = table_->nrow() ;548 602 TableIterator iter5( t4, "STATE_ID" ) ; 549 603 while( !iter5.pastEnd() ) { … … 560 614 561 615 Int nrow = t5.nrow() ; 562 Int prevnr = oldnr ;563 Int addednr = 0 ;564 Int nloop = 0 ;565 616 time1 = gettimeofday_sec() ; 566 617 os_ << "end 5th iteration init: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 567 568 // SPECTRA, FLAG 569 time0 = gettimeofday_sec() ; 570 os_ << "start fill SPECTRA and FLAG: " << time0 << LogIO::POST ; 571 ROArrayColumn<Bool> mFlagCol( t5, "FLAG" ) ; 572 //Cube<Bool> mFlagArr = mFlagCol.getColumn() ; 618 619 Cube<Float> spArr ; 620 Cube<Bool> flArr ; 573 621 if ( isFloatData_ ) { 574 //os_ << "FLOAT_DATA exists" << LogIO::POST;622 ROArrayColumn<Bool> mFlagCol( t5, "FLAG" ) ; 575 623 ROArrayColumn<Float> mFloatDataCol( t5, "FLOAT_DATA" ) ; 576 //Cube<Float> mFloatDataArr = mFloatDataCol.getColumn() ; 577 addednr = nrow*npol ; 578 oldnr += addednr ; 579 stab.addRow( addednr ) ; 580 nloop = npol ; 581 for ( Int irow = 0 ; irow < nrow ; irow++ ) { 582 Matrix<Float> sp = mFloatDataCol( irow ) ; 583 //Matrix<Float> sp = mFloatDataArr.xyPlane( irow ) ; 584 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 585 spCol->put( prevnr+ipol*nrow+irow, sp.row(ipol) ) ; 586 } 587 } 588 for ( Int irow = 0 ; irow < nrow ; irow++ ) { 589 Matrix<Bool> flb = mFlagCol( irow ) ; 590 //Matrix<Bool> flb = mFlagArr.xyPlane( irow ) ; 591 Matrix<uChar> fl( flb.shape() ) ; 592 convertArray( fl, flb ) ; 593 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 594 flCol->put( prevnr+ipol*nrow+irow, fl.row(ipol) ) ; 595 } 596 } 624 spArr = mFloatDataCol.getColumn() ; 625 flArr = mFlagCol.getColumn() ; 597 626 if ( sdh.fluxunit == "" ) { 598 const TableRecord dataColKeys = mFloatDataCol.keywordSet() ;627 const TableRecord &dataColKeys = mFloatDataCol.keywordSet() ; 599 628 if ( dataColKeys.isDefined( "UNIT" ) ) 600 629 sdh.fluxunit = dataColKeys.asString( "UNIT" ) ; … … 602 631 } 603 632 else if ( isData_ ) { 604 //os_ << "DATA exists" << LogIO::POST ; 633 spArr.resize( npol, nchan, nrow ) ; 634 ROArrayColumn<Bool> mFlagCol( t5, "FLAG" ) ; 605 635 ROArrayColumn<Complex> mDataCol( t5, "DATA" ) ; 606 //Cube<Complex> mDataArr = mDataCol.getColumn() ;607 addednr = nrow*npol ;608 oldnr += addednr ;609 stab.addRow( addednr ) ;610 nloop = npol ;611 636 for ( Int irow = 0 ; irow < nrow ; irow++ ) { 612 637 Bool crossOK = False ; 613 638 Matrix<Complex> sp = mDataCol( irow ) ; 614 //Matrix<Complex> sp = mDataArr.xyPlane( irow ) ;639 Matrix<Bool> fl = mFlagCol( irow ) ; 615 640 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 616 641 if ( corrtype[ipol] == Stokes::XY || corrtype[ipol] == Stokes::YX 617 642 || corrtype[ipol] == Stokes::RL || corrtype[ipol] == Stokes::LR ) { 618 643 if ( !crossOK ) { 619 crossOK = True ; 620 Int pidx = prevnr + ipol * nrow + irow ; 621 spCol->put( pidx, real( sp.row(ipol) ) ) ; 622 if ( corrtype[ipol] == Stokes::XY || corrtype[ipol] == Stokes::RL ) 623 spCol->put( pidx+nrow, imag( sp.row(ipol) ) ) ; 624 else 625 spCol->put( pidx+nrow, imag( conj(sp.row(ipol)) ) ) ; 644 spArr( Slicer(Slice(ipol,1),Slice(0,nchan),Slice(irow,1)) ) = real( sp.row( ipol ) ) ; 645 flArr( Slicer(Slice(ipol,1),Slice(0,nchan),Slice(irow,1)) ) = fl.row( ipol ) ; 646 if ( corrtype[ipol] == Stokes::XY || corrtype[ipol] == Stokes::RL ) { 647 spArr( Slicer(Slice(ipol+1,1),Slice(0,nchan),Slice(irow,1)) ) = imag( sp.row( ipol ) ) ; 648 flArr( Slicer(Slice(ipol+1,1),Slice(0,nchan),Slice(irow,1)) ) = fl.row( ipol ) ; 649 } 650 else { 651 spArr( Slicer(Slice(ipol+1,1),Slice(0,nchan),Slice(irow,1)) ) = imag( conj( sp.row( ipol ) ) ) ; 652 flArr( Slicer(Slice(ipol+1,1),Slice(0,nchan),Slice(irow,1)) ) = fl.row( ipol ) ; 653 } 626 654 } 627 655 } 628 656 else { 629 spCol->put( prevnr+ipol*nrow+irow, real( sp.row(ipol) ) ) ; 630 } 631 } 632 } 633 for ( Int irow = 0 ; irow < nrow ; irow++ ) { 634 Bool crossOK = False ; 635 Matrix<Bool> flb = mFlagCol( irow ) ; 636 //Matrix<Bool> flb = mFlagArr.xyPlane( irow ) ; 637 Matrix<uChar> fl( flb.shape() ) ; 638 convertArray( fl, flb ) ; 639 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 640 if ( corrtype[ipol] == Stokes::XY || corrtype[ipol] == Stokes::YX 641 || corrtype[ipol] == Stokes::RL || corrtype[ipol] == Stokes::LR ) { 642 if ( !crossOK ) { 643 crossOK = True ; 644 Int pidx = prevnr + ipol * nrow + irow ; 645 flCol->put( pidx, fl.row(ipol) ) ; 646 flCol->put( pidx+nrow, fl.row(ipol+1) ) ; 647 } 648 } 649 else { 650 flCol->put( prevnr+ipol*nrow+irow, fl.row(ipol) ) ; 657 spArr( Slicer(Slice(ipol,1),Slice(0,nchan),Slice(irow,1)) ) = real( sp.row( ipol ) ) ; 658 flArr( Slicer(Slice(ipol,1),Slice(0,nchan),Slice(irow,1)) ) = fl.row( ipol ) ; 651 659 } 652 660 } 653 661 } 654 662 if ( sdh.fluxunit == "" ) { 655 const TableRecord dataColKeys = mDataCol.keywordSet() ;663 const TableRecord &dataColKeys = mDataCol.keywordSet() ; 656 664 if ( dataColKeys.isDefined( "UNIT" ) ) 657 sdh.fluxunit = dataColKeys.asString( "UNIT" ) ; 658 } 665 sdh.fluxunit = dataColKeys.asString( "UNIT" ) ; 666 } 659 667 } 660 time1 = gettimeofday_sec() ;661 os_ << "end fill SPECTRA and FLAGTRA: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ;662 663 // number of rows added in this cycle664 //os_ << "prevnr = " << prevnr << LogIO::POST ;665 //os_ << "addednr = " << addednr << LogIO::POST ;666 RefRows rows( prevnr, prevnr+addednr-1 ) ;667 668 // TIME669 time0 = gettimeofday_sec() ;670 os_ << "start fill TIME: " << time0 << LogIO::POST ;671 668 ROScalarMeasColumn<MEpoch> *mTimeCol = new ROScalarMeasColumn<MEpoch>( t5, "TIME" ) ; 672 Int tidx = prevnr ; 673 for ( Int i = 0 ; i < nloop ; i++ ) { 674 for ( Int j = 0 ; j < nrow ; j++ ) { 675 timeCol->put( tidx++, (*mTimeCol)( j ) ) ; 676 } 669 Block<MEpoch> mTimeB( nrow ) ; 670 for ( Int irow = 0 ; irow < nrow ; irow++ ) 671 mTimeB[irow] = (*mTimeCol)( irow ) ; 672 ROTableColumn *mIntervalCol = tpoolr->construct( t5, "INTERVAL" ) ; 673 ROTableColumn *mFlagRowCol = tpoolr->construct( t5, "FLAG_ROW" ) ; 674 //Vector<Double> sysCalTime( nrow, -1.0 ) ; 675 Block<Double> sysCalTime( nrow, -1.0 ) ; 676 Block<Int> sysCalIdx( nrow, -1 ) ; 677 if ( isSysCal_ ) { 678 //getSysCalTime( caltabsel, mTimeB, sysCalTime, sysCalIdx ) ; 679 getSysCalTime( scTime, scInterval, mTimeB, sysCalTime, sysCalIdx ) ; 677 680 } 678 time1 = gettimeofday_sec() ; 679 os_ << "end fill TIME: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 680 681 // TSYS 682 time0 = gettimeofday_sec() ; 683 os_ << "start fill TSYS: " << time0 << LogIO::POST ; 684 Vector<Double> sysCalTime( nrow, -1.0 ) ; 685 //Vector<Double> sysCalTime ; 686 if ( isSysCal_ ) { 687 //sysCalTime = getSysCalTime( caltabsel, *mTimeCol ) ; 688 getSysCalTime( caltabsel, *mTimeCol, sysCalTime ) ; 689 tidx = prevnr ; 690 uInt calidx = 0 ; 691 for ( Int i = 0 ; i < nrow ; i++ ) { 692 Matrix<Float> tsys ; 693 calidx = getTsys( calidx, tsys, caltabsel, sysCalTime(i) ) ; 694 //os_ << "tsys = " << tsys << LogIO::POST ; 695 uInt ncol = tsys.ncolumn() ; 696 if ( ncol == 0 ) { 697 IPosition cShape = IPosition( 2, npol, 1 ) ; 698 tsys.resize( cShape ) ; 699 tsys = 1.0 ; 700 } 701 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 702 //floatArrCol->put( prevnr+i+nrow*ipol, tsys.row( ipol ) ) ; 703 tsysCol->put( prevnr+i+nrow*ipol, tsys.row( ipol ) ) ; 704 } 705 } 706 } 707 else { 708 Vector<Float> tsys( 1, 1.0 ) ; 709 for ( Int i = prevnr ; i < prevnr+addednr ; i++ ) 710 tsysCol->put( i, tsys ) ; 711 } 712 time1 = gettimeofday_sec() ; 713 os_ << "end fill TSYS: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 714 715 716 // INTERVAL 717 time0 = gettimeofday_sec() ; 718 os_ << "start fill INTERVAL: " << time0 << LogIO::POST ; 719 tcolr = tpoolr->construct( t5, "INTERVAL" ) ; 720 //Vector<Double> integ = mIntervalCol->getColumn() ; 721 for ( int i = 0 ; i < nloop ; i++ ) { 722 //Int startrow = prevnr + i ; 723 //Int endrow = startrow + nrow - 1 ; 724 //RefRows prows( startrow, endrow ) ; 725 //intervalCol->putColumnCells( prows, integ ) ; 726 for ( int j = 0 ; j < nrow ; j++ ) { 727 intervalCol->putScalar( prevnr+i*nrow+j, (Double)(tcolr->asdouble( j )) ) ; 728 } 729 } 730 tpoolr->destroy( tcolr ) ; 731 time1 = gettimeofday_sec() ; 732 os_ << "end fill INTERVAL: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 733 734 // SRCTYPE 735 time0 = gettimeofday_sec() ; 736 os_ << "start fill SRCTYPE: " << time0 << LogIO::POST ; 681 delete mTimeCol ; 682 Matrix<Float> defaulttsys( npol, 1, 1.0 ) ; 683 //uInt calidx = 0 ; 737 684 Int srcType = getSrcType( stateId, tpoolr ) ; 738 for ( int i = 0 ; i < addednr ; i++ ) {739 srctypeCol->putScalar( prevnr+i, srcType ) ;740 }741 //Vector<Int> *srcType = new Vector<Int>( addednr, getSrcType( stateId ) ) ;742 //srcTypeCol->putColumnCells( rows, *srcType ) ;743 time1 = gettimeofday_sec() ;744 os_ << "end fill SRCTYPE: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ;745 746 // DIRECTION, AZIMUTH, ELEVATION, SCANRATE747 time0 = gettimeofday_sec() ;748 os_ << "start fill DIRECTION, AZIMUTH, ELEVATION, SCANRATE: " << time0 << LogIO::POST ;749 685 Vector<Double> defaultScanrate( 2, 0.0 ) ; 750 686 uInt diridx = 0 ; 751 687 MDirection::Types dirType ; 752 if ( getPt_ ) { 753 for ( Int i = 0 ; i < nrow ; i++ ) { 688 uInt wid = 0 ; 689 Int pidx = 0 ; 690 Bool crossOK = False ; 691 Block<uInt> polnos( npol, 99 ) ; 692 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 693 Block<uInt> p = getPolNo( corrtype[ipol] ) ; 694 if ( p.size() > 1 ) { 695 if ( crossOK ) continue ; 696 else { 697 polnos[pidx] = p[0] ; 698 pidx++ ; 699 polnos[pidx] = p[1] ; 700 pidx++ ; 701 crossOK = True ; 702 } 703 } 704 else { 705 polnos[pidx] = p[0] ; 706 pidx++ ; 707 } 708 } 709 710 // SRCTYPE 711 intRF.attachToRecord( trec, "SRCTYPE" ) ; 712 *intRF = srcType ; 713 714 for ( Int irow = 0 ; irow < nrow ; irow++ ) { 715 // CYCLENO 716 *cycleRF = cycle ; 717 718 // FLAGROW 719 *flrRF = (uInt)mFlagRowCol->asBool( irow ) ; 720 721 // SPECTRA, FLAG 722 //Matrix<Float> sp = mFloatDataCol( irow ) ; 723 //Matrix<Bool> flb = mFlagCol( irow ) ; 724 Matrix<Float> sp = spArr.xyPlane( irow ) ; 725 Matrix<Bool> flb = flArr.xyPlane( irow ) ; 726 Matrix<uChar> fl( flb.shape() ) ; 727 convertArray( fl, flb ) ; 728 729 // TIME 730 *timeRF = mTimeB[irow].get("s").getValue() ; 731 732 // INTERVAL 733 *intervalRF = (Double)(mIntervalCol->asdouble( irow )) ; 734 735 // TSYS 736 Matrix<Float> tsys ; 737 //calidx = getTsys( calidx, tsys, caltabsel, sysCalTime[irow] ) ; 738 // calidx = getTsys( calidx, tsys, scTsysCol, scTime, sysCalTime[irow] ) ; 739 // if ( tsys.nelements() == 0 ) 740 // tsys = defaulttsys ; 741 if ( sysCalIdx[irow] != -1 ) 742 tsys = scTsysCol( irow ) ; 743 else 744 tsys = defaulttsys ; 745 746 // TCAL_ID 747 Block<uInt> tcalids = getTcalId( feedId, spwId, sysCalTime[irow] ) ; 748 if ( tcalids.size() == 0 ) { 749 tcalids.resize( npol ) ; 750 tcalids = 0 ; 751 } 752 753 // WEATHER_ID 754 if ( isWeather_ ) 755 wid = getWeatherId( wid, mTimeB[irow].get("s").getValue() ) ; 756 *widRF = wid ; 757 758 759 // DIRECTION, AZEL, SCANRATE 760 if ( getPt_ ) { 754 761 Vector<Double> dir ; 755 762 Vector<Double> scanrate ; 756 763 String refString ; 757 diridx = getDirection( diridx, dir, scanrate, refString, pointtab, (*mTimeCol)(i).get("s").getValue() ) ; 758 //os_ << "diridx = " << diridx << " dmTimeCol(" << i << ") = " << mTimeCol(i).get("s").getValue()-mTimeCol(0).get("s").getValue() << LogIO::POST ; 759 //os_ << "dir = " << dir << LogIO::POST ; 760 //os_ << "scanrate = " << scanrate << LogIO::POST ; 761 //os_ << "refString = " << refString << LogIO::POST ; 764 diridx = getDirection( diridx, dir, scanrate, refString, pointtab, mTimeB[irow].get("s").getValue() ) ; 762 765 MDirection::getType( dirType, refString ) ; 763 //os_ << "dirType = " << dirType << LogIO::POST ; 764 mf.resetEpoch( (*mTimeCol)(i) ) ; 766 mf.resetEpoch( mTimeB[irow] ) ; 765 767 mf.resetDirection( MDirection( MVDirection(dir), dirType ) ) ; 766 768 if ( refString == "J2000" ) { 767 //os_ << "J2000" << LogIO::POST ; 768 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 769 dirCol->put( prevnr+i+nrow*ipol, dir ) ; 770 } 769 *dirRF = dir ; 771 770 MDirection::Convert toazel( dirType, MDirection::Ref( MDirection::AZEL, mf ) ) ; 772 771 Vector<Double> azel = toazel( dir ).getAngle("rad").getValue() ; 773 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 774 azCol->putScalar( prevnr+i+nrow*ipol, (Float)azel(0) ) ; 775 elCol->putScalar( prevnr+i+nrow*ipol, (Float)azel(1) ) ; 776 } 772 *azRF = (Float)azel(0) ; 773 *elRF = (Float)azel(1) ; 777 774 } 778 775 else if ( refString(0,4) == "AZEL" ) { 779 //os_ << "AZEL" << LogIO::POST ; 780 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 781 azCol->putScalar( prevnr+i+nrow*ipol, (Float)dir(0) ) ; 782 elCol->putScalar( prevnr+i+nrow*ipol, (Float)dir(1) ) ; 783 } 776 *azRF = (Float)dir(0) ; 777 *elRF = (Float)dir(1) ; 784 778 MDirection::Convert toj2000( dirType, MDirection::Ref( MDirection::J2000, mf ) ) ; 785 779 Vector<Double> newdir = toj2000( dir ).getAngle("rad").getValue() ; 786 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 787 dirCol->put( prevnr+i+nrow*ipol, newdir ) ; 788 } 780 *dirRF = newdir ; 789 781 } 790 782 else { 791 //os_ << "OTHER: " << refString << LogIO::POST ;792 783 MDirection::Convert toazel( dirType, MDirection::Ref( MDirection::AZEL, mf ) ) ; 793 784 Vector<Double> azel = toazel( dir ).getAngle("rad").getValue() ; 794 785 MDirection::Convert toj2000( dirType, MDirection::Ref( MDirection::J2000, mf ) ) ; 795 786 Vector<Double> newdir = toj2000( dir ).getAngle("rad").getValue() ; 796 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 797 dirCol->put( prevnr+i+nrow*ipol, newdir ) ; 798 azCol->putScalar( prevnr+i+nrow*ipol, (Float)dir(0) ) ; 799 elCol->putScalar( prevnr+i+nrow*ipol, (Float)dir(1) ) ; 800 } 787 *dirRF = newdir ; 788 *azRF = (Float)azel(0) ; 789 *elRF = (Float)azel(1) ; 801 790 } 802 791 if ( scanrate.size() != 0 ) { 803 //os_ << "scanrate.size() = " << scanrate.size() << LogIO::POST ; 804 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 805 scanrateCol->put( prevnr+i+nrow*ipol, scanrate ) ; 806 } 792 *scrRF = scanrate ; 807 793 } 808 794 else { 809 //os_ << "scanrate.size() = " << scanrate.size() << LogIO::POST ; 810 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 811 scanrateCol->put( prevnr+i+nrow*ipol, defaultScanrate ) ; 812 } 795 *scrRF = defaultScanrate ; 813 796 } 814 797 } 815 } 816 else { 817 // All directions are set to source direction 818 //ROArrayMeasColumn<MDirection> dmcol( pointtab, "DIRECTION" ) ; 819 //ROArrayColumn<Double> dcol( pointtab, "DIRECTION" ) ; 820 //IPosition ip( dmcol(0).shape().nelements(), 0 ) ; 821 //IPosition outp( 1, 2 ) ; 822 //String ref = dmcol(0)(ip).getRefString() ; 823 String ref = md.getRefString() ; 824 //Slice ds( 0, 2, 1 ) ; 825 //Slice ds0( 0, 1, 1 ) ; 826 //Slicer dslice0( ds, ds0 ) ; 827 //Vector<Double> defaultDir = dcol(0)(dslice0).reform(outp) ; 828 Vector<Double> defaultDir = srcDir ; 829 MDirection::getType( dirType, "J2000" ) ; 830 //mf.resetDirection( MDirection( MVDirection(srcDir), dirType ) ) ; 831 if ( ref != "J2000" ) { 832 ROScalarMeasColumn<MEpoch> tmCol( pointtab, "TIME" ) ; 833 mf.resetEpoch( tmCol( 0 ) ) ; 834 MDirection::Convert toj2000( dirType, MDirection::Ref( MDirection::J2000, mf ) ) ; 835 defaultDir = toj2000( defaultDir ).getAngle("rad").getValue() ; 798 else { 799 String ref = md.getRefString() ; 800 Vector<Double> defaultDir = srcDir ; 801 MDirection::getType( dirType, "J2000" ) ; 802 if ( ref != "J2000" ) { 803 ROScalarMeasColumn<MEpoch> tmCol( pointtab, "TIME" ) ; 804 mf.resetEpoch( tmCol( 0 ) ) ; 805 MDirection::Convert toj2000( dirType, MDirection::Ref( MDirection::J2000, mf ) ) ; 806 defaultDir = toj2000( defaultDir ).getAngle("rad").getValue() ; 807 } 808 mf.resetEpoch( mTimeB[irow] ) ; 809 MDirection::Convert toazel( dirType, MDirection::Ref( MDirection::AZEL, mf ) ) ; 810 Vector<Double> azel = toazel( defaultDir ).getAngle("rad").getValue() ; 811 *azRF = (Float)azel(0) ; 812 *elRF = (Float)azel(1) ; 813 *dirRF = defaultDir ; 814 *scrRF = defaultScanrate ; 836 815 } 837 for ( Int i = 0 ; i < nrow ; i++ ) { 838 mf.resetEpoch( (*mTimeCol)(i) ) ; 839 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 840 Int localidx = prevnr+i+nrow*ipol ; 841 MDirection::Convert toazel( dirType, MDirection::Ref( MDirection::AZEL, mf ) ) ; 842 Vector<Double> azel = toazel( defaultDir ).getAngle("rad").getValue() ; 843 azCol->putScalar( localidx, (Float)azel(0) ) ; 844 elCol->putScalar( localidx, (Float)azel(1) ) ; 845 dirCol->put( localidx, defaultDir ) ; 846 scanrateCol->put( localidx, defaultScanrate ) ; 847 } 816 817 // Polarization dependent things 818 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 819 // POLNO 820 *polnoRF = polnos[ipol] ; 821 822 *spRF = sp.row( ipol ) ; 823 *ucarrRF = fl.row( ipol ) ; 824 *tsysRF = tsys.row( ipol ) ; 825 *tcalidRF = tcalids[ipol] ; 826 827 // Commit row 828 stab.addRow() ; 829 row.put( stab.nrow()-1 ) ; 848 830 } 849 831 } 832 //delete mTimeCol ; 833 tpoolr->destroy( mIntervalCol ) ; 834 tpoolr->destroy( mFlagRowCol ) ; 835 850 836 time1 = gettimeofday_sec() ; 851 os_ << "end fill DIRECTION, AZIMUTH, ELEVATION, SCANRATE: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 852 853 // CYCLENO 854 time0 = gettimeofday_sec() ; 855 os_ << "start fill CYCLENO: " << time0 << LogIO::POST ; 856 for ( int i = 0 ; i < nloop ; i++ ) { 857 for ( int j = 0 ; j < nrow ; j++ ) { 858 cyclenoCol->putScalar( prevnr+nrow*i+j, cycle+j ) ; 859 } 860 } 861 cycle += nrow ; 862 time1 = gettimeofday_sec() ; 863 os_ << "end fill CYCLENO: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 864 865 // POLNO 866 time0 = gettimeofday_sec() ; 867 os_ << "start fill POLNO: " << time0 << LogIO::POST ; 868 Int pidx = 0 ; 869 Bool crossOK = False ; 870 for ( int i = 0 ; i < npol ; i++ ) { 871 Vector<uInt> polnos = getPolNo( corrtype[i] ) ; 872 if ( polnos.size() > 1 ) { 873 if ( crossOK ) continue ; 874 else crossOK = True ; 875 } 876 for ( uInt j = 0 ; j < polnos.size() ; j++ ) { 877 for ( Int irow = 0 ; irow < nrow ; irow++ ) { 878 polnoCol->putScalar( prevnr+pidx*nrow+irow, polnos[j] ) ; 879 } 880 pidx++ ; 881 } 882 } 883 time1 = gettimeofday_sec() ; 884 os_ << "end fill POLNO: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 885 886 // FLAGROW 887 time0 = gettimeofday_sec() ; 888 os_ << "start fill FLAGROW: " << time0 << LogIO::POST ; 889 tcolr = tpoolr->construct( t5, "FLAG_ROW" ) ; 890 for ( int i = 0 ; i < nloop ; i++ ) { 891 for ( int j = 0 ; j < nrow ; j++ ) { 892 flagrowCol->putScalar( prevnr+nrow*i+j, (uInt)(tcolr->asBool( j )) ) ; 893 } 894 } 895 tpoolr->destroy( tcolr ) ; 896 time1 = gettimeofday_sec() ; 897 os_ << "end fill FLAGROW: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 898 899 // TCAL_ID 900 time0 = gettimeofday_sec() ; 901 os_ << "start fill TCAL_ID: " << time0 << LogIO::POST ; 902 if ( isSysCal_ ) { 903 for( Int irow = 0 ; irow < nrow ; irow++ ) { 904 Vector<uInt> tcalids = getTcalId( feedId, spwId, sysCalTime[irow] ) ; 905 if ( tcalids.size() == 0 ) { 906 tcalids.resize( npol ) ; 907 tcalids = 0 ; 908 } 909 for ( Int ipol = 0 ; ipol < npol ; ipol++ ) { 910 tcalidCol->putScalar( prevnr+irow+nrow*ipol, tcalids[ipol] ) ; 911 } 912 } 913 } 914 else { 915 //Vector<uInt> tcalid( addednr, 0 ) ; 916 //uIntCol->putColumnCells( rows, tcalid ) ; 917 uInt tcalid = 0 ; 918 for ( Int irow = 0 ; irow < addednr ; irow++ ) 919 tcalidCol->putScalar( prevnr+irow, tcalid ) ; 920 } 921 time1 = gettimeofday_sec() ; 922 os_ << "end fill TCAL_ID: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 923 924 // WEATHER_ID 925 time0 = gettimeofday_sec() ; 926 os_ << "start fill WEATHER_ID: " << time0 << LogIO::POST ; 927 if ( isWeather_ ) { 928 uInt wid = 0 ; 929 for ( Int irow = 0 ; irow < nrow ; irow++ ) { 930 wid = getWeatherId( wid, (*mTimeCol)(irow).get("s").getValue() ) ; 931 for ( Int ipol = 0 ; ipol < nloop ; ipol++ ) { 932 weatheridCol->putScalar( prevnr+ipol*nrow+irow, wid ) ; 933 } 934 } 935 } 936 time1 = gettimeofday_sec() ; 937 os_ << "end fill WEATHER_ID: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 938 939 delete mTimeCol ; 940 941 //os_ << "field: " << fieldId << " scan: " << scanNum << " obs: " << obsId << " state: " << stateId << " ddid: " << ddId << endl ; 942 os_ << "addednr = " << addednr << endl ; 943 added5 += addednr ; 837 os_ << "end 5th iteration: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 838 944 839 iter5.next() ; 945 840 } 946 841 947 // SCANNO 948 // MS: 1-base 949 // Scantable: 0-base 950 time0 = gettimeofday_sec() ; 951 os_ << "start fill SCANNO: " << time0 << LogIO::POST ; 952 Int dest5 = current5 + added5 ; 953 scanNum -= 1 ; 954 for ( Int irow = current5 ; irow < dest5 ; irow++ ) 955 scannoCol->putScalar( irow, (uInt)scanNum ) ; 956 time1 = gettimeofday_sec() ; 957 os_ << "end fill SCANNO: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 958 959 os_ << "added5 = " << added5 << endl ; 960 added4 += added5 ; 842 961 843 iter4.next() ; 962 844 } 963 845 964 // IFNO 965 time0 = gettimeofday_sec() ; 966 os_ << "start fill IFNO: " << time0 << LogIO::POST ; 967 Int dest4 = current4 + added4 ; 968 for ( Int irow = current4 ; irow < dest4 ; irow++ ) 969 ifnoCol->putScalar( irow, (uInt)spwId ) ; 970 time1 = gettimeofday_sec() ; 971 os_ << "end fill IFNO: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 972 973 // FREQ_ID 974 time0 = gettimeofday_sec() ; 975 os_ << "start fill FREQ_ID: " << time0 << LogIO::POST ; 976 uInt fid = ifmap[spwId] ; 977 for ( Int irow = current4 ; irow < dest4 ; irow++ ) 978 freqidCol->putScalar( irow, fid ) ; 979 time1 = gettimeofday_sec() ; 980 os_ << "end fill FREQ_ID: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 981 982 // MOLECULE_ID 983 time0 = gettimeofday_sec() ; 984 os_ << "start fill MOLECULE_ID: " << time0 << LogIO::POST ; 985 for ( Int irow = current4 ; irow < dest4 ; irow++ ) 986 molidCol->putScalar( irow, molId ) ; 987 time1 = gettimeofday_sec() ; 988 os_ << "end fill MOLECULE_ID: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 989 990 // SRCNAME 991 time0 = gettimeofday_sec() ; 992 os_ << "start fill SRCNAME: " << time0 << LogIO::POST ; 993 for ( Int irow = current4 ; irow < dest4 ; irow++ ) 994 srcnameCol->putScalar( irow, srcName ) ; 995 time1 = gettimeofday_sec() ; 996 os_ << "end fill SRCNAME: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 997 998 // SRCVELOCITY, SRCPROPERMOTION and SRCDIRECTION 999 // no reference conversion for direction at the moment (assume J2000) 1000 // no reference conversion for velocity at the moment (assume LSRK) 1001 time0 = gettimeofday_sec() ; 1002 os_ << "start fill SRCPROPERMOTION: " << time0 << LogIO::POST ; 1003 for ( Int irow = current4 ; irow < dest4 ; irow++ ) 1004 srcpmCol->put( irow, srcPM ) ; 1005 time1 = gettimeofday_sec() ; 1006 os_ << "end fill SRCPROPERMOTION: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 1007 time0 = gettimeofday_sec() ; 1008 os_ << "start fill SRCDIRECTION: " << time0 << LogIO::POST ; 1009 for ( Int irow = current4 ; irow < dest4 ; irow++ ) 1010 srcdirCol->put( irow, srcDir ) ; 1011 time1 = gettimeofday_sec() ; 1012 os_ << "end fill SRCDIRECTION: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 1013 time0 = gettimeofday_sec() ; 1014 os_ << "start fill SRCVELOCITY: " << time0 << LogIO::POST ; 1015 for ( Int irow = current4 ; irow < dest4 ; irow++ ) 1016 srcvelCol->putScalar( irow, sysVel ) ; 1017 time1 = gettimeofday_sec() ; 1018 os_ << "end fill SRCVELOCITY: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 1019 1020 os_ << "added4 = " << added4 << endl ; 1021 added3 += added4 ; 846 1022 847 iter3.next() ; 1023 848 } 1024 849 1025 // FIELDNAME1026 time0 = gettimeofday_sec() ;1027 os_ << "start fill FIELDNAME: " << time0 << LogIO::POST ;1028 Int dest3 = current3 + added3 ;1029 for ( Int irow = current3 ; irow < dest3 ; irow++ )1030 fieldnameCol->putScalar( irow, fieldName ) ;1031 time1 = gettimeofday_sec() ;1032 os_ << "end fill FIELDNAME: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ;1033 850 1034 os_ << "added3 = " << added3 << endl ;1035 added2 += added3 ;1036 851 iter2.next() ; 1037 852 } 1038 853 1039 // BEAMNO 1040 time0 = gettimeofday_sec() ; 1041 os_ << "start fill BEAMNO: " << time0 << LogIO::POST ; 1042 Int dest2 = current2 + added2 ; 1043 for ( Int irow = current2 ; irow < dest2 ; irow++ ) 1044 beamnoCol->putScalar( irow, (uInt)feedId ) ; 1045 time1 = gettimeofday_sec() ; 1046 os_ << "end fill BEAMNO: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 1047 1048 // FOCUS_ID 1049 // tentative 1050 time0 = gettimeofday_sec() ; 1051 os_ << "start fill FOCUS_ID: " << time0 << LogIO::POST ; 1052 uInt focusId = 0 ; 1053 for ( Int irow = current2 ; irow < dest2 ; irow++ ) 1054 focusidCol->putScalar( irow, focusId ) ; 1055 time1 = gettimeofday_sec() ; 1056 os_ << "end fill FOCUS_ID: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ; 1057 1058 os_ << "added2 = " << added2 << endl ; 1059 added1 += added2 ; 854 1060 855 iter1.next() ; 1061 856 } 1062 857 if ( sdh.nbeam < nbeam ) sdh.nbeam = nbeam ; 1063 858 1064 os_ << "added1 = " << added1 << endl ;1065 added0 += added1 ;1066 859 iter0.next() ; 1067 860 } 1068 861 1069 os_ << "added0 = " << added0 << endl ;1070 1071 // REFBEAMNO1072 // set 0 at the moment1073 time0 = gettimeofday_sec() ;1074 os_ << "start fill REFBEAMNO: " << time0 << LogIO::POST ;1075 tcolw = tpoolw->construct( stab, "REFBEAMNO" ) ;1076 for ( Int irow = current0 ; irow < added0 ; irow++ )1077 tcolw->putScalar( irow, 0 ) ;1078 tpoolw->destroy( tcolw ) ;1079 time1 = gettimeofday_sec() ;1080 os_ << "end fill REFBEAMNO: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ;1081 1082 // FIT_ID1083 // nothing to do1084 time0 = gettimeofday_sec() ;1085 os_ << "start fill FIT_ID: " << time0 << LogIO::POST ;1086 tcolw = tpoolw->construct( stab, "FIT_ID" ) ;1087 for ( Int irow = current0 ; irow < added0 ; irow++ )1088 tcolw->putScalar( irow, -1 ) ;1089 tpoolw->destroy( tcolw ) ;1090 time1 = gettimeofday_sec() ;1091 os_ << "end fill FIT_ID: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ;1092 1093 // OPACITY1094 // not used?1095 time0 = gettimeofday_sec() ;1096 os_ << "start fill OPACITY: " << time0 << LogIO::POST ;1097 tcolw = tpoolw->construct( stab, "OPACITY" ) ;1098 for ( Int irow = current0 ; irow < added0 ; irow++ )1099 tcolw->putScalar( irow, 0.0 ) ;1100 tpoolw->destroy( tcolw ) ;1101 time1 = gettimeofday_sec() ;1102 os_ << "end fill OPACITY: " << time1 << " (" << time1-time0 << "sec)" << LogIO::POST ;1103 1104 // delete Scantable columns1105 tpoolw->destroy( scannoCol ) ;1106 tpoolw->destroy( cyclenoCol ) ;1107 tpoolw->destroy( beamnoCol ) ;1108 tpoolw->destroy( ifnoCol ) ;1109 tpoolw->destroy( polnoCol ) ;1110 tpoolw->destroy( freqidCol ) ;1111 tpoolw->destroy( molidCol ) ;1112 tpoolw->destroy( flagrowCol ) ;1113 delete timeCol ;1114 tpoolw->destroy( intervalCol ) ;1115 tpoolw->destroy( srcnameCol ) ;1116 tpoolw->destroy( srctypeCol ) ;1117 tpoolw->destroy( fieldnameCol ) ;1118 delete spCol ;1119 delete flCol ;1120 delete tsysCol ;1121 delete dirCol ;1122 tpoolw->destroy( azCol ) ;1123 tpoolw->destroy( elCol ) ;1124 tpoolw->destroy( tcalidCol ) ;1125 tpoolw->destroy( focusidCol ) ;1126 tpoolw->destroy( weatheridCol ) ;1127 delete srcpmCol ;1128 delete srcdirCol ;1129 tpoolw->destroy( srcvelCol ) ;1130 delete scanrateCol ;1131 862 1132 863 delete tpoolr ; 1133 delete tpoolw ;1134 864 1135 865 … … 1306 1036 } 1307 1037 1308 Vector<uInt> MSFiller::getPolNo( Int corrType ) 1038 //Vector<uInt> MSFiller::getPolNo( Int corrType ) 1039 Block<uInt> MSFiller::getPolNo( Int corrType ) 1309 1040 { 1310 1041 double startSec = gettimeofday_sec() ; 1311 1042 os_ << "start MSFiller::getPolNo() startSec=" << startSec << LogIO::POST ; 1312 Vector<uInt> polno( 1 ) ;1043 Block<uInt> polno( 1 ) ; 1313 1044 1314 1045 if ( corrType == Stokes::I || corrType == Stokes::RR || corrType == Stokes::XX ) { … … 1460 1191 } 1461 1192 1462 void MSFiller::fillTcal( boost::object_pool<ROTableColumn> *tpoolr, 1463 boost::object_pool<TableColumn> *tpoolw ) 1193 void MSFiller::fillTcal( boost::object_pool<ROTableColumn> *tpoolr ) 1464 1194 { 1465 1195 double startSec = gettimeofday_sec() ; 1466 1196 os_ << "start MSFiller::fillTcal() startSec=" << startSec << LogIO::POST ; 1467 1197 1468 //MSSysCal sctab = mstable_.sysCal() ;1469 1198 Table sctab = mstable_.sysCal() ; 1470 1199 if ( sctab.nrow() == 0 ) { … … 1482 1211 //os_ << "fillTcal(): npol = " << npol << LogIO::POST ; 1483 1212 Table tab = table_->tcal().table() ; 1484 TableColumn *idCol = tpoolw->construct( tab, "ID" ) ;1485 TableColumn *timeCol = tpoolw->construct( tab, "TIME" ) ;1486 1213 ArrayColumn<Float> tcalCol( tab, "TCAL" ) ; 1487 1214 ROTableColumn *sharedCol ; 1488 ROArrayColumn<Float> scTcalCol ;1489 1215 uInt oldnr = 0 ; 1490 1216 uInt newnr = 0 ; 1217 TableRow row( tab ) ; 1218 TableRecord &trec = row.record() ; 1219 RecordFieldPtr<uInt> idRF( trec, "ID" ) ; 1220 RecordFieldPtr<String> timeRF( trec, "TIME" ) ; 1221 RecordFieldPtr< Array<Float> > tcalRF( trec, "TCAL" ) ; 1491 1222 TableIterator iter0( sctabsel, "FEED_ID" ) ; 1492 1223 while( !iter0.pastEnd() ) { … … 1495 1226 Int feedId = sharedCol->asInt( 0 ) ; 1496 1227 tpoolr->destroy( sharedCol ) ; 1497 //String ffield = "FEED" + String::toString( feedId ) ;1498 //Record rec0 ;1499 1228 TableIterator iter1( t0, "SPECTRAL_WINDOW_ID" ) ; 1500 1229 while( !iter1.pastEnd() ) { … … 1503 1232 Int spwId = sharedCol->asInt( 0 ) ; 1504 1233 tpoolr->destroy( sharedCol ) ; 1505 //String spwfield = "SPW" + String::toString( spwId ) ; 1506 //Record rec1 ; 1507 TableIterator iter2( t1, "TIME" ) ; 1508 while( !iter2.pastEnd() ) { 1509 Table t2 = iter2.table() ; 1510 uInt nrow = t2.nrow() ; // must be 1 1511 //os_ << "fillTcal::t2.nrow = " << nrow << LogIO::POST ; 1512 ROScalarQuantColumn<Double> scTimeCol( t2, "TIME" ) ; 1513 scTcalCol.attach( t2, colTcal_ ) ; 1514 tab.addRow( nrow*npol ) ; 1515 newnr += nrow*npol ; 1516 String sTime = MVTime( scTimeCol(0) ).string( MVTime::YMD ) ; 1234 tmpTcalCol = new ROArrayColumn<Float>( t1, colTcal_ ) ; 1235 ROScalarQuantColumn<Double> scTimeCol( t1, "TIME" ) ; 1236 Vector<uInt> idminmax( 2, oldnr ) ; 1237 for ( uInt irow = 0 ; irow < t1.nrow() ; irow++ ) { 1238 String sTime = MVTime( scTimeCol(irow) ).string( MVTime::YMD ) ; 1239 *timeRF = sTime ; 1240 uInt idx = oldnr ; 1241 Matrix<Float> subtcal = (*tmpTcalCol)( irow ) ; 1517 1242 for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) { 1518 timeCol->putScalar( oldnr+ipol, sTime ) ; 1243 *idRF = idx++ ; 1244 *tcalRF = subtcal.row( ipol ) ; 1245 1246 // commit row 1247 tab.addRow() ; 1248 row.put( tab.nrow()-1 ) ; 1249 1250 newnr++ ; 1519 1251 } 1520 uInt idx = oldnr ; 1521 for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) { 1522 idCol->putScalar( oldnr+ipol, idx++ ) ; 1523 } 1524 Vector<uInt> idminmax( 2, oldnr ) ; 1525 Matrix<Float> subtcal = scTcalCol( 0 ) ; 1526 for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) { 1527 tcalCol.put( oldnr+ipol, subtcal.row( ipol ) ) ; 1528 } 1252 idminmax[0] = oldnr ; 1529 1253 idminmax[1] = newnr - 1 ; 1530 1254 oldnr = newnr ; 1531 1255 1532 //String key = ffield+":"+spwfield+":"+sTime ;1533 1256 String key = keyTcal( feedId, spwId, sTime ) ; 1534 1257 tcalrec_.define( key, idminmax ) ; 1535 //rec1.define( sTime, idminmax ) ;1536 iter2++ ;1537 1258 } 1538 //rec0.defineRecord( spwfield, rec1 );1259 delete tmpTcalCol ; 1539 1260 iter1++ ; 1540 1261 } 1541 //tcalrec_.defineRecord( ffield, rec0 ) ;1542 1262 iter0++ ; 1543 1263 } 1544 1545 tpoolw->destroy( idCol ) ;1546 tpoolw->destroy( timeCol ) ;1547 1264 1548 1265 //tcalrec_.print( std::cout ) ; … … 1588 1305 } 1589 1306 1590 //Vector<Double> MSFiller::getSysCalTime( MSSysCal &tab, MEpoch::ROScalarColumn &tcol ) 1591 void MSFiller::getSysCalTime( MSSysCal &tab, MEpoch::ROScalarColumn &tcol, Vector<Double> &tstr ) 1307 void MSFiller::getSysCalTime( Block<MEpoch> &scTime, Vector<Double> &scInterval, Block<MEpoch> &tcol, Block<Double> &tstr, Block<Int> &tidx ) 1592 1308 { 1593 1309 double startSec = gettimeofday_sec() ; 1594 1310 os_ << "start MSFiller::getSysCalTime() startSec=" << startSec << LogIO::POST ; 1595 //uInt nrow = tcol.table().nrow() ;1596 1311 uInt nrow = tstr.nelements() ; 1597 //Vector<Double> tstr( nrow, -1.0 ) ; 1598 if ( tab.nrow() == 0 ) 1599 //return tstr ; 1312 if ( scTime.nelements() == 0 ) 1600 1313 return ; 1601 uInt scnrow = tab.nrow() ; 1602 ROScalarMeasColumn<MEpoch> scTimeCol( tab, "TIME" ) ; 1603 ROScalarQuantColumn<Double> scIntervalCol( tab, "INTERVAL" ) ; 1314 uInt scnrow = scTime.nelements() ; 1604 1315 uInt idx = 0 ; 1605 1316 const Double half = 0.5e0 ; 1606 1317 for ( uInt i = 0 ; i < nrow ; i++ ) { 1607 Double t = tcol( i ).get( "s" ).getValue() ; 1318 Double t = tcol[i].get( "s" ).getValue() ; 1319 os_ << "getSysCalTime() irow = " << i << " idx = " << idx << LogIO::POST ; 1608 1320 for ( uInt j = idx ; j < scnrow-1 ; j++ ) { 1609 Double tsc1 = scTime Col( j ).get( "s" ).getValue() ;1610 Double dt1 = scInterval Col( j ).getValue("s");1611 Double tsc2 = scTime Col( j+1 ).get( "s" ).getValue() ;1612 Double dt2 = scInterval Col( j+1 ).getValue("s");1321 Double tsc1 = scTime[j].get( "s" ).getValue() ; 1322 Double dt1 = scInterval[j] ; 1323 Double tsc2 = scTime[j+1].get( "s" ).getValue() ; 1324 Double dt2 = scInterval[j+1] ; 1613 1325 if ( t > tsc1-half*dt1 && t <= tsc2-half*dt2 ) { 1614 1326 tstr[i] = tsc1 ; 1327 tidx[i] = j ; 1615 1328 idx = j ; 1616 1329 break ; … … 1618 1331 } 1619 1332 if ( tstr[i] == -1.0 ) { 1620 Double tsc = scTime Col( scnrow-1 ).get( "s" ).getValue() ;1621 Double dt = scInterval Col( scnrow-1 ).getValue( "s" );1622 if ( t <= tsc+0.5*dt ) 1333 Double tsc = scTime[scnrow-1].get( "s" ).getValue() ; 1334 Double dt = scInterval[scnrow-1] ; 1335 if ( t <= tsc+0.5*dt ) { 1623 1336 tstr[i] = tsc ; 1337 tidx[i] = scnrow-1 ; 1338 } 1624 1339 } 1625 1340 } 1626 1341 double endSec = gettimeofday_sec() ; 1627 os_ << "end MSFiller::getSysCalTime() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ; 1628 //return tstr ; 1342 os_ << "end MSFiller::getSysCalTime() endSec=" << endSec << " (" << endSec-startSec << "sec) scnrow = " << scnrow << "tcol.nelements() = " << tcol.nelements() << LogIO::POST ; 1629 1343 return ; 1630 1344 } 1631 1345 1632 uInt MSFiller::getTsys( uInt idx, Matrix<Float> &tsys, MSSysCal &tab, Double t ) 1633 { 1634 double startSec = gettimeofday_sec() ; 1635 os_ << "start MSFiller::getTsys() startSec=" << startSec << LogIO::POST ; 1636 uInt nrow = tab.nrow() ; 1637 if ( nrow == 0 ) { 1638 os_ << "No SysCal rows" << LogIO::POST ; 1639 tsys.resize( IPosition(0) ) ; 1640 return 0 ; 1641 } 1642 ROScalarMeasColumn<MEpoch> scTimeCol( tab, "TIME" ) ; 1643 ROArrayColumn<Float> mTsysCol ; 1644 mTsysCol.attach( tab, colTsys_ ) ; 1645 for ( uInt i = idx ; i < nrow ; i++ ) { 1646 Double tref = scTimeCol( i ).get( "s" ).getValue() ; 1647 if ( t == tref ) { 1648 tsys.reference( mTsysCol( i ) ) ; 1649 idx = i ; 1650 break ; 1651 } 1652 } 1653 //os_ << "MSFiller::getTsys() idx = " << idx << " tsys = " << tsys << LogIO::POST ; 1654 double endSec = gettimeofday_sec() ; 1655 os_ << "end MSFiller::getTsys() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ; 1656 return idx ; 1657 } 1658 1659 Vector<uInt> MSFiller::getTcalId( Int fid, Int spwid, Double t ) 1346 // uInt MSFiller::getTsys( uInt idx, Matrix<Float> &tsys, ROArrayColumn<Float> &mTsysCol, Block<MEpoch> &scTimeCol, Double t ) 1347 // { 1348 // double startSec = gettimeofday_sec() ; 1349 // os_ << "start MSFiller::getTsys() startSec=" << startSec << LogIO::POST ; 1350 // uInt nrow = mTsysCol.nrow() ; 1351 // if ( nrow == 0 ) { 1352 // os_ << "No SysCal rows" << LogIO::POST ; 1353 // tsys.resize( IPosition(0) ) ; 1354 // return 0 ; 1355 // } 1356 // for ( uInt i = idx ; i < nrow ; i++ ) { 1357 // Double tref = scTimeCol[i].get("s").getValue() ; 1358 // if ( t == tref ) { 1359 // tsys.reference( mTsysCol( i ) ) ; 1360 // idx = i ; 1361 // break ; 1362 // } 1363 // } 1364 // //os_ << "MSFiller::getTsys() idx = " << idx << " tsys = " << tsys << LogIO::POST ; 1365 // double endSec = gettimeofday_sec() ; 1366 // os_ << "end MSFiller::getTsys() endSec=" << endSec << " (" << endSec-startSec << "sec)" << LogIO::POST ; 1367 // return idx ; 1368 // } 1369 1370 Block<uInt> MSFiller::getTcalId( Int fid, Int spwid, Double t ) 1660 1371 { 1661 1372 double startSec = gettimeofday_sec() ; … … 1663 1374 if ( table_->tcal().table().nrow() == 0 ) { 1664 1375 os_ << "No TCAL rows" << LogIO::POST ; 1665 Vector<uInt> tcalids( 0 ) ;1376 Block<uInt> tcalids( 0 ) ; 1666 1377 return tcalids ; 1667 1378 } 1668 //String feed = "FEED" + String::toString(fid) ;1669 //String spw = "SPW" + String::toString(spwid) ;1670 1379 String sctime = MVTime( Quantum<Double>(t,"s") ).string(MVTime::YMD) ; 1671 //String key = feed + ":" + spw + ":" + sctime ;1672 1380 String key = keyTcal( fid, spwid, sctime ) ; 1673 //Record rec = tcalrec_.asRecord(feed).asRecord(spw) ;1674 //if ( !rec.isDefined( sctime ) ) {1675 1381 if ( !tcalrec_.isDefined( key ) ) { 1676 1382 os_ << "No TCAL rows" << LogIO::POST ; 1677 Vector<uInt> tcalids( 0 ) ;1383 Block<uInt> tcalids( 0 ) ; 1678 1384 return tcalids ; 1679 1385 } 1680 //Vector<uInt> ids = rec.asArrayuInt(sctime) ;1681 1386 Vector<uInt> ids = tcalrec_.asArrayuInt( key ) ; 1682 1387 uInt npol = ids[1] - ids[0] + 1 ; 1683 Vector<uInt> tcalids( npol ) ;1388 Block<uInt> tcalids( npol ) ; 1684 1389 tcalids[0] = ids[0] ; 1685 1390 tcalids[1] = ids[1] ; -
trunk/src/MSFiller.h
r1993 r2002 29 29 30 30 #include <casa/Containers/Record.h> 31 #include <casa/Containers/Block.h> 31 32 32 33 //#include <tables/Tables/TableColumn.h> … … 66 67 //void fillHistory() ; 67 68 //void fillFit() ; 68 void fillTcal( boost::object_pool<casa::ROTableColumn> *poolr, 69 boost::object_pool<casa::TableColumn> *poolw ) ; 69 void fillTcal( boost::object_pool<casa::ROTableColumn> *poolr ) ; 70 70 71 71 // get SRCTYPE from STATE_ID 72 72 casa::Int getSrcType( casa::Int stateId, boost::object_pool<casa::ROTableColumn> *pool ) ; 73 //casa::Int getSrcType( casa::Int stateId ) ;74 73 75 74 // get POLNO from CORR_TYPE 76 casa:: Vector<casa::uInt> getPolNo( casa::Int corrType ) ;75 casa::Block<casa::uInt> getPolNo( casa::Int corrType ) ; 77 76 78 77 // get poltype from CORR_TYPE … … 85 84 // assume that tab is selected by ANTENNA_ID, FEED_ID, SPECTRAL_WINDOW_ID 86 85 // and sorted by TIME 87 //casa::Vector<casa::Double> getSysCalTime( casa::MSSysCal &tab, casa::MEpoch::ROScalarColumn &tcol ) ; 88 void getSysCalTime( casa::MSSysCal &tab, casa::MEpoch::ROScalarColumn &tcal, casa::Vector<casa::Double> &scTime ) ; 86 void getSysCalTime( casa::Block<casa::MEpoch> &scTimeIn, casa::Vector<casa::Double> &scInterval, casa::Block<casa::MEpoch> &tcol, casa::Block<casa::Double> &scTimeOut, casa::Block<casa::Int> &tidx ) ; 89 87 90 88 // get tsys by time stamp 91 89 // assume that tab is selected by ANTENNA_ID, FEED_ID, SPECTRAL_WINDOW_ID 92 90 // and sorted by TIME 93 casa::uInt getTsys( casa::uInt idx, casa::Matrix<casa::Float> &tsys, casa::MSSysCal &tab, casa::Double t ) ;91 // casa::uInt getTsys( casa::uInt idx, casa::Matrix<casa::Float> &tsys, casa::ROArrayColumn<casa::Float> &tsysCol, casa::Block<casa::MEpoch> &scTime, casa::Double t ) ; 94 92 95 93 // get TCAL_ID 96 casa:: Vector<casa::uInt> getTcalId( casa::Int feedId, casa::Int spwId, casa::Double t ) ;94 casa::Block<casa::uInt> getTcalId( casa::Int feedId, casa::Int spwId, casa::Double t ) ; 97 95 98 96 // get direction for DIRECTION, AZIMUTH, and ELEVATION columns
Note:
See TracChangeset
for help on using the changeset viewer.