Changeset 2012 for trunk/src/STLineFinder.cpp

Timestamp:

02/25/11 16:51:50 (13 years ago)

Author:

WataruKawasaki

Message:

New Development: Yes

JIRA Issue: Yes (CAS-2373, CAS-2620)

Ready for Test: Yes

Interface Changes: Yes

What Interface Changed: For Scantable::polyBaseline(), parameters and return value have been changed.

Test Programs:

Put in Release Notes: Yes

Module(s): sdbaseline, sd.linefinder

Description: (1) CAS-2373-related:

(1.1) Modified Scantable::polyBaseline() to have the row-based loop inside.

Now it fits and subtracts baseline for all rows and also output info
about the fitting result to logger and text file, while in the
previous version this method just did baseline fitting/subtraction
for one row only and had to be put inside a row-based loop at the
python side ("poly_baseline()" in scantable.py) and result output had
also to be controlled at the python side. Using a test data from NRO
45m telescope (348,000 rows, 512 channels), the processing time of
scantable.poly_baseline() has reduced from 130 minutes to 5-6 minutes.

(1.2) For accelerating the another polynomial fitting function, namely

scantable.auto_poly_baseline(), added a method
Scantable::autoPolyBaseline(). This basically does the same thing
with Scantable::polyBaseline(), but uses linefinfer also to
automatically flag the line regions.

(1.3) To make linefinder usable in Scantable class, added a method

linefinder.setdata(). This makes it possible to apply linefinder
for a float-list data given as a parameter, without setting scantable,
while it was indispensable to set scantable to use linefinger previously.

(2) CAS-2620-related:

Added Scantable::cubicSplineBaseline() and autoCubicSplineBaseline() for
fit baseline using the cubic spline function. Parameters include npiece
(number of polynomial pieces), clipthresh (clipping threshold), and
clipniter (maximum iteration number).

File:

• trunk/src/STLineFinder.cpp (modified) (8 diffs)

trunk/src/STLineFinder.cpp

@@ -925 +925 @@
   itsNoiseBox = in_noise_box;
   itsUseMedian = in_median;
+
+  useScantable = true;
 }
 
@@ -934 +936 @@
   scan=in_scan.getCP();
   AlwaysAssert(!scan.null(),AipsError);
-
+}
+
+// set spectrum data to work with. this is a method to allow linefinder work 
+// without setting scantable for the purpose of using linefinder inside some 
+// method in scantable class. (Dec 22, 2010 by W.Kawasaki)
+void STLineFinder::setData(const std::vector<float> &in_spectrum)
+{
+  spectrum = Vector<Float>(in_spectrum);
+  useScantable = false;
 }
 
@@ -948 +958 @@
                 const casa::uInt &whichRow) throw(casa::AipsError)
 {
-  if (scan.null())
+  if (useScantable && scan.null())
       throw AipsError("STLineFinder::findLines - a scan should be set first,"
                       " use set_scan");
 
-  uInt nchan = scan->nchan(scan->getIF(whichRow));
+  uInt nchan = useScantable ? scan->nchan(scan->getIF(whichRow)) : spectrum.nelements();
   // set up mask and edge rejection
   // no mask given...
@@ -962 +972 @@
   }
   if (mask.nelements()!=nchan)
-      throw AipsError("STLineFinder::findLines - in_scan and in_mask have different"
-            "number of spectral channels.");
+      throw AipsError("STLineFinder::findLines - in_scan and in_mask, or in_spectrum "
+                      "and in_mask have different number of spectral channels.");
 
   // taking flagged channels into account
