From 4fe5ff2a766531372e1925ed4d6570038767bb07 Mon Sep 17 00:00:00 2001 From: Paolo Cignoni cignoni Date: Fri, 24 Jul 2009 15:28:36 +0000 Subject: [PATCH] refactoring: ransac loop has been divided in more loops, each one performing a specific task, such as base selection, matching, ranking, short consensus and full consensus. timer in millisec have been added to each loop. many log infos have been added in DEBUG mode. Random generator for base selection has been initialized with iteration, instead of clock() that is dangerous 'cause iterations are really close one another. --- .../feature_alignment.h | 172 ++++++++---------- .../filter_feature_alignment.cpp | 46 +++-- 2 files changed, 106 insertions(+), 112 deletions(-) diff --git a/src/meshlabplugins/filter_feature_alignment/feature_alignment.h b/src/meshlabplugins/filter_feature_alignment/feature_alignment.h index 7c8bd13cd..59414eaa2 100644 --- a/src/meshlabplugins/filter_feature_alignment/feature_alignment.h +++ b/src/meshlabplugins/filter_feature_alignment/feature_alignment.h @@ -375,13 +375,14 @@ template class FeatureAlignment int errorCode; Matrix44Type tr; float bestConsensus; - int numSkippedIter; + int numBasesFound; int numWonShortCons; int numWonFullCons; int numMatches; int totalTime; int baseSelectionTime; int matchingTime; + int rankTime; int shortConsTime; int fullConsTime; int initTime; @@ -392,7 +393,7 @@ template class FeatureAlignment errorCode = -1; tr = Matrix44Type::Identity(); bestConsensus = 0.0f; - numSkippedIter = 0; + numBasesFound = 0; numWonShortCons = 0; numWonFullCons = 0; numMatches = 0; @@ -402,6 +403,7 @@ template class FeatureAlignment shortConsTime = 0; fullConsTime = 0; baseSelectionTime = 0; + rankTime = 0; errorMsg = "An unkown error occurred."; } }; @@ -414,6 +416,7 @@ template class FeatureAlignment FeatureType** matchPtr; FeatureType** basePtr; + CandidateType(FeatureType** base, FeatureType** match, int cons = -1, float spd = numeric_limits::max()):basePtr(base),matchPtr(match),shortCons(cons),summedPointDist(spd){} static bool SortByDistance(CandidateType i,CandidateType j) { return (i.summedPointDistj.shortCons); } }; @@ -472,6 +475,8 @@ template class FeatureAlignment QTime tot_timer, timer; //timers tot_timer.start(); + assert(vecFMov); assert(vecFFix); + if(param.nBase>int(vecFMov->size())){ setError(2,res); return res; } //not enough features to pick a base if(param.k>int(vecFFix->size())){ setError(3, res); return res; } //not enough features to pick k neighboors @@ -496,102 +501,99 @@ template class FeatureAlignment //auxiliary vectors needed inside the loop vector* baseVec = new vector(); - vector* matchesVec = new vector(); vector* candidates = new vector(); //variables needed for progress bar callback - float progBar = 0.0f; - float offset = 30.0f/param.ransacIter; + float progBar = 1.0f; + if(cb) cb(int(progBar),"Selecting bases..."); //loop of ransac iterations to find all possible matches at once + timer.start(); for(int i = 0; isize()-1], *matchesVec, param); - - assert(baseVec->size()>=1); - int currBase = baseVec->size()-1; - - //scan all the new matches found and computes datas for each one - for(unsigned int j=candidates->size(); jsize(); j++) - { - CandidateType data; - data.shortCons = -1; //this is needed to get structures well sorted later - data.matchPtr = (*matchesVec)[j]; //set the pointer to the match - data.basePtr = (*baseVec)[currBase]; //set the pointer to the relative base - - //computes the rigid transformation matrix that overlaps the two points sets - Matrix44Type tr; - errCode = FeatureAlignment::FindRigidTransformation(mFix, mMov, (*baseVec)[currBase], (*matchesVec)[j], param.nBase, tr); - if(errCode) { res.numSkippedIter++; continue; } //can't find a rigid transformation, skip this iteration - - data.tr = tr; //store transformation - - Matrix44Type oldTr = ApplyTransformation(mMov, tr); //apply transformation - data.summedPointDist = SummedPointDistances(mFix, mMov, (*baseVec)[currBase], (*matchesVec)[j], param.nBase); //compute and store the sum of points distances - ResetTransformation(mMov, oldTr); //restore old tranformation - - candidates->push_back(data); - } - assert(candidates->size()==matchesVec->size()); + int errCode = FeatureAlignment::SelectBase(*vecFMov,*baseVec, param,i); + if(errCode) continue; //can't find a base, skip this iteration } - if(param.log) param.log(3,"%i bases found.", param.ransacIter-res.numSkippedIter); + res.numBasesFound=baseVec->size(); + res.baseSelectionTime = timer.elapsed(); + if(param.log) param.log(3,"%i bases found in %i msec.", baseVec->size(), res.baseSelectionTime); - - //res.matchingTime = int(difftime(end_loop,start_loop)); + timer.start(); + for(int i = 0; isize(); i++) + { + int errCode = FeatureAlignment::Matching(*vecFFix, *vecFMov, fkdTree, (*baseVec)[i], *candidates, param); + if(cb){ progBar+=(20.0f/baseVec->size()); cb(int(progBar),"Matching..."); } + } + res.matchingTime = timer.elapsed(); res.numMatches = candidates->size(); + if(param.log) param.log(3,"%i candidates found in %i msec.", candidates->size(), res.matchingTime); + + timer.start(); + for(unsigned int j=0; jsize(); j++) + { + CandidateType& currCandidate = (*candidates)[j]; + + //computes the rigid transformation matrix that overlaps the two points sets + Matrix44Type tr; + int errCode = FeatureAlignment::FindRigidTransformation(mFix, mMov, currCandidate.basePtr, currCandidate.matchPtr, param.nBase, tr); + if(errCode) continue; //can't find a rigid transformation, skip this iteration + + currCandidate.tr = tr; //store transformation + + Matrix44Type oldTr = ApplyTransformation(mMov, tr); //apply transformation + currCandidate.summedPointDist = SummedPointDistances(mFix, mMov, currCandidate.basePtr, currCandidate.matchPtr, param.nBase); //compute and store the sum of points distance + ResetTransformation(mMov, oldTr); //restore old tranformation + + if(cb){ progBar+=(20.0f/candidates->size()); cb(int(progBar),"Ranking candidates..."); } + } + res.rankTime = timer.elapsed(); + if(param.log) param.log(3,"Ranking performed in %i msec.", res.rankTime); //sort candidates by summed point distances sort(candidates->begin(), candidates->end(), CandidateType::SortByDistance); - +for(int i=0; i<50;i++)if(param.log) param.log(3,"%i spd %.2f.", i,(*candidates)[i].summedPointDist); //variable needed for progress bar callback - offset = (40.0f/math::Min(param.maxNumShortConsensus,int(candidates->size()))); + float offset = (20.0f/math::Min(param.maxNumShortConsensus,int(candidates->size()))); - //time(&start_loop); //evaluetes at most first g candidates + timer.start(); for(unsigned int j=0; jsize() && j= short_cons_succ) res.numWonShortCons++; //count how many won, and use this as bound later - } + if(currCandidate.shortCons >= short_cons_succ) res.numWonShortCons++; //count how many won, and use this as bound later - //res.shortConsTime = int(difftime(end_loop,start_loop)); + if(cb){ progBar+=offset; cb(int(progBar),"Short consensus..."); } + } + res.shortConsTime = timer.elapsed(); + if(param.log) param.log(3,"Short consensus performed in %i msec.", res.shortConsTime); //sort candidates by short consensus sort(candidates->begin(), candidates->end(), CandidateType::SortByScore); //variables needed for progress bar callback - offset = (40.0f/math::Min(param.maxNumFullConsensus,int(candidates->size()))); + offset = (20.0f/math::Min(param.maxNumFullConsensus,res.numWonShortCons)); - //time(&start_loop); //evaluetes at most maxNumFullConsensus candidates beetween those who won short consensus + timer.start(); for(int j=0; j= cons_succ) res.numWonFullCons++; @@ -604,11 +606,13 @@ template class FeatureAlignment res.bestConsensus = 100.0f*(float(bestConsensus)/param.fullConsensusSamples); bestConsIdx = j; if(consensus >= ransac_succ) break; //very good alignment, no more iterations are done - } - } + } - //res.fullConsTime = int(difftime(end_loop,start_loop)); + if(cb){ progBar+=offset; cb(int(progBar),"Full consensus..."); } + } + res.fullConsTime = timer.elapsed(); res.totalTime = tot_timer.elapsed(); + if(param.log) param.log(3,"Full consensus performed in %i msec.", res.fullConsTime); //if flag 'points' is checked, clear old picked points and save the new points if(param.pickPoints){ @@ -619,9 +623,9 @@ template class FeatureAlignment } //Clean structures... - FeatureAlignment::CleanTuplesVector(baseVec, true); - FeatureAlignment::CleanTuplesVector(matchesVec, true); - candidates->clear(); if(candidates) delete candidates; candidates = NULL; + for(int i=0; isize();i++){ delete[] (*baseVec)[i];} + for(int i=0; isize();i++){ delete[] (*candidates)[i].basePtr; delete[] (*candidates)[i].matchPtr;} + delete baseVec; delete candidates; return res; } @@ -653,7 +657,7 @@ template class FeatureAlignment res.errorMsg = FeatureAlignment::getErrorMsg(errorCode); } - static FEATURE_TYPE** FeatureUniform(vector& vecF, int* sampleNum) + static FEATURE_TYPE** FeatureUniform(vector& vecF, int* sampleNum, int seed = clock()) { typedef MESH_TYPE MeshType; typedef FEATURE_TYPE FeatureType; @@ -666,7 +670,7 @@ template class FeatureAlignment assert(vec->size()==vecF.size()); - tri::SurfaceSampling >::SamplingRandomGenerator().initialize(clock()); + tri::SurfaceSampling >::SamplingRandomGenerator().initialize(seed); unsigned int (*p_myrandom)(unsigned int) = tri::SurfaceSampling >::RandomInt; std::random_shuffle(vec->begin(),vec->end(), p_myrandom); @@ -765,32 +769,18 @@ template class FeatureAlignment if(close) annClose(); //done with ANN; clean memory shared among all the kdTrees } - static void CleanTuplesVector(vector* tuplesVec, bool deleteVec) - { - //clean vectors contents... - for(unsigned int j=0; jsize(); j++){ - delete [] (*tuplesVec)[j]; (*tuplesVec)[j] = NULL; - } - tuplesVec->clear(); - - //if deleteVec is set to true, deallocates the vector too... - if(deleteVec){ - delete tuplesVec; tuplesVec = NULL; - } - } - - static int SelectBase(vector& vecFMov, vector& baseVec, Parameters& param) + static int SelectBase(vector& vecFMov, vector& baseVec, Parameters& param,int seed = clock()) { assert(param.nBase<=int(vecFMov.size())); //not enough features to pick a base float baseDist = param.sparseBaseDist*(param.mMovBBoxDiag/100.0f); //compute needed params - FeatureType** base = FeatureAlignment::FeatureUniform(vecFMov, &(param.nBase)); //randomly chooses a base of features from vecFFix + FeatureType** base = FeatureAlignment::FeatureUniform(vecFMov, &(param.nBase), seed); //randomly chooses a base of features from vecFFix if(!VerifyBaseDistances(base, param.nBase, baseDist)) return 4; //if base point are not enough sparse, skip baseVec.push_back(base); return 0; } - static int Matching(vector& vecFFix, vector& vecFMov, ANNkd_tree* kdTree, FEATURE_TYPE** base, vector& matchesVec, Parameters& param, CallBackPos *cb=NULL) + static int Matching(vector& vecFFix, vector& vecFMov, ANNkd_tree* kdTree, FEATURE_TYPE** base, vector& candidates, Parameters& param, CallBackPos *cb=NULL) { float pBar = 0, offset = 100.0f/param.nBase; //used for progresss bar callback if(cb) cb(0, "Matching..."); @@ -832,7 +822,7 @@ template class FeatureAlignment //branch and bound int* curSolution = new int[param.nBase]; for(int i=0; i class FeatureAlignment static int FindRigidTransformation(MESH_TYPE& mFix, MESH_TYPE& mMov, FEATURE_TYPE* fixF[], FEATURE_TYPE* movF[], int nBase, Matrix44& tr, CallBackPos *cb=NULL) { - typedef MESH_TYPE MeshType; - typedef FEATURE_TYPE FeatureType; - typedef typename MESH_TYPE::ScalarType ScalarType; - if(cb) cb(0,"Computing rigid transformation..."); //computes the rigid transformation matrix that overlaps the two points sets @@ -954,7 +940,7 @@ template class FeatureAlignment } } - static void Match(FEATURE_TYPE** base, vector* >& matchedVec, int nBase, int level, int curSolution[], vector& solutionsVec, float errDist, CallBackPos *cb = NULL) + static void Match(FEATURE_TYPE** base, vector* >& matchedVec, int nBase, int level, int curSolution[], vector& candidates, float errDist, CallBackPos *cb = NULL) { assert(level class FeatureAlignment for(int h=0; h* baseVec = new vector(); - vector* matchesVec = new vector(); + vector* candidates = new vector(); AlignerType aligner; ResultType res = aligner.init(mFix.cm, mMov.cm, param); @@ -558,14 +559,15 @@ typename ALIGNER_TYPE::Result FilterFeatureAlignment::MatchingOperation(MeshMode //execute matching procedure with requested parameters; int errCode = AlignerType::SelectBase(*(aligner.vecFMov),*baseVec,param); if(errCode){ AlignerType::setError(errCode, res); return res; } - errCode = AlignerType::Matching(*(aligner.vecFFix), *(aligner.vecFMov), aligner.fkdTree, (*baseVec)[baseVec->size()-1], *matchesVec, param, cb); + errCode = AlignerType::Matching(*(aligner.vecFFix), *(aligner.vecFMov), aligner.fkdTree, (*baseVec)[baseVec->size()-1], *candidates, param, cb); if(errCode){ AlignerType::setError(errCode, res); return res; } - res.numMatches = matchesVec->size(); //store the numeber of matches found + res.numMatches = candidates->size(); //store the numeber of matches found - //cleaning baseVec and matchesVec... - AlignerType::CleanTuplesVector(baseVec, true); - AlignerType::CleanTuplesVector(matchesVec, true); + //Clean structures... + for(int i=0; isize();i++){ delete[] (*baseVec)[i];} + for(int i=0; isize();i++){ delete[] (*candidates)[i].basePtr; delete[] (*candidates)[i].matchPtr;} + delete baseVec; delete candidates; res.totalTime = timer.elapsed(); @@ -581,13 +583,14 @@ typename ALIGNER_TYPE::Result FilterFeatureAlignment::RigidTransformationOperati typedef typename MeshType::ScalarType ScalarType; typedef typename AlignerType::Result ResultType; typedef Matrix44 Matrix44Type; + typedef typename AlignerType::CandidateType CandidateType; QTime timer; timer.start(); //start timer //create vectors to hold tuples of bases and matches vector* baseVec = new vector(); - vector* matchesVec = new vector(); + vector* candidates = new vector(); AlignerType aligner; ResultType res = aligner.init(mFix.cm, mMov.cm, param); @@ -596,25 +599,27 @@ typename ALIGNER_TYPE::Result FilterFeatureAlignment::RigidTransformationOperati //execute matching procedure with requested parameters; int errCode = AlignerType::SelectBase(*(aligner.vecFMov),*baseVec,param); if(errCode){ AlignerType::setError(errCode, res); return res; } - errCode = AlignerType::Matching(*(aligner.vecFFix), *(aligner.vecFMov), aligner.fkdTree, (*baseVec)[baseVec->size()-1], *matchesVec, param, cb); + errCode = AlignerType::Matching(*(aligner.vecFFix), *(aligner.vecFMov), aligner.fkdTree, (*baseVec)[baseVec->size()-1], *candidates, param, cb); if(errCode){ AlignerType::setError(errCode, res); return res; } - res.numMatches = matchesVec->size(); //store the numeber of matches found + res.numMatches = candidates->size(); //store the numeber of matches found assert(baseVec->size()==1); //now baseVec must hold exactly one base of features - for(unsigned int j=0; jsize(); j++) + for(unsigned int j=0; jsize(); j++) { + CandidateType currCandidate = (*candidates)[j]; Matrix44Type tr; - errCode = AlignerType::FindRigidTransformation(mFix.cm, mMov.cm, (*baseVec)[0], (*matchesVec)[j], param.nBase, tr, cb); + errCode = AlignerType::FindRigidTransformation(mFix.cm, mMov.cm, currCandidate.basePtr, currCandidate.matchPtr, param.nBase, tr, cb); if(errCode) { AlignerType::setError(errCode,res); return res; } res.numWonFullCons++; //use this variable to increase num of found transform res.exitCode = ResultType::ALIGNED; //this means at least one matrix found } - //cleaning baseVec and matchesVec... - AlignerType::CleanTuplesVector(baseVec, true); - AlignerType::CleanTuplesVector(matchesVec, true); + //Clean structures... + for(int i=0; isize();i++){ delete[] (*baseVec)[i];} + for(int i=0; isize();i++){ delete[] (*candidates)[i].basePtr; delete[] (*candidates)[i].matchPtr;} + delete baseVec; delete candidates; res.totalTime = timer.elapsed(); @@ -695,12 +700,14 @@ typename ALIGNER_TYPE::Result FilterFeatureAlignment::RansacDiagramOperation(Mes fprintf(file,"Fix Mesh#%s\nMove Mesh#%s\nFeature#%s\nNum. of vertices of Fix Mesh#%i\nNum. of vertices of Move Mesh#%i\nOverlap#%.2f%%\nShort consensus threshold#%.2f%% overlap#%.2f%% of Move Mesh#\nFull consensus threshold#%.2f%% overlap#%.2f%% of Move Mesh#\nTrials#%i\n",mFix.fileName.c_str(),mMov.fileName.c_str(),FeatureType::getName(),mFix.cm.VertexNumber(),mMov.cm.VertexNumber(),param.overlap,param.shortConsOffset,(param.shortConsOffset*param.overlap/100.0f),param.consOffset,(param.consOffset*param.overlap/100.0f),trials); fflush(file); for(int h=0; h<4; h++) { - fprintf(file,"Iterazioni#Tempo di inizializzazione#Tempo di esecuzione#Prob. Succ.#Prob. Fall. per Sec\n0#0#0#0\n"); fflush(file); + fprintf(file,"Iterazioni#Tempo di inizializzazione#Tempo di esecuzione#Prob. Succ.#Prob. Fall. per Sec#Num. medio basi\n0#0#0#0#0\n"); fflush(file); AlignerType aligner; float probSucc = 0.0f, meanTime = 0.0f, meanInitTime = 0.0f, failPerSec = -1.0f; - int numWon = 0, trialsTotTime = 0, trialsInitTotTime = 0; + int numWon = 0, trialsTotTime = 0, trialsInitTotTime = 0, numBases = 0; param.ransacIter = from; + //param.log = &mylogger; //this is the way to assign a pointer to log function + //move res here while(param.ransacIter<=to) { @@ -714,7 +721,8 @@ for(int h=0; h<4; h++) res = aligner.align(mFix.cm, mMov.cm, param); trialsTotTime+=res.totalTime; - if(res.exitCode==ResultType::ALIGNED)numWon++; + numBases+=res.numBasesFound; + if(res.exitCode==ResultType::ALIGNED) numWon++; if(res.exitCode==ResultType::FAILED) return res; //failure: stop everything and return error } @@ -722,8 +730,8 @@ for(int h=0; h<4; h++) meanTime = trialsTotTime/float(trials); //t=sec elapsed to perform N ransac iterations meanInitTime = trialsInitTotTime/float(trials); failPerSec = std::pow(1-probSucc,1.0f/(meanTime/1000)); //fail rate per sec is: (1-k)^(1/t) - fprintf(file,"%i#%.2f#%.2f#%.2f#%.2f\n", param.ransacIter, meanInitTime/1000, meanTime/1000, probSucc, failPerSec); fflush(file); - numWon = 0; trialsTotTime = 0; trialsInitTotTime=0; progBar=0.0f; + fprintf(file,"%i#%.2f#%.2f#%.2f#%.2f#%i\n", param.ransacIter, meanInitTime/1000, meanTime/1000, probSucc, failPerSec,numBases/trials); fflush(file); + numWon = 0; trialsTotTime = 0; trialsInitTotTime=0; progBar=0.0f; numBases=0; param.ransacIter+=step; } }