PLearn 0.1
|
00001 // -*- C++ -*- 00002 00003 // ChemicalICP.cc 00004 // 00005 // Copyright (C) 2006 Pascal Lamblin 00006 // 00007 // Redistribution and use in source and binary forms, with or without 00008 // modification, are permitted provided that the following conditions are met: 00009 // 00010 // 1. Redistributions of source code must retain the above copyright 00011 // notice, this list of conditions and the following disclaimer. 00012 // 00013 // 2. Redistributions in binary form must reproduce the above copyright 00014 // notice, this list of conditions and the following disclaimer in the 00015 // documentation and/or other materials provided with the distribution. 00016 // 00017 // 3. The name of the authors may not be used to endorse or promote 00018 // products derived from this software without specific prior written 00019 // permission. 00020 // 00021 // THIS SOFTWARE IS PROVIDED BY THE AUTHORS ``AS IS'' AND ANY EXPRESS OR 00022 // IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES 00023 // OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN 00024 // NO EVENT SHALL THE AUTHORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, 00025 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED 00026 // TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR 00027 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF 00028 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING 00029 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS 00030 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 00031 // 00032 // This file is part of the PLearn library. For more information on the PLearn 00033 // library, go to the PLearn Web site at www.plearn.org 00034 00035 // Authors: Pascal Lamblin 00036 00039 #include "ChemicalICP.h" 00040 #include <plearn/base/stringutils.h> 00041 #include <plearn/io/openFile.h> 00042 #include <plearn/var/UnaryVariable.h> 00043 #include <plearn/var/VarColumnsVariable.h> 00044 #include <plearn/io/pl_log.h> 00045 #include "geometry.h" 00046 00047 namespace PLearn { 00048 using namespace std; 00049 00050 PLEARN_IMPLEMENT_OBJECT( 00051 ChemicalICP, 00052 "Iterative Closest Points (ICP) aligner.", 00053 "" 00054 ); 00055 00056 ChemicalICP::ChemicalICP(): 00057 mol_feat_indices( new UnaryVariable() ), 00058 template_feat_indices( new UnaryVariable() ), 00059 matching_neighbors( new UnaryVariable() ), 00060 weighting_method( "features_sigmoid" ), 00061 weighting_params( Vec(2,1) ), 00062 matching_method( "exhaustive" ), 00063 memory( "none" ), 00064 initial_angles_step( 0 ), 00065 max_iter( 50 ), 00066 error_t( 0 ), 00067 angle_t( 0.5 ), 00068 trans_t( 0 ), 00069 rotation( 3, 3 ), 00070 translation( 3 ), 00071 all_mol_features( 0 ), 00072 all_template_features( 0 ), 00073 all_template_feat_dev( 0 ), 00074 mol_coordinates( 0 ), 00075 used_mol_features( new VarColumnsVariable( all_mol_features, 00076 mol_feat_indices ) ), 00077 template_coordinates( 0 ), 00078 template_geom_dev( 0 ), 00079 used_template_features( new VarColumnsVariable(all_template_features, 00080 template_feat_indices) ), 00081 used_template_feat_dev( new VarColumnsVariable(all_template_feat_dev, 00082 template_feat_indices) ) 00083 {} 00084 00085 ChemicalICP::ChemicalICP( const MolTemplate& the_template, 00086 const Mol& the_molecule, 00087 const TVec<string>& the_feature_names, 00088 string the_weighting_method, 00089 const Var& the_weighting_params, 00090 string the_matching_method ): 00091 mol_feat_indices( new UnaryVariable() ), 00092 template_feat_indices( new UnaryVariable() ), 00093 matching_neighbors( new UnaryVariable() ), 00094 mol_template(the_template), 00095 molecule(the_molecule), 00096 feature_names(the_feature_names), 00097 weighting_method(the_weighting_method), 00098 weighting_params(the_weighting_params), 00099 matching_method( the_matching_method ), 00100 memory( "none" ), 00101 initial_angles_step( 0 ), 00102 max_iter( 50 ), 00103 error_t( 0 ), 00104 angle_t( 0.5 ), 00105 trans_t( 0 ), 00106 rotation( 3, 3 ), 00107 translation( 3 ), 00108 all_mol_features( 0 ), // calls SourceVariable() 00109 all_template_features( 0 ), 00110 all_template_feat_dev( 0 ), 00111 mol_coordinates( 0 ), 00112 used_mol_features( new VarColumnsVariable( all_mol_features, 00113 mol_feat_indices ) ), 00114 template_coordinates( 0 ), 00115 template_geom_dev( 0 ), 00116 used_template_features( new VarColumnsVariable(all_template_features, 00117 template_feat_indices) ), 00118 used_template_feat_dev( new VarColumnsVariable(all_template_feat_dev, 00119 template_feat_indices) ) 00120 { 00121 //pout << "Big constructor called" << endl; 00122 build(); 00123 } 00124 00125 void ChemicalICP::setMolecule( const Mol& the_molecule ) 00126 { 00127 molecule = the_molecule; 00128 computeUsedFeatures(); 00129 computeVariables(); 00130 00131 if( matching_method == "exhaustive" ) 00132 cacheFeatureDistances(); 00133 } 00134 00135 void ChemicalICP::computeUsedFeatures() 00136 { 00137 // compute intersection between mol_template->feature_names and 00138 // feature_names 00139 TVec<string> template_names = mol_template->feature_names; 00140 TVec<string> common_names; 00141 TVec<int> common_indices; 00142 00143 if( feature_names.length() == 0 ) 00144 { 00145 common_names = template_names.copy(); 00146 common_indices = TVec<int>( 0, template_names.length() - 1, 1 ); 00147 } 00148 else if( feature_names[0] == "none" ) 00149 { 00150 used_feat_names.resize( 0 ); 00151 mol_feat_indices->resize( 0, 0 ); 00152 template_feat_indices->resize( 0, 0 ); 00153 return; 00154 } 00155 else 00156 { 00157 common_names.resize(0, feature_names.length()); // to have some space 00158 common_indices.resize(0, feature_names.length()); 00159 TVec<int> indices = feature_names.find( template_names ); 00160 00161 for( int i=0 ; i<indices.length() ; i++ ) 00162 { 00163 if( indices[i] >= 0 ) 00164 { 00165 // template_names[i] is present in feature_names 00166 common_names.append( template_names[i] ); 00167 common_indices.append( i ); 00168 } 00169 } 00170 } 00171 00172 // then compute intersection with molecule->feature_names 00173 TVec<string> mol_names = molecule->feature_names; 00174 used_feat_names.resize(0, common_names.length()); 00175 00176 // storage for indices 00177 Vec template_feat_indices_; 00178 Vec mol_feat_indices_; 00179 00180 TVec<int> indices = mol_names.find( common_names ); 00181 for( int i=0 ; i<indices.length() ; i++ ) 00182 { 00183 if( indices[i] >= 0 ) 00184 { 00185 // common_names[i] is present in mol_names at position indices[i] 00186 used_feat_names.append( common_names[i] ); 00187 template_feat_indices_.append( common_indices[i] ); 00188 mol_feat_indices_.append( indices[i] ); 00189 } 00190 } 00191 00192 // update Var's value 00193 template_feat_indices->defineValueLocation( 00194 template_feat_indices_.toMat( template_feat_indices_.size(), 1 ) ); 00195 00196 mol_feat_indices->defineValueLocation( 00197 mol_feat_indices_.toMat( mol_feat_indices_.size(), 1 ) ); 00198 } 00199 00200 void ChemicalICP::computeVariables() 00201 { 00202 // make mol_coordinates point to the new molecule's coordinates 00203 mol_coordinates->defineValueLocation( molecule->coordinates ); 00204 00205 // make mol_features point to the new molecule's features 00206 all_mol_features->defineValueLocation( molecule->features ); 00207 00208 // update used_... Var's so they have the right size 00209 used_mol_features->sizefprop(); 00210 used_template_features->sizefprop(); 00211 used_template_feat_dev->sizefprop(); 00212 } 00213 00214 void ChemicalICP::cacheFeatureDistances() 00215 { 00216 //pout << "cacheFeatureDistances()" << endl; 00217 Mat t_features = used_template_features->matValue; 00218 Mat t_feat_dev = used_template_feat_dev->matValue; 00219 Mat m_features = used_mol_features->matValue; 00220 00221 int n_template_points = t_features.length(); 00222 int n_mol_points = m_features.length(); 00223 int n_features = m_features.width(); 00224 00225 feat_distances2.resize( n_template_points, n_mol_points ); 00226 feat_distances2.fill(0); 00227 00228 for( int i=0 ; i<n_template_points ; i++ ) 00229 { 00230 for( int j=0 ; j<n_mol_points ; j++ ) 00231 { 00232 for( int k=0 ; k<n_features ; k++ ) 00233 { 00234 real diff = (t_features(i, k) - m_features(j, k)) 00235 / t_feat_dev(i, k); 00236 feat_distances2(i, j) += diff * diff; 00237 } 00238 } 00239 } 00240 /*pout << "t_features(0) = " << t_features(0) << endl; 00241 pout << "t_feat_dev(0) = " << t_feat_dev(0) << endl; 00242 pout << "feat_distances2(0) = " << feat_distances2(0) << endl;*/ 00243 } 00244 00246 // forgetMemorizedAlignments // 00248 void ChemicalICP::forgetMemorizedAlignments() 00249 { 00250 mem_aligned_molecules.resize(0); 00251 mem_error.resize(0); 00252 mem_rotation.resize(0); 00253 mem_translation.resize(0); 00254 mem_matching.resize(0); 00255 mem_count.resize(0); 00256 } 00257 00259 // run // 00261 void ChemicalICP::run() 00262 { 00263 if( initial_angles_list.length() > 0 ) 00264 { 00265 if( initial_angles_step > 0 ) 00266 { 00267 PLWARNING( "ChemicalICP::run - both 'initial_angles_step' and" 00268 " 'initial_angles_list'\n" 00269 "are provided. Setting 'initial_angles_step' to 0.\n" ); 00270 initial_angles_step = 0; 00271 } 00272 } 00273 else 00274 { 00275 if( fast_is_equal( initial_angles_step, 0. ) ) 00276 initial_angles_step = 360; 00277 00278 initial_angles_list.resize( 0, 3 ); 00279 for( real rx=0. ; rx<360. ; rx += initial_angles_step ) 00280 for( real ry=0. ; ry<360. ; ry += initial_angles_step ) 00281 for( real rz=0. ; rz<180. ; rz += initial_angles_step ) 00282 { 00283 Vec angles( 3 ); 00284 angles[0] = rx; 00285 angles[1] = ry; 00286 angles[2] = rz; 00287 initial_angles_list.appendRow( angles ); 00288 } 00289 } 00290 EXTREME_NAMED_LOG("ICP+NNet") << "begin run()" << endl; 00291 00292 // Only perform the alignment if we do not just obtain the result from a 00293 // previously memorized alignment. 00294 bool realign = true; 00295 int find_mol = -1; 00296 if (memory != "none") { 00297 find_mol = mem_aligned_molecules.find(molecule); 00298 if (find_mol == -1) { 00299 // This molecule had not been previously aligned. 00300 find_mol = mem_aligned_molecules.length(); 00301 mem_aligned_molecules.append(molecule); 00302 // Dummy additions just to resize the vectors correctly. 00303 mem_error.append(0); 00304 mem_rotation.append(Mat()); 00305 mem_translation.append(Vec()); 00306 mem_matching.append(TVec<int>()); 00307 mem_count.append(0); 00308 } 00309 } 00310 00311 if (every_x > 0) { 00312 PLASSERT( find_mol >= 0 ); 00313 if (mem_count[find_mol] % every_x != 0) 00314 realign = false; 00315 } 00316 00317 if (realign) { 00318 00319 int n_points = mol_template->n_points(); 00320 if( n_points < 3 ) 00321 PLERROR( "ChemicalICP::run() - not enough points in template (%d).\n", 00322 n_points ); 00323 00324 Mat best_rotation(3, 3); 00325 Vec best_translation(3); 00326 TVec<int> best_matching( n_points ); 00327 Vec best_weights( n_points ); 00328 real best_error = REAL_MAX; 00329 00330 // transformed template coordinates 00331 Mat tr_template_coords( n_points, 3 ); 00332 00333 // coordinates of molecule points matched to template ones 00334 Mat matched_mol_coords( n_points, 3 ); 00335 00336 int n_initial_angles = initial_angles_list.length(); 00337 EXTREME_NAMED_LOG("ICP+NNet") << n_initial_angles << " initial angles to try" << endl; 00338 for( int i=0 ; i<n_initial_angles ; i++ ) 00339 { 00340 EXTREME_NAMED_LOG("ICP+NNet") << "global iteration number " << i << endl; 00341 // initialization 00342 rotation = rotationMatrixFromAngles( initial_angles_list(i) ); 00343 translation.fill(0); 00344 matching = TVec<int>( n_points, -1 ); 00345 weights = Vec( n_points ); 00346 00347 int n_iter = 0; 00348 real delta_rot_length = REAL_MAX; 00349 real delta_trans_length = REAL_MAX; 00350 00351 applyGeomTransformation( rotation, translation, 00352 template_coordinates->matValue, 00353 tr_template_coords ); 00354 EXTREME_NAMED_LOG("ICPDEBUG") << "beginning ICP:" << endl 00355 << "rotation = " << endl << rotation << endl 00356 << "translation = " << translation << endl 00357 << "template_coordinates->matValue = " << endl 00358 << template_coordinates->matValue << endl 00359 << "tr_template_coords = " << endl << tr_template_coords << endl 00360 << endl; 00361 // main loop 00362 do 00363 { 00364 EXTREME_NAMED_LOG("ICPDEBUG") << " beginning of main loop" << endl; 00365 EXTREME_NAMED_LOG("ICPDEBUG") << " matchNearestNeighbors()" << endl; 00366 matchNearestNeighbors( tr_template_coords, matched_mol_coords ); 00367 EXTREME_NAMED_LOG("ICPDEBUG") << " minimizeWeightedDistance()" << endl; 00368 minimizeWeightedDistance( tr_template_coords, matched_mol_coords, 00369 delta_rot_length, delta_trans_length ); 00370 EXTREME_NAMED_LOG("ICPDEBUG") << " applyGeomTransformation()" << endl; 00371 applyGeomTransformation( rotation, translation, 00372 template_coordinates->matValue, 00373 tr_template_coords ); 00374 EXTREME_NAMED_LOG("ICPDEBUG") << "tr_template_coords = " << endl << tr_template_coords << endl; 00375 error = computeWeightedDistance( tr_template_coords, 00376 matched_mol_coords ); 00377 n_iter++; 00378 EXTREME_NAMED_LOG("ICPDEBUG") << "end of main loop" << endl; 00379 EXTREME_NAMED_LOG("ICPDEBUG") << " iteration = " << n_iter << " / " << max_iter << endl 00380 << " error = " << error << " / " << error_t << endl 00381 << " delta_rot_length = " << delta_rot_length << " / " 00382 << angle_t << endl 00383 << " delta_trans_length = " << delta_trans_length << " / " 00384 << trans_t << endl 00385 << endl; 00386 } 00387 while( n_iter < max_iter && 00388 error > error_t && 00389 delta_rot_length > angle_t && 00390 delta_trans_length > trans_t ); 00391 00392 // keep the best one 00393 if( error < best_error ) 00394 { 00395 best_error = error; 00396 best_rotation << rotation; 00397 best_translation << translation; 00398 best_matching << matching; 00399 } 00400 EXTREME_NAMED_LOG("ICP+NNet") << "end global iteration number " << i << endl; 00401 } 00402 00403 // get best parameters 00404 error = best_error; 00405 rotation = best_rotation; 00406 translation = best_translation; 00407 matching = best_matching; 00408 00409 if (memory != "none") { 00410 // Remember this alignment. 00411 PLASSERT( find_mol >= 0 ); 00412 mem_error[find_mol] = error; 00413 mem_rotation[find_mol].resize(rotation.length(), rotation.width()); 00414 mem_rotation[find_mol] << rotation; 00415 mem_translation[find_mol].resize(translation.length()); 00416 mem_translation[find_mol] << translation; 00417 mem_matching[find_mol].resize(matching.length()); 00418 mem_matching[find_mol] << matching; 00419 } 00420 00421 } else { 00422 PLASSERT( find_mol >= 0 ); 00423 error = mem_error[find_mol]; 00424 // Copies make the code simpler, though it would be slightly more 00425 // efficient to use resize and the << operator. 00426 rotation = mem_rotation[find_mol].copy(); 00427 translation = mem_translation[find_mol].copy(); 00428 matching = mem_matching[find_mol].copy(); 00429 } 00430 00431 if (find_mol >= 0) 00432 // Increment counter of alignments for this molecule. 00433 mem_count[find_mol]++; 00434 00435 // Update the 'matching_neighbors' variable. 00436 matching_neighbors->resize(matching.length(), 1); 00437 for (int i = 0; i < matching.length(); i++) 00438 matching_neighbors->value[i] = matching[i]; 00439 00440 if( !fast_is_equal( initial_angles_step, 0. ) ) 00441 initial_angles_list.resize( 0, 3 ); 00442 } 00443 00445 // matchNearestNeighbors // 00447 void ChemicalICP::matchNearestNeighbors( const Mat& tr_template_coords, 00448 const Mat& matched_mol_coords ) 00449 { 00450 //pout << "matchNearestNeighbors()" << endl; 00451 Mat mol_coords = mol_coordinates->matValue; 00452 int n_template_points = tr_template_coords.length(); 00453 int n_mol_points = mol_coords.length(); 00454 00455 if( matching_method == "exhaustive" ) 00456 { 00457 // bruteforce searche 00458 for( int i=0 ; i<n_template_points ; i++ ) 00459 { 00460 Vec t_point = tr_template_coords( i ); 00461 //pout << "tr_template_coords(" << i << ") = " << t_point << endl; 00462 Vec dists( n_template_points ); 00463 real closest_dist2 = REAL_MAX; 00464 //pout << "closest_dist2 = " << closest_dist2 << endl; 00465 00466 for( int j=0 ; j<n_mol_points ; j++ ) 00467 { 00468 // compute distance 00469 Vec m_point = mol_coords( j ); 00470 real dist2 = powdistance( t_point, m_point, 2 ) 00471 + feat_distances2(i, j); 00472 /*pout << "m_point = " << m_point << endl; 00473 pout << "powdistance( t_point, m_point, 2 ) = " 00474 << powdistance( t_point, m_point, 2 ) << endl; 00475 pout << "feat_distances2(" << i << "," << j << ") = " << 00476 feat_distances2(i, j) << endl; 00477 pout << "dist2 = " << dist2 << endl;*/ 00478 00479 // keep the smallest 00480 if( dist2 < closest_dist2 ) 00481 { 00482 //pout << "beep " << endl; 00483 //pout << " matching[" << i << "] = " << j << endl; 00484 closest_dist2 = dist2; 00485 matching[i] = j; 00486 } 00487 } 00488 matched_mol_coords( i ) << mol_coords( matching[i] ); 00489 } 00490 } 00491 } 00492 00493 void ChemicalICP::minimizeWeightedDistance( const Mat& tr_template_coords, 00494 const Mat& matched_mol_coords, 00495 real& delta_rot_length, 00496 real& delta_trans_length ) 00497 { 00498 Mat delta_rot( 3, 3 ); 00499 Vec delta_trans( 3 ); 00500 real err = REAL_MAX; 00501 00502 //pout << "before computeWeights()" << endl 00503 // << tr_template_coords(0) << endl; 00504 computeWeights( tr_template_coords, matched_mol_coords ); 00505 //pout << "after computeWeights() / before transformationFromWeightedMatchedPoints" << endl 00506 // << tr_template_coords(0) << endl; 00507 transformationFromWeightedMatchedPoints( tr_template_coords, 00508 matched_mol_coords, 00509 weights, 00510 delta_rot, delta_trans, 00511 err ); 00512 /*pout << "after transformationFromWeightedMatchedPoints" << endl 00513 << tr_template_coords(0) << endl; 00514 pout << "weights = " << weights << endl 00515 << "delta_rot = " << endl << delta_rot << endl 00516 << "delta_trans = " << delta_trans << endl;*/ 00517 00518 delta_rot_length = norm( anglesFromRotationMatrix( delta_rot ), 2 ); 00519 delta_trans_length = norm( delta_trans ); 00520 00521 // accumulate transformation ensuring normalization 00522 Vec angles = anglesFromRotationMatrix( product( rotation, delta_rot ) ); 00523 rotation << rotationMatrixFromAngles( angles ); 00524 00525 // translation = delta_trans + delta_rot * translation 00526 productAcc( delta_trans, delta_rot, translation ); 00527 translation << delta_trans; 00528 //pout << "rotation = " << endl << rotation << endl 00529 // << "translation = " << translation << endl; 00530 } 00531 00532 real ChemicalICP::computeWeightedDistance( const Mat& tr_template_coords, 00533 const Mat& matched_mol_coords ) 00534 { 00535 //pout << "computeWeightedDistance()" << endl; 00536 real err = 0; 00537 int n_points = tr_template_coords.length(); 00538 if( matching_method == "exhaustive" ) 00539 { 00540 for( int i=0 ; i<n_points ; i++ ) 00541 { 00542 real diff2 = feat_distances2( i, matching[i] ) + 00543 powdistance( tr_template_coords(i), matched_mol_coords(i), 2 ); 00544 err += weights[i] * sqrt( diff2 ); 00545 /*pout << "i = " << i << endl 00546 << "tr_template_coords = " << tr_template_coords(i) << endl 00547 << "matched_mol_coords = " << matched_mol_coords(i) << endl 00548 << "weights = " << weights[i] << endl 00549 << "delta_err = " << weights[i] * sqrt( diff2 ) << endl 00550 << endl;*/ 00551 } 00552 } 00553 return err; 00554 } 00555 00556 void ChemicalICP::computeWeights( const Mat& tr_template_coords, 00557 const Mat& matched_mol_coords ) 00558 { 00559 //pout << "computeWeights()" << endl; 00560 int n_points = tr_template_coords.length(); 00561 //pout << "n_points = " << n_points << endl; 00562 00563 if( weighting_method == "none" ) 00564 { 00565 //pout << "pouet!" << endl; 00566 weights.fill( 1. / real( n_points ) ); 00567 //pout << "1. / real( n_points ) = " << 1. / real( n_points ) << endl; 00568 //pout << "weights[0] = " << weights[0] << endl; 00569 } 00570 else if( weighting_method == "features_sigmoid" ) 00571 { 00572 //pout << "plouf!" << endl; 00573 real total_weight = 0; 00574 real mid = weighting_params->value[0]; 00575 real slope = weighting_params->value[1]; 00576 00577 if( matching_method == "exhaustive" ) 00578 { 00579 for( int i=0 ; i<n_points ; i++ ) 00580 { 00581 real diff2 = feat_distances2( i, matching[i] ); 00582 weights[i] = sigmoid( slope * (mid - sqrt(diff2)) ); 00583 total_weight += weights[i]; 00584 /*pout << "i = " << i << endl 00585 << "diff2 = " << diff2 << endl 00586 << "slope = " << slope << endl 00587 << "mid = " << mid << endl 00588 << "weights = " << weights[i] << endl 00589 << endl;*/ 00590 } 00591 weights /= total_weight; 00592 //pout << "total_weight = " << total_weight << endl; 00593 } 00594 } 00595 } 00596 00597 void ChemicalICP::saveMatch( const PPath& filename ) 00598 { 00599 PStream file = openFile( filename, PStream::plearn_ascii, "w" ); 00600 file<< "rotation =" << endl 00601 << rotation << endl 00602 << "translation =" << endl 00603 << translation << endl 00604 << "matching =" << endl 00605 << matching << endl 00606 << "error = " << error << endl; 00607 } 00608 00609 // ### Nothing to add here, simply calls build_ 00610 void ChemicalICP::build() 00611 { 00612 inherited::build(); 00613 build_(); 00614 } 00615 00617 // makeDeepCopyFromShallowCopy // 00619 void ChemicalICP::makeDeepCopyFromShallowCopy(CopiesMap& copies) 00620 { 00621 inherited::makeDeepCopyFromShallowCopy(copies); 00622 00623 varDeepCopyField(mol_feat_indices, copies); 00624 varDeepCopyField(template_feat_indices, copies); 00625 varDeepCopyField(matching_neighbors, copies); 00626 deepCopyField(mol_template, copies); 00627 deepCopyField(molecule, copies); 00628 deepCopyField(feature_names, copies); 00629 varDeepCopyField(weighting_params, copies); 00630 deepCopyField(initial_angles_list, copies); 00631 deepCopyField(rotation, copies); 00632 deepCopyField(translation, copies); 00633 deepCopyField(matching, copies); 00634 deepCopyField(weights, copies); 00635 deepCopyField(used_properties, copies); 00636 deepCopyField(other_base_properties, copies); 00637 00638 deepCopyField(mem_aligned_molecules, copies); 00639 deepCopyField(mem_error, copies); 00640 deepCopyField(mem_rotation, copies); 00641 deepCopyField(mem_translation, copies); 00642 deepCopyField(mem_matching, copies); 00643 deepCopyField(mem_count, copies); 00644 00645 deepCopyField(used_feat_names, copies); 00646 deepCopyField(feat_distances2, copies); 00647 00648 varDeepCopyField(all_mol_features, copies); 00649 varDeepCopyField(all_template_features, copies); 00650 varDeepCopyField(all_template_feat_dev, copies); 00651 varDeepCopyField(mol_coordinates, copies); 00652 varDeepCopyField(used_mol_features, copies); 00653 varDeepCopyField(template_coordinates, copies); 00654 varDeepCopyField(template_geom_dev, copies); 00655 varDeepCopyField(used_template_features, copies); 00656 varDeepCopyField(used_template_feat_dev, copies); 00657 } 00658 00660 // declareOptions // 00662 void ChemicalICP::declareOptions(OptionList& ol) 00663 { 00664 // declareOption(ol, "myoption", &ChemicalICP::myoption, OptionBase::buildoption, 00665 // "Help text describing this option"); 00666 // ... 00667 00668 declareOption(ol, "mol_template", &ChemicalICP::mol_template, 00669 OptionBase::buildoption, 00670 "The template we try to align on the molecule"); 00671 00672 declareOption(ol, "molecule", &ChemicalICP::molecule, 00673 OptionBase::buildoption, 00674 "The molecule"); 00675 00676 declareOption(ol, "feature_names", &ChemicalICP::feature_names, 00677 OptionBase::buildoption, 00678 "Names of features to use during alignment.\n" 00679 "Empty TVec means 'use all available features.\n" 00680 "Use '[ \"none\" ]' if you don't want to use any feature.\n" 00681 ); 00682 00683 declareOption(ol, "weighting_method", &ChemicalICP::weighting_method, 00684 OptionBase::buildoption, 00685 "Method used to compute the weight of a pair of point." 00686 " One of:\n" 00687 " - \"features_sigmoid\": sigmoid of feature distance,\n" 00688 " - \"none\": same weight for each pair.\n" 00689 ); 00690 00691 declareOption(ol, "weighting_params", &ChemicalICP::weighting_params, 00692 OptionBase::buildoption, 00693 "Var containing parameters used during weighting.\n" 00694 "Size and meaning depends on value of 'weighting_method'.\n" 00695 ); 00696 00697 declareOption(ol, "matching_method", &ChemicalICP::matching_method, 00698 OptionBase::buildoption, 00699 "Method used to find the nearest neighbors. For the moment," 00700 " only one:\n" 00701 " - \"exhaustive\": exhaustive search (caching feature" 00702 " distances).\n" 00703 ); 00704 00705 declareOption(ol, "memory", &ChemicalICP::memory, OptionBase::buildoption, 00706 "The kind of memory allowing ICP to re-use previous alignments:\n" 00707 "- 'none' : no memory\n" 00708 "- 'every-X': the alignment for a given molecule is re-computed only\n" 00709 " once every X times it is asked\n" 00710 ); 00711 00712 declareOption(ol, "initial_angles_step", &ChemicalICP::initial_angles_step, 00713 OptionBase::buildoption, 00714 "Tries initial rotations every \"initial_angles_step\"" 00715 " degrees"); 00716 00717 declareOption(ol, "initial_angles_list", &ChemicalICP::initial_angles_list, 00718 OptionBase::buildoption, 00719 "Explicit list of initial rotations angles"); 00720 00721 declareOption(ol, "max_iter", &ChemicalICP::max_iter, 00722 OptionBase::buildoption, 00723 "Maximum number of iterations to perform during alignment"); 00724 00725 declareOption(ol, "error_t", &ChemicalICP::error_t, 00726 OptionBase::buildoption, 00727 "Stop alignment if error falls below this threshold"); 00728 00729 declareOption(ol, "angle_t", &ChemicalICP::angle_t, 00730 OptionBase::buildoption, 00731 "Stop alignment if angles falls below this threshold"); 00732 00733 declareOption(ol, "trans_t", &ChemicalICP::trans_t, 00734 OptionBase::buildoption, 00735 "Stop alignment if translation falls below this threshold"); 00736 00737 declareOption(ol, "rotation", &ChemicalICP::rotation, 00738 OptionBase::learntoption, 00739 "Learned rotation matrix"); 00740 00741 declareOption(ol, "translation", &ChemicalICP::translation, 00742 OptionBase::learntoption, 00743 "Learned translation vector"); 00744 00745 declareOption(ol, "matching", &ChemicalICP::matching, 00746 OptionBase::learntoption, 00747 "matching[i] is the index of the molecule point being\n" 00748 "the nearest neighbor of template point i.\n"); 00749 00750 declareOption(ol, "weights", &ChemicalICP::weights, 00751 OptionBase::learntoption, 00752 "Weight of the pair of points (i, matching[i])"); 00753 00754 declareOption(ol, "error", &ChemicalICP::error, 00755 OptionBase::learntoption, 00756 "Weigted error of the alignment"); 00757 00758 // Now call the parent class' declareOptions 00759 inherited::declareOptions(ol); 00760 } 00761 00763 // build_ // 00765 void ChemicalICP::build_() 00766 { 00767 00768 #ifdef BOUNDCHECK 00769 // Variable names for debugging. 00770 used_mol_features->setName("used_mol_features"); 00771 used_template_features->setName("used_template_features"); 00772 used_template_feat_dev->setName("used_template_feat_dev"); 00773 #endif 00774 00775 if (memory == "none") { 00776 // A negative value for 'every_x' means it is not used. 00777 every_x = -1; 00778 } else if (memory.find("every-") == 0) { 00779 string x = memory.substr(6); 00780 every_x = atoi(x.c_str()); 00781 if (every_x <= 0) 00782 PLERROR("In ChemicalICP::build_ - The X value in 'every-X' (option" 00783 " 'memory') must be positive"); 00784 } else 00785 PLERROR("In ChemicalICP::build_ - Invalid value for option 'memory'"); 00786 00787 if( feature_names.size() > 0 && 00788 lowerstring( feature_names[0] ) == "none" ) 00789 { 00790 // no feature will be used during the alignment nor score computation 00791 feature_names[0] = "none"; 00792 if( feature_names.size() > 1 ) 00793 { 00794 PLWARNING("First element of 'feature_names' is 'none', but" 00795 " other features are present.\n" 00796 "Resizing 'feature_names' to 1.\n"); 00797 feature_names.resize( 1 ); 00798 } 00799 } 00800 00801 if( lowerstring( weighting_method ) == "none" || weighting_method == "" ) 00802 { 00803 weighting_method = "none"; 00804 weighting_params->resize(0, 0); 00805 } 00806 else if( lowerstring( weighting_method ) == "features_sigmoid" ) 00807 { 00808 weighting_method = "features_sigmoid"; 00809 weighting_params->resize(2, 1); 00810 } 00811 else 00812 PLERROR( "ChemicalICP::build_ - weighting_method '%s' is unknown.\n", 00813 weighting_method.c_str() ); 00814 00815 00816 if( lowerstring( matching_method ) == "exhaustive" ) 00817 matching_method = "exhaustive"; 00818 else 00819 PLERROR( "ChemicalICP::build_ - matching_method '%s' is unknown.\n", 00820 matching_method.c_str() ); 00821 00822 if( mol_template ) 00823 { 00824 // make the Var's relative to the template have the right storage 00825 template_coordinates->defineValueLocation( mol_template->coordinates ); 00826 template_geom_dev->defineValueLocation( 00827 mol_template->geom_dev.toMat( mol_template->n_points(), 1 ) ); 00828 all_template_features->defineValueLocation( mol_template->features ); 00829 all_template_feat_dev->defineValueLocation( mol_template->feat_dev ); 00830 00831 if( molecule ) 00832 { 00833 // make as if 'setMolecule' were called 00834 computeUsedFeatures(); 00835 computeVariables(); 00836 00837 if( matching_method == "exhaustive" ) 00838 cacheFeatureDistances(); 00839 } 00840 00841 // Resize 'matching_neighbors' so that it has the correct size (i.e. 00842 // the number of points in the template). 00843 matching_neighbors->resize(mol_template->n_points(), 1); 00844 } 00845 00846 // build VarArray 00847 used_properties = mol_coordinates 00848 & template_coordinates & template_geom_dev 00849 & used_mol_features & used_template_features & used_template_feat_dev; 00850 00851 other_base_properties = all_mol_features 00852 & all_template_features & all_template_feat_dev; 00853 //pout << "end build_()" << endl; 00854 00855 } 00856 00857 00858 } // end of namespace PLearn 00859 00860 00861 /* 00862 Local Variables: 00863 mode:c++ 00864 c-basic-offset:4 00865 c-file-style:"stroustrup" 00866 c-file-offsets:((innamespace . 0)(inline-open . 0)) 00867 indent-tabs-mode:nil 00868 fill-column:79 00869 End: 00870 */ 00871 // vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:encoding=utf-8:textwidth=79 :