-  vector<bool> flaggedChannels = scan->getMask(whichRow);
-  if (flaggedChannels.size()) {
+  if (useScantable) {
+    vector<bool> flaggedChannels = scan->getMask(whichRow);
+    if (flaggedChannels.size()) {
       // there is a mask set for this row
       if (flaggedChannels.size() != mask.nelements()) {
-          throw AipsError("STLineFinder::findLines - internal inconsistency: number of mask elements do not match the number of channels");
+          throw AipsError("STLineFinder::findLines - internal inconsistency: number of "
+                          "mask elements do not match the number of channels");
       }
       for (size_t ch = 0; ch<mask.nelements(); ++ch) {
            mask[ch] &= flaggedChannels[ch];
       }
+    }
   }
 
@@ -1015 +1028 @@
       throw AipsError("STLineFinder::findLines - noise_box is supposed to be at least 2 elements");
 
-  spectrum.resize();
-  spectrum = Vector<Float>(scan->getSpectrum(whichRow));
+  if (useScantable) {
+    spectrum.resize();
+    spectrum = Vector<Float>(scan->getSpectrum(whichRow));
+  }
 
   lines.resize(0); // search from the scratch
@@ -1141 +1156 @@
 {
   try {
-       if (scan.null())
-           throw AipsError("STLineFinder::getMask - a scan should be set first,"
-                      " use set_scan followed by find_lines");
-       DebugAssert(mask.nelements()==scan->getChannels(last_row_used), AipsError);
-       /*
-       if (!lines.size())
-           throw AipsError("STLineFinder::getMask - one have to search for "
+    if (useScantable) {
+      if (scan.null())
+        throw AipsError("STLineFinder::getMask - a scan should be set first,"
+                        " use set_scan followed by find_lines");
+      DebugAssert(mask.nelements()==scan->getChannels(last_row_used), AipsError);
+    }
+    /*
+    if (!lines.size())
+       throw AipsError("STLineFinder::getMask - one have to search for "
                            "lines first, use find_lines");
-       */
-       std::vector<bool> res_mask(mask.nelements());
-       // iterator through lines
-       std::list<std::pair<int,int> >::const_iterator cli=lines.begin();
-       for (int ch=0;ch<int(res_mask.size());++ch) {
-            if (ch<edge.first || ch>=edge.second) res_mask[ch]=false;
-            else if (!mask[ch]) res_mask[ch]=false;
-            else {
-                    res_mask[ch]=!invert; // no line by default
-                    if (cli!=lines.end())
-                        if (ch>=cli->first && ch<cli->second)
-                             res_mask[ch]=invert; // this is a line
-            }
-            if (cli!=lines.end())
-                if (ch>=cli->second) {
-                    ++cli; // next line in the list
-                }
-       }
-       return res_mask;
+    */
+    std::vector<bool> res_mask(mask.nelements());
+    // iterator through lines
+    std::list<std::pair<int,int> >::const_iterator cli=lines.begin();
+    for (int ch=0;ch<int(res_mask.size());++ch) {
+      if (ch<edge.first || ch>=edge.second) res_mask[ch]=false;
+      else if (!mask[ch]) res_mask[ch]=false;
+      else {
+        res_mask[ch]=!invert; // no line by default
+        if (cli!=lines.end())
+          if (ch>=cli->first && ch<cli->second)
+            res_mask[ch]=invert; // this is a line
+      }
+      if (cli!=lines.end())
+        if (ch>=cli->second)
+          ++cli; // next line in the list
+    }
+    return res_mask;
   }
   catch (const AipsError &ae) {
-      throw;
+    throw;
   }
   catch (const exception &ex) {
-      throw AipsError(String("STLineFinder::getMask - STL error: ")+ex.what());
+    throw AipsError(String("STLineFinder::getMask - STL error: ")+ex.what());
   }
 }
@@ -1182 +1198 @@
                              const throw(casa::AipsError)
 {
-  // convert to required abscissa units
-  std::vector<double> vel=scan->getAbcissa(last_row_used);
+  std::vector<double> vel;
+  if (useScantable) {
+    // convert to required abscissa units
+    vel = scan->getAbcissa(last_row_used);
+  } else {
+    for (int i = 0; i < spectrum.nelements(); ++i)
+      vel.push_back((double)i);
+  }
   std::vector<int> ranges=getLineRangesInChannels();
   std::vector<double> res(ranges.size());
@@ -1202 +1224 @@
 {
   try {
-       if (scan.null())
-           throw AipsError("STLineFinder::getLineRangesInChannels - a scan should be set first,"
-                      " use set_scan followed by find_lines");
-       DebugAssert(mask.nelements()==scan->getChannels(last_row_used), AipsError);
-
-       if (!lines.size())
-           throw AipsError("STLineFinder::getLineRangesInChannels - one have to search for "
-                           "lines first, use find_lines");
-
-       std::vector<int> res(2*lines.size());
-       // iterator through lines & result
-       std::list<std::pair<int,int> >::const_iterator cli=lines.begin();
-       std::vector<int>::iterator ri=res.begin();
-       for (;cli!=lines.end() && ri!=res.end();++cli,++ri) {
-            *ri=cli->first;
-            if (++ri!=res.end())
-                *ri=cli->second-1;
-       }
-       return res;
-  }
-  catch (const AipsError &ae) {
-      throw;
-  }
-  catch (const exception &ex) {
-      throw AipsError(String("STLineFinder::getLineRanges - STL error: ")+ex.what());
+    if (useScantable) {
+      if (scan.null())
+        throw AipsError("STLineFinder::getLineRangesInChannels - a scan should be set first,"
+                        " use set_scan followed by find_lines");
+      DebugAssert(mask.nelements()==scan->getChannels(last_row_used), AipsError);
+    }
+
+    if (!lines.size())
+      throw AipsError("STLineFinder::getLineRangesInChannels - one have to search for "
+                      "lines first, use find_lines");
+
+    std::vector<int> res(2*lines.size());
+    // iterator through lines & result
+    std::list<std::pair<int,int> >::const_iterator cli = lines.begin();
+    std::vector<int>::iterator ri = res.begin();
+    for (; cli != lines.end() && ri != res.end(); ++cli,++ri) {
+      *ri = cli->first;
+      if (++ri != res.end())
+        *ri = cli->second - 1;
+    }
+    return res;
+  } catch (const AipsError &ae) {
+    throw;
+  } catch (const exception &ex) {
+    throw AipsError(String("STLineFinder::getLineRanges - STL error: ") + ex.what());
   }
 }