10
$\begingroup$

MorphologicalComponents can break an image into components.

How can we measure the pairwise distance between these components efficiently? Performance is important.

Here's an example (that I took from here and cleaned up).

enter image description here

It has 112 components:

im = Import["https://i.stack.imgur.com/U7myR.png"];

MorphologicalComponents[im] // Max
(* 112 *)

Thus we need to find Binomial[112, 2] == 6216 distances efficiently.

Computing the pairwise distances of components would allow us to determine which are adjacent. This image is just one example. I am looking for a general approach that would work on any map.

I suspect that there must be something built-in that I did not yet find. I am aware of DistanceTransform, but I am unable to make it perform sufficiently well. I solved this particular map efficiently using ImageMesh, but that approach is not generalizable: the mesh is not usually accurate. It is also difficult to map back to the original component indices.

$\endgroup$
  • $\begingroup$ If i understand correctly can you use this to find neighbour indices: ComponentMeasurements[ Dilation[MorphologicalComponents@img, 1], "Neighbors"] and once we know the neighbours we can compute the pairwise distances? $\endgroup$ – Ali Hashmi May 1 '18 at 13:45
  • $\begingroup$ @AliHashmi Yes, that works in this particular case for finding neighbours. It would not work if the distance between components were not 1 or 2 pixels only. But the question is how to find the pairwise distances, not just how to find neighbours. $\endgroup$ – Szabolcs May 1 '18 at 13:48
  • 2
    $\begingroup$ It's not clear what you mean by the "distance between components." Do you mean the Hausdorf distance? The L2 distance from the centroids of components? $\endgroup$ – bill s May 1 '18 at 14:14
  • $\begingroup$ @bills The smallest distance between any two points $p_1$ and $p_2$ such that $p_1$ is in the first component and $p_2$ is in the second one. $\endgroup$ – Szabolcs May 1 '18 at 14:27
  • $\begingroup$ @Szabolcs does MorphologicalGraph help in general case? if you use img//ColorNegate//MorphologicalGraph ? $\endgroup$ – Ali Hashmi May 1 '18 at 17:58
7
$\begingroup$

How is this? It first groups the points of each component together and applies Nearest and to find minimal distances in the sense of

$$d(A,B) = \inf_{a \in A} \inf_{b\in B} |a-b|.$$

This extracts the boundary pixels of each component (also of the background component):

im = Import["https://i.stack.imgur.com/U7myR.png"];
comp = SparseArray@With[{
       perimeter = 
        SparseArray[ImageData[MorphologicalPerimeter[im]]],
       comp = MorphologicalComponents[im]},
      Plus[
       Times[UnitStep[-comp], 
        MorphologicalTransform[perimeter, Max]],
       Times[comp + 1, perimeter]
       ]
      ]; // AbsoluteTiming // First
Colorize@comp

0.029138

enter image description here

This groups the boundary points' coordinates according to their components:

comppts = Developer`ToPackedArray /@ KeySort@Association@Reap[
         MapThread[
          Sow[#1, #2] &, 
          {comp["NonzeroPositions"],comp["NonzeroValues"]}],
         _, Rule
         ][[2]]; // AbsoluteTiming // First

0.041497

Now we use Nearest to compute pairwise distances:

n = Length[comppts];
data = ParallelTable[
    With[{nf = Nearest[comppts[i] -> "Distance"]},
     {i, j} -> Min[nf[comppts[j]]]
     ],
    {i, 1, n}, {j, i + 1, n}
    ]; // AbsoluteTiming // First

A = Normal[SparseArray[Join @@ data, {n, n}]];
A += Transpose[A];

MatrixPlot[A]

0.424216

enter image description here

So, everything is done in roughly half a second.

If you are not interested in the background component, you can use

comppts = 
    Association@ComponentMeasurements[im, "PerimeterPositions"]; // 
  AbsoluteTiming // First
n = Length[comppts];
data = ParallelTable[
     With[{nf = Nearest[comppts[i][[1]] -> "Distance"]},
      {i, j} -> Min[nf[comppts[j][[1]]]]],
     {i, 1, n},
     {j, i + 1, n}
     ]; // AbsoluteTiming // First
A = Normal[SparseArray[Join @@ data, {n, n}]];
A += Transpose[A];
MatrixPlot[A]

0.006328

0.384601

$\endgroup$
2
$\begingroup$

I am not sure that I would call this "efficient", but I think it does what you are asking: given a set of component masks and a pair of component labels, it should return the minimum distance between the components with the given labels.

Get the image:

maze = Import["https://i.stack.imgur.com/U7myR.png"];

Get an association of labels to component masks, and assemble the pairs of component labels:

compMasks = Association @@ ComponentMeasurements[maze, "Mask"];
compLabels = Keys[compMasks];
labelPairs = Tuples[compLabels, 2];

Define a function meant to take the sparse array mask of two regions as separate arguments and return the minimum Euclidean distance between them:

compDistance[mask1_, mask2_] := Block[
    {m1Back, genDist, m2Im},
    m1Back = ColorNegate[Image[mask1]];
    genDist = DistanceTransform[m1Back];
    m2Im = Image[mask2];
    ImageMeasurements[m2Im genDist, "Min", Masking -> m2Im]
]

Define another function that calls that one to compute distance between labelled components:

labelledDistance[maskAssoc_, {label1_, label2_}] := compDistance[maskAssoc[label1], maskAssoc[label2]]

See the results and timing on the first 10 component pairs:

AbsoluteTiming[{#, labelledDistance[compMasks, #]} & /@ labelPairs[[1 ;; 10]]]

returns:

{0.193874, {{{1, 1}, 0.}, {{1, 2}, 3.}, {{1, 3}, 3.}, {{1, 4}, 3.}, {{1, 5}, 3.}, {{1, 6}, 3.}, {{1, 7}, 23.1948}, {{1, 8}, 63.0714}, {{1, 9}, 103.044}, {{1, 10}, 143.031}}}

To compute the distances for all pairs takes ~4m.

 AbsoluteTiming[{#, labelledDistance[compMasks, #]} & /@ labelPairs;]
 {242.331, Null}
$\endgroup$
2
$\begingroup$

not a general answer (it is specific to the maze above or when the pixel distance is very small). Assuming distance to be Euclidean, the code below finds distance between neighbours (ignores any pair that is not a neighbour):

comp = MorphologicalComponents@img // Dilation[#, 1] &;
pairs = DeleteDuplicates[Sort/@ Flatten@Map[Thread, ComponentMeasurements[comp,"Neighbors"]]];
neigh = GroupBy[pairs, Keys -> Values];
cent = <|ComponentMeasurements[comp, "Centroid"]|>;

distance[{p_, q___}] := EuclideanDistance[p, #] & /@ {q};

dist = KeyValueMap[distance@FlattenAt[{Lookup[cent, #1], Lookup[cent, #2]}, {2}] &,neigh];

res = Thread[(pairs /. Rule -> List) -> Flatten@dist];

SparseArray[Join[MapAt[Reverse, res, {All, 1}], res]]//MatrixPlot

enter image description here

Absolute timing: 0.0315194

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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