Changeset 2761 for trunk/external-alma/atnf/PKSIO

Timestamp:

02/06/13 13:13:20 (12 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...

Added new option, freqref, to NRO filler. Posible values are:
1) 'rest' to import frequency in REST frame, which results in an exactly
same frequency label as NEWSTAR, and 2) 'vref' to import frequency
in the frame that source velocity refers, which results in the same
velocity label as NEWSTAR. The option must be given to scantable
constructor.

Location:

trunk/external-alma/atnf/PKSIO

Files:

• NROReader.cc (modified) (10 diffs)
• NROReader.h (modified) (10 diffs)

trunk/external-alma/atnf/PKSIO/NROReader.cc

@@ -46 +46 @@
 #include <measures/Measures/MeasConvert.h>
 
+#include <casa/Logging/LogIO.h>
 #include <casa/IO/RegularFileIO.h>
 #include <casa/OS/File.h>
@@ -185 +186 @@
   : dataset_( NULL ),
     srcdir_( 0 ),
-    msrcdir_( 0 )
+    msrcdir_( 0 ),
+    freqRefFromVREF_( false ) // import frequency as REST
 {
   // initialization
@@ -199 +201 @@
     dataset_ = NULL ;
   }
+}
+
+// set frequency reference frame
+void NROReader::setFreqRefFromVREF( bool fromVREF )
+{
+  os_.origin( LogOrigin( "NROReader", "setFreqRefFromVREF", WHERE ) ) ;
+  os_ << ((fromVREF) ? "Take frequency reference frame from VREF" : 
+          "Use frequency reference frame REST") << LogIO::POST ;
+
+  freqRefFromVREF_ = fromVREF ; 
 }
 
@@ -387 +399 @@
         for ( int ip = 0 ; ip < 3 ; ip++ )
           qantpos[ip] = Quantity( antpos[ip], "m" ) ;
-        MPosition mp( MVPosition( qantpos ), MPosition::ITRF ) ;
-        mf_->set( mp ) ;
+        mp_ = MPosition( MVPosition( qantpos ), MPosition::ITRF ) ;
+        mf_->set( mp_ ) ;
       }
       converter_ = new MDirection::Convert( MDirection::AZEL,
@@ -397 +409 @@
 
   if ( coord_ == 2 ) {
-    MEpoch me( Quantity( t, "d" ), MEpoch::UTC ) ;
-    mf_->set( me ) ;
+    me_ = MEpoch( Quantity( t, "d" ), MEpoch::UTC ) ;
+    mf_->set( me_ ) ;
   }
 }
@@ -458 +470 @@
     }
   }
+  else if ( vref.compare( 0, 3, "GAL" ) == 0 ) {
+    vref = "GALACTO" ;
+  }
+  else if (vref.compare( 0, 3, "HEL" ) == 0 ) {
+    os_.origin( LogOrigin( "NROReader", "getHeaderInfo", WHERE ) ) ;
+    os_ << LogIO::WARN << "Heliocentric frame is not supported. Use Barycentric frame instead." << LogIO::POST ;
+    vref = "BARY" ;
+  }
   //freqref = vref ;
   //freqref = "LSRK" ;
-  freqref = "REST" ;
+  //freqref = "REST" ;
+  freqref = freqRefFromVREF_ ? vref : "REST" ; 
   //cout << "freqref = " << freqref << endl ;
   NRODataRecord *record = dataset_->getRecord( 0 ) ;
-  reffreq = record->FREQ0 ;
+  reffreq = record->FREQ0 ; // rest frequency
+
   //cout << "reffreq = " << reffreq << endl ;
   bw = dataset_->getBEBW()[0] ;
@@ -538 +560 @@
                             Float &winddir,      
                             Double &srcvel,
-                            Vector<Double> &propermotion,
+                            Vector<Double> &/*propermotion*/,
                             Vector<Double> &srcdir,
