Changeset 2262 for branches/parallel/external-alma

Timestamp:

08/04/11 12:58:21 (13 years ago)

Author:

Takeshi Nakazato

Message:

New Development: No

JIRA Issue: No

Ready for Test: No

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

Preparation for check in tuning on NRO filler.
Merge changes in trunk.

r2154
r2156
r2158
r2198-2203
r2212
r2261

Location:

branches/parallel

Files:

• . (modified) (1 prop)
• external-alma (modified) (1 prop)
• external-alma/atnf/PKSIO/ASTEReader.cc (modified) (1 diff)
• external-alma/atnf/PKSIO/ASTEReader.h (modified) (1 diff)
• external-alma/atnf/PKSIO/NRO45Reader.cc (modified) (1 diff)
• external-alma/atnf/PKSIO/NRO45Reader.h (modified) (1 diff)
• external-alma/atnf/PKSIO/NRODataset.cc (modified) (6 diffs)
• external-alma/atnf/PKSIO/NRODataset.h (modified) (3 diffs)
• external-alma/atnf/PKSIO/NROFITSDataset.cc (modified) (2 diffs)
• external-alma/atnf/PKSIO/NROFITSDataset.h (modified) (1 diff)
• external-alma/atnf/PKSIO/NROReader.cc (modified) (10 diffs)
• external-alma/atnf/PKSIO/NROReader.h (modified) (2 diffs)
• external-alma/atnf/pks/pks_maths.cc (modified) (1 prop)

branches/parallel/external-alma/atnf/PKSIO/ASTEReader.cc

@@ -118 +118 @@
 }
 
+Int ASTEReader::getNumIF() 
+{
+  vector<Bool> v ;
+  vector<int> arry = dataset_->getARRY() ;
+  for ( uInt i = 0 ; i < arry.size() ; i++ ) {
+    if ( arry[i] != 0 ) {
+      v.push_back( True ) ;
+    }
+  }
+  Int nif = v.size() ;
+  return nif ;
+}
+
+Int ASTEReader::getNumBeam()
+{
+  // ASTE doesn't have array receiver so far
+  return 1 ;
+}

branches/parallel/external-alma/atnf/PKSIO/ASTEReader.h

@@ -81 +81 @@
   virtual int read() ;
 
+  // get number of IF
+  virtual Int getNumIF() ;
+
+  // get number of Beams
+  virtual Int getNumBeam() ;
+
  protected:
   // Get Antenna Position in ITRF coordinate

branches/parallel/external-alma/atnf/PKSIO/NRO45Reader.cc

@@ -124 +124 @@
 }
 
+Int NRO45Reader::getNumIF() 
+{
+  Int nif = 1 ;
+  if ( (dataset_->getRX()[0]).find("MULT2-") == string::npos ) {
+    vector<Bool> v ;
+    vector<int> arry = dataset_->getARRY() ;
+    for ( uInt i = 0 ; i < arry.size() ; i++ ) {
+      if ( arry[i] != 0 ) {
+        v.push_back( True ) ;
+      }
+    }
+    nif = v.size() ;
+  }
+  return nif ;
+}
+
+Int NRO45Reader::getNumBeam()
+{
+  Int nbeam = 25 ;
+  if ( (dataset_->getRX()[0]).find("MULT2-") == string::npos ) {
+    nbeam = 1 ;
+  }
+  return nbeam ;
+}

branches/parallel/external-alma/atnf/PKSIO/NRO45Reader.h

@@ -83 +83 @@
   virtual int read() ;
 
+  // get number of IF
+  virtual Int getNumIF() ;
+
+  // get number of Beams
+  virtual Int getNumBeam() ;
+
  protected:
   // Get Antenna Position in ITRF coordinate

branches/parallel/external-alma/atnf/PKSIO/NRODataset.cc

@@ -83 +83 @@
   }
   same_ = -1 ;
+
+  // Record for frequency setting
+  frec_ = Record() ;
 }
 
@@ -602 +605 @@
 }
 
