29
$\begingroup$

After happily using the v9 image assistant to crop elliptically an image, and then the drawing tools to put a white disk in the middle, I turned this image:

into the following one:

Mathematica graphics

that can be imported with

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

The overall objective is to programatically count the number of those radial lines. Due to the lighting, there are parts of the image in which those lines are darker than its surroundings and others where it is lighter.

So far I haven't found a good way worth posting, so any pointer to a good solution would be appreciated. I have the feeling the image processing people will see better ways of dealing with this and I will be grateful to learn something. Thanks a lot

EDIT

A first approach with @RahulNarain's suggestion would be

int = ListInterpolation[
   ImageData[ColorConvert[im, "GrayLevel"]], {{-1, 1}, {-1, 1}}];
polInt = Function[t, int[0.95 Cos[2 \[Pi] t], 0.95 Sin[2 \[Pi] t]]];

Plot[polInt[t], {t, 0, 1}, AspectRatio -> 0.2, ImageSize -> Large, 
 PlotRange -> Full]

Mathematica graphics

Now,

ListLinePlot[Abs@Fourier[polInt@Range[0, 1 - 0.001, 0.001]], 
 PlotRange -> {{0, 500}, {0, 50}}]

Mathematica graphics

A better zoom shows the maximum at 98

However, manual counting (could be wrong) gave me 96, and nikie's approach is suggesting 97. 1 or 2 off count could be due to the light changes making the real ridge be the local minimum at some places and local maximum at others

$\endgroup$
9
  • $\begingroup$ I'm by no means proficient in image manipulation, but Sharpen may be what you're looking for - nest it 10 times. $\endgroup$
    – VF1
    Dec 8, 2012 at 6:10
  • $\begingroup$ Maybe: build an interpolating function from the image, sample it over a circular path something like $(0.95\cos\theta, 0.95\sin\theta)$, and then analyze the variation of intensity as a function of $\theta$. That reduces it to a one-dimensional signal processing problem. $\endgroup$
    – user484
    Dec 8, 2012 at 7:40
  • $\begingroup$ @RahulNarain nice idea, I'll give it a shot now :) $\endgroup$
    – Rojo
    Dec 8, 2012 at 7:41
  • $\begingroup$ The big version image is rather tiny... $\endgroup$
    – Yves Klett
    Dec 8, 2012 at 7:58
  • $\begingroup$ @YvesKlett haha, I thought it wasn't appropriate to upload the big version, but I can do it if you want $\endgroup$
    – Rojo
    Dec 8, 2012 at 7:59

2 Answers 2

24
$\begingroup$

A simple first try would be to rotate the image, then measure the distance to the original image:

img = ColorConvert[Import["https://i.stack.imgur.com/NNzNM.png"], 
   "Grayscale"];

Monitor[t = 
   Table[{i, 
     ImageDistance[img, ImageRotate[img, 360°/i, Full, Background -> White]]}, {i, 10, 
     200}], i];

ListLinePlot[t]

Mathematica graphics

Obviously, if the angle is exactly 360° / [number of radial lines], the match should be lowest, so the estimated count would be:

count = Extract[t, Position[t[[All, 2]], Min[t[[All, 2]]]]][[1, 1]]

which is 97

If I overlay 97 radial lines over your image, it seems as if the count wasn't too far of. I can't tell if it's exact, though:

center = 0.5 ImageDimensions[img];
Show[img,
 Graphics[
  {Red, Table[
    Line[{center + 80 {Cos[\[Phi]], Sin[\[Phi]]}, 
      center + 250 {Cos[\[Phi]], Sin[\[Phi]]}}], {\[Phi], 0, 2 \[Pi], 
     2 \[Pi]/count}]}]]

Mathematica graphics


EDIT: I've been playing with this some more, especially with the FFT idea. First, I polar-transform the image:

polar = ImageTransformation[
  img, #[[2]]*{Cos[#[[1]]], Sin[#[[1]]]} &, {500, 20}, 
  PlotRange -> {{0, 2 \[Pi]}, {0.9, 1.0}}, 
  DataRange -> {{-1, 1}, {-1, 1}}]

Mathematica graphics

Then I've applied a windowed Fourier transform to the mean of that signal:

mean = Mean[ImageData[polar]];
window = Array[HannWindow, Length[mean], {-1.5, 1.5}];
stft = Table[(Abs[Fourier[mean*RotateLeft[window, i]]][[
       80 ;; 120]]^2) // #/Max[#] &, {i, 0, 500}];
ArrayPlot[stft\[Transpose], ColorFunction -> GrayLevel, 
 DataRange -> {{1, 500}, {80, 120}}, FrameTicks -> True]

Mathematica graphics

The windowed Fourier transform looks as if the frequency isn't constant over the whole area. Which would make sense, if the center isn't perfect or if there's an affine/perspective transformation. Sadly, I'm not sure what to do with this, but I thought I'd post it, in case it gives somebody else an idea.

$\endgroup$
4
  • $\begingroup$ Very smart approach! Thanks a lot! What I'm trying using Rahul's suggestion so far is giving me 98 or 99, but counting them manually gives me 96. Haven't triple checked manually, you can imagine why $\endgroup$
    – Rojo
    Dec 8, 2012 at 8:48
  • $\begingroup$ If you want 96, you could use DistanceFunction -> "MutualInformationVariation". That's not very scientific, though. ;-) $\endgroup$ Dec 8, 2012 at 8:49
  • $\begingroup$ After your polar transform, do a column sum to make it 1d, and then MaxDetect would seem suited for the problem? $\endgroup$
    – JxB
    Dec 8, 2012 at 14:23
  • $\begingroup$ Nice edit. Perhaps you could have used SpectrogramArray for the windowed fourier transform $\endgroup$
    – Rojo
    Dec 8, 2012 at 18:13
10
$\begingroup$

I think the big bright spots in the lower half of the image are confusing the Fourier transform, because they lie right in between where two bright lines should be, and so they are exactly out of phase with the rest of the signal. How about I just throw them away?*

polInt = Function[t, 
   If[(t > 0.11 && t < 0.16) || (t > 0.89), 0.2, 
    int[0.95 Cos[2 \[Pi] t], 0.95 Sin[2 \[Pi] t]]]];
Plot[polInt[t], {t, 0, 1}, AspectRatio -> 0.2, ImageSize -> Large, 
 PlotRange -> Full]

enter image description here

Then the Fourier transform has a maximum at 96, as you want.

enter image description here

enter image description here


* This procedure is entirely unscientific. Chopping data at hard boundaries is not recommended by signal processing experts. Do not try this at home. Void where prohibited.

$\endgroup$
1
  • 2
    $\begingroup$ Actually... Does the output of Fourier have frequency 0 (the DC component) at position 1? Then the peak here is at 95 cycles, which might be a bad thing. $\endgroup$
    – user484
    Dec 8, 2012 at 11:16

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.