Changeset 2711 for trunk

Timestamp:

12/27/12 17:37:43 (12 years ago)

Author:

Kana Sugimoto

Message:

New Development: No

JIRA Issue: Yes (CAS-4141/CSV-1526)

Ready for Test: Yes

Interface Changes: Yes

What Interface Changed: Added asap.sbseparator.set_lo1root() = SBSeparator::set_lo1root(),

and STSideBandSep::setLO1Root()

Test Programs:

Put in Release Notes: No

Module(s): sbseparator and STSideBandSep

Description:

Added a new python method asap.sbseparator.set_lo1root() to set MS
name to search for LO1 frequency.
Underlying c++ functions are SBSeparator::set_lo1root (interface) and STSideBandSep::setLO1Root().
Also developed a capability to get relevant table names from scantable header and
fill frequencies of image side band.

Location:

trunk

Files:

• python/sbseparator.py (modified) (5 diffs)
• src/STSideBandSep.cpp (modified) (9 diffs)
• src/STSideBandSep.h (modified) (3 diffs)
• src/python_PlotHelper.cpp (modified) (1 diff)
• src/python_STSideBandSep.cpp (modified) (1 diff)

trunk/python/sbseparator.py

@@ -49 +49 @@
         self.solveother = False
         self.dirtol = [1.e-5, 1.e-5] # direction tolerance in rad (2 arcsec)
-        self.lo1 = 0.
+        #self.lo1 = 0.
 
         self.tables = []
@@ -57 +57 @@
         self.set_data(infiles)
         
-        #self.separator = sidebandsep()
+        self.separator = SBSeparator()
 
     @asaplog_post_dec