-                            Vector<Double> &scanrate )
+                            Vector<Double> &/*scanrate*/ )
 {
   static const IPosition oneByOne( 1, 1 );
@@ -580 +602 @@
   //cout << "freqs = [" << freqs[0] << ", " << freqs[1] << ", " << freqs[2] << "]" << endl ;
 
+  if ( freqRefFromVREF_ ) {
+    freqs = shiftFrequency( freqs, 
+                            dataset_->getURVEL(), 
+                            dataset_->getVDEF() ) ;
+  }
+
   // restfreq (for MOLECULE_ID)
   restfreq.resize( oneByOne ) ;
@@ -609 +637 @@
   vector<double> spec = dataset_->getSpectrum( irow ) ;
   spectra.resize( spec.size() ) ;
-  int index = 0 ;
+  //int index = 0 ;
   Bool b ;
   Float *fp = spectra.getStorage( b ) ;
@@ -708 +736 @@
   return dataset_->getRowNum() ;
 }
+
+vector<double> NROReader::shiftFrequency( const vector<double> &f, 
+                                          const double &v, 
+                                          const string &vdef )
+{
+  vector<double> r( f ) ;
+  double factor = v / 2.99792458e8 ;
+  if ( vdef.compare( 0, 3, "RAD" ) == 0 ) {
+    factor = 1.0 / ( 1.0 + factor ) ;
+    r[1] *= factor ;
+    r[2] *= factor ;
+  }
+  else if ( vdef.compare( 0, 3, "OPT" ) == 0 ) {
+    factor += 1.0 ;
+    r[1] *= factor ;
+    r[2] *= factor ;
+  }
+  else {
+    cout << "vdef=" << vdef << " is not supported." << endl;
+  }
+  return r ;
+}

trunk/external-alma/atnf/PKSIO/NROReader.h

@@ -42 +42 @@
 #include <casa/BasicSL/String.h>
 #include <measures/Measures/MPosition.h>
+#include <measures/Measures/MEpoch.h>
 #include <measures/Measures/MCPosition.h>
 #include <measures/Measures/MDirection.h>
@@ -56 +57 @@
 #include <atnf/PKSIO/NRODataset.h>
 #include <atnf/PKSIO/NRODataRecord.h>
-
-using namespace std ;
 
 // <summary>
@@ -103 +102 @@
   // Destructor.
   virtual ~NROReader() ;
+
+  // determine whether to import frequency as REST (frequency is same as 
+  // NEWSTAR/NOSTAR) or as is (velocity is same as NEWSTAR/NOSTAR)
+  void setFreqRefFromVREF( bool fromVREF ) ;
 
   // Read data header
@@ -133 +136 @@
                            uInt &beamno,
                            uInt &polno,
-                           vector<double> &freqs,   
+                           std::vector<double> &freqs,   
                            Vector<Double> &restfreq,  
                            uInt &refbeamno,
@@ -171 +174 @@
 
   // Get IF settings
-  virtual vector<Bool> getIFs() ;
+  virtual std::vector<Bool> getIFs() ;
 
   // Get Number of IFs
@@ -177 +180 @@
 
   // Get Beam settings
-  virtual vector<Bool> getBeams() ;
+  virtual std::vector<Bool> getBeams() ;
 
   // Get Number of Beams
@@ -188 +191 @@
 
   // Get spectrum
-  virtual vector< vector<double> > getSpectrum() ;
+  virtual std::vector< std::vector<double> > getSpectrum() ;
 
   // Get number of polarization
@@ -196 +199 @@
   virtual double getStartTime() ;
   virtual double getEndTime() ;
-  virtual vector<double> getStartIntTime() ;
+  virtual std::vector<double> getStartIntTime() ;
   //virtual double getStartIntTime( int i ) ;
 
   // Get Antenna Position in ITRF coordinate
-  virtual vector<double> getAntennaPosition() = 0 ;
+  virtual std::vector<double> getAntennaPosition() = 0 ;
 
   // Get SRCDIRECTION in RADEC(J2000)
@@ -209 +212 @@
   virtual void initConvert( int icoord, double t, char *epoch ) ;
 
+  // Shift frequency by given velocity with respect to specified 
+  // velocity reference
+  std::vector<double> shiftFrequency( const std::vector<double> &f, const double &v, const string &vref ) ; 
+
   // filename 
   string filename_ ;
@@ -222 +229 @@
   CountedPtr<MDirection::Convert> converter_ ;
   CountedPtr<MeasFrame> mf_ ;
+  MEpoch me_ ;
+  MPosition mp_ ;
   int coord_ ;
 
+  bool freqRefFromVREF_ ;
+
   // Logger
-  //LogIO os ;
+  LogIO os_ ;
 };