contrib/oxl/mvl/ImageMetric.cxx
Go to the documentation of this file.
00001 // This is oxl/mvl/ImageMetric.cxx
00002 #ifdef VCL_NEEDS_PRAGMA_INTERFACE
00003 #pragma implementation
00004 #endif
00005 //:
00006 // \file
00007 
00008 #include "ImageMetric.h"
00009 
00010 #include <vcl_iostream.h>
00011 #include <vcl_cassert.h>
00012 #include <vnl/vnl_double_2.h>
00013 #include <vnl/vnl_double_3x3.h>
00014 #include <vnl/vnl_transpose.h>
00015 #include <vnl/vnl_math.h>
00016 #include <vgl/vgl_point_2d.h>
00017 #include <vgl/vgl_homg_point_2d.h>
00018 #include <vgl/vgl_homg_line_2d.h>
00019 #include <vgl/algo/vgl_homg_operators_2d.h>
00020 #include <vgl/vgl_line_segment_2d.h>
00021 
00022 #include <mvl/HomgPoint2D.h>
00023 #include <mvl/HomgLineSeg2D.h>
00024 #include <mvl/HomgOperator2D.h>
00025 
00026 static vcl_ostream& warning(char const * fn)
00027 {
00028   return vcl_cerr << fn << " WARNING: ";
00029 }
00030 
00031 // TRANSFORMATIONS
00032 
00033 //: Condition the 2D point p.
00034 //  Default implementation is simply to return p in homogeneous coordinates
00035 vgl_homg_point_2d<double> ImageMetric::image_to_homg(vgl_point_2d<double> const& p) const
00036 {
00037   assert(!"ImageMetric::image_to_homg should be implemented for efficiency");
00038   return vgl_homg_point_2d<double>(p.x(), p.y(), 1.0);
00039 }
00040 
00041 //: Condition the 2D point p.
00042 //  Default implementation is simply to return p in homogeneous coordinates
00043 HomgPoint2D ImageMetric::image_to_homg(const vnl_double_2& p) const
00044 {
00045   assert(!"ImageMetric::image_to_homg should be implemented for efficiency");
00046   return HomgPoint2D(p[0], p[1], 1.0);
00047 }
00048 
00049 //: Condition 2D point (x,y)
00050 HomgPoint2D ImageMetric::image_to_homg(double x, double y) const
00051 {
00052   assert(!"ImageMetric::image_to_homg should be implemented for efficiency");
00053   return HomgPoint2D(x, y, 1.0);
00054 }
00055 
00056 //: Convert conditioned point p to image coordinates
00057 vgl_point_2d<double> ImageMetric::homg_to_image(vgl_homg_point_2d<double> const& p) const
00058 {
00059   assert(!"ImageMetric::homg_to_image should be implemented for efficiency");
00060   return p;
00061 }
00062 
00063 vnl_double_2 ImageMetric::homg_to_image(const HomgPoint2D& p) const
00064 {
00065   assert(!"ImageMetric::homg_to_image should be implemented for efficiency");
00066   double x,y;
00067   p.get_nonhomogeneous(x, y);
00068   return vnl_double_2(x, y);
00069 }
00070 
00071 //: Convert homogeneous point in image coordinates to one in conditioned coordinates
00072 vgl_homg_point_2d<double> ImageMetric::imagehomg_to_homg(vgl_homg_point_2d<double> const& x) const
00073 {
00074   return x;
00075 }
00076 
00077 HomgPoint2D ImageMetric::imagehomg_to_homg(const HomgPoint2D& x) const
00078 {
00079   return x;
00080 }
00081 
00082 //: Convert homogeneous point in conditioned coordinates to one in image coordinates
00083 vgl_homg_point_2d<double> ImageMetric::homg_to_imagehomg(vgl_homg_point_2d<double> const& x) const
00084 {
00085   return x;
00086 }
00087 
00088 HomgPoint2D ImageMetric::homg_to_imagehomg(const HomgPoint2D& x) const
00089 {
00090   return x;
00091 }
00092 
00093 //: Convert homogeneous line in conditioned coordinates to one in image coordinates
00094 vgl_homg_line_2d<double> ImageMetric::homg_to_image_line(vgl_homg_line_2d<double> const& l) const
00095 {
00096   if (is_linear())
00097     return get_C_inverse().transpose() * l;
00098 
00099   // get points, decondition, and rejoin
00100   vgl_homg_point_2d<double> p1, p2;
00101   l.get_two_points(p1,p2);
00102   vgl_homg_point_2d<double> i1 = homg_to_imagehomg(p1);
00103   vgl_homg_point_2d<double> i2 = homg_to_imagehomg(p2);
00104   return vgl_homg_line_2d<double>(i1,i2);
00105 }
00106 
00107 HomgLine2D ImageMetric::homg_to_image_line(const HomgLine2D& l) const
00108 {
00109   if (is_linear())
00110     return HomgLine2D(vnl_transpose(get_C_inverse()) * l.get_vector().as_ref());
00111 
00112   // get points, decondition, and rejoin
00113   HomgPoint2D p1, p2;
00114   l.get_2_points_on_line(&p1, &p2);
00115   HomgPoint2D i1 = homg_to_imagehomg(p1);
00116   HomgPoint2D i2 = homg_to_imagehomg(p2);
00117   return HomgOperator2D::join(i1, i2);
00118 }
00119 
00120 vgl_homg_line_2d<double> ImageMetric::image_to_homg_line(const vgl_homg_line_2d<double>& l) const
00121 {
00122   if (is_linear())
00123     return get_C().transpose() * l;
00124 
00125   // get points, condition, and rejoin
00126   vgl_homg_point_2d<double> p1, p2;
00127   l.get_two_points(p1,p2);
00128   vgl_homg_point_2d<double> i1 = imagehomg_to_homg(p1);
00129   vgl_homg_point_2d<double> i2 = imagehomg_to_homg(p2);
00130   return vgl_homg_line_2d<double>(i1,i2);
00131 }
00132 
00133 HomgLine2D ImageMetric::image_to_homg_line(const HomgLine2D& l) const
00134 {
00135   if (is_linear())
00136     return HomgLine2D(vnl_transpose(get_C()) * l.get_vector().as_ref());
00137 
00138   // get points, condition, and rejoin
00139   HomgPoint2D p1, p2;
00140   l.get_2_points_on_line(&p1, &p2);
00141   HomgPoint2D i1 = imagehomg_to_homg(p1);
00142   HomgPoint2D i2 = imagehomg_to_homg(p2);
00143   return HomgOperator2D::join(i1, i2);
00144 }
00145 
00146 
00147 //: Convert homogeneous line segment in conditioned coordinates to one in image coordinates
00148 vgl_line_segment_2d<double> ImageMetric::homg_line_to_image(vgl_line_segment_2d<double> const& l) const
00149 {
00150   // get points, decondition, and rejoin
00151   vgl_homg_point_2d<double> i1 = homg_to_imagehomg(vgl_homg_point_2d<double>(l.point1()));
00152   vgl_homg_point_2d<double> i2 = homg_to_imagehomg(vgl_homg_point_2d<double>(l.point2()));
00153   return vgl_line_segment_2d<double>(i1,i2);
00154 }
00155 
00156 HomgLineSeg2D ImageMetric::homg_line_to_image(const HomgLineSeg2D& l) const
00157 {
00158   // get points, decondition, and rejoin
00159   HomgPoint2D i1 = homg_to_imagehomg(l.get_point1());
00160   HomgPoint2D i2 = homg_to_imagehomg(l.get_point2());
00161   return HomgLineSeg2D(i1, i2);
00162 }
00163 
00164 vgl_line_segment_2d<double> ImageMetric::image_to_homg_line(vgl_line_segment_2d<double> const& l) const
00165 {
00166   vgl_homg_point_2d<double> i1 = imagehomg_to_homg(vgl_homg_point_2d<double>(l.point1()));
00167   vgl_homg_point_2d<double> i2 = imagehomg_to_homg(vgl_homg_point_2d<double>(l.point2()));
00168   return vgl_line_segment_2d<double>(i1,i2);
00169 }
00170 
00171 HomgLineSeg2D ImageMetric::image_to_homg_line(const HomgLineSeg2D& l) const
00172 {
00173   HomgPoint2D i1 = imagehomg_to_homg(l.get_point1());
00174   HomgPoint2D i2 = imagehomg_to_homg(l.get_point2());
00175   return HomgLineSeg2D(i1, i2);
00176 }
00177 
00178 // == MEASUREMENTS ==
00179 
00180 //: Compute perpendicular distance in image coordinates between point p and line l, expressed in conditioned coordinates.
00181 double ImageMetric::perp_dist_squared(vgl_homg_point_2d<double> const& p,
00182                                       vgl_homg_line_2d<double> const& l) const
00183 {
00184   return vgl_homg_operators_2d<double>::perp_dist_squared(homg_to_imagehomg(p), homg_to_image_line(l));
00185 }
00186 
00187 double ImageMetric::perp_dist_squared(const HomgPoint2D& p, const HomgLine2D& l) const
00188 {
00189   return HomgOperator2D::perp_dist_squared(homg_to_imagehomg(p), homg_to_image_line(l));
00190 }
00191 
00192 //: Project point onto line.
00193 vgl_homg_point_2d<double> ImageMetric::perp_projection(vgl_homg_line_2d<double> const& l,
00194                                                        vgl_homg_point_2d<double> const& p) const
00195 {
00196   if (l.ideal())
00197     vcl_cerr << "ImageMetric::perp_projection -- line at infinity\n";
00198 
00199   return vgl_homg_operators_2d<double>::perp_projection(homg_to_image_line(l), homg_to_imagehomg(p));
00200 }
00201 
00202 HomgPoint2D ImageMetric::perp_projection(const HomgLine2D & l, const HomgPoint2D & p) const
00203 {
00204   if (p.ideal()) {
00205     vcl_cerr << "ImageMetric::perp_projection -- point at infinity\n";
00206   }
00207 
00208   if (l.ideal()) {
00209     vcl_cerr << "ImageMetric::perp_projection -- line at infinity\n";
00210   }
00211 
00212   return HomgOperator2D::perp_projection(homg_to_image_line(l), homg_to_imagehomg(p));
00213 }
00214 
00215 //: Get perpendicular distance in image.
00216 double ImageMetric::distance_squared(vgl_homg_point_2d<double> const& p1,
00217                                      vgl_homg_point_2d<double> const& p2) const
00218 {
00219   return vgl_homg_operators_2d<double>::distance_squared(homg_to_imagehomg(p1), homg_to_imagehomg(p2));
00220 }
00221 
00222 double ImageMetric::distance_squared(const HomgPoint2D & p1, const HomgPoint2D & p2) const
00223 {
00224   return HomgOperator2D::distance_squared(homg_to_imagehomg(p1), homg_to_imagehomg(p2));
00225 }
00226 
00227 //: Get distance between a line segment and an infinite line.
00228 //  The metric used is the maximum of the two endpoint perp distances.
00229 double ImageMetric::distance_squared(vgl_line_segment_2d<double> const& segment,
00230                                      vgl_homg_line_2d<double> const& line) const
00231 {
00232   // ca_distance_squared_lineseg_to_line
00233   return vnl_math_max(this->perp_dist_squared(vgl_homg_point_2d<double>(segment.point1()), line),
00234                       this->perp_dist_squared(vgl_homg_point_2d<double>(segment.point2()), line));
00235 }
00236 
00237 double ImageMetric::distance_squared(const HomgLineSeg2D& segment, const HomgLine2D& line) const
00238 {
00239   // ca_distance_squared_lineseg_to_line
00240   return vnl_math_max(this->perp_dist_squared(segment.get_point1(), line),
00241                       this->perp_dist_squared(segment.get_point2(), line));
00242 }
00243 
00244 bool ImageMetric::is_within_distance(vgl_homg_point_2d<double> const& p1,
00245                                      vgl_homg_point_2d<double> const& p2,
00246                                      double distance) const
00247 {
00248   return distance_squared(p1, p2) < distance*distance;
00249 }
00250 
00251 bool ImageMetric::is_within_distance(const HomgPoint2D& p1, const HomgPoint2D& p2, double distance) const
00252 {
00253   return distance_squared(p1, p2) < distance*distance;
00254 }
00255 
00256 //: Convert a distance in image coordinates to one in conditioned coordinates.
00257 // This is only possible for similarity transformations, but where it does make
00258 // sense it can mean significant increases in speed.
00259 double ImageMetric::homg_to_image_distance(double image_distance) const
00260 {
00261   warning("ImageMetric::invert_distance()") << "returning 0\n";
00262   return image_distance;
00263 }
00264 
00265 //: Convert a distance in image coordinates to one in conditioned coordinates.
00266 // This is only possible for similarity transformations, but where it does make
00267 // sense it can mean significant increases in speed.
00268 double ImageMetric::image_to_homg_distance(double image_distance) const
00269 {
00270   warning("ImageMetric::invert_distance()") << "returning 0\n";
00271   return image_distance;
00272 }
00273 
00274 //: Return true if the invert_distance function makes sense.
00275 bool ImageMetric::can_invert_distance() const
00276 {
00277   return false;
00278 }
00279 
00280 // == MATRIX REPRESENTATION ==
00281 
00282 //: Return true if the action of the conditioner can be represented as a planar homography.
00283 bool ImageMetric::is_linear() const
00284 {
00285   return false;
00286 }
00287 
00288 #include <vnl/vnl_identity_3x3.h>
00289 
00290 //: Return conditioning matrix C that converts homogeneous image points to homogeneous conditioned points.
00291 //  If the ImageMetric used is nonlinear, then we'll have to make other arrangements...
00292 vnl_double_3x3 ImageMetric::get_C() const
00293 {
00294   static vnl_identity_3x3 I;
00295   warning("ImageMetric::get_C()") << "returning identity\n";
00296   return I;
00297 }
00298 
00299 //: Return conditioning matrix C that converts homogeneous conditioned points to image coordinates.
00300 vnl_double_3x3 ImageMetric::get_C_inverse() const
00301 {
00302   static vnl_identity_3x3 I;
00303   warning("ImageMetric::get_C_inverse()") << "returning identity\n";
00304   return I;
00305 }
00306 
00307 #include <mvl/FMatrix.h>
00308 
00309 FMatrix ImageMetric::decondition(const FMatrix& F, const ImageMetric* c1, const ImageMetric* c2)
00310 {
00311   if (!c1->is_linear())
00312     warning("ImageMetric::decondition(FMatrix...)") << "ImageMetric for image 1 is nonlinear\n";
00313 
00314   if (!c2->is_linear())
00315     warning("ImageMetric::decondition(FMatrix...)") << "ImageMetric for image 2 is nonlinear\n";
00316 
00317   vnl_double_3x3 C1inv = c1->get_C_inverse();
00318   vnl_double_3x3 C2inv = c2->get_C_inverse();
00319   return FMatrix( vnl_transpose(C1inv) * F.get_matrix() * C2inv );
00320 }
00321 
00322 vcl_ostream & ImageMetric::print(vcl_ostream& s) const
00323 {
00324   return s << "Empty ImageMetric";
00325 }