Changeset 343 for trunk


Ignore:
Timestamp:
02/01/05 13:21:27 (20 years ago)
Author:
vor010
Message:

SDLineFinder - now algorithm does many iterations
to ensure that weak lines, which are close to the strong ones, are also found

Location:
trunk/src
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/SDLineFinder.cc

    r331 r343  
    3535
    3636// STL
     37#include <functional>
     38#include <algorithm>
    3739#include <iostream>
    3840
     
    4143using namespace std;
    4244using namespace boost::python;
     45
     46///////////////////////////////////////////////////////////////////////////////
     47//
     48// LFRunningMean - a running mean algorithm for line detection
     49//
     50//
    4351
    4452namespace asap {
     
    117125} // namespace asap
    118126
     127//
     128///////////////////////////////////////////////////////////////////////////////
     129
     130
    119131///////////////////////////////////////////////////////////////////////////////
    120132//
     
    280292///////////////////////////////////////////////////////////////////////////////
    281293
    282 
     294///////////////////////////////////////////////////////////////////////////////
     295//
     296// SDLineFinder::IntersectsWith  -  An auxiliary object function to test
     297// whether two lines have a non-void intersection
     298//
     299
     300
     301// line1 - range of the first line: start channel and stop+1
     302SDLineFinder::IntersectsWith::IntersectsWith(const std::pair<int,int> &in_line1) :
     303                          line1(in_line1) {}
     304
     305
     306// return true if line2 intersects with line1 with at least one
     307// common channel, and false otherwise
     308// line2 - range of the second line: start channel and stop+1
     309bool SDLineFinder::IntersectsWith::operator()(const std::pair<int,int> &line2)
     310                          const throw()
     311{
     312  if (line2.second<line1.first) return false; // line2 is at lower channels
     313  if (line2.first>line1.second) return false; // line2 is at upper channels
     314  return true; // line2 has an intersection or is adjacent to line1
     315}
     316
     317//
     318///////////////////////////////////////////////////////////////////////////////
     319
     320///////////////////////////////////////////////////////////////////////////////
     321//
     322// SDLineFinder::BuildUnion - An auxiliary object function to build a union
     323// of several lines to account for a possibility of merging the nearby lines
     324//
     325
     326// set an initial line (can be a first line in the sequence)
     327SDLineFinder::BuildUnion::BuildUnion(const std::pair<int,int> &line1) :
     328                             temp_line(line1) {}
     329
     330// update temp_line with a union of temp_line and new_line
     331// provided there is no gap between the lines
     332void SDLineFinder::BuildUnion::operator()(const std::pair<int,int> &new_line)
     333                                   throw()
     334{
     335  if (new_line.first<temp_line.first) temp_line.first=new_line.first;
     336  if (new_line.second>temp_line.second) temp_line.second=new_line.second;
     337}
     338
     339// return the result (temp_line)
     340const std::pair<int,int>& SDLineFinder::BuildUnion::result() const throw()
     341{
     342  return temp_line;
     343}
     344
     345//
     346///////////////////////////////////////////////////////////////////////////////
     347
     348///////////////////////////////////////////////////////////////////////////////
     349//
     350// SDLineFinder::LaterThan - An auxiliary object function to test whether a
     351// specified line is at lower spectral channels (to preserve the order in
     352// the line list)
     353//
     354
     355// setup the line to compare with
     356SDLineFinder::LaterThan::LaterThan(const std::pair<int,int> &in_line1) :
     357                         line1(in_line1) {}
     358
     359// return true if line2 should be placed later than line1
     360// in the ordered list (so, it is at greater channel numbers)
     361bool SDLineFinder::LaterThan::operator()(const std::pair<int,int> &line2)
     362                          const throw()
     363{
     364  if (line2.second<line1.first) return false; // line2 is at lower channels
     365  if (line2.first>line1.second) return true; // line2 is at upper channels
     366 
     367  // line2 intersects with line1. We should have no such situation in
     368  // practice
     369  return line2.first>line1.first;
     370}
     371
     372//
     373///////////////////////////////////////////////////////////////////////////////
     374
     375
     376///////////////////////////////////////////////////////////////////////////////
     377//
    283378// SDLineFinder  -  a class for automated spectral line search
     379//
     380//
    284381
    285382SDLineFinder::SDLineFinder() throw() : edge(0,0)
     
    371468  lines.resize(0); // search from the scratch
    372469  Vector<Bool> temp_mask(mask);
    373   size_t cursz;
    374   do {
    375      cursz=lines.size();
     470  Bool need2iterate=True;
     471  while (need2iterate) {
    376472     // line find algorithm
    377      LFRunningMean lfalg(spectrum,temp_mask,edge,max_box_nchan,min_nchan,threshold);
    378      lfalg.findLines(lines);
     473     LFRunningMean lfalg(spectrum,temp_mask,edge,max_box_nchan,min_nchan,
     474                         threshold);
     475     std::list<pair<int,int> > new_lines;
     476     lfalg.findLines(new_lines);
     477     if (!new_lines.size()) need2iterate=False;
    379478     temp_mask=getMask();
    380   } while (cursz!=lines.size());
     479     // update the list (lines) merging intervals, if necessary
     480     addNewSearchResult(new_lines);
     481  }
    381482  return int(lines.size());
    382483}
     
    474575      for (std::list<pair<int,int> >::const_iterator cli=newlines.begin();
    475576           cli!=newlines.end();++cli) {
    476            // search for a right place for the new line
    477            //TODO
     577           
     578           // the first item, which has a non-void intersection or touches
     579           // the new line
     580           std::list<pair<int,int> >::iterator pos_beg=find_if(lines.begin(),
     581                          lines.end(), IntersectsWith(*cli));           
     582           // the last such item         
     583           std::list<pair<int,int> >::iterator pos_end=find_if(pos_beg,
     584                          lines.end(), not1(IntersectsWith(*cli)));
     585
     586           // extract all lines which intersect or touch a new one into
     587           // a temporary buffer. This may invalidate the iterators
     588           // line_buffer may be empty, if no lines intersects with a new
     589           // one.
     590           std::list<pair<int,int> > lines_buffer;
     591           lines_buffer.splice(lines_buffer.end(),lines, pos_beg, pos_end);
     592
     593           // build a union of all intersecting lines
     594           pair<int,int> union_line=for_each(lines_buffer.begin(),
     595                   lines_buffer.end(),BuildUnion(*cli)).result();
     596           
     597           // search for a right place for the new line (union_line) and add
     598           std::list<pair<int,int> >::iterator pos2insert=find_if(lines.begin(),
     599                          lines.end(), LaterThan(union_line));
     600           lines.insert(pos2insert,union_line);
    478601      }
    479602  }
     
    484607      throw AipsError(String("SDLineFinder::addNewSearchResult - STL error: ")+ex.what());
    485608  }
    486 
    487 }
     609}
  • trunk/src/SDLineFinder.h

    r331 r343  
    8787protected:
    8888   // concatenate two lists preserving the order. If two lines appear to
    89    // be adjacent, they are joined into the new one
     89   // be adjacent or have a non-void intersection, they are joined into
     90   // the new line
    9091   void addNewSearchResult(const std::list<std::pair<int, int> > &newlines)
    9192                           throw(casa::AipsError);
     93                           
     94   // An auxiliary object function to test whether two lines have a non-void
     95   // intersection
     96   class IntersectsWith : public std::unary_function<pair<int,int>, bool> {
     97       std::pair<int,int> line1;           // range of the first line
     98                                           // start channel and stop+1
     99   public:
     100        IntersectsWith(const std::pair<int,int> &in_line1);
     101        // return true if line2 intersects with line1 with at least one
     102        // common channel, and false otherwise
     103        bool operator()(const std::pair<int,int> &line2) const throw();
     104   };
     105
     106   // An auxiliary object function to build a union of several lines
     107   // to account for a possibility of merging the nearby lines
     108   class BuildUnion {
     109       std::pair<int,int> temp_line;       // range of the first line
     110                                           // start channel and stop+1
     111   public:
     112        BuildUnion(const std::pair<int,int> &line1);
     113        // update temp_line with a union of temp_line and new_line
     114        // provided there is no gap between the lines
     115        void operator()(const std::pair<int,int> &new_line) throw();
     116        // return the result (temp_line)
     117        const std::pair<int,int>& result() const throw();
     118   };
     119   
     120   // An auxiliary object function to test whether a specified line
     121   // is at lower spectral channels (to preserve the order in the line list)
     122   class LaterThan : public std::unary_function<pair<int,int>, bool> {
     123       std::pair<int,int> line1;           // range of the first line
     124                                           // start channel and stop+1
     125   public:
     126        LaterThan(const std::pair<int,int> &in_line1);
     127
     128        // return true if line2 should be placed later than line1
     129        // in the ordered list (so, it is at greater channel numbers)
     130        bool operator()(const std::pair<int,int> &line2) const throw();
     131   };
     132   
    92133private:
    93134   casa::CountedConstPtr<SDMemTable> scan; // the scan to work with
Note: See TracChangeset for help on using the changeset viewer.