+double NRODataset::getScanTime( int i ) 
+{
+  double startTime = getStartIntTime( i ) ;
+  double interval = getIPTIM() / 86400.0 ; // [sec]->[day]
+  return startTime+0.5*interval ;
+}
+
 vector<bool> NRODataset::getIFs()
 {
@@ -657 +667 @@
   NRODataRecord *record = getRecord( i ) ;
   string arryt = string( record->ARRYT ) ;
-  //string sbeamno = arryt.substr( 1, arryt.size()-1 ) ;
-  //uInt ib = atoi( sbeamno.c_str() ) - 1 ; 
   uInt ib = getArrayId( arryt ) ;
-
-  int ia = -1 ;
+  string rxname = getRX()[0] ;
+  string key = arryt ;
+  if ( rxname.find("MULT2") != string::npos )
+    key = "BEARS" ;
+
+  if ( frec_.isDefined( key ) ) {
+    // frequency for the array is already calculated
+    Vector<Double> f =  frec_.asArrayDouble( key ) ;
+    Double *f_p = f.data() ;
+    for ( int i = 0 ; i < 3 ; i++ )
+      v[i] = (double)(f_p[i]) ;
+    return v ;
+  }
+
+  //int ia = -1 ;
   bool isAOS = false ;
   //cout << "NRODataset::getFrequencies()  record->ARRYT=" << record->ARRYT << endl ;
   //cout << "NRODataset::getFrequencies()  ib = " << ib << endl ;
 
-  if ( strncmp( record->ARRYT, "X", 1 ) == 0 ) {
-    // FX
-    if ( strncmp( (record->ARRYT)+1, "1", 1 ) == 0 
-         || strncmp( (record->ARRYT)+1, "3", 1 ) ) {
-      // FX1, 3
-      ia = 2 ;
-    }
-    else {
-      // FX2, 4
-      ia = 1 ;
-    }
-  }
-  else if ( strncmp( record->ARRYT, "A", 1 ) == 0 )
-    ia = 2 ;  // AC
-  else if ( strncmp( record->ARRYT, "W", 1 ) == 0 ) {
-    // AOS-W    
-    ia = 2 ;  
+  if ( arryt[0] == 'W' || arryt[0] == 'U' || arryt[0] == 'H' )
     isAOS = true ;
-  }
-  else if ( strncmp( record->ARRYT, "U", 1 ) == 0 ) {
-    // AOS-U
-    ia = 2 ;  
-    isAOS = true ;
-  }
-  else if ( strncmp( record->ARRYT, "H", 1 ) == 0 ) {
-    // AOS-H
-    isAOS = true ;
-    //cout << record->ARRYT << " " <<  strlen(record->ARRYT) << endl ;
-    //cout << (record->ARRYT)+1 << endl ;
-    if ( strncmp( (record->ARRYT)+2, " ", 1 ) == 0 ) {
-      // H1-9
-      if ( strncmp( (record->ARRYT)+1, "9", 1 ) == 0 ) {
-        // H9
-        ia = 2 ;
-      }
-      else {
-        // H1-8
-        ia = 1 ;
-      }
-    }
-    else {
-      // H10-16
-      ia = 2 ;
-    }
-  }
-
-  int iu ;
+
+  Bool isUSB ;
   if ( record->FQIF1 > 0 )
-    iu = 1 ;  // USB
+    isUSB = True ;  // USB
   else 
-    iu = 2 ;  // LSB
+    isUSB = False ;  // LSB
 
   int ivdef = -1 ;
-  //if ( strncmp( (dataset_->getVDEF()).c_str(), "RAD", 3 ) == 0 )
   if ( (getVDEF()).compare( 0, 3, "RAD" ) == 0 )
     ivdef = 0 ;
-  //else if ( strncmp( dataset_->getVDEF(), "OPT", 3 ) == 0 )
   else if ( (getVDEF()).compare( 0, 3, "OPT" ) == 0 )
     ivdef = 1 ;
@@ -734 +710 @@
 
   int ncal = getNFCAL()[ib] ;
-  vector<double> freqs( ncal ) ;
   double cw = 0.0 ;
   vector<double> fqcal = getFQCAL()[ib] ;
   vector<double> chcal = getCHCAL()[ib] ;
-
+  double f0cal = getF0CAL()[ib] ;
+  Vector<Double> freqs( ncal, fq0-f0cal ) ;
+
+  double factor = vel / cvel ;
+  if ( ivdef == 0 )
+    factor = 1.0 / ( 1.0 - factor ) ;
   for ( int ii = 0 ; ii < ncal ; ii++ ) {
-    freqs[ii] = fqcal[ii] ;
-    freqs[ii] -= getF0CAL()[ib] ;
-    if ( ia == 1 ) {
-      if ( iu == 1 ) {
-        freqs[ii] = fq0 + freqs[ii] ;
-      }
-      else if ( iu == 2 ) {
-        freqs[ii] = fq0 - freqs[ii] ;
-      }
-    }
-    else if ( ia == 2 ) {
-      if ( iu == 1 ) {
-        freqs[ii] = fq0 - freqs[ii] ;
-      }
-      else if ( iu == 2 ) {
-        freqs[ii] = fq0 + freqs[ii] ;
-      }
-    }     
-//       if ( ivdef == 0 ) {
-//         double factor = 1.0 / ( 1.0 - vel / cvel ) ;
-//         freqs[ii] = freqs[ii] * factor - record->FQTRK * ( factor - 1.0 ) ;
-//       }
-//       else if ( ivdef == 1 ) {
-//         double factor = vel / cvel ;
-//         freqs[ii] = freqs[ii] * ( 1.0 + factor ) - record->FQTRK * factor ;
-//       }
+    freqs[ii] += fqcal[ii] ;
+    if ( ivdef == 0 ) {
+      freqs[ii] = freqs[ii] * factor + record->FQTRK * ( 1.0 - factor ) ;
+    }
+    else if ( ivdef == 1 ) {
+      freqs[ii] = freqs[ii] * ( 1.0 + factor ) - record->FQTRK * factor ;
+    }
+
+    //ofstream ofs("freqs.txt",ios::out|ios::app) ;
+    //ofs << setprecision(16) ;
+    //ofs << i << " " << record->ARRYT << " " << chcal[ii] << " " << freqs[ii] << " " << record->ISCAN << " " << fqcal[ii] << " " << f0cal << " " << record->FQTRK << " " << vel << endl ; 
+    //ofs.close() ;
+
   }
 
@@ -790 +757 @@
   }
   else {
-    cw = getBERES()[ib] ;
+
+    cw = ( freqs[1] - freqs[0] ) 
+      / ( chcal[1] - chcal[0] ) ;    
+
+    if ( isUSB ) {
+      // channel frequency inversion 
+      cw *= -1.0 ;
+      Double tmp = freqs[1] ;
+      freqs[1] = freqs[0] ; 
+      freqs[0] = tmp ;
+    }
     
-    if ( iu == 2 )
-      cw *= -1.0 ;
-
-    if ( cw == 0.0 ) {
-      cw = ( freqs[1] - freqs[0] ) 
-        / ( chcal[1] - chcal[0] ) ;
-//           if ( cw < 0.0 ) 
-//             cw = - cw ;
-    }
     v[0] = chcal[0] - 1 ; // 0-base
     v[1] = freqs[0] ;
@@ -806 +774 @@
   }
 
