contrib/mul/mfpf/mfpf_edge_finder.cxx
Go to the documentation of this file.
00001 #include "mfpf_edge_finder.h"
00002 //:
00003 // \file
00004 // \brief Locates strongest edge along a profile
00005 // \author Tim Cootes
00006 
00007 #include <vsl/vsl_binary_loader.h>
00008 #include <vnl/vnl_vector.h>
00009 #include <vgl/vgl_point_2d.h>
00010 #include <vgl/vgl_vector_2d.h>
00011 #include <vcl_cmath.h>
00012 
00013 #include <vimt/vimt_bilin_interp.h>
00014 #include <vimt/vimt_sample_profile_bilin.h>
00015 
00016 //=======================================================================
00017 // Dflt ctor
00018 //=======================================================================
00019 
00020 mfpf_edge_finder::mfpf_edge_finder()
00021 {
00022 }
00023 
00024 //=======================================================================
00025 // Destructor
00026 //=======================================================================
00027 
00028 mfpf_edge_finder::~mfpf_edge_finder()
00029 {
00030 }
00031 
00032 //: Radius of circle containing modelled region
00033 double mfpf_edge_finder::radius() const
00034 {
00035   return 1.0;
00036 }
00037 
00038 //: Generate points in ref frame that represent boundary
00039 //  Points of a closed contour around the shape.
00040 //  Used for display purposes.
00041 void mfpf_edge_finder::get_outline(vcl_vector<vgl_point_2d<double> >& pts) const
00042 {
00043   pts.resize(2);
00044   pts[0]=vgl_point_2d<double>(-0.5,0);
00045   pts[1]=vgl_point_2d<double>( 0.5,0);
00046 }
00047 
00048 
00049 //: Evaluate match at p, using u to define scale and orientation
00050 // Returns -1*edge strength at p along direction u
00051 double mfpf_edge_finder::evaluate(const vimt_image_2d_of<float>& image,
00052                                   const vgl_point_2d<double>& p,
00053                                   const vgl_vector_2d<double>& u)
00054 {
00055   double v1 = vimt_bilin_interp_safe(image,p+0.5*step_size_*u);
00056   double v2 = vimt_bilin_interp_safe(image,p-0.5*step_size_*u);
00057   return -1.0*vcl_fabs(v1-v2);
00058 }
00059 
00060    //: Evaluate match at in a region around p
00061    // Returns a quality of fit at a set of positions.
00062    // response image (whose size and transform is set inside the
00063    // function), indicates the points at which the function was
00064    // evaluated.  response(i,j) is the fit at the point
00065 // response.world2im().inverse()(i,j).  The world2im() transformation
00066 // may be affine.
00067 void mfpf_edge_finder::evaluate_region(
00068                         const vimt_image_2d_of<float>& image,
00069                         const vgl_point_2d<double>& p,
00070                         const vgl_vector_2d<double>& u,
00071                         vimt_image_2d_of<double>& response)
00072 {
00073   int n=1+2*search_ni_;
00074   vnl_vector<double> v(n+1);
00075   vgl_vector_2d<double> u1=step_size_*u;
00076   const vgl_point_2d<double> p0 = p-(search_ni_+0.5)*u1;
00077   vimt_sample_profile_bilin(v,image,p0,u1,n+1);
00078   response.image().set_size(n,1);
00079   double* r = response.image().top_left_ptr();
00080   for (int i=0;i<n;++i,++r)
00081   {
00082     *r = -1*vcl_fabs(v[i+1]-v[i]);
00083   }
00084 
00085   // Set up transformation parameters
00086 
00087   // Point (i,j) in dest corresponds to p1+i.u+j.v,
00088   // an affine transformation for image to world
00089   const vgl_point_2d<double> p1 = p-search_ni_*u1;
00090 
00091   vimt_transform_2d i2w;
00092   i2w.set_similarity(vgl_point_2d<double>(u1.x(),u1.y()),p1);
00093   response.set_world2im(i2w.inverse());
00094 }
00095 
00096    //: Search given image around p, using u to define scale and orientation
00097    //  On exit, new_p and new_u define position, scale and orientation of
00098    //  the best nearby match.  Returns a quality of fit measure at that
00099    //  point (the smaller the better).
00100 double mfpf_edge_finder::search_one_pose(
00101                                 const vimt_image_2d_of<float>& image,
00102                                 const vgl_point_2d<double>& p,
00103                                 const vgl_vector_2d<double>& u,
00104                                 vgl_point_2d<double>& new_p)
00105 {
00106   int n=1+2*search_ni_;
00107   vnl_vector<double> v(n+1);
00108   vgl_vector_2d<double> u1=step_size_*u;
00109   const vgl_point_2d<double> p0 = p-(search_ni_+0.5)*u1;
00110   vimt_sample_profile_bilin(v,image,p0,u1,n+1);
00111   int best_i=0;
00112   double best_e = vcl_fabs(v[1]-v[0]);
00113   for (int i=1;i<n;++i)
00114   {
00115     double e = vcl_fabs(v[i+1]-v[i]);
00116     if (e>best_e) { best_e=e; best_i=i; }
00117   }
00118   new_p = p+(best_i-search_ni_)*u1;
00119   return -1.0 * best_e;
00120 }
00121 
00122 //=======================================================================
00123 // Method: is_a
00124 //=======================================================================
00125 
00126 vcl_string mfpf_edge_finder::is_a() const
00127 {
00128   return vcl_string("mfpf_edge_finder");
00129 }
00130 
00131 //: Create a copy on the heap and return base class pointer
00132 mfpf_point_finder* mfpf_edge_finder::clone() const
00133 {
00134   return new mfpf_edge_finder(*this);
00135 }
00136 
00137 //=======================================================================
00138 // Method: print
00139 //=======================================================================
00140 
00141 void mfpf_edge_finder::print_summary(vcl_ostream& os) const
00142 {
00143   os<<"{ ";
00144   mfpf_point_finder::print_summary(os);
00145   os<<" }";
00146 }
00147 
00148 short mfpf_edge_finder::version_no() const
00149 {
00150   return 1;
00151 }
00152 
00153 
00154 void mfpf_edge_finder::b_write(vsl_b_ostream& bfs) const
00155 {
00156   vsl_b_write(bfs,version_no());
00157   mfpf_point_finder::b_write(bfs);  // Save baseclass
00158 }
00159 
00160 //=======================================================================
00161 // Method: load
00162 //=======================================================================
00163 
00164 void mfpf_edge_finder::b_read(vsl_b_istream& bfs)
00165 {
00166   if (!bfs) return;
00167   short version;
00168   vsl_b_read(bfs,version);
00169   switch (version)
00170   {
00171     case 1:
00172       mfpf_point_finder::b_read(bfs);  // Load in baseclass
00173       break;
00174     default:
00175       vcl_cerr << "I/O ERROR: vsl_b_read(vsl_b_istream&)\n"
00176                << "           Unknown version number "<< version << vcl_endl;
00177       bfs.is().clear(vcl_ios::badbit); // Set an unrecoverable IO error on stream
00178       return;
00179   }
00180 }