contrib/mul/mbl/mbl_thin_plate_spline_weights_3d.h
Go to the documentation of this file.
00001 // This is mul/mbl/mbl_thin_plate_spline_weights_3d.h
00002 #ifndef mbl_thin_plate_spline_weights_3d_h_
00003 #define mbl_thin_plate_spline_weights_3d_h_
00004 //:
00005 // \file
00006 // \brief Construct thin plate spline to map 3D to 3D
00007 // \author Tim Cootes
00008 
00009 #include <vgl/vgl_point_3d.h>
00010 #include <vcl_vector.h>
00011 #include <vnl/vnl_vector.h>
00012 #include <vnl/vnl_matrix.h>
00013 #include <vsl/vsl_binary_io.h>
00014 
00015 //=======================================================================
00016 //: Construct thin plate spline to map 3D to 3D.
00017 // I.e. does some mapping (x',y',z') = f(x,y,z). (See Booksteins work, e.g. IPMI 1993)
00018 // The warp is `guided' by a set of
00019 // landmarks p(0) .. p(n-1) in the source plane which are to be
00020 // mapped to a (possibly deformed) set q(0)..q(n-1) in the destination.
00021 // Thus the mapping is constrained so that f(p(i)) = q(i) for i = 0..n-1.
00022 // The points are given to the build() function to set up the object.
00023 //
00024 // If one wishes to map a set of source points to multiple target points,
00025 // use set_source_pts(src_pts);  then build(target_pts); for each target set.
00026 class mbl_thin_plate_spline_weights_3d
00027 {
00028   vnl_vector<double> Wx_,Wy_,Wz_;
00029   double Ax0_, AxX_, AxY_, AxZ_;
00030   double Ay0_, AyX_, AyY_, AyZ_;
00031   double Az0_, AzX_, AzY_, AzZ_;
00032   double energy_x_,energy_y_,energy_z_;
00033 
00034   vcl_vector<vgl_point_3d<double> > src_pts_, pt_wts_;
00035 
00036     //: Used to estimate weights in set_source_points()
00037   vnl_matrix<double> L_inv_;
00038 
00039     //: Build from small number of points
00040   void build_pure_affine(const vcl_vector<vgl_point_3d<double> >& source_pts,
00041                          const vcl_vector<vgl_point_3d<double> >& dest_pts);
00042 
00043    //: Set parameters from vectors
00044   void set_params(const vnl_vector<double>& W1,
00045                   const vnl_vector<double>& W2,
00046                   const vnl_vector<double>& W3);
00047 
00048   void set_up_rhs(vnl_vector<double>& Bx,
00049                   vnl_vector<double>& By,
00050                   vnl_vector<double>& Bz,
00051                   const vcl_vector<vgl_point_3d<double> >& dest_pts);
00052 
00053    //: Compute spline-bending energy
00054   void compute_energy(vnl_vector<double>& W1,
00055                   vnl_vector<double>& W2,
00056                   vnl_vector<double>& W3,
00057                   const vnl_matrix<double>& L);
00058 
00059  public:
00060 
00061     //: Dflt ctor
00062   mbl_thin_plate_spline_weights_3d();
00063 
00064     //: Destructor
00065   virtual ~mbl_thin_plate_spline_weights_3d();
00066 
00067     //: Sets up internal transformation to map source_pts onto dest_pts
00068   void build(const vcl_vector<vgl_point_3d<double> >& source_pts,
00069         const vcl_vector<vgl_point_3d<double> >& dest_pts,
00070                    bool compute_the_energy=false);
00071 
00072     //: Define source point positions
00073     //  Performs pre-computations so that build(dest_points) can be
00074     //  called multiple times efficiently
00075   void set_source_pts(const vcl_vector<vgl_point_3d<double> >& source_pts);
00076 
00077     //: Define source point weights
00078     //  Sets x, y, z weights for each of the source points
00079   void set_pt_wts( const vcl_vector<vgl_point_3d<double> >& pt_wts );
00080 
00081     //: Sets up internal transformation to map source_pts onto dest_pts
00082   void build(const vcl_vector<vgl_point_3d<double> >& dest_pts);
00083 
00084        //: Return transformed version of (x,y,z)
00085   vgl_point_3d<double>  operator()(double x, double y, double z) const;
00086 
00087        //: Return transformed version of (x,y,z)
00088   vgl_point_3d<double>  operator()(const vgl_point_3d<double>&  p) const
00089   { return operator()(p.x(),p.y(),p.z()); }
00090 
00091     //: Bending energy of X component (zero for pure affine)
00092     //  A measure of total amount of non-linear deformation
00093   double bendingEnergyX() const { return energy_x_; }
00094 
00095     //: Bending energy of Y component (zero for pure affine)
00096     //  A measure of total amount of non-linear deformation
00097   double bendingEnergyY() const { return energy_y_; }
00098 
00099     //: Bending energy of Z component (zero for pure affine)
00100     //  A measure of total amount of non-linear deformation
00101   double bendingEnergyZ() const { return energy_z_; }
00102 
00103     //: Version number for I/O
00104   short version_no() const;
00105 
00106     //: Print class to os
00107   void print_summary(vcl_ostream& os) const;
00108 
00109     //: Save class to binary file stream
00110   void b_write(vsl_b_ostream& bfs) const;
00111 
00112     //: Load class from binary file stream
00113   void b_read(vsl_b_istream& bfs);
00114 
00115     //: Comparison operator
00116   bool operator==(const mbl_thin_plate_spline_weights_3d& tps) const;
00117 };
00118 
00119   //: Binary file stream output operator for class reference
00120 void vsl_b_write(vsl_b_ostream& bfs, const mbl_thin_plate_spline_weights_3d& b);
00121 
00122   //: Binary file stream input operator for class reference
00123 void vsl_b_read(vsl_b_istream& bfs, mbl_thin_plate_spline_weights_3d& b);
00124 
00125   //: Stream output operator for class reference
00126 vcl_ostream& operator<<(vcl_ostream& os,const mbl_thin_plate_spline_weights_3d& b);
00127 
00128 #endif
00129 
00130