mrgloom 242 Жалоба Опубликовано August 27, 2012 http://en.wikipedia.org/wiki/Iterative_Closest_Point ICP(Iterative Closest Point) http://www.morethantechnical.com/2010/06/06/iterative-closest-point-icp-with-opencv-w-code/ для 2д, для 3д есть в PCL. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 27, 2012 библиотека посерьёзней и более приятная. (нужен буст, но собралась без проблем, только в дебаг версии 1 функция почему то не была определена check_query_in_bound) http://www.cvlibs.net/software/libicp.html прямая ссылка http://www.cvlibs.net/downloads/libicp.zip пока правда не понял какие ограничения в демо, только поворот и сдвиг. 1 Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 27, 2012 попробовал, почему то не всегда работает нормально при малом кол-ве точек. и опять же только rigid регистрация(т.е. только поворот и сдвиг) и не всегда нормально сходится, наверно надо как то давать приближение по другому. пример работы. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 27, 2012 http://code.google.com/p/dbicp/ dbicp показано сравнение с ICP Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 28, 2012 я вот подумал, а что если ransac приспособить под это дело. т.е. берем рендомно из 1 сета точек 4 точки, так же из второго, считаем матрицу перехода, затем трансформируем точки первого сета, считаем метрику(какую не уверен ну например сумма мин расстояний для каждой точки) повторяем всё это k раз. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 28, 2012 попробовал прикрутить ранзак для этих дел. но почеу то работает неправильно, видимо выбрана неправильная метрика. я для каждой точки нахожу ближайшую точку, потом суммирую и делю на кол-во точек. хотя возможно надо просто какое то ограничение на матрицу трансформации давать. еще есть идея, что если мы не имеем пар точек от дескрипторов, то можно попытаться найти пары точек, определяя дескриптор как shape context. начальное состояние расстрел из 1000 попыток Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 29, 2012 похоже я понял почему ничего не получилось, т.к. когда я беру рэндомно по 4 точки и считаю матрицу перехода, то если мы имеем лишь похожие фигуры(т.е. ни одна пара не образует нормального преобразования), то перебрав даже все пары мы не получим того преобразования которое требуется. т.е. получается, что надо не точки выбирать рэндомно, а саму матрицу, а это никаких ограничений, хотя возможно можно вывести примерные ограничения, например взять ограничивающие прямоугольники и из этого вывести ограничение на матрицу. http://www.cse.psu.edu/~rcollins/CSE486/lecture15_6pp.pdf Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 30, 2012 всё таки непонятно по обычному ICP (где только поворот и сдвиг) там написано для каждой точки первого сета найти ближайшую точку из второго сета и на основе этого сделать трансформацию. http://www.cs.technion.ac.il/~cs236329/tutorials/ICP.pdf так вот как делать эту трансформацию? вроде как можно как СЛАУ через МНК как я пробовал тут http://www.compvision.ru/forum/index.php?showtopic=999 но может надо подругому? Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 30, 2012 vector<Point> vec_pair; for(int i=0;i<vec_M.size();++i) { double min_dist=INT_MAX; int id=-1; for(int j=0;j<vec_T.size();++j) { double metric= sqrt(double(vec_T[j].x-vec_M[i].x)*(vec_T[j].x-vec_M[i].x)+ (vec_T[j].y-vec_M[i].y)*(vec_T[j].y-vec_M[i].y)); if(min_dist>metric) { min_dist=metric; id=j; } } line(img,vec_M[i],vec_T[id],cvScalar(0,0,255)); vec_pair.push_back(vec_T[id]); } Mat m= estimateRigidTransform(vec_M,vec_pair,1); cout<<m; попробовал так, но не дает нормальный результат, возможно надо отсекать какие либо пары? Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
Smorodov 579 Жалоба Опубликовано August 30, 2012 А здесь не то что Вам нужно?: http://code.google.com/p/gmmreg/ Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 31, 2012 да я это уже видел,и даже оно завелось, но там алгоритмы сложные и не очень понятно как управлять параметрами. но похоже придётся всё равно к этому вернутся. там экзешник, где все алгоритмы требует, какую то довольно большую стороннюю либу.(какую не помню давно собирал) впрочем там экзешник есть, можно не пересобирать. а чтобы завести под питон мне пришлось ставить еще 4 сторонние библиотеки для питона. причем scipy почему то через pip ставится отказался. на матлабе тоже пробовал, но тоже давно. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 31, 2012 вообщем попробовал на тестовом примере там где рыбу из двух частей хотят составить. EM_TPS EM_GRBF TPS_L2 GRBF_L2 TPS_KC GRBF_KC rigid в итоге если не указывать ctrl_pts. (хотя я не понял, что это такое то ли начальное приближение, то ли их можно как то сформировать из имеющихся точек) # 'ctrl_pts' serve as the control points when thin plate splines (TPS) or # Gaussian radial basis functions (GRBF) are used in the nonrigid registration # if 'ctrl_pts' is not provided, model will be used as ctrl_pts # the program 'gmmreg_aux' can be used to generate ctrl pts from regular grid pts сработали нормально только TPS_L2,TPS_KC. но на моем сете 6к и 8к точек, алгоритм выделил слишком много памяти и благополучно вылетел. вроде память выделяется когда как раз не указаны ctrl_pts. вот так работает rigid - видимо обычный ICP Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 31, 2012 забавно, но эти оказались просто сеткой с каким то шагом, причем если их использовать память не кушается, странно почему алгоритм не может сам генерировать эту сетку, а надо ее делать отдельно. но все равно ничего хорошего не получилось,т.к. надо видимо примерное приближение всё таки. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано September 4, 2012 http://en.wikipedia.org/wiki/Weierstrass_transform начал смотреть rigid версию в gmmreg, похоже код там поддерживает еще и скейл. втретил там некий GaussTransform и похоже это не размытие, как то связано с методом сопряженного градиента(conjugate gradient method) или каким то еще итеративным методом, но причём тут гауссиан? double energy = GaussTransform(gmmreg->transformed_model, gmmreg->scene, scale, gradient); вообщем похоже все хвосты ведут на какой то минимайзер, у которого виртуальные функции определенные в самом классе rigid gmmreg. только всё таки непонятно, что это зам метод. //: An object that represents a function from R^n -> R. // It is commonly used to express the // interface of a minimizer. class vnl_cost_function : public vnl_unary_function<double, vnl_vector<double> > { public: //! Default constructor vnl_cost_function():dim(0) {} //! Construct with a specified number of unknowns vnl_cost_function(int number_of_unknowns):dim(number_of_unknowns) {} virtual ~vnl_cost_function() {} //: The main function. Given the parameter vector x, compute the value of f(x). virtual double f(vnl_vector<double> const& x); //: Calculate the gradient of f at parameter vector x. virtual void gradf(vnl_vector<double> const& x, vnl_vector<double>& gradient); //: Compute one or both of f and g. // Normally implemented in terms of the above two, but may be faster if specialized. f != 0 => compute f virtual void compute(vnl_vector<double> const& x, double *f, vnl_vector<double>* g); //: Return the number of unknowns int get_number_of_unknowns() const { return dim; } //: Compute finite-difference gradient void fdgradf(vnl_vector<double> const& x, vnl_vector<double>& gradient, double stepsize = 1e-5); //: Called when error is printed for user. virtual double reported_error(double f_value) { return f_value; } //: Conveniences for printing grad, fdgrad vnl_vector<double> gradf(vnl_vector<double> const& x); vnl_vector<double> fdgradf(vnl_vector<double> const& x); protected: //! Set number of unknowns. void set_number_of_unknowns(int number_of_unknowns) { dim=number_of_unknowns; } public: int dim; }; Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
Smorodov 579 Жалоба Опубликовано September 4, 2012 Контрольные точки, насколько я понял, это точки по которым считается сплайн-деформация. Например для лица это будут глаза, уголки губ и т.д. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
Smorodov 579 Жалоба Опубликовано September 4, 2012 Есть еще ресурс с OpenCV-шным кодом: http://www.morethantechnical.com/2010/06/06/iterative-closest-point-icp-with-opencv-w-code/ Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
kilop 7 Жалоба Опубликовано September 5, 2012 а можете пояснить для какого рода задач вы используете алгоритм? Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано September 5, 2012 Контрольные точки, насколько я понял, это точки по которым считается сплайн-деформация. Например для лица это будут глаза, уголки губ и т.д. какие контрольные точки? последний код отсюда http://www.morethantechnical.com/2010/06/06/iterative-closest-point-icp-with-opencv-w-code/ какой то мутный я его пробовал, вроде даже он не заработал на моих примерах. вообщем что то там мне не понравилось, во всяком случае чтобы его исправить надо понимать что там делается. по поводу conjugate gradient method там это помоему задача energy minimization типа есть transformed_model, а есть scene и видимо итеративно смотрят как хорошо они наложились при каком либо преобразовании. написал еще китайцу который делал gmmreg, но видимо это бесполезно.а код довольно жёсткий, минимум комментариев, иногда вставлены куски-аналоги из матлаба. а можете пояснить для какого рода задач вы используете алгоритм? регистрация 2-х наборов точек. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано September 5, 2012 http://www.morethantechnical.com/2012/08/12/registering-point-clouds-rigidly-with-scale-using-pcl-wcode/#more-1129 а тут вот еще rigid через PCL и как то скейл определяют непонятно через PCA double s = sqrt(ev_B[0])/sqrt(ev_A[0]); через первые собственные значения. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано September 5, 2012 вообщем похоже там по коду используется этот алгоритм minimization with constraints http://en.wikipedia.org/wiki/Limited-memory_BFGS //: Limited memory Broyden Fletcher Goldfarb Shannon minimization with constraints. // Lower and upper bounds may be specified for the variables to be optimized. // The algorithm miminizes a nonlinear function f(x) of n variables // subject to simple bound constraints of l <= x <= u. осталось понять какая задаётся функция для минимизации. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано September 5, 2012 вообщем вся соль похоже в этих функциях, но почему они так считаются не понятно. возможно это связан ос темой "robust point set registration using mixture of Gaussians" //непонятно назначение, по сути просто меняет знак для энергии и градиента double gmmreg_rigid_func::eval(double& f, vnl_matrix<double>& g) { /* L2 version and KC version are equivalent in the rigid case */ g = -g; return -f; } //основная функция которая вычисляет функцию от 3-х переменных(х,у смещение+ поворот) //непонятно каким боком тут GaussTransform(вроде бы используется градиент с предыдущей итерации) double gmmreg_rigid_func::f(const vnl_vector<double>& x) { // std::cout << "f begin" << std::endl; gmmreg->perform_transform(x); //x это вроде параметры которые мы пытаемся определить // std::cout << "transform done" << std::endl; double energy = GaussTransform(gmmreg->transformed_model, gmmreg->scene, scale, gradient); // std::cout << "gauss transform done" << std::endl; return eval(energy, gradient); // std::cout << "f end" << std::endl; } //вычисление градиента void gmmreg_rigid_func::gradf(const vnl_vector<double>& x, vnl_vector<double>& g) { //код из матлаба как называется этот метод? //похоже как то связано с "robust point set registration using mixture of Gaussians" т.к. rigid_costfunc указывает на это /* case 'rigid2d' * [f, grad] = rigid_costfunc(transformed_model, scene, scale); * grad = grad'; * g(1) = sum(grad(1,); * g(2) = sum(grad(2,); * grad = grad*model; * theta = param(3); * r = [-sin(theta) -cos(theta); * cos(theta) -sin(theta)]; * g(3) = sum(sum(grad.*r)); * */ int i = 0; vnl_matrix<double> r; vnl_matrix<double> gm; if (d == 2) { //rigid2d g[0] = gradient.get_column(0).sum(); g[1] = gradient.get_column(1).sum(); r.set_size(2,2); double theta = x[2]; r[0][0] = -sin(theta); r[0][1] = -cos(theta); r[1][0] = cos(theta); r[1][1] = -sin(theta); gm = gradient.transpose()*gmmreg->model; //что дает перемножение градиента на модель? g[2] = 0; for (i=0;i<2;++i) { for (int j=0;j<2;++j) { g[2] += r[i][j]*gm[i][j]; } } } else if (d == 3) { //rigid3d /* * case 'rigid3d' * [f,grad] = rigid_costfunc(transformed_model, scene, scale); * [r,gq] = quaternion2rotation(param(1:4)); * grad = grad'; * gm = grad*model; * g(1) = sum(sum(gm.*gq{1})); * g(2) = sum(sum(gm.*gq{2})); * g(3) = sum(sum(gm.*gq{3})); * g(4) = sum(sum(gm.*gq{4})); * g(5) = sum(grad(1,); * g(6) = sum(grad(2,); * g(7) = sum(grad(3,); */ g[4] = gradient.get_column(0).sum(); g[5] = gradient.get_column(1).sum(); g[6] = gradient.get_column(2).sum(); vnl_vector<double> q; q.set_size(4); for (i=0;i<4;++i) q[i] = x[i]; r.set_size(3,3); vnl_matrix<double> g1,g2,g3,g4; g1.set_size(3,3); g2.set_size(3,3); g3.set_size(3,3); g4.set_size(3,3); quaternion2rotation(q,r,g1,g2,g3,g4); gm = gradient.transpose()*gmmreg->model; g[0] = 0; for (i=0;i<3;++i) { for (int j=0;j<3;++j) { g[0] += g1[i][j]*gm[i][j]; } } g[1] = 0; for (i=0;i<3;++i) { for (int j=0;j<3;++j) { g[1] += g2[i][j]*gm[i][j]; } } g[2] = 0; for (i=0;i<3;++i) { for (int j=0;j<3;++j) { g[2] += g3[i][j]*gm[i][j]; } } g[3] = 0; for (i=0;i<3;++i) { for (int j=0;j<3;++j) { g[3] += g4[i][j]*gm[i][j]; } } } } [/code] Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
Smorodov 579 Жалоба Опубликовано September 5, 2012 GaussTransform - похоже нужен для получения более менее гладкого градиента, чтобы минимизатор работал лучше. На кусочно-гладкой функции производная будет "вылетать" в бесконечность на стыках гладких кусков. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано September 5, 2012 там используется метод Robust Point Set Registration Using Gaussian Mixture Models т.е. сначала наборы точек представляются как Gaussian Mixture Models ,а потом уже ка кто совмещаются. но как это относится к коду пока непонятно. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
Smorodov 579 Жалоба Опубликовано September 5, 2012 При совмещении распределений обычно минимизируется дивергенция Кульбака — Лейблера. http://ru.wikipedia.org/wiki/%D0%A0%D0%B0%D1%81%D1%81%D1%82%D0%BE%D1%8F%D0%BD%D0%B8%D0%B5_%D0%9A%D1%83%D0%BB%D1%8C%D0%B1%D0%B0%D0%BA%D0%B0_%E2%80%94_%D0%9B%D0%B5%D0%B9%D0%B1%D0%BB%D0%B5%D1%80%D0%B0 Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано September 5, 2012 они там предлагают в пейпэре минимизировать L2 distance между микстурами гауссиан. и там это формулируется как inner product of gaussian mixture, только непонятно как выглядят сами микстуры гауссиан, как векторы? и опять же я не очень понимаю, ну допустим мы сформировали положительно определенную cost function и пытаемся каким то там методом найти её минимум, допустим она зависит от 3-х параметров, как мы определяем куда двигать эти параметры чтобы сошлось? (хотя можно минимизатор рассматривать как черный ящик, а задать только целевую функцию и ограничения) тут вот нашел разные методы\метрики для того чтобы посчитать разницу между микстурами гауссиан. (правда тулбокс для аудио) http://kom.aau.dk/project/isound/isp_gmmdistance.html что еще означает affinity matrix которая часто фигурирует. там аж 3 метода с этим связанных во введении описывается Instead of assuming the one-to-one correspondence based on the nearest neighbor criterion, one-to-many relaxations have been proposed to allow for fuzzy correspondences. One of the most notable contributions along these lines is the work in Chui and Rangarajan [2], wherein a soft assignment technique and the deterministic annealing are combined to alternatively estimate the transformation and update the correspondence. The key idea of their work is to assume that each model point corresponds to a weighted sum of the scene points, instead of the closest scene point alone. These weights are taken from an affinity matrix whose entries are proportional to a Gaussian function of the pairwise distances between the moving model and the fixed scene. In other words, this affinity matrix assigns each model point to all points in the scene. The closer the scene point is to the model point, the larger is the assigned weight. Clearly, the closest point correspondence in ICP can be viewed as a binary version of this graduated assignment scheme. In Chui and Rangarajan [2], the registration problem is expressed as an joint optimization over the transformation parameter and correspondence matrix. However, in practice, it is implemented as an alternating update strategy: at each iteration, a transformation is estimated from the virtual correspondence between model points and the linear combinations of scene points, then the correspondence matrix gets modified as the model points move towards the scene and the scale parameter in the Gaussian function for computing the weights decreases. Gaussian affinity matrix Recently Myronenko et al. [6] proposed yet another robust non-rigid point set registration algorithm. Like in Chui and Rangarajan [2], Myronenko et al. [6] maintain the same Gaussian affinity matrix and also adopt a similar alternating update strategy interpreted in an expectation maximization (EM) [20] framework. The major difference between these two algorithms is that in Chui and Rangarajan [2] the non-rigid deformation is modeled by thin-plate splines while in Myronenko et al. [6] they use Gaussian radial basis functions. Another trick is to model each of the two point sets by a probability distribution and then a distance measure between the two distributions is minimized over the transformation space to yield the desired transformation. For instance, Tsin and Kanade [4] proposed a kernel correlation (KC) based point set registration approach where the cost function is proportional to the correlation of two kernel density estimates. Glaunes et al. [24] formulated the problem of aligning two unlabeled point sets by first modeling the point sets as weighted sums of Dirac measures and then finding the optimal diffeomorphic transformation between the two discrete distributions. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах