6
$\begingroup$

overlay- fiber network, in red - skeletonize fibers, blue - position of junctions

I work on a structural analysis of SEM image of a gel. I extracted the center of fibers by binarizing ane seletonizing the picture in ImageJ (in red). Then I extracted the positions of junctions where fibers cross using Ultimate Points also in ImageJ(in blue). I want to extract the information on how many fibers cross at a given junction (e.g. histogram). I thought of making a matrix of size 5 or 7 with the position of junction in the middle and then summing up all integers at the sides of the matrix, but I have no idea how to do it. Any idea how to proceed?

Files with positions of points for fibers and junctions below: https://www.dropbox.com/s/93onef13rrkdxrz/fibers.txt?dl=0 https://www.dropbox.com/s/1pjbf2nu8ww17pe/junctions.txt?dl=0

$\endgroup$
4
$\begingroup$

I try to show how you might use Nearest for detecting neighboting white pixels and Graph for connectedness analysis.

This is not perfect as it counts only the neigboring pixels of a junction: For me, a junction is a pixel with more than two direct neighbors. So there might be many false positives.

<< Developer`;
SetDirectory[FileNameJoin[{$HomeDirectory, "Downloads"}]]
fibers0 = ToPackedArray@Import["fibers.txt", "Table"];
(*cutting out a reactangular piece of the data*)

A = SparseArray[fibers0 -> 1];
i0 = 1; i1 = 80;
j0 = 1; j1 = 50;
fibers = A[[i0 ;; i1, j0 ;; j1]]["NonzeroPositions"];
edges = DeleteDuplicates[Sort /@ Join @@ Map[
      x \[Function] Thread[UndirectedEdge[First[x], Rest[x]]],
      Nearest[fibers -> Automatic, fibers, {\[Infinity], 1}, 
       DistanceFunction -> ChessboardDistance]
      ]
   ];
G = Graph[Range[Length[fibers]], edges,
   VertexCoordinates -> fibers,
   VertexSize -> .5,
   VertexLabels -> Placed["Name", Tooltip]];
degrees = AssociationThread[VertexList[G], VertexDegree[G]];
degreegroups = KeySort@GroupBy[Normal[degrees], Last -> First];
junctions = 
 Flatten[Values[KeyDrop[degreegroups, {0, 1, 2}]]]; counts = 
 Length /@ degreegroups;
BarChart[counts, PlotRange -> All, ChartLabels -> Keys[counts]]

enter image description here

Here is the colorization of the graph with respect to junction degree:

HighlightGraph[G, degreegroups /@ Keys[KeyDrop[degreegroups, {1, 2}]],
  ImageSize -> Large]

enter image description here

As you can see, there are many triangles and some quadrangles that ought to be collapsed into single junctions. At least the triangles can be handled by the Mathematica package IGraphM by Szabolcs.

$\endgroup$
  • $\begingroup$ A typo in the argument of Sort: Sort[junctions0] == Sort[junctions] ? $\endgroup$ – Anjan Kumar Feb 26 '18 at 7:53
  • $\begingroup$ @AnjanKumar Good point! No wonder that it evaluates to True. That's a bit embarassing, you know. =/ $\endgroup$ – Henrik Schumacher Feb 26 '18 at 9:28
  • $\begingroup$ very interesting approach, but definitely biased with your understanding of junction (ultimate Points gave 1200 points whereas here there are more than 20 000). And a question, why there are no three way connections there then? Nevertheless, thumbs up! $\endgroup$ – dziakku Feb 26 '18 at 9:48
  • $\begingroup$ Thanks. I correct the a bug with axis labels. It turns out: most "junctions" have valence three. Could tell us how the term junction is defined? what are "Ultimate Points"? I am still quite sure that a combination of Nearest and Graph might help you. $\endgroup$ – Henrik Schumacher Feb 26 '18 at 10:01
  • $\begingroup$ Now, that it shows that three is the most abundant number of fibers at a junction it makes so much more sense. Ultimate Points are maxima of the Euclidian distance map (you erode until you end up with a point). Most likely I will use your approach, but redefined the junctions using this ImageJ Ultimate Points to reduce number of those false positives. Nevertheless, once more, insightful approach! $\endgroup$ – dziakku Feb 26 '18 at 10:34
0
$\begingroup$

here I how I approached to problem (code is not clean, but works)

skeleton

junctions

Once more, I sued ImageJ to skeletonize image and find the Ultimate Points to get junctions (Pics above) First, I import files, read in the position of pixels and visualize it:

Import["skeleton.jpg"] // ImageRotate[#, -90 Degree] & // Binarize // 
    ImageData // (mask = #) & // 
  Position[#, 0] & // (points = #) & // Point //   
Graphics[{Red, #}] & // (sk = #) &;
Import["junctions.jpg"] // ImageRotate[#, -90 Degree] & // Binarize //
    ImageData // Position[#, 0] & // (pointsUlti = #) & // Point //
Graphics[{PointSize[0.01], Green, #}] & // (ultiPoints = #) &;
Show[sk, ultiPoints]

Overlay

Then I delete the points for junctions that are too close to edges:

DeleteCases[pointsUlti, 
x_ /; x[[1]] < 4 || x[[1]] > 836 || x[[2]] < 4 || 
x[[2]] > 836] // (newPoints = #) &;

Next, because Mathematica is interpreting data differently I reverse the matrix with data and define how many neighours there are in 3x3 matrix with the center positioned where the juncion is:

(mask - 1) // Abs // (list = #) &;
i = 1;
l = 1;
threeFibers3 = {};
fourFibers3 = {};
other3 = {};
threex3 = {};
Table[newPoints[[i]] // 
   list[[#[[1]] - l ;; #[[1]] + l, #[[2]] - l ;; #[[2]] + l]] & //
   Total // Total // (value = # - 1) & // AppendTo[threex3, #] &;
If[value == 3, AppendTo[threeFibers3, newPoints[[i]]], 
If[value == 4, AppendTo[fourFibers3, newPoints[[i]]], 
 AppendTo[other3, newPoints[[i]]]]], {i, 1, (newPoints // Length)}];
threeFibers3 // Point // 
Graphics[{PointSize[0.008], Green, #}] & // (ttt = #) &;
fourFibers3 // Point // 
Graphics[{PointSize[0.008], Blue, #}] & // (tttt = #) &;
other3 // Point // 
Graphics[{PointSize[0.01], Black, #}] & // (tt = #) &;
Show[sk, tt, ttt, tttt]
Histogram[threex3, {1}, "Probability"]

Then you get junctions with 3 fibers in green, with 4 blue, and all other with black 3x3matrix

Next, I check if the results with a bigger matrix 5x5 give a more realistic results. I use 5x5 atrix centered at junction position - 3x3 matrix just to get values from edges:

 i = 1;
 l = 2;
 threeFibers5 = {};
 fourFibers5 = {};
 other5 = {};
 fivex5 = {};
 Table[newPoints[[i]] // 
    list[[#[[1]] - l ;; #[[1]] + l, #[[2]] - l ;; #[[2]] + l]] & //
    Total // Total // (value = # - (threex3[[i]] + 1)) & // 
 AppendTo[fivex5, #] &;;
If[value == 3, AppendTo[threeFibers5, newPoints[[i]]],
 If[value == 4, AppendTo[fourFibers5, newPoints[[i]]],
 AppendTo[other5, newPoints[[i]]]]]
, {i, 1, (newPoints // Length)}];
threeFibers5 // Point // 
Graphics[{PointSize[0.008], Darker[Green], #}] & // (ttt5 = #) &;
fourFibers5 // Point // 
Graphics[{PointSize[0.008], Blue, #}] & // (tttt5 = #) &;
other5 // Point // 
Graphics[{PointSize[0.01], Black, #}] & // (tt5 = #) &;
Show[sk, tt5, ttt5, tttt5]
Histogram[
Table[3, (threeFibers5 // Length)]~Join~
Table[4, (fourFibers5 // Length)]~Join~
Table[5, (other5 // Length)], {1}, "Probability"]

This gives you a more probable results if you look at them, when you compare the approach with 3x3 matrix and 5x5 matrix 5x5matrix

Once more, sorry for messy code, but that's current working condition. This way you have a realistic view on the junctions ~80% have 3 fibers,~17% 4 and rest is "other"

$\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.