Changeset 2383


Ignore:
Timestamp:
12/22/11 17:05:22 (13 years ago)
Author:
Takeshi Nakazato
Message:

New Development: No

JIRA Issue: Yes CAS-2816

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

Re-defined row-by-row gridding such that each polarization component
is processed separately. After this change, sorting table rows is
no more necessary.

Location:
trunk/src
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/STGrid.cpp

    r2382 r2383  
    214214}
    215215
     216Bool STGrid::pastEnd()
     217{
     218  LogIO os( LogOrigin("STGrid","pastEnd",WHERE) ) ;
     219  Bool b = nprocessed_ >= nrow_ ;
     220  return b ;
     221}
     222
     223void STGrid::grid()
     224{
     225  // data selection
     226  selectData() ;
     227  setupArray() ;
     228
     229  if ( npol_ > 1 ) {
     230    LogIO os( LogOrigin("STGrid", "grid", WHERE) ) ;
     231    os << LogIO::WARN << "STGrid doesn't support assigning polarization-dependent weight. Use averaged weight over polarization." << LogIO::POST ;
     232  }
     233 
     234  if ( wtype_.compare("UNIFORM") != 0 &&
     235       wtype_.compare("TINT") != 0 &&
     236       wtype_.compare("TSYS") != 0 &&
     237       wtype_.compare("TINTSYS") != 0 ) {
     238    LogIO os( LogOrigin("STGrid", "grid", WHERE) ) ;
     239    os << LogIO::WARN << "Unsupported weight type '" << wtype_ << "', apply UNIFORM weight" << LogIO::POST ;
     240    wtype_ = "UNIFORM" ;
     241  }
     242
     243  Bool doAll = examine() ;
     244
     245  if ( doAll )
     246    gridPerPol() ;
     247  else
     248    gridPerRow() ;
     249}
     250
     251Bool STGrid::examine()
     252{
     253  // TODO: nchunk_ must be determined from nchan_, npol_, and (nx_,ny_)
     254  //       by considering data size to be allocated for ggridsd input/output
     255  nchunk_ = 100 ;
     256  Bool b = nchunk_ >= nrow_ ;
     257  nchunk_ = min( nchunk_, nrow_ ) ;
     258  return b ;
     259}
     260
    216261void STGrid::gridPerRow()
    217262{
     
    262307  IPosition dshape( 2, 2, nchunk_ ) ;
    263308  IPosition wshape( 2, nchan_, nchunk_ ) ;
    264   Array<Complex> spectra( cshape ) ;
    265   Array<Double> direction( dshape ) ;
    266   Array<Int> flagtra( cshape ) ;
     309  Array<Complex> spectra( wshape ) ;
     310  Array<Int> flagtra( wshape ) ;
    267311  Array<Int> rflag( vshape ) ;
    268312  Array<Float> weight( wshape ) ;
     313  Array<Double> direction( dshape ) ;
    269314  Array<Double> xypos( dshape ) ;
    270315
    271   spectraF_.resize( cshape ) ;
    272   flagtraUC_.resize( cshape ) ;
     316  spectraF_.resize( wshape ) ;
     317  flagtraUC_.resize( wshape ) ;
    273318  rflagUI_.resize( vshape ) ;
    274319
     
    282327    }
    283328  }
    284   Int *polMap = new Int[npol_] ;
    285   {
    286     Int *work_p = polMap ;
    287     for ( Int i = 0 ; i < npol_ ; i++ ) {
    288       *work_p = i ;
    289       work_p++ ;
    290     }
    291   }
     329  Int *polMap = new Int[1] ;
     330  Int nvispol = 1 ;
    292331
    293332  // for performance check
     
    297336  double eToInt = 0.0 ;
    298337
    299   // prepare pointer
    300   Float *conv_p = convFunc.getStorage( deleteConv ) ;
    301   Complex *gdata_p = gdataArrC.getStorage( deleteDataG ) ;
    302   Float *wdata_p = gwgtArr.getStorage( deleteWgtG ) ;
    303   const Complex *data_p = spectra.getStorage( deleteData ) ;
    304   const Float *wgt_p = weight.getStorage( deleteWgt ) ;
    305   const Int *flag_p = flagtra.getStorage( deleteFlag ) ;
    306   const Int *rflag_p = rflag.getStorage( deleteFlagR ) ;
    307   Double *xypos_p = xypos.getStorage( deletePos ) ;
    308 
    309   while( !pastEnd() ) {
    310     // retrieve data
    311     t0 = mathutil::gettimeofday_sec() ;
    312     Int nrow = getDataChunk( spectra, direction, flagtra, rflag, weight ) ;
    313     //os << "nrow = " << nrow << LogIO::POST ;
    314     t1 = mathutil::gettimeofday_sec() ;
    315     eGetDataChunk += t1-t0 ;
    316     eToInt += subtime_ ;
     338  for ( Int ipol = 0 ; ipol < npol_ ; ipol++ ) {
     339    initPol( ipol ) ;
    317340   
    318     // world -> pixel
    319     t0 = mathutil::gettimeofday_sec() ;
    320     toPixel( direction, xypos ) ;
    321     t1 = mathutil::gettimeofday_sec() ;
    322     eToPixel += t1-t0 ;
     341    polMap[0] = ipol ;
     342
     343    while( !pastEnd() ) {
     344      // retrieve data
     345      t0 = mathutil::gettimeofday_sec() ;
     346      Int nrow = getDataChunk( spectra, direction, flagtra, rflag, weight ) ;
     347      //os << "nrow = " << nrow << LogIO::POST ;
     348      t1 = mathutil::gettimeofday_sec() ;
     349      eGetDataChunk += t1-t0 ;
     350      eToInt += subtime_ ;
     351     
     352      // world -> pixel
     353      t0 = mathutil::gettimeofday_sec() ;
     354      toPixel( direction, xypos ) ;
     355      t1 = mathutil::gettimeofday_sec() ;
     356      eToPixel += t1-t0 ;
    323357   
    324     // call ggridsd
    325     irow = -1 ;
    326     t0 = mathutil::gettimeofday_sec() ;
    327     call_ggridsd( xypos_p,
    328                   data_p,
    329                   &npol_,
    330                   &nchan_,
    331                   flag_p,
    332                   rflag_p,
    333                   wgt_p,
    334                   &nrow,
    335                   &irow,
    336                   gdata_p,
    337                   wdata_p,
    338                   &gnx,
    339                   &gny,
    340                   &npol_,
    341                   &nchan_,
    342                   &convSupport_,
    343                   &convSampling_,
    344                   conv_p,
    345                   chanMap,
    346                   polMap ) ;
    347     t1 = mathutil::gettimeofday_sec() ;
    348     eGGridSD += t1-t0 ;
    349 
     358      // prepare pointer
     359      Float *conv_p = convFunc.getStorage( deleteConv ) ;
     360      Complex *gdata_p = gdataArrC.getStorage( deleteDataG ) ;
     361      Float *wdata_p = gwgtArr.getStorage( deleteWgtG ) ;
     362      const Complex *data_p = spectra.getStorage( deleteData ) ;
     363      const Float *wgt_p = weight.getStorage( deleteWgt ) ;
     364      const Int *flag_p = flagtra.getStorage( deleteFlag ) ;
     365      const Int *rflag_p = rflag.getStorage( deleteFlagR ) ;
     366      Double *xypos_p = xypos.getStorage( deletePos ) ;
     367
     368      // call ggridsd
     369      irow = -1 ;
     370      t0 = mathutil::gettimeofday_sec() ;
     371      call_ggridsd( xypos_p,
     372                    data_p,
     373                    //&npol_,
     374                    &nvispol,
     375                    &nchan_,
     376                    flag_p,
     377                    rflag_p,
     378                    wgt_p,
     379                    &nrow,
     380                    &irow,
     381                    gdata_p,
     382                    wdata_p,
     383                    &gnx,
     384                    &gny,
     385                    &npol_,
     386                    &nchan_,
     387                    &convSupport_,
     388                    &convSampling_,
     389                    conv_p,
     390                    chanMap,
     391                    polMap ) ;
     392      t1 = mathutil::gettimeofday_sec() ;
     393      eGGridSD += t1-t0 ;
     394
     395      // finalization
     396      convFunc.putStorage( conv_p, deleteConv ) ;
     397      gdataArrC.putStorage( gdata_p, deleteDataG ) ;
     398      gwgtArr.putStorage( wdata_p, deleteWgtG ) ;
     399      xypos.putStorage( xypos_p, deletePos ) ;
     400      spectra.freeStorage( data_p, deleteData ) ;
     401      weight.freeStorage( wgt_p, deleteWgt ) ;
     402      flagtra.freeStorage( flag_p, deleteFlag ) ;
     403      rflag.freeStorage( rflag_p, deleteFlagR ) ;
     404
     405    }
     406
     407    nprocessed_ = 0 ;
    350408  }
    351409  os << "getDataChunk: elapsed time is " << eGetDataChunk << " sec." << LogIO::POST ;
     
    353411  os << "ggridsd: elapsed time is " << eGGridSD << " sec." << LogIO::POST ;
    354412  os << "toInt: elapsed time is " << eToInt << " sec." << LogIO::POST ;
    355 
    356   // finalization
    357   convFunc.putStorage( conv_p, deleteConv ) ;
    358   gdataArrC.putStorage( gdata_p, deleteDataG ) ;
    359   gwgtArr.putStorage( wdata_p, deleteWgtG ) ;
    360   xypos.putStorage( xypos_p, deletePos ) ;
    361   spectra.freeStorage( data_p, deleteData ) ;
    362   weight.freeStorage( wgt_p, deleteWgt ) ;
    363   flagtra.freeStorage( flag_p, deleteFlag ) ;
    364   rflag.freeStorage( rflag_p, deleteFlagR ) ;
     413 
    365414  delete polMap ;
    366415  delete chanMap ;
     
    369418  setData( gdataArrC, gwgtArr ) ;
    370419
    371 }
    372 
    373 Bool STGrid::pastEnd()
    374 {
    375   LogIO os( LogOrigin("STGrid","pastEnd",WHERE) ) ;
    376   Bool b = nprocessed_ >= nrow_ ;
    377   return b ;
    378 }
    379 
    380 void STGrid::grid()
    381 {
    382   // data selection
    383   selectData() ;
    384   setupArray() ;
    385 
    386   if ( npol_ > 1 ) {
    387     LogIO os( LogOrigin("STGrid", "grid", WHERE) ) ;
    388     os << LogIO::WARN << "STGrid doesn't support assigning polarization-dependent weight. Use averaged weight over polarization." << LogIO::POST ;
    389   }
    390  
    391   if ( wtype_.compare("UNIFORM") != 0 &&
    392        wtype_.compare("TINT") != 0 &&
    393        wtype_.compare("TSYS") != 0 &&
    394        wtype_.compare("TINTSYS") != 0 ) {
    395     LogIO os( LogOrigin("STGrid", "grid", WHERE) ) ;
    396     os << LogIO::WARN << "Unsupported weight type '" << wtype_ << "', apply UNIFORM weight" << LogIO::POST ;
    397     wtype_ = "UNIFORM" ;
    398   }
    399 
    400   Bool doAll = examine() ;
    401 
    402   if ( doAll )
    403     gridPerPol() ;
    404   else {
    405     sortData() ;
    406     gridPerRow() ;
    407   }
    408 }
    409 
    410 Bool STGrid::examine()
    411 {
    412   // TODO: nchunk_ must be determined from nchan_, npol_, and (nx_,ny_)
    413   //       by considering data size to be allocated for ggridsd input/output
    414   nchunk_ = 10000 ;
    415   Bool b = nchunk_ >= nrow_ ;
    416   nchunk_ = min( nchunk_, nrow_ ) ;
    417   return b ;
    418420}
    419421
     
    459461  // to read data from the table
    460462  IPosition mshape( 2, nchan_, nrow_ ) ;
     463  IPosition dshape( 2, 2, nrow_ ) ;
    461464  Array<Complex> spectra( mshape ) ;
    462   Array<Double> direction( IPosition(2,2,nrow_) ) ;
     465  Array<Double> direction( dshape ) ;
    463466  Array<Int> flagtra( mshape ) ;
    464467  Array<Int> rflag( IPosition(1,nrow_) ) ;
    465468  Array<Float> weight( mshape ) ;
     469  Array<Double> xypos( dshape ) ;
    466470
    467471  // maps
     
    477481  Int nvispol = 1 ;
    478482 
     483  // for performance check
     484  double eInitPol = 0.0 ;
     485  double eGetData = 0.0 ;
     486  double eToPixel = 0.0 ;
     487  double eGGridSD = 0.0 ;
     488  double eToInt = 0.0 ;
     489
    479490  for ( Int ipol = 0 ; ipol < npol_ ; ipol++ ) {
    480491    t0 = mathutil::gettimeofday_sec() ;
    481492    initPol( ipol ) ;
    482493    t1 = mathutil::gettimeofday_sec() ;
    483     os << "initPol: elapsed time is " << t1-t0 << " sec." << LogIO::POST ;
     494    eInitPol += t1-t0 ;
    484495
    485496    // retrieve data
     
    487498    getDataPerPol( spectra, direction, flagtra, rflag, weight ) ;
    488499    t1 = mathutil::gettimeofday_sec() ;
    489     os << "getData: elapsed time is " << t1-t0 << " sec." << LogIO::POST ;
     500    eGetData += t1-t0 ;
     501   
    490502    IPosition sshape = spectra.shape() ;
    491503   
    492504    // world -> pixel
    493     Array<Double> xypos( direction.shape(), 0.0 ) ;
    494505    t0 = mathutil::gettimeofday_sec() ;
    495506    toPixel( direction, xypos ) ; 
     
    543554  }
    544555
     556  os << "initPol: elapsed time is " << eInitPol << " sec." << LogIO::POST ;
     557  os << "getData: elapsed time is " << eGetData << " sec." << LogIO::POST ;
    545558  // delete maps
    546559  delete polMap ;
     
    734747    os << LogIO::EXCEPTION ;
    735748  }
    736   attach( tab_ ) ;
     749  //attach( tab_ ) ;
    737750}
    738751
     
    770783{
    771784  LogIO os( LogOrigin("STGrid","getDataPerPol",WHERE) ) ;
    772   Vector<uInt> rflagVec( rflag ) ;
     785
    773786  // 2011/12/22 TN
    774787  // weight and tsys shares its storage
    775788  Array<Float> tsys( weight ) ;
    776789  Array<Double> tint( rflag.shape() ) ;
     790
     791  Vector<uInt> rflagVec( rflag ) ;
    777792  Vector<Double> tintVec( tint ) ;
    778793
     
    790805  }
    791806
    792   getWeightPerPol( weight, tsys, tint ) ;
     807  getWeight( weight, tsys, tint ) ;
    793808}
    794809
     
    812827  toInt( rflagUI, rflag ) ;
    813828  t1 = mathutil::gettimeofday_sec() ;
    814   os << "toInt: elapsed time is " << t1-t0 << " sec." << LogIO::POST ;
     829  //os << "toInt: elapsed time is " << t1-t0 << " sec." << LogIO::POST ;
    815830}
    816831
     
    822837{
    823838  Int nrow = getDataChunk( spectraF_, direction, flagtraUC_, rflagUI_, weight ) ;
     839  if ( nrow < nchunk_ ) {
     840    spectra.resize( spectraF_.shape() ) ;
     841    flagtra.resize( flagtraUC_.shape() ) ;
     842    rflag.resize( rflagUI_.shape() ) ;
     843  }
    824844  double t0, t1 ;
    825845  t0 = mathutil::gettimeofday_sec() ;
     
    840860{
    841861  LogIO os( LogOrigin("STGrid","getDataChunk",WHERE) ) ;
    842   Int nrow = min( spectra.shape()[2], nrow_-nprocessed_ ) ;
    843   Array<Float> tsys( spectra.shape() ) ;
    844   Array<Double> tint( IPosition(2,spectra.shape()[0],spectra.shape()[2]) ) ;
    845   IPosition m( 2, 0, 2 ) ;
    846   IPosition v( 1, 0 ) ;
    847   ArrayIterator<Float> spi( spectra, m, False ) ;
    848   ArrayIterator<uChar> fli( flagtra, m, False ) ;
    849   ArrayIterator<Float> tsi( tsys, m, False ) ;
    850   ArrayIterator<Double> di( direction, v ) ;
    851   Bool bfr, bti ;
    852   uInt *fr_p = rflag.getStorage( bfr ) ;
    853   Double *ti_p = tint.getStorage( bti ) ;
    854   uInt *wfr_p = fr_p ;
    855   Double *wti_p = ti_p ;
    856   Int offset = nprocessed_ * npol_ ;
    857   Int start = offset ;
    858   Int end = start + nrow * npol_ ;
    859   for ( Int irow = 0 ; irow < npol_ * nrow ; irow++ ) {
    860     Array<Float> spSlice = spi.array() ;
    861     Array<uChar> flSlice = fli.array() ;
    862     Array<Float> tsSlice = tsi.array() ;
    863 
    864     uInt idx = rows_[offset+irow] ;
    865     spectraCol_.get( idx, spSlice ) ;
    866 
    867     flagtraCol_.get( idx, flSlice ) ;
    868     Vector<Float> tmpF = tsysCol_( idx ) ;
    869     if ( tmpF.nelements() == (uInt)nchan_ ) {
    870       tsSlice = tmpF ;
    871     }
    872     else {
    873       tsSlice = tmpF[0] ;
    874     }
    875     *wti_p = intervalCol_( idx ) ;
    876 
    877     spi.next() ;
    878     fli.next() ;
    879     tsi.next() ;
    880     wti_p++ ;
    881   }
    882   tint.putStorage( ti_p, bti ) ;
    883   for ( Int irow = 0 ; irow < nrow ; irow += npol_ ) {
    884     uInt idx = rows_[offset+irow] ;
    885     directionCol_.get( idx, di.array() ) ;
    886     di.next() ;
    887 
    888     *wfr_p = flagRowCol_( idx ) ;
    889     for ( Int ipol = 1 ; ipol < npol_ ; ipol++ ) {
    890       *wfr_p = max( *wfr_p, flagRowCol_( rows_[offset+irow+ipol] ) ) ;
    891     }
    892     wfr_p++ ;
    893   }
    894   rflag.putStorage( fr_p, bfr ) ;
    895 //   Bool bti, bfr ;
    896 //   Double *ti_p = tint.getStorage( bti ) ;
    897 //   Double *wti_p = ti_p ;
    898 //   uInt *fr_p = rflag.getStorage( bfr ) ;
    899 //   uInt *wfr_p = fr_p - 1 ;
    900 //   Int start = nprocessed_ * npol_ ;
    901 //   Int end = start + nrow * npol_ ;
    902 //   for ( Int irow = start ; irow < end ; irow++ ) {
    903 //     uInt idx = rows_[irow] ;
    904 //     row_.get( idx ) ;
    905 //     // SPECTRA
    906 // //     os << "spi.array().shape()=" << spi.array().shape() << LogIO::POST ;
    907 // //     os << "(*spectraRF_).shape()=" << (*spectraRF_).shape() << LogIO::POST ;
    908 //     spi.array() = *spectraRF_ ;
    909    
    910 //     // FLAGTRA
    911 //     fli.array() = *flagtraRF_ ;
    912 
    913 //     // INTERVAL
    914 //     *wti_p = *intervalRF_ ;
    915 
    916 //     // TSYS
    917 // //     os << "tsi.array().shape()=" << tsi.array().shape() << LogIO::POST ;
    918 // //     os << "(*tsysRF_).shape()=" << (*tsysRF_).shape() << LogIO::POST ;
    919 // //     os << "(*tsysRF_).nelements()=" << (*tsysRF_).nelements() << LogIO::POST ;
    920 //     if ( (*tsysRF_).nelements() == (uInt)nchan_ )
    921 //       tsi.array() = *tsysRF_ ;
    922 //     else
    923 //       tsi.array() = *((*tsysRF_).data()) ;
    924 
    925 //     // DIRECTION and FLAGROW
    926 //     Int mod = irow % npol_ ;
    927 //     if ( mod == 0 ) {
    928 // //       os << "di.array().shape()=" << di.array().shape() << LogIO::POST ;
    929 // //       os << "(*directionRF_).shape()=" << (*directionRF_).shape() << LogIO::POST ;
    930 //       di.array() = *directionRF_ ;
    931 //       wfr_p++ ;
    932 //       *wfr_p = *flagRowRF_ ;
    933 //     }
    934 //     else {
    935 //       *wfr_p += *flagRowRF_ ;
    936 //     }
    937    
    938 //     // increment
    939 //     spi.next() ;
    940 //     fli.next() ;
    941 //     tsi.next() ;
    942 //     di.next() ;
    943 //     *wti_p++ ;
    944 //   }
    945 //   tint.putStorage( ti_p, bti ) ;
    946 //   rflag.putStorage( fr_p, bfr ) ;
    947 
     862  Int nrow = spectra.shape()[1] ;
     863  Int remainingRow = nrow_ - nprocessed_ ;
     864  if ( remainingRow < nrow ) {
     865    nrow = remainingRow ;
     866    IPosition mshape( 2, nchan_, nrow ) ;
     867    IPosition vshape( 1, nrow ) ;
     868    spectra.resize( mshape ) ;
     869    flagtra.resize( mshape ) ;
     870    direction.resize( IPosition(2,2,nrow) ) ;
     871    rflag.resize( vshape ) ;
     872    weight.resize( mshape ) ;
     873  }
     874  // 2011/12/22 TN
     875  // tsys shares its storage with weight
     876  Array<Float> tsys( weight ) ;
     877  Array<Double> tint( rflag.shape() ) ;
     878
     879  Vector<uInt> rflagVec( rflag ) ;
     880  Vector<Double> tintVec( tint ) ;
     881
     882  RefRows rows( nprocessed_, nprocessed_+nrow-1, 1 ) ;
     883  os << "rows.nrows()=" << rows.nrows() << LogIO::POST ;
     884  spectraCol_.getColumnCells( rows, spectra ) ;
     885  flagtraCol_.getColumnCells( rows, flagtra ) ;
     886  directionCol_.getColumnCells( rows, direction ) ;
     887  flagRowCol_.getColumnCells( rows, rflagVec ) ;
     888  intervalCol_.getColumnCells( rows, tintVec ) ;
     889  Vector<Float> tsysTemp = tsysCol_( nprocessed_ ) ;
     890  if ( tsysTemp.nelements() == (uInt)nchan_ )
     891    tsysCol_.getColumnCells( rows, tsys ) ;
     892  else
     893    tsys = tsysTemp[0] ;
     894 
    948895  getWeight( weight, tsys, tint ) ;
    949896
     
    998945
    999946void STGrid::getWeight( Array<Float> &w,
    1000                         Array<Float> &tsys,
    1001                         Array<Double> &tint )
    1002 {
    1003   LogIO os( LogOrigin("STGrid","getWeight",WHERE) ) ;
    1004 //   double t0, t1 ;
    1005 //   t0 = mathutil::gettimeofday_sec() ;
    1006   // resize
    1007   IPosition refShape = tsys.shape() ;
    1008   Int nchan = refShape[1] ;
    1009   Int nrow = refShape[2] ;
    1010 
    1011   // set weight
    1012   if ( wtype_.compare( "UNIFORM" ) == 0 ) {
    1013     w = 1.0 ;
    1014   }
    1015   else if ( wtype_.compare( "TINT" ) == 0 ) {
    1016     Bool b0, b1 ;
    1017     Float *w_p = w.getStorage( b0 ) ;
    1018     Float *w0_p = w_p ;
    1019     const Double *ti_p = tint.getStorage( b1 ) ;
    1020     const Double *w1_p = ti_p ;
    1021     for ( Int irow = 0 ; irow < nrow ; irow++ ) {
    1022       Float val = (Float)(polMean( w1_p )) ;
    1023       for ( Int ichan = 0 ; ichan < nchan ; ichan++ ) {
    1024         *w0_p = val ;
    1025         w0_p++ ;
    1026       }
    1027     }
    1028     w.putStorage( w_p, b0 ) ;
    1029     tint.freeStorage( ti_p, b1 ) ;
    1030   }
    1031   else if ( wtype_.compare( "TSYS" ) == 0 ) {
    1032     Bool b0, b1 ;
    1033     Float *w_p = w.getStorage( b0 ) ;
    1034     Float *w0_p = w_p ;
    1035     const Float *ts_p = tsys.getStorage( b1 ) ;
    1036     const Float *w1_p = ts_p ;
    1037     for ( Int irow = 0 ; irow < nrow ; irow++ ) {
    1038       for ( Int ichan = 0 ; ichan < nchan ; ichan++ ) {
    1039         Float val = polMean( w1_p ) ;
    1040         *w0_p = 1.0 / ( val * val ) ;
    1041         w0_p++ ;
    1042       }
    1043     }
    1044     w.putStorage( w_p, b0 ) ;
    1045     tsys.freeStorage( ts_p, b1 ) ;
    1046   }
    1047   else if ( wtype_.compare( "TINTSYS" ) == 0 ) {
    1048     Bool b0, b1, b2 ;
    1049     Float *w_p = w.getStorage( b0 ) ;
    1050     Float *w0_p = w_p ;
    1051     const Double *ti_p = tint.getStorage( b1 ) ;
    1052     const Double *w1_p = ti_p ;
    1053     const Float *ts_p = tsys.getStorage( b2 ) ;
    1054     const Float *w2_p = ts_p ;
    1055     for ( Int irow = 0 ; irow < nrow ; irow++ ) {
    1056       Float interval = (Float)(polMean( w1_p )) ;
    1057       for ( Int ichan = 0 ; ichan < nchan ; ichan++ ) {
    1058         Float temp = polMean( w2_p ) ;
    1059         *w0_p = interval / ( temp * temp ) ;
    1060         w0_p++ ;
    1061       }
    1062     }
    1063     w.putStorage( w_p, b0 ) ;
    1064     tint.freeStorage( ti_p, b1 ) ;
    1065     tsys.freeStorage( ts_p, b2 ) ;
    1066   }
    1067   else {
    1068     //LogIO os( LogOrigin("STGrid", "getWeight", WHERE) ) ;
    1069     //os << LogIO::WARN << "Unsupported weight type '" << wtype_ << "', apply UNIFORM weight" << LogIO::POST ;
    1070     w = 1.0 ;
    1071   }
    1072 
    1073 //   t1 = mathutil::gettimeofday_sec() ;
    1074 //   os << "getWeight: elapsed time is " << t1-t0 << " sec" << LogIO::POST ;
    1075 }
    1076 
    1077 void STGrid::getWeightPerPol( Array<Float> &w,
    1078947                              Array<Float> &tsys,
    1079948                              Array<Double> &tint )
    1080949{
    1081   LogIO os( LogOrigin("STGrid","getWeightPerPol",WHERE) ) ;
    1082   //os << "start getWeightPerPol" << LogIO::POST ;
     950  LogIO os( LogOrigin("STGrid","getWeight",WHERE) ) ;
     951  //os << "start getWeight" << LogIO::POST ;
    1083952  double t0, t1 ;
    1084953  t0 = mathutil::gettimeofday_sec() ;
     
    11511020
    11521021  t1 = mathutil::gettimeofday_sec() ;
    1153   os << "getWeight: elapsed time is " << t1-t0 << " sec" << LogIO::POST ;
    1154 }
    1155 
    1156 Float STGrid::polMean( const Float *p )
    1157 {
    1158   Float v = 0.0 ;
    1159   for ( Int i = 0 ; i < npol_ ; i++ ) {
    1160     v += *p ;
    1161     p++ ;
    1162   }
    1163   v /= npol_ ;
    1164   return v ;
    1165 }
    1166 
    1167 Double STGrid::polMean( const Double *p )
    1168 {
    1169   Double v = 0.0 ;
    1170   for ( Int i = 0 ; i < npol_ ; i++ ) {
    1171     v += *p ;
    1172     p++ ;
    1173   }
    1174   v /= npol_ ;
    1175   return v ;
     1022  //os << "getWeight: elapsed time is " << t1-t0 << " sec" << LogIO::POST ;
     1023  subtime_ += t1-t0 ;
    11761024}
    11771025
     
    14091257}
    14101258
    1411 void STGrid::sortData()
    1412 {
    1413   LogIO os( LogOrigin("STGRid","sortData",WHERE) ) ;
    1414   double t0, t1 ;
    1415   t0 = mathutil::gettimeofday_sec() ;
    1416   rows_.resize( npol_*nrow_ ) ;
    1417   Table tab = tab_( tab_.col("POLNO").in(pollist_) ) ;
    1418   Block<String> cols( 2 ) ;
    1419   cols[0] = "TIME" ;
    1420   cols[1] = "BEAMNO" ;
    1421   TableIterator iter( tab, cols ) ;
    1422   uInt idx = 0 ;
    1423   while( !iter.pastEnd() ) {
    1424     Table t = iter.table().sort( "POLNO" ) ;
    1425     Vector<uInt> rows = t.rowNumbers() ;
    1426     //os << "rows=" << rows << LogIO::POST ;
    1427     for ( uInt i = 0 ; i < rows.nelements() ; i++ ) {
    1428       rows_[idx] = rows[i] ;
    1429       idx++ ;
    1430     }
    1431     iter.next() ;
    1432   }
    1433   t1 = mathutil::gettimeofday_sec() ;
    1434   os << "sortData: elapsed time is " << t1-t0 << " sec" << LogIO::POST ;
    1435   //os << "rows_ = " << rows_ << LogIO::POST ;
    1436 }
    1437 }
     1259}
  • trunk/src/STGrid.h

    r2382 r2383  
    109109                  Array<Float> &tsys,
    110110                  Array<Double> &tint ) ;
    111   void getWeightPerPol( Array<Float> &w,
    112                         Array<Float> &tsys,
    113                         Array<Double> &tint ) ;
    114111 
    115112  void toInt( Array<uChar> &u, Array<Int> &v ) ;
     
    124121  void setConvFunc( Vector<Float> &convFunc ) ;
    125122
    126   Float polMean( const Float *p ) ;
    127   Double polMean( const Double *p ) ;
    128 
    129123  void prepareTable( Table &tab, String &name ) ;
    130124
     
    133127  void selectData() ;
    134128  void setupArray() ;
    135   void sortData() ;
    136129
    137130  Bool examine() ;
     
    189182  ROArrayColumn<Float> tsysCol_ ;
    190183  ROScalarColumn<Double> intervalCol_ ;
    191 //   ROTableRow row_ ;
    192 //   RORecordFieldPtr< Array<Float> > spectraRF_ ;
    193 //   RORecordFieldPtr< Array<uChar> > flagtraRF_ ;
    194 //   RORecordFieldPtr< Array<Double> > directionRF_ ;
    195 //   RORecordFieldPtr<uInt> flagRowRF_ ;
    196 //   RORecordFieldPtr< Array<Float> > tsysRF_ ;
    197 //   RORecordFieldPtr<Double> intervalRF_ ;
    198184  Int nprocessed_ ;
    199185  Vector<uInt> rows_ ;
Note: See TracChangeset for help on using the changeset viewer.