Find matching points from two coordinate systems
显示 更早的评论
Hello,
I have 10 points with x and y coordinate (cell locations imaged with microscope). I have used two different methods to image these cells, so I have two sets of coordinates that could be rotated, sheared, shifted, mirrored to each other. I would like to assign each point from one coordinate system to the matching point in the other coordinate system. How can I do that? In theory, some kind of transformation optimizing minimum distance between the points would be desirable.
In the end I would like to assign properties I obtained from the cells in one method to properties obtained with the other method.
Thank you.
7 个评论
darova
2019-10-7
Can you please attach some sample data? How final iamges should look like?
Yangfan Peng
2019-10-8
编辑：Yangfan Peng
2019-10-8
darova
2019-10-8
And what do you want to do with these (x,y)?
Yangfan Peng
2019-10-8
coordinate_a = [376 455;421 489;465 537;353 512;355 535;329 571;377 593;417 598;482 575;355 634]
coordinate_b = [168 88;107 138;69 194;126 229;163 232;199 267;272 239;228 210;235 155;215 68]
This is the code for the two different set of coordinates. They represent the same cells in different microscope images. However, the points are mixed up, so point 1 in coordinate_a is not point 1 in coordinate_b (see image below). I would like to match the points to know which points belong to which in the different coordinate systems, e.g. point 1 in coordinate_a is point 7 in coordinate_b
Generally speaking, it is an impossible problem. For example, suppose there were 3 points in each image with the A-points forming an equilateral triangle and the B-points an equilateral triangle in a different position. Then there is no unique matching between the points that can be determined purely from their relative spacing - the symmetry among the points is too perfect. There would have to be some a priori known asymmetry among the points in order to match them.
darova
2019-10-8
What if just find closest pairs of points? Using pdist2() for example?
Yangfan Peng
2019-10-8
采纳的回答
If the shear component of the deformation isn't too strong, the matchpoints function defined below might work. It will sort the rows of B to correspond with A and also return the corresponding permutation indices. It relies on absor, mentioned by Bruno, which you will have to download
function [permIndices,Bsorted]=matchpoints(A,B)
%
%IN:
%
% A: an Nx2 matrix of points
% B: an Nx2 matrix of points from A transformed and unordered
%
%OUT:
%
% permIndices: permutation indices of rows of B thought to match A
% Bsorted: the Nx2 permuted version of B
La=landmarks(A);
Lb=landmarks(B);
B3=B.'; B3(3,:)=0;
reg=absor( Lb,La,'doScale',1);
C3=(reg.s*reg.R)*B3+reg.t;
C=C3(1:2,:).';
[~,permIndices]=pdist2(C,A,'euclidean','Smallest',1);
if nargout>1
Bsorted=B(permIndices,:);
end
function L=landmarks(P)
G=pdist2(P,P); G(~G)=nan;
[i,j]=find( G==min(G(:)) ,1);
I=P(i,:);
J=P(j,:);
K=mean(P,1);
if norm(I-K)<norm(J-K)
[I,J]=deal(J,I);
end
L=[I;J;K].';
L(3,:)=0;
end
end
Applying it to your example data, I obtain,
A = [376 455;421 489;465 537;353 512;355 535;329 571;377 593;417 598;482 575;355 634];
B = [168 88;107 138;69 194;126 229;163 232;199 267;272 239;228 210;235 155;215 68];
[permIndices,Bsorted]=matchpoints(A,B)
plot(graph(1:10,1:10),'EdgeColor','none','XData',A(:,1),'YData',A(:,2));
hold on
plot(graph(1:10,1:10),'EdgeColor','none','XData',Bsorted(:,1),'YData',Bsorted(:,2));
hold off
8 个评论
darova
2019-10-8
Why so many times edited? Do you have some problems?
Matt J
2019-10-10
Just making improvements.
Yangfan Peng
2020-4-29
Dear Matt,
thanks you very much for the code and the help! I am incredibly sorry that I am just getting back to you now. I have been distracted and busy with other projects, but am continuing on this problem now. Your function works very well for "easy" cases (great help and thanks for that!!). However, if two points get very close together, these two points seem to overlap and match to one coordinate in Bsorted (see attached image). Is there a way to restrict this overlap? Can we add some margin of minimal distance between the Bsorted points while minimizing the distances?
I know that this answer comes late, but I would be grateful for your help.
Thanks. Best, Yangfan
Matt J
2020-4-30
I cannot copy/paste the A and B data for this example.
Yangfan Peng
2020-4-30
Matt J
2020-4-30
Do you have the Optimization Toolbox?
Yangfan Peng
2020-5-1
Here is a modification of the code that works on your example. You have a very large outlier in your coordinates, however, so it can't always be guaranteed that such cases would be robustly handled. Note that the code uses the FEX file minL1intlin.m (Download).
function [permIndices,Bsorted]=matchpoints(A,B)
%
%IN:
%
% A: an Nx2 matrix of points
% B: an Nx2 matrix of points from A transformed and unordered
%
%OUT:
%
% permIndices: permutation indices of rows of B thought to match A
% Bsorted: the Nx2 permuted version of B
La=landmarks(A);
Lb=landmarks(B);
B3=B.'; B3(3,:)=0;
reg=absor( Lb,La,'doScale',1);
C3=(reg.s*reg.R)*B3+reg.t;
C=C3(1:2,:).';
N=size(A,1);
e=ones(1,N);
E=speye(N);
Aeq=[ kron(E,e) ; kron(e,E) ]; beq=[e,e].';
Q=kron(E,C.');
d=reshape(A.',[],1);
lb=zeros(1,N^2);
ub=ones(1,N^2);
intcon=1:N^2;
P=minL1intlin(Q,d,intcon,[],[],Aeq,beq,lb,ub);
P=round(reshape(P,N,[]));
permIndices=(1:N)*P;
if nargout>1
Bsorted=B(permIndices,:);
end
function L=landmarks(P)
G=pdist2(P,P); G(~G)=nan;
[i,j]=find( G==min(G(:)) ,1);
I=P(i,:);
J=P(j,:);
K=mean(P,1);
if norm(I-K)<norm(J-K)
[I,J]=deal(J,I);
end
L=[I;J;K].';
L(3,:)=0;
end
end
更多回答（1 个）
Bruno Luong
2019-10-8
0 个投票
Look for the literature of image registration.
For simple rotation/scaling/translation you can use Matt J's submission
For more complex deformation, you need to apply spline deformation
There are a bunch of intermediate method for camera which take into account for camera cushion distortion or higher order. Pick one that is suitable for your need.
类别
在 帮助中心 和 File Exchange 中查找有关 Geometric Transformation and Image Registration 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
