Changeset 2550 for branches


Ignore:
Timestamp:
05/29/12 11:24:15 (13 years ago)
Author:
Takeshi Nakazato
Message:

New Development: No

JIRA Issue: No

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

Clean up code.


Location:
branches/hpc33/src
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • branches/hpc33/src/MSFiller.cpp

    r2549 r2550  
    304304    // set dummy epoch
    305305    mf.set( currentTime ) ;
    306 
    307     // initialize dirType
    308     //dirType = MDirection::N_Types ;
    309306
    310307    //
     
    957954                    pointingDirection.xyPlane( idx-1 ), pointingDirection.xyPlane( idx ) ) ;
    958955    }
    959     //mf.resetEpoch( currentTime ) ;
    960956    mf.set( currentTime ) ;
    961957    Quantum< Vector<Double> > tmp( d.column( 0 ), Unit( "rad" ) ) ;
    962958    if ( dirType != MDirection::J2000 ) {
    963       //MDirection::Convert toj2000( dirType, MDirection::Ref( MDirection::J2000, mf ) ) ;
    964959      dir = toj2000( tmp ).getAngle( "rad" ).getValue() ;
    965960    }
     
    968963    }
    969964    if ( dirType != MDirection::AZELGEO ) {
    970       //MDirection::Convert toazel( dirType, MDirection::Ref( MDirection::AZELGEO, mf ) ) ;
    971       //MDirection::Convert toazel( dirType, MDirection::Ref( MDirection::AZEL, mf ) ) ;
    972965      azel = toazel( tmp ).getAngle( "rad" ).getValue() ;
    973966    }
     
    981974  {
    982975    dir = sourceDir.getAngle( "rad" ).getValue() ;
    983     mf.resetEpoch( currentTime ) ;
    984     //MDirection::Convert toazel( dirType, MDirection::Ref( MDirection::AZELGEO, mf ) ) ;
     976    mf.set( currentTime ) ;
    985977    azel = toazel( Quantum< Vector<Double> >( dir, Unit("rad") ) ).getAngle( "rad" ).getValue() ;
    986978    if ( dirType != MDirection::J2000 ) {
    987       //MDirection::Convert toj2000( dirType, MDirection::Ref( MDirection::J2000, mf ) ) ;
    988979      dir = toj2000( Quantum< Vector<Double> >( dir, Unit("rad") ) ).getAngle( "rad" ).getValue() ;
    989980    }
  • branches/hpc33/src/STIdxIter.cpp

    r2548 r2550  
    1616  prod_m.resize( nfield_m ) ;
    1717  idx_m.resize( nfield_m ) ;
    18 //   prod_m[0] = 1 ;
    19 //   idx_m[0] = 0 ;
    20 //   for ( uInt i = 1 ; i < nfield_m ; i++ ) {
    21 //     prod_m[i] = shape[i-1] * prod_m[i-1] ;
    22 //     idx_m[i] = 0 ;
    23 //   }
    24 //   maxiter_m = prod_m[nfield_m-1] * shape[nfield_m-1] ;
    2518  prod_m[nfield_m-1] = 1 ;
    2619  idx_m[nfield_m-1] = 0 ;
     
    4235  niter_m++ ;
    4336  uInt x = niter_m ;
    44 //   for ( Int i = nfield_m-1 ; i >= 0 ; i-- ) {
    4537  for ( Int i = 0 ; i < nfield_m ; i++ ) {
    4638    idx_m[i] = x / prod_m[i] ;
     
    136128  // ...
    137129  storage_m.resize( arr_m.nelements()+nrow_m ) ;
    138 //   uInt *p = storage_m.storage() + arr_m.nelements() ;
     130  // initialize constant temporary storage
    139131  uInt *p = storage_m.storage() ;
    140132  for ( uInt i = 0 ; i < nrow_m ; i++ ) {
     
    142134    p++ ;
    143135  }
     136  // len_m layout
     137  // len[0]: length of temporary storage that incidates whole rows (constant)
     138  // len[1]: number of matched row for column 0 selection
     139  // len[2]: number of matched row for column 1 selection
     140  // ...
    144141  len_m.resize( ncol_m+1 ) ;
    145142  p = len_m.storage() ;
     
    149146//     cout << "len[" << i << "]=" << len_m[i] << endl ;
    150147  }
     148  // skip_m layout
     149  // skip_m[0]: True if column 0 is filled by unique value
     150  // skip_m[1]: True if column 1 is filled by unique value
     151  // ...
     152  skip_m.resize( ncol_m ) ;
     153  for ( uInt i = 0 ; i < ncol_m ; i++ ) {
     154    skip_m[i] = (Bool)(idxlist_m[i].size()==1) ;
     155//     cout << "skip_m[" << i << "]=" << skip_m[i] << endl ;
     156  }
    151157  prev_m = iter_m->current() ;
    152158  p = prev_m.storage() ;
     
    156162  }
    157163//   cout << "prev_m=" << Vector<uInt>(IPosition(1,ncol_m),prev_m.storage()) << endl ;
    158   skip_m.resize( ncol_m ) ;
    159   for ( uInt i = 0 ; i < ncol_m ; i++ ) {
    160     skip_m[i] = (Bool)(idxlist_m[i].size()==1) ;
    161 //     cout << "skip_m[" << i << "]=" << skip_m[i] << endl ;
    162   }
    163164}
    164165
     
    175176//   cout << "v=" << Vector<uInt>(IPosition(1,v.nelements()),v.storage(),SHARE) << endl ;
    176177//   cout << "c=" << c << endl ;
    177 //   if ( c < 0 ) {
    178178  if ( c > ncol_m-1 ) {
    179 //     pos_m[0] = len_m[0] ;
    180 //     uInt offset = 0 ;
    181 //     while( offset < ncol_m && skip_m[offset] )
    182 //       offset++ ;
    183179    pos_m[0] = len_m[ncol_m] ;
    184180    Int offset = ncol_m - 1 ;
     
    189185    return Vector<uInt>( pos_m, storage_m.storage()+offset*nrow_m, policy ) ;
    190186  }
    191 //   Int offset = c + 1 ;
    192 //   while( offset < ncol_m && skip_m[offset] )
    193 //     offset++ ;
    194187  Int offset = c - 1 ;
    195188  while( offset >= 0 && skip_m[offset] )
     
    200193//   cout << "len_m[c+1]=" << len_m[c+1] << endl ;
    201194//   cout << "base=" << Vector<uInt>(IPosition(1,len_m[c+1]),base,SHARE) << endl ;
    202 //   for ( Int i = c ; i >= 0 ; i-- ) {
    203 //     base = updateStorage( i, base, idxlist_m[i][v[i]] ) ;
    204   for ( Int i = c+1 ; i <= ncol_m ; i++ ) {
    205     base = updateStorage( i, base, idxlist_m[i-1][v[i-1]] ) ;
     195  for ( Int i = c ; i < ncol_m ; i++ ) {
     196    base = updateStorage( i, base, idxlist_m[i][v[i]] ) ;
    206197//     cout << "len_m[" << i << "]=" << len_m[i] << endl ;
    207198//     cout << "base=" << Vector<uInt>(IPosition(1,len_m[i]),base,SHARE) << endl ;
    208199  }
    209 //   pos_m[0] = len_m[0] ;
    210200  pos_m[0] = len_m[ncol_m] ;
    211201//   cout << "pos_m=" << pos_m << endl ;
     
    216206Int ArrayIndexIteratorAcc::isChanged( Block<uInt> &idx )
    217207{
    218 //   Int i = ncol_m - 1 ;
    219 //   while( i >= 0 && idx[i] == prev_m[i] )
    220 //     i-- ;
    221208  Int i = 0 ;
    222209  while( i < ncol_m && idx[i] == prev_m[i] )
     
    229216                                            uInt &v )
    230217{
     218  // CAUTION:
     219  // indexes for storage_m and len_m differ from index for skip_m by 1
     220  // (skip_m[0] corresponds to storage_m[1] and len_m[1]) since there
     221  // is additional temporary storage at first segment in storage_m
    231222  uInt *p ;
    232 //   if ( skip_m[icol] ) {
    233   if ( skip_m[icol-1] ) {
     223  if ( skip_m[icol] ) {
    234224    // skip update, just update len_m[icol] and pass appropriate pointer
    235225//     cout << "skip " << icol << endl ;
    236226    p = base ;
    237 //     len_m[icol] = len_m[icol+1] ;
    238     len_m[icol] = len_m[icol-1] ;
     227    len_m[icol+1] = len_m[icol] ;
    239228  }
    240229  else {
    241230//     cout << "update " << icol << endl ;
    242231    uInt len = 0 ;
    243     p = storage_m.storage() + icol * nrow_m ;
     232    p = storage_m.storage() + (icol+1) * nrow_m ;
    244233    uInt *work = p ;
    245234    Bool b ;
     
    247236    long offset = 0 ;
    248237    // warr_p points a first element of (icol)-th column
    249 //     const uInt *warr_p = arr_p + icol * nrow_m ;
    250     const uInt *warr_p = arr_p + (icol-1) * nrow_m ;
    251 //     for ( uInt i = 0 ; i < len_m[icol+1] ; i++ ) {
    252     for ( uInt i = 0 ; i < len_m[icol-1] ; i++ ) {
     238    const uInt *warr_p = arr_p + icol * nrow_m ;
     239    for ( uInt i = 0 ; i < len_m[icol] ; i++ ) {
    253240      // increment warr_p by (*(base)-*(base-1))
    254241      warr_p += *base - offset ;
     
    267254    }
    268255    arr_m.freeStorage( arr_p, b ) ;
    269     len_m[icol] = len ;
     256    len_m[icol+1] = len ;
    270257  }
    271258  return p ;
     
    464451{
    465452  convertArray( out, in ) ;
    466   //uInt len = in.nelements() ;
    467   //Vector<Int> tmp = in.copy() ;
    468   //uInt n = genSort( tmp, Sort::Ascending, Sort::QuickSort|Sort::NoDuplicates ) ;
    469   //val.resize( n ) ;
    470   //for ( uInt i = 0 ; i < n ; i++ ) {
    471   //  val[i] = tmp[i] ;
    472 //     cout << "val[" << i << "]=" << val[i] << endl ;
    473   //}
    474   //out.resize( len ) ;
    475   //uInt idx = -1 ;
    476   //cout << "idx=" << idx << endl ;
    477   //int idx2 = (Int)idx ;
    478   //cout << "idx2=" << idx2 << endl ;
    479   //idx = (uInt) idx2 ;
    480   //cout << "idx=" << idx << endl ;
    481   //if ( n == 1 ) {
    482   //  cout << "n=1" << endl ;
    483   //  out = 0 ;
    484   //}
    485   //else if ( n == 2 ) {
    486   //  cout << "n=2" << endl ;
    487   //  for ( uInt i = 0 ; i < len ; i++ ) {
    488   //    out[i] = (in[i] == val[0]) ? 0 : 1 ;
    489   //  }
    490   //}
    491   //else {
    492   //  cout << "n=" << n << endl ;
    493   //  map<Int,uInt> m ;
    494   // for( uInt i = 0 ; i < n ; i++ )
    495   //    m[val[i]] = i ;
    496   // for ( uInt i = 0 ; i < len ; i++ ) {
    497      //for ( uInt j = 0 ; j < n ; j++ ) {
    498       //  if ( in[i] == val[j] ) {
    499       //    out[i] = j ;
    500       //    break ;
    501       //  }
    502       //}
    503   //    out[i] = m[in[i]] ;
    504   //  }
    505   // }
    506453}
    507454
     
    518465//     cout << "val[" << i << "]=" << val[i] << endl ;
    519466  }
    520   //out.resize( len ) ;
    521467  if ( n == 1 ) {
    522468    //cout << "n=1" << endl ;
     
    535481      m[val[i]] = i ;
    536482    for ( uInt i = 0 ; i < len ; i++ ) {
    537       //for ( uInt j = 0 ; j < n ; j++ ) {
    538       //  if ( in[i] == val[j] ) {
    539       //    out[i] = j ;
    540       //    break ;
    541       //  }
    542       //}
    543483      out[i] = m[in[i]] ;
    544484    }
     
    548488Int STIdxIterExAcc::getSrcType()
    549489{
    550   //if ( srctype_m.nelements() > 0 )
    551490  if ( srctypeid_m >= 0 )
    552     //return srctype_m[iter_m->current()[srctypeid_m]] ;
    553491    return (Int)(iter_m->current()[srctypeid_m]) ;
    554492  else
  • branches/hpc33/src/STMath.cpp

    r2546 r2550  
    184184    tout.addRow();
    185185//     t4 = mathutil::gettimeofday_sec() ;
    186     //TableCopy::copyRows(tout, subt, outrowCount, 0, 1);
    187     //TableCopy::copyRows(tout, baset, outrowCount, rows[0], 1);
    188186    // skip to copy SPECTRA, FLAGTRA, and TSYS since those heavy columns are
    189187    // overwritten in the following process
     
    235233//     outrowCount += rowNum ;
    236234
    237 // merge loop
     235    // merge loop
    238236    uInt i = outrowCount ;
    239 //  for (uInt i=0; i < tout.nrow(); ++i) {
    240 
    241     // in[0] is already selected by TableItertor
     237    // in[0] is already selected by iterator
    242238    specCol.attach(baset,"SPECTRA");
    243239    flagCol.attach(baset,"FLAGTRA");
     
    248244    Vector<uChar> flag;
    249245    Double inter,time;
    250 //     for (uInt k = 0; k < subt.nrow(); ++k ) {
     246
    251247    for (uInt l = 0; l < rows.nelements(); ++l ) {
    252248      uInt k = rows[l] ;
     
    273269    // in[0] is already selected by TableIterator so that index is
    274270    // started from 1
    275     //for ( int j=0; j < int(in.size()); ++j ) {
    276271    for ( int j=1; j < int(in.size()); ++j ) {
    277272      const Table& tin = in[j]->table();
     
    362357    const Vector<Bool>& msk = acc.getMask();
    363358    if ( allEQ(msk, False) ) {
    364       //uint n = rowstodelete.nelements();
    365       //rowstodelete.resize(n+1, True);
    366       //rowstodelete[n] = i;
    367359      rowstodelB[nrowdel] = i ;
    368360      nrowdel++ ;
     
    408400  }
    409401
    410   //if (rowstodelete.nelements() > 0) {
    411402  if ( nrowdel > 0 ) {
    412403    Vector<uInt> rowstodelete( IPosition(1,nrowdel), rowstodelB.storage(), SHARE ) ;
     
    53155306}
    53165307
    5317 void STMath::calibrateALMA( CountedPtr<Scantable>& on,
    5318                             CountedPtr<Scantable>& off )
    5319 {
    5320   if ( on->nrow() == 0 )
    5321     return ;
    5322   //   string reftime = on->getTime( index ) ;
    5323   ROScalarColumn<Double> timeCol( off->table(), "TIME" ) ;
    5324   Vector<Double> timeVec = timeCol.getColumn() ;
    5325   //ROTableColumn timeCol( on->table(), "TIME" ) ;
    5326   timeCol.attach( on->table(), "TIME" ) ;
    5327   ArrayColumn<Float> tsysCol( on->table(), "TSYS" ) ;
    5328   vector<float> sp( on->nchan( on->getIF(0) ) ) ;
    5329   unsigned int spsize = sp.size() ;
    5330   for ( int index = 0 ; index < on->nrow() ; index++ ) {
    5331     double reftime = timeCol.asdouble(index) ;
    5332     Vector<Float> spoff = getSpectrumFromTime( reftime, timeVec, off, "linear" ) ;
    5333     vector<float> spec = on->getSpectrum( index ) ;
    5334     Vector<Float> tsys = tsysCol( index ) ;
    5335     // ALMA Calibration
    5336     //
    5337     // Ta* = Tsys * ( ON - OFF ) / OFF
    5338     //
    5339     // 2010/01/07 Takeshi Nakazato
    5340     unsigned int tsyssize = tsys.nelements() ;
    5341     for ( unsigned int j = 0 ; j < sp.size() ; j++ ) {
    5342       float tscale = 0.0 ;
    5343       if ( tsyssize == spsize )
    5344         tscale = tsys[j] ;
    5345       else
    5346         tscale = tsys[0] ;
    5347       float v = tscale * ( spec[j] - spoff[j] ) / spoff[j] ;
    5348       sp[j] = v ;
    5349     }
    5350     on->setSpectrum( sp, index ) ;
    5351   }
    5352 }
    5353 
    5354 void STMath::calibrateALMA( CountedPtr<Scantable>& on,
    5355                             CountedPtr<Scantable>& off,
    5356                             Vector<uInt>& rows )
    5357 {
    5358   // if rows is empty, just return
    5359   if ( rows.nelements() == 0 )
    5360     return ;
    5361 //   string reftime = on->getTime( index ) ;
    5362   ROScalarColumn<Double> timeCol( off->table(), "TIME" ) ;
    5363   Vector<Double> timeVec = timeCol.getColumn() ;
    5364   timeCol.attach( on->table(), "TIME" ) ;
    5365   ArrayColumn<Float> tsysCol( on->table(), "TSYS" ) ;
    5366   vector<float> sp( on->nchan( on->getIF(rows[0]) ) ) ;
    5367   unsigned int spsize = sp.size() ;
    5368   // I know that the data is contiguous
    5369   const uInt *p = rows.data() ;
    5370   for ( int irow = 0 ; irow < rows.nelements() ; irow++ ) {
    5371     double reftime = timeCol.asdouble(*p) ;
    5372     Vector<Float> spoff = getSpectrumFromTime( reftime, timeVec, off, "linear" ) ;
    5373     vector<float> spec = on->getSpectrum( *p ) ;
    5374     Vector<Float> tsys = tsysCol( *p ) ;
    5375     // ALMA Calibration
    5376     //
    5377     // Ta* = Tsys * ( ON - OFF ) / OFF
    5378     //
    5379     // 2010/01/07 Takeshi Nakazato
    5380     unsigned int tsyssize = tsys.nelements() ;
    5381     for ( unsigned int j = 0 ; j < sp.size() ; j++ ) {
    5382       float tscale = 0.0 ;
    5383       if ( tsyssize == spsize )
    5384         tscale = tsys[j] ;
    5385       else
    5386         tscale = tsys[0] ;
    5387       float v = tscale * ( spec[j] - spoff[j] ) / spoff[j] ;
    5388       sp[j] = v ;
    5389     }
    5390     on->setSpectrum( sp, *p ) ;
    5391     p++ ;
    5392   }
    5393 }
    5394 
    53955308void STMath::calibrateALMA( CountedPtr<Scantable>& out,
    53965309                            const CountedPtr<Scantable>& on,
     
    54095322  ROArrayColumn<Float> tsysCol( on->table(), "TSYS" ) ;
    54105323  unsigned int spsize = on->nchan( on->getIF(rows[0]) ) ;
    5411   //Vector<Float> spoff( spsize ) ;
    5412   //Vector<Float> spec( spsize ) ;
    54135324  // I know that the data is contiguous
    54145325  const uInt *p = rows.data() ;
     
    54175328    Vector<Float> spoff = getSpectrumFromTime( reftime, timeVec, off, "linear" ) ;
    54185329    Vector<Float> spec = on->specCol_( *p ) ;
    5419     //getSpectrumFromTime2( spoff, reftime, timeVec, off, "linear" ) ;
    5420     //on->specCol_.get( *p, spec ) ;
    54215330    Vector<Float> tsys = tsysCol( *p ) ;
    54225331    // ALMA Calibration
  • branches/hpc33/src/STMath.h

    r2546 r2550  
    422422                                      casa::CountedPtr<Scantable>& off,
    423423                                      int index ) ;
    424   void calibrateALMA( casa::CountedPtr<Scantable>& on,
    425                       casa::CountedPtr<Scantable>& off ) ;
    426   void calibrateALMA( casa::CountedPtr<Scantable>& on,
    427                       casa::CountedPtr<Scantable>& off,
    428                       casa::Vector<casa::uInt>& rows ) ;
    429424  void calibrateALMA( casa::CountedPtr<Scantable>& out,
    430425                      const casa::CountedPtr<Scantable>& on,
Note: See TracChangeset for help on using the changeset viewer.