@@ -221 +221 @@
         lo1 : LO1 frequency in float
         """
+        lo1val = -1.
         if isinstance(lo1, dict) and lo1["unit"] == "Hz":
-            self.lo1 = lo1["value"]
-        else:
-            self.lo1 = float(lo1)
-        if self.lo1 <= 0.:
+            lo1val = lo1["value"]
+        else:
+            lo1val = float(lo1)
+        if lo1val <= 0.:
             asaplog.push("Got negative LO1 frequency. It will be ignored.")
             asaplog.post("WARN")
+        else:
+            self.separator.set_lo1(lo1val)
+
+
+    def set_lo1root(self, name):
+        """
+        Set MS name which stores LO1 frequency of signal side band.
+        It is used to calculate frequency of image sideband.
+
+        name : MS name which contains 'ASDM_SPECTRALWINDOW' and
+               'ASDM_RECEIVER' tables.
+        """
+        self.separator.set_lo1root(name)
 
 
@@ -289 +303 @@
             asaplog.post()
             asaplog.push("Saving IMAGE sideband.")
-            asaplog.push("Note, frequency information of IMAGE sideband cannot be properly filled so far. (future development)")
-            asaplog.push("Storing frequency setting of SIGNAL sideband in FREQUENCIES table for now.")
-            asaplog.post("WARN")
+            #asaplog.push("Note, frequency information of IMAGE sideband cannot be properly filled so far. (future development)")
+            #asaplog.push("Storing frequency setting of SIGNAL sideband in FREQUENCIES table for now.")
+            #asaplog.post("WARN")
 
             imagename = outname + ".imageband"
@@ -300 +314 @@
                     shutil.rmtree(imagename)
             # Update frequency information
-            sbsep = SBSeparator()
-            sbsep.set_imgtable(imagetab)
-            if self.lo1 > 0.:
-                sbsep.set_lo1(self.lo1)
-            sbsep.solve_imgfreq()
+            self.separator.set_imgtable(imagetab)
+            self.separator.solve_imgfreq()
             imagetab.save(imagename)
-            del sbsep
 
 

trunk/src/STSideBandSep.cpp

@@ -21 +21 @@
 #include <measures/Measures/MCFrequency.h>
 
+#include <tables/Tables/TableRow.h>
 #include <tables/Tables/TableRecord.h>
 #include <tables/Tables/TableVector.h>
@@ -98 +99 @@
 };
 
-void STSideBandSep::setLO1Asdm(string asdmname)
-{
-  // Check for existance of the asdm.
-  if (!checkFile(asdmname)) {
+void STSideBandSep::setLO1Root(string name)
+{
+   LogIO os(LogOrigin("STSideBandSep","setLO1Root()", WHERE));
+   os << "Searching for '" << name << "'..." << LogIO::POST;
+  // Check for existance of the file
+  if (!checkFile(name)) {
      throw(AipsError("File does not exist"));
   }
-  if (asdmname[(asdmname.size()-1)] == '/')
-    asdmname = asdmname.substr(0,(asdmname.size()-2));
-
-  asdmName_ = asdmname;
-
-#ifdef KS_DEBUG
-  cout << "STSideBandSep::setLO1Asdm" << endl;
-  cout << "asdm name = " << asdmName_ << endl;
+  if (name[(name.size()-1)] == '/')
+    name = name.substr(0,(name.size()-2));
+
+  if (checkFile(name+"/Receiver.xml", "file") &&
+      checkFile(name+"/SpectralWindow.xml", "file")){
+    os << "Found '" << name << "/Receiver.xml' ... got an ASDM name." << LogIO::POST;
+    asdmName_ = name;
+  } else if (checkFile(name+"/ASDM_RECEIVER") &&
+             checkFile(name+"/ASDM_SPECTRALWINDOW")){
+    os << "Found '" << name << "/ASDM_RECEIVER' ... got a Table name." << LogIO::POST;
+    asisName_ = name;
+  } else {
+    throw(AipsError("Invalid file name. Set an MS or ASDM name."));
+  }
+
+#ifdef KS_DEBUG
+  cout << "STSideBandSep::setLO1Root" << endl;
+  if (!asdmName_.empty())
+    cout << "asdm name = " << asdmName_ << endl;
+  if (!asisName_.empty())
+    cout << "MS name = " << asisName_ << endl;
 #endif
 };
@@ -126 +142 @@
     throw AipsError("STSideBandSep::solveImageFreqency - an image side band scantable should be set first");
 
-  // Check for the availability of LO1
-  os << "Looking for LO1. lo1Freq_ = " << lo1Freq_ << LogIO::POST; 
-  if (lo1Freq_ > 0.) {
-    os << "Using user defined LO1 frequency." << LogIO::POST;
-  } else if (!asdmName_.empty() > 0) {
-      // ASDM name is set.
-    os << "Using user defined ASDM: " << asdmName_ << LogIO::POST;
-    if (!getLo1FromAsdm(asdmName_)) {
-      throw AipsError("Failed to get LO1 frequency from ASDM");
-    }
-  } else {
-    // Try getting ASDM name from scantable header
-    os << "Try getting information from scantable header" << LogIO::POST;
-    if (!getLo1FromTab(imgTab_p)) {
-      throw AipsError("Failed to get LO1 frequency from asis table");
-    }
-  }
-  // LO1 should now be ready.
-  if (lo1Freq_ < 0.)
-    throw(AipsError("Got negative LO1 Frequency"));
-
+  // Convert frequency REFVAL to the value in frame of LO.
   // The code assumes that imgTab_p has only an IF and only a FREQ_ID
   // is associated to an IFNO
@@ -160 +156 @@
   // assuming single freqID per IFNO
   uInt freqid = freqIdVec(0);
+  int nChan = imgTab_p->nchan(imgTab_p->getIF(0));
   double refpix, refval, increment ;
   freqTab_.getEntry(refpix, refval, increment, freqid);
@@ -203 +200 @@
     md = MDirection(ra.getValue("rad"), dec.getValue("rad"),epoch);
   }
-  // Summary (LO1)
+  MeasFrame mframe( me, mp, md );
+  MFrequency::Convert tobframe(loFrame_, MFrequency::Ref(tabframe, mframe));
+  MFrequency::Convert toloframe(tabframe, MFrequency::Ref(loFrame_, mframe));
+  // Convert refval to loFrame_
+  double sigrefval;
+  if (tabframe == loFrame_)
+    sigrefval = refval;
+  else
+    sigrefval = toloframe(Quantum<Double>(refval, "Hz")).get("Hz").getValue();
+
+
+  // Check for the availability of LO1
+  if (lo1Freq_ > 0.) {
+    os << "Using user defined LO1 frequency." << LogIO::POST;
+  } else if (!asisName_.empty()) {
+      // MS name is set.
+    os << "Using user defined MS (asis): " << asisName_ << LogIO::POST;
+    if (!getLo1FromAsisTab(asisName_, sigrefval, refpix, increment, nChan)) {
+      throw AipsError("Failed to get LO1 frequency from MS");
+    }
+  } else if (!asdmName_.empty()) {
+      // ASDM name is set.
+    os << "Using user defined ASDM: " << asdmName_ << LogIO::POST;
+    if (!getLo1FromAsdm(asdmName_, sigrefval, refpix, increment, nChan)) {
+      throw AipsError("Failed to get LO1 frequency from ASDM");
+    }
+  } else {
+    // Try getting ASDM name from scantable header
+    os << "Try getting information from scantable header" << LogIO::POST;
+    if (!getLo1FromScanTab(imgTab_p, sigrefval, refpix, increment, nChan)) {
+      //throw AipsError("Failed to get LO1 frequency from asis table");
+      os << LogIO::WARN << "Failed to get LO1 frequency using information in scantable." << LogIO::POST;
+      os << LogIO::WARN << "Could not fill frequency information of IMAGE sideband properly." << LogIO::POST;
+      os << LogIO::WARN << "Storing values of SIGNAL sideband in FREQUENCIES table" << LogIO::POST;
+      return;
+    }
+  }
+  // LO1 should now be ready.
+  if (lo1Freq_ < 0.)
+    throw(AipsError("Got negative LO1 Frequency"));
+
+  // Print summary (LO1)
   Vector<Double> dirvec = md.getAngle(Unit(String("rad"))).getValue();
   os << "[LO1 settings]" << LogIO::POST;
@@ -212 +250 @@
   os << "- Reference direction: [" << dirvec[0] << ", " << dirvec[1]
      << "] (" << md.getRefString() << ") " << LogIO::POST;
-  //
-  MeasFrame mframe( me, mp, md );
-  MFrequency::Convert tobframe(loFrame_, MFrequency::Ref(tabframe, mframe));
-  MFrequency::Convert toloframe(tabframe, MFrequency::Ref(loFrame_, mframe));
-  // convert refval to loFrame_
-  double sigrefval;
-  if (tabframe == loFrame_)
-    sigrefval = refval;
-  else
-    sigrefval = toloframe(Quantum<Double>(refval, "Hz")).get("Hz").getValue();
-
-  //Summary (signal)
+
+  //Print summary (signal)
   os << "[Signal side band]" << LogIO::POST;
   os << "- IFNO: " << imgTab_p->getIF(0) << " (FREQ_ID = " << freqid << ")"
@@ -233 +261 @@
   os << "- Reference pixel: " << refpix  << LogIO::POST;
   os << "- Increment: " << increment << " [Hz]" << LogIO::POST;
-  // calculate image band incr and refval in loFrame_
+
+  // Calculate image band incr and refval in loFrame_
   Double imgincr = -increment;
   Double imgrefval = 2 * lo1Freq_ - sigrefval;
   Double imgrefval_tab = imgrefval;
-  // convert imgrefval back to table base frame
+  // Convert imgrefval back to table base frame
   if (tabframe != loFrame_)
     imgrefval = tobframe(Quantum<Double>(imgrefval, "Hz")).get("Hz").getValue();
-  // set new frequencies table
+  // Set new frequencies table
   uInt fIDnew = freqTab_.addEntry(refpix, imgrefval, imgincr);
-  // update freq_ID in table.
+  // Update FREQ_ID in table.
   freqIdVec = fIDnew;
-  // Summary (Image side band)
+
+  // Print summary (Image side band)
   os << "[Image side band]" << LogIO::POST;
   os << "- IFNO: " << imgTab_p->getIF(0) << " (FREQ_ID = " << freqIdVec(0)
@@ -278 +308 @@
 };
 
-bool STSideBandSep::getLo1FromAsdm(const string asdmname)
+bool STSideBandSep::getLo1FromAsdm(const string asdmname,
+                                   const double refval,
+                                   const double refpix,
+                                   const double increment,
+                                   const int nChan)
 {
   // Check for relevant tables.
@@ -292 +326 @@
 
 
-bool STSideBandSep::getLo1FromTab(CountedPtr< Scantable > &scantab)
-{
+bool STSideBandSep::getLo1FromScanTab(CountedPtr< Scantable > &scantab, 
+                                      const double refval,
+                                      const double refpix,
+                                      const double increment,
+                                      const int nChan)
+{
+  LogIO os(LogOrigin("STSideBandSep","getLo1FromScanTab()", WHERE));
   Table& tab = scantab->table();
   // Check for relevant tables.
-  const String spwname = tab.keywordSet().asString("ASDM_SPECTRALWINDOW");
-  const String recname = tab.keywordSet().asString("ASDM_RECEIVER");
+  String spwname, recname;
+  try {
+    spwname = tab.keywordSet().asString("ASDM_SPECTRALWINDOW");
+    recname = tab.keywordSet().asString("ASDM_RECEIVER");
+  } catch (...) {
+    // keywords are not there
+    os << LogIO::WARN
+       << "Could not find necessary table names in scantable header."
+       << LogIO::POST;
+    return false;
+  }
   if (!checkFile(spwname,"directory") || !checkFile(recname,"directory")) {
-    throw(AipsError("Could not find subtables in MS"));
-  }
-  // get
-
+    throw(AipsError("Could not find relevant subtables in MS"));
+  }
+  // Get root MS name
+  string msname;
+  const String recsuff = "/ASDM_RECEIVER";
+  String::size_type pos;
+  pos = recname.size()-recsuff.size()-1;
+  if (recname.substr(pos) == recsuff)
+    msname = recname.substr(0, pos-1);
+  else
+    throw(AipsError("Internal error in parsing table name from a scantable keyword."));
+
+  if (!checkFile(msname))
+    throw(AipsError("Internal error in parsing MS name from a scantable keyword."));
+
+  return getLo1FromAsisTab(msname, refval, refpix, increment, nChan);
+
+};
+
+bool STSideBandSep::getLo1FromAsisTab(const string msname,
+                                      const double refval,
+                                      const double refpix,
+                                      const double increment,
+                                      const int nChan)
+{
+  LogIO os(LogOrigin("STSideBandSep","getLo1FromAsisTab()", WHERE));
+  os << "Searching an LO1 frequency in '" << msname << "'" << LogIO::POST;
+  // Check for relevant tables.
+  const string spwname = msname + "/ASDM_SPECTRALWINDOW";
+  const string recname = msname + "/ASDM_RECEIVER";
+  if (!checkFile(spwname,"directory") || !checkFile(recname,"directory")) {
+    throw(AipsError("Could not find relevant tables in MS"));
+  }
+
+  Table spwtab_ = Table(spwname);
+  String asdmSpw;
+  ROTableRow spwrow(spwtab_);
+  const Double rtol = 0.01;
+  for (uInt idx = 0; idx < spwtab_.nrow(); idx++){
+    const TableRecord& rec = spwrow.get(idx);
+    // Compare nchan
+    if (rec.asInt("numChan") != (Int) nChan)
+      continue;
+    // Compare increment
+    Double asdminc;
+    Array<Double> incarr = rec.asArrayDouble("chanWidthArray");
+    if (incarr.size() > 0)
+      asdminc = incarr(IPosition(1, (uInt) refpix));
+    else
+      asdminc = rec.asDouble("chanWidth");
+    if (abs(asdminc - abs(increment)) > rtol * abs(increment))
+      continue;
+    // Compare refval
+    Double asdmrefv;
+    Array<Double> refvarr = rec.asArrayDouble("chanFreqArray");
+    if (refvarr.size() > 0){
+      const uInt iref = (uInt) refpix;
+      const Double ratio = refpix - (Double) iref;
+      asdmrefv = refvarr(IPosition(1, iref))*(1.-ratio)
+        + refvarr(IPosition(1,iref+1))*ratio;
+    }
+    else {
+      const Double ch0 = rec.asDouble("chanFreqStart");
+      const Double chstep = rec.asDouble("chanFreqStep");
+      asdmrefv = ch0 + chstep * refpix;
+    }
+    if (abs(asdmrefv - refval) < 0.5*abs(asdminc)){
+      asdmSpw = rec.asString("spectralWindowId");
+      break;
+    }
+  }
+
+  if (asdmSpw.empty()){
+    os << LogIO::WARN << "Could not find relevant SPW ID in " << spwname << LogIO::POST;
+    return false;
+  }
+  else {
+    os << asdmSpw << " in " << spwname
+       << " matches the freqeuncies of signal side band." << LogIO::POST;
+  }
+
+  Table rectab_ = Table(recname);
+  ROTableRow recrow(rectab_);
+  for (uInt idx = 0; idx < rectab_.nrow(); idx++){
+    const TableRecord& rec = recrow.get(idx);
+    if (rec.asString("spectralWindowId") == asdmSpw){
+      const Array<Double> loarr = rec.asArrayDouble("freqLO");
+      lo1Freq_ = loarr(IPosition(1,0));
+      os << "Found LO1 Frequency in " << recname << ": "
+         << lo1Freq_ << " [Hz]" << LogIO::POST;
+      return true;
+    }
+  }
+  os << LogIO::WARN << "Could not find " << asdmSpw << " in " << recname
+     << LogIO::POST;
   return false;
-
-};
-
+};
+
+// String STSideBandSep::()
+// {
+
+// };
 
 } //namespace asap

trunk/src/STSideBandSep.h

@@ -54 +54 @@
    **/
   void setLO1(const double lo1, string frame="TOPO", double reftime=-1, string refdir="");
-  void setLO1Asdm(const string asdmname);
+  void setLO1Root(const string name);
   /**
    * Actual calculation of frequencies of image sideband
@@ -62 +62 @@
 private:
   Bool checkFile(const string name, string type="");
-  bool getLo1FromAsdm(string asdmname);
-  bool getLo1FromTab(casa::CountedPtr< Scantable > &scantab);
+  bool getLo1FromAsdm(const string asdmname,
+                      const double refval, const double refpix,
+                      const double increment, const int nChan);
+  bool getLo1FromAsisTab(const string msname,
+                         const double refval, const double refpix,
+                         const double increment, const int nChan);
+  bool getLo1FromScanTab(casa::CountedPtr< Scantable > &scantab,
+                         const double refval, const double refpix,
+                         const double increment, const int nChan);
 
   unsigned int sigIfno_;
@@ -71 +78 @@
   double loTime_;
   string loDir_;
-  string asdmName_;
+  string asdmName_, asisName_;
 
   CountedPtr< Scantable > imgTab_p, sigTab_p;

trunk/src/python_PlotHelper.cpp

@@ -1 +1 @@
 //#---------------------------------------------------------------------------
-//# python_PlotHelper.cc: python exposure of c++ PlotHelper class
+//# python_PlotHelper.cpp: python exposure of c++ PlotHelper class
+//#---------------------------------------------------------------------------
+//# Author: Kanako Sugimoto, (C) 2012
+//#
+//# Copyright: See COPYING file that comes with this distribution
+//#
 //#---------------------------------------------------------------------------
 #include <boost/python.hpp>
 #include <boost/python/args.hpp>
 
-#include <string>
-#include <vector>
-
 #include "PlotHelper.h"
-#include "ScantableWrapper.h"
+//#include "ScantableWrapper.h"
 
 using namespace boost::python;

trunk/src/python_STSideBandSep.cpp

@@ -27 +27 @@
            boost::python::arg("reftime")=-1,
            boost::python::arg("refdir")="") )
-    //.def( "set_asdm", &STSideBandSep::setLO1Asdm )
+    .def( "set_lo1root", &STSideBandSep::setLO1Root )
     // temporal methods
     .def( "set_imgtable", &STSideBandSep::setImageTable )