-  // conversion from TOPO to LSRK
-  v[1] = toLSR( v[1], getStartIntTime( i ), record_->SCX, record_->SCY ) ;
-
-  if ( refFreq_[ib] != 0.0 ) {
-    v[1] = refFreq_[ib] ;
-  }
-  else {
+  if ( refFreq_[ib] == 0.0 ) 
     refFreq_[ib] = v[1] ;
-  }
+
+  // register frequency setting to Record
+  Vector<Double> f( v ) ;
+  frec_.define( key, f ) ;
 
   return v ;

branches/parallel/external-alma/atnf/PKSIO/NRODataset.h

@@ -41 +41 @@
 //#include <casa/aips.h>
 #include <casa/Logging/LogIO.h>
+#include <casa/Containers/Record.h>
 #include <atnf/PKSIO/NRODataRecord.h>
 #include <casa/namespace.h>
@@ -194 +195 @@
   virtual vector<double> getStartIntTime() ;
   virtual double getStartIntTime( int i ) ;
+  virtual double getScanTime( int i ) ;
   virtual double getMJD( char *time ) ;
   virtual vector<bool> getIFs() ;
@@ -512 +514 @@
   // reference frequency for each array
   vector<double> refFreq_ ;
+
+  // record to store REFPIX, REFVAL, INCREMENT pair for each array
+  Record frec_ ;
 } ;
 

branches/parallel/external-alma/atnf/PKSIO/NROFITSDataset.cc

@@ -622 +622 @@
 //   for ( int i = 0 ; i < ARYNM ; i++ ) {
 //     if ( i == 0 ) {
-//       //cout << "RX      " ;
+//       cout << "RX      " ;
 //     }
 //     else if ( ( i % 5 ) == 0 ) {
-//       //cout << endl << "        " ;
-//     }
-//     //cout << "\'" << RX[i] << "\' " ;
+//       cout << endl << "        " ;
+//     }
+//     cout << "\'" << RX[i] << "\' " ;
 //   }
-//   //cout << endl ;
+//   cout << endl ;
   //
   if ( readColumn( HPBW, "HPBW", sameEndian ) != 0 ) {
@@ -2734 +2734 @@
 }
 
+double NROFITSDataset::getScanTime( int i ) 
+{
+  double startTime = getStartIntTime( i ) ;
+  double interval = getIPTIM() ;
+  interval /= 86400.0 ;
+  return startTime+0.5*interval ;
+}
+
 // double NROFITSDataset::toLSR( double v, double t, double x, double y ) 
 // {

branches/parallel/external-alma/atnf/PKSIO/NROFITSDataset.h

@@ -96 +96 @@
   virtual uInt getArrayId( string type ) ;
   virtual double getStartIntTime( int i ) ;
+  virtual double getScanTime( int i ) ;
 
  protected:

branches/parallel/external-alma/atnf/PKSIO/NROReader.cc

@@ -39 +39 @@
 #include <atnf/PKSIO/ASTEDataset.h>
 
+#include <measures/Measures/MEpoch.h>
+#include <measures/Measures/MPosition.h>
 #include <measures/Measures/MDirection.h>
 #include <measures/Measures/MCDirection.h>
+#include <measures/Measures/MeasFrame.h>
 #include <measures/Measures/MeasConvert.h>
 
@@ -314 +317 @@
     //<< srcra << "," << srcdec << ") B1950 to (" 
     //<< v( 0 ) << ","<< v( 1 ) << ") J2000" << endl ;
-    os << LogIO::NORMAL << "SRCDIRECTION convert from (" 
-       << srcra << "," << srcdec << ") B1950 to (" 
-       << v( 0 ) << ","<< v( 1 ) << ") J2000" << LogIO::POST ;
+    //os << LogIO::NORMAL << "SRCDIRECTION convert from (" 
+    //   << srcra << "," << srcdec << ") B1950 to (" 
+    //   << v( 0 ) << ","<< v( 1 ) << ") J2000" << LogIO::POST ;
   }
   else if ( strncmp( epoch, "J2000", 5 ) == 0 ) {
@@ -353 +356 @@
     //<< dirx << "," << diry << ") LB to (" 
     //<< v( 0 ) << ","<< v( 1 ) << ") RADEC" << endl ;
-    os << LogIO::NORMAL << "DIRECTION convert from (" 
-       << dirx << "," << diry << ") LB to (" 
-       << v( 0 ) << ","<< v( 1 ) << ") RADEC" << LogIO::POST ;
+    //os << LogIO::NORMAL << "DIRECTION convert from (" 
+    //   << dirx << "," << diry << ") LB to (" 
+    //   << v( 0 ) << ","<< v( 1 ) << ") RADEC" << LogIO::POST ;
   }
   else if ( icoord == 2 ) {
     // convert from AZEL to RADEC
+    vector<double> antpos = getAntennaPosition() ;
+    Vector<Quantity> qantpos( 3 ) ;
+    for ( int ip = 0 ; ip < 3 ; ip++ )
+      qantpos[ip] = Quantity( antpos[ip], "m" ) ;
+    Double scantime = Double( dataset_->getScanTime( i ) ) ;
+    MEpoch me( Quantity( scantime, "d" ), MEpoch::UTC ) ;
+    MPosition mp( MVPosition( qantpos ), MPosition::ITRF ) ;
+    MeasFrame mf( me, mp ) ;
     MDirection result = 
       MDirection::Convert( MDirection( Quantity( dirx, "rad" ), 
                                        Quantity( diry, "rad" ),
                                        MDirection::Ref( MDirection::AZEL ) ),
-                           MDirection::Ref( MDirection::J2000 ) ) () ;
+                           MDirection::Ref( MDirection::J2000, mf ) ) () ;
     v = result.getAngle().getValue() ;
     //cout << "NROReader::getDirection()  DIRECTION convert from (" 
     //<< dirx << "," << diry << ") AZEL to (" 
     //<< v( 0 ) << ","<< v( 1 ) << ") RADEC" << endl ;
-    os << LogIO::NORMAL << "DIRECTION convert from (" 
-       << dirx << "," << diry << ") AZEL to (" 
-       << v( 0 ) << ","<< v( 1 ) << ") RADEC" << LogIO::POST ;
+    //os << LogIO::NORMAL << "DIRECTION convert from (" 
+    //   << dirx << "," << diry << ") AZEL to (" 
+    //   << v( 0 ) << ","<< v( 1 ) << ") RADEC" << LogIO::POST ;
   }
   else if ( icoord == 0 ) {
@@ -387 +398 @@
       //<< dirx << "," << diry << ") B1950 to (" 
       //<< v( 0 ) << ","<< v( 1 ) << ") J2000" << endl ;
-      os << LogIO::NORMAL << "DIRECTION convert from (" 
-         << dirx << "," << diry << ") B1950 to (" 
-         << v( 0 ) << ","<< v( 1 ) << ") J2000" << LogIO::POST ;
+      //os << LogIO::NORMAL << "DIRECTION convert from (" 
+      //   << dirx << "," << diry << ") B1950 to (" 
+      //   << v( 0 ) << ","<< v( 1 ) << ") J2000" << LogIO::POST ;
     }
     else if ( strncmp( epoch, "J2000", 5 ) == 0 ) {
@@ -441 +452 @@
   //cout << "antpos = " << antpos << endl ;
   string eq = dataset_->getEPOCH() ;
-  if ( eq.compare( 0, 5, "B1950" ) == 0 )
-    equinox = 1950.0 ;
-  else if ( eq.compare( 0, 5, "J2000" ) == 0 ) 
-    equinox = 2000.0 ;
+//   if ( eq.compare( 0, 5, "B1950" ) == 0 )
+//     equinox = 1950.0 ;
+//   else if ( eq.compare( 0, 5, "J2000" ) == 0 ) 
+//     equinox = 2000.0 ;
+  // equinox is always 2000.0
+  equinox = 2000.0 ;
   //cout << "equinox = " << equinox << endl ;
   string vref = dataset_->getVREF() ;
@@ -456 +469 @@
   }
   //freqref = vref ;
-  freqref = "LSRK" ;
+  //freqref = "LSRK" ;
+  freqref = "REST" ;
   //cout << "freqref = " << freqref << endl ;
   NRODataRecord *record = dataset_->getRecord( 0 ) ;
@@ -483 +497 @@
   //cout << "poltype = " << poltype << endl ;
 
-  vector<Bool> ifs = getIFs() ;
-  nif = ifs.size() ;
+  //vector<Bool> ifs = getIFs() ;
+  //nif = ifs.size() ;
+  nif = getNumIF() ;
   //cout << "nif = " << nif << endl ;
 
-  vector<Bool> beams = getBeams() ;
-  nbeam = beams.size() ;
+  //vector<Bool> beams = getBeams() ;
+  //nbeam = beams.size() ;
+  nbeam = getNumBeam() ;
   //cout << "nbeam = " << nbeam << endl ;
 
@@ -505 +521 @@
                             uInt &scanno,
                             uInt &cycleno,
+                            uInt &ifno,
                             uInt &beamno,
                             uInt &polno,
@@ -548 +565 @@
   //cout << "cycleno = " << cycleno << endl ;
 
-  // beamno
-  string arryt = string( record->ARRYT ) ;
-  string sbeamno = arryt.substr( 1, arryt.size()-1 ) ;
-  uInt ibeamno = atoi( sbeamno.c_str() ) ; 
-  beamno = ibeamno - 1 ;
+  // beamno and ifno
+  string rxname = dataset_->getRX()[0] ;
+  if ( rxname.find("MULT2") != string::npos ) {
+    string arryt = string( record->ARRYT ) ;
+    beamno = dataset_->getArrayId( arryt ) ;
+    ifno = 0 ;
+  }
+  else {
+    beamno = 0 ;
+    string arryt = string( record->ARRYT ) ;
+    ifno = dataset_->getArrayId( arryt ) ;
+  }
   //cout << "beamno = " << beamno << endl ;
 
@@ -575 +599 @@
 
   // scantime
-  scantime = Double( dataset_->getStartIntTime( irow ) ) ;
+  //scantime = Double( dataset_->getStartIntTime( irow ) ) ;
+  scantime = Double( dataset_->getScanTime( irow ) ) ;
   //cout << "scantime = " << scantime << endl ;
 

branches/parallel/external-alma/atnf/PKSIO/NROReader.h

@@ -129 +129 @@
                            uInt &scanno,
                            uInt &cycleno,
+                           uInt &ifno,
                            uInt &beamno,
                            uInt &polno,
@@ -171 +172 @@
   virtual vector<Bool> getIFs() ;
 
+  // Get Number of IFs
+  virtual Int getNumIF() = 0 ;
+
   // Get Beam settings
   virtual vector<Bool> getBeams() ;
+
+  // Get Number of Beams
+  virtual Int getNumBeam() = 0 ;
 
  protected: