I'm interested in Centroidal Voronoi tessellation. Here are two examples of unit disk tessellations:

enter image description here

They are the same with respect to rotation.
Recently, I've found FindGeometricTransform. For a square it works well, but for disk I see a problem.
So here are data of mentioned tessellations (we need to find transformation only for seeds):

    seeds1 = {{0.775108, -0.241578}, {-0.792184, -0.056252}, {-0.101257, \
-0.787851}, {-0.182002, 0.791042}, {-0.150633, 
  0.393956}, {0.380565, -0.180766}, {0.713349, 
  0.293661}, {-0.340764, -0.0153295}, {0.208997, 
  0.192923}, {0.412246, -0.662962}, {-0.0408312, -0.337878}, \
{-0.63377, 0.455064}, {-0.556108, -0.526441}, {0.350942, 0.68688}}
    seeds2 = {{-0.802537, 0.109056}, {-0.227327, -0.761387}, {0.761821, 
  0.164738}, {-0.423125, -0.0215327}, {0.291033, 
  0.300411}, {0.705383, -0.364998}, {0.300652, -0.705072}, {0.488612, 
  0.648416}, {-0.0351534, 0.769755}, {-0.119228, 
  0.255533}, {0.306248, -0.150524}, {-0.0965345, -0.329819}, \
{-0.523693, 0.565791}, {-0.651289, -0.434905}}

Define distance between sets as:

arrDist[arr1_, arr2_] := Max[SquaredEuclideanDistance[
       First@Nearest[arr1, #], #] & /@ arr2];

and start to rotate one of set. We get such picture:

data = Table[arrDist[seeds1, RotationTransform[t][seeds2]], 
{t, 0, 2 \[Pi], 0.01}]; ListPlot[data, Joined -> True]

enter image description here

Obviously there are many local minima and one global. And if value of global minimum less than some tolerance (which is derived from region discretization properties), sets treated as the same.
But how to find it?!

FindGeometricTransform failed:

FindGeometricTransform[seeds1, seeds2]

enter image description here
FindMinimum failed 'cause we need global minimum:

goal[arr1_,arr2_,t_]:= arrDist[arr1, RotationTransform[t][arr2]] /; 
  MatrixQ[arr1, NumericQ] && MatrixQ[arr2, NumericQ] && 

FindMinimum[{goal[seeds1,seeds2,t], 0<=t<=2 \[Pi]},t]
{0.062312, {t -> 1.275423}}

NMinimize is better, but it's not global minimum:

NMinimize[{goal[seeds1, seeds2, t], 0 <= t <= 2 \[Pi]}, t,
 Method -> "DifferentialEvolution"]
{0.0207009, {t -> 0.027393}}

Manually (using above picture) we easily get

goal[seeds1, seeds2, 5.06]

What does all this mean and how do I end up comparing sets in terms of rotations?!


3 Answers 3


I consider it to be a bug or at least inefficient algorithm.

When you provide permutation of seeds2 that is consistent with order of seeds1 then FindGeometricTransform is able to find the correct transformation.

I found the permutation of seeds2 manually.

On the picture you can see in red points of seeds1, in green points of seeds2 and in black points that are transformation of seeds2.

pseeds2 = seeds2[[{5, 7, 2, 9, 3, 8, 6, 1, 13, 10, 11, 12, 4}]];

FindGeometricTransform[seeds1, pseeds2, TransformationClass -> "Rigid"]

Graphics[{{PointSize[0.04], Red, Point[seeds1]}, {PointSize[0.04], 
   Green, Point[seeds2]}, {PointSize[0.02], Black, 
   GeometricTransformation[Point[seeds2], %[[2]]]}}]

enter image description here


I used DeleteDuplicates on seeds in case there are duplicates, then sorted both seeds by angle (with error 0.1) and by norm.

Then I found the minimum of total squared distance for all 14 rotations of sseeds2 and sseeds1 and choose the best one.

It is at rotation angle fi->-1.22166 which is consistent with yours t->5.06 because -1.22166+2*π=5.06153. Total squared distance is 0.000866629.



sseeds1=SortBy[seeds1,{Round[Mod[ArcTan[#[[1]],#[[2]]],2 π],0.1]&,Norm[#]&}];
sseeds2=SortBy[seeds2,{Round[Mod[ArcTan[#[[1]],#[[2]]],2 π],0.1]&,Norm[#]&}];






Well, for hours of exploration looks like next combination works. Firstly, distance function should be cumulative, not maximise:

arrDist[arr1_, arr2_] := Total[SquaredEuclideanDistance[
       First@Nearest[arr1, #], #] & /@ arr2];

This makes goal function smoother:
enter image description here

Secondly, Method "RandomSearch" better for this task. So NMinimize with such options find global minimum in all cases that I've explored and may be good alternative for FindGeometricTransform in case of continuous rotations.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.