contrib/rpl/rgrl/rgrl_trans_homo2d_proj_rad.cxx
Go to the documentation of this file.
00001 #include "rgrl_trans_homo2d_proj_rad.h"
00002 //:
00003 // \file
00004 #include <vcl_cassert.h>
00005 
00006 #include <vnl/vnl_vector_fixed.h>
00007 #include <vnl/vnl_transpose.h>
00008 #include <vnl/vnl_double_2.h>
00009 #include <vnl/vnl_double_2x2.h>
00010 #include <rgrl/rgrl_util.h>
00011 
00012 rgrl_trans_homo2d_proj_rad::
00013 rgrl_trans_homo2d_proj_rad()
00014   : rgrl_est_proj_rad_func<2,2>( 0, true )
00015   , H_( 0.0 )
00016   , rad_k_(0)
00017 {
00018   const vnl_vector_fixed<double, 2> zeroc( 0, 0 );
00019   this->set_centres( zeroc, zeroc, zeroc );
00020 }
00021 
00022 
00023 rgrl_trans_homo2d_proj_rad::
00024 rgrl_trans_homo2d_proj_rad( vnl_matrix_fixed<double, 3, 3> const& H,
00025                             vcl_vector<double>             const& k,
00026                             vnl_vector_fixed< double, 2 >  const& image_centre)
00027   : rgrl_est_proj_rad_func<2,2>( k.size(), true ),
00028     H_(H),
00029     rad_k_(k)
00030 {
00031   const vnl_vector_fixed<double, 2> zeroc( 0, 0 );
00032   this->set_centres( zeroc, zeroc, image_centre );
00033   set_up_rad_k(k);
00034 }
00035 
00036 rgrl_trans_homo2d_proj_rad::
00037 rgrl_trans_homo2d_proj_rad( vnl_matrix_fixed<double, 3, 3> const& H,
00038                             vcl_vector<double>             const& k,
00039                             vnl_vector_fixed< double, 2 >  const& image_centre,
00040                             vnl_matrix<double> const& covar,
00041                             vnl_vector<double> const& from_centre,
00042                             vnl_vector<double> const& to_centre )
00043   : rgrl_transformation( covar ),
00044     rgrl_est_proj_rad_func<2,2>( k.size(), true ),
00045     H_(H),
00046     rad_k_(k)
00047 {
00048   this->set_centres( from_centre, to_centre, image_centre );
00049   set_up_rad_k( k );
00050 }
00051 
00052 rgrl_trans_homo2d_proj_rad::
00053 rgrl_trans_homo2d_proj_rad( vcl_vector<double>             const& norm_k,
00054                             vnl_matrix_fixed<double, 3, 3> const& H,
00055                             vnl_vector_fixed< double, 2 >  const& image_centre )
00056   : rgrl_est_proj_rad_func<2,2>( norm_k.size(), true ),
00057     H_(H),
00058     rad_k_(norm_k)
00059 {
00060   const vnl_vector_fixed<double, 2> zeroc( 0, 0 );
00061   this->set_centres( zeroc, zeroc, image_centre );
00062 }
00063 
00064 rgrl_trans_homo2d_proj_rad::
00065 rgrl_trans_homo2d_proj_rad( vcl_vector<double>             const& norm_k,
00066                             vnl_matrix_fixed<double, 3, 3> const& H,
00067                             vnl_vector_fixed< double, 2 >  const& image_centre,
00068                             vnl_matrix<double> const& covar,
00069                             vnl_vector<double> const& from_centre,
00070                             vnl_vector<double> const& to_centre )
00071   : rgrl_transformation( covar ),
00072     rgrl_est_proj_rad_func<2,2>( norm_k.size(), true ),
00073     H_(H),
00074     rad_k_(norm_k)
00075 {
00076   this->set_centres( from_centre, to_centre, image_centre );
00077 }
00078 
00079 void
00080 rgrl_trans_homo2d_proj_rad::
00081 set_up_rad_k(vcl_vector<double> const & rad_k)
00082 {
00083   rad_k_.resize(rad_k.size());
00084   const double centre_mag_norm = centre_mag_norm_const();
00085 
00086   double base = 1;
00087   for ( unsigned i=0; i<rad_k.size(); ++i ) {
00088     base *= centre_mag_norm;
00089     rad_k_[i] = rad_k[i] * base;
00090   }
00091 }
00092 
00093 void
00094 rgrl_trans_homo2d_proj_rad::
00095 map_loc( vnl_vector<double> const& from,
00096          vnl_vector<double>      & to  ) const
00097 {
00098   // use vnl_double_2 to reduce memory allocation
00099   vnl_double_2 pt = from;
00100   vnl_double_2 mapped;
00101   rgrl_est_proj_rad_func<2,2>::map_loc( mapped, H_, rad_k_, pt );
00102   to = mapped.as_ref();
00103 }
00104 
00105 void
00106 rgrl_trans_homo2d_proj_rad::
00107 map_dir( vnl_vector<double> const& from_loc,
00108          vnl_vector<double> const& from_dir,
00109          vnl_vector<double>      & to_dir  ) const
00110 {
00111   assert ( from_loc.size() == 2 );
00112   assert ( from_dir.size() == 2 );
00113 
00114   const vnl_double_2 from_begin( from_loc );
00115   vnl_double_2 from_end( from_loc );
00116   from_end += from_dir;
00117 
00118   vnl_double_2 to_loc_begin, to_loc_end;
00119   rgrl_est_proj_rad_func<2,2>::map_loc( to_loc_begin, H_, rad_k_, from_begin );
00120   rgrl_est_proj_rad_func<2,2>::map_loc( to_loc_end,   H_, rad_k_, from_end );
00121 
00122   to_dir = (to_loc_end - to_loc_begin).as_ref();
00123   to_dir.normalize();
00124 }
00125 
00126 rgrl_transformation_sptr
00127 rgrl_trans_homo2d_proj_rad::
00128 scale_by( double scale ) const
00129 {
00130   vnl_matrix_fixed<double,3,3> new_H( H_ );
00131 
00132   // scale
00133   new_H(0,2) *= scale;
00134   new_H(1,2) *= scale;
00135 
00136   // move the scale on the fixed coordinate,
00137   // and divide the 3rd row by this scale
00138   new_H(2,0) /= scale;
00139   new_H(2,1) /= scale;
00140 
00141   // normalize
00142   new_H /= new_H.fro_norm();
00143 
00144   // centers
00145   const vnl_vector_fixed<double,2> from = from_centre_ * scale;
00146   const vnl_vector_fixed<double,2> to = to_centre_ * scale;
00147   const vnl_vector_fixed<double,2> ic = image_centre_ * scale;
00148 
00149   // radial terms
00150   vcl_vector<double> radk = radial_params();
00151   const double sq_scale = scale*scale;
00152 
00153   double base = 1;
00154   for ( unsigned int i=0; i<radk.size(); ++i ) {
00155     base *= sq_scale;
00156     radk[i] /= base;
00157   }
00158 
00159   rgrl_transformation_sptr xform
00160     =  new rgrl_trans_homo2d_proj_rad( new_H,
00161                                        radk,
00162                                        ic,
00163                                        vnl_matrix<double>(),
00164                                        from.as_ref(), to.as_ref() );
00165   xform->set_scaling_factors( this->scaling_factors() );
00166   return xform;
00167 }
00168 
00169 
00170 vnl_matrix_fixed<double, 3, 3>
00171 rgrl_trans_homo2d_proj_rad::
00172 H( ) const
00173 {
00174   return uncentre_proj( H_ );
00175 }
00176 
00177 vnl_matrix<double>
00178 rgrl_trans_homo2d_proj_rad::
00179 transfer_error_covar( vnl_vector<double> const& p ) const
00180 {
00181   assert ( p.size() == 2 );
00182   vnl_matrix<double> jac;
00183   full_proj_rad_jacobian( jac, H_, rad_k_, p );
00184   return jac*covar_*vnl_transpose(jac);
00185 }
00186 
00187 //: Return the jacobian of the transform.
00188 void
00189 rgrl_trans_homo2d_proj_rad::
00190 jacobian_wrt_loc( vnl_matrix<double>& jacobian, vnl_vector<double> const& from_loc ) const
00191 {
00192   vnl_double_2x2 jac_loc;
00193   proj_jac_wrt_loc( jac_loc, H_, rad_k_, from_loc );
00194   jacobian = jac_loc.as_ref();
00195 }
00196 
00197 // for output CENTERED transformation
00198 void
00199 rgrl_trans_homo2d_proj_rad::
00200 write(vcl_ostream& os ) const
00201 {
00202   // tag
00203   os << "HOMOGRAPHY2D+RADIAL\n"
00204   // parameters
00205      << H_ << vcl_endl;
00206 
00207   // radial terms
00208   os << rad_k_.size() << "  ";
00209   for ( unsigned i=0; i<rad_k_.size(); ++i )
00210     os << rad_k_[i] << ' ';
00211 
00212   os << '\n' << from_centre_ << "   " << to_centre_
00213      << '\n' << image_centre_ << "   " << centre_mag_norm_const_
00214      << vcl_endl;
00215 
00216   // parent
00217   rgrl_transformation::write( os );
00218 }
00219 
00220 // for input
00221 bool
00222 rgrl_trans_homo2d_proj_rad::
00223 read(vcl_istream& is )
00224 {
00225   // skip empty lines
00226   rgrl_util_skip_empty_lines( is );
00227 
00228   vcl_string str;
00229   vcl_getline( is, str );
00230 
00231   // The token should appear at the beginning of line
00232   if ( str.find( "HOMOGRAPHY2D+RADIAL" ) != 0 ) {
00233    return false;
00234   }
00235 
00236   // input global xform
00237   is >> H_;
00238 
00239   // input radial terms
00240   {
00241     int size = -1;
00242     is >> size;
00243     if ( size < 0 || !is.good() )
00244       return false;
00245 
00246     rad_k_.resize( size );
00247     for ( int i=0; i<size; ++i )
00248       is >> rad_k_[i];
00249 
00250     // set H_dof_
00251     camera_dof_ = size;
00252   }
00253 
00254   is >> from_centre_ >> to_centre_;
00255   is >> image_centre_;
00256   is >> centre_mag_norm_const_;
00257 
00258   // parent
00259   return is.good() && rgrl_transformation::read( is );
00260 }
00261 
00262 //: make a clone copy
00263 rgrl_transformation_sptr
00264 rgrl_trans_homo2d_proj_rad::
00265 clone() const
00266 {
00267   return new rgrl_trans_homo2d_proj_rad( *this );
00268 }
00269 
00270 //: return radial parameters
00271 vcl_vector<double>
00272 rgrl_trans_homo2d_proj_rad::
00273 radial_params() const
00274 {
00275   vcl_vector<double> ori_radk( rad_k_ );
00276 
00277   const double centre_mag_norm = centre_mag_norm_const();
00278 
00279   double base = 1;
00280   for ( unsigned i=0; i<ori_radk.size(); ++i ) {
00281     base *= centre_mag_norm;
00282     ori_radk[i] /= base;
00283   }
00284   return ori_radk;
00285 }