20

I start with a photograph of a shape (physically made by the flow of a liquid into another), of which I can extract the border, manually or using Mathematica's feature detect feature :

Example

Using the method described here, I extract some points of the shape, in a standard list form ({{x,y},{x,y}, ...}). Here is the example plot of such data :

A plot of the obtained data

How to :

  1. Obtain the shape's area (area of the enclosed zone)
  2. How to get an algebraic fit of the shape (and/or part of it)?
  3. How to compare one shape against another, on their respective "jaggedness".

Question 3 is more of a mathematical question, but I'm just searching for an approximate comparison tool, more in the spirit of the following example than something absolute:

this example

I expected CornerFilter to work, but it seems to give no result whatsoever. As for 2, I can fit small part of the curve using Fit[], but the general shape has multiple point with the same x, which forbid this.

William Briand
  • 303
  • 2
  • 6
  • 3
    For the first: here's a short routine for computing the signed area, assuming your points have already been ordered either clockwise or anticlockwise: PolygonSignedArea[pts_?MatrixQ] := Total[Det /@ Partition[pts, 2, 1, {1, 1}]]/2. If the points aren't already sorted, Sjoerd says to look at ListCurvePathPlot[] (and the related function FindCurvePath[]). – J. M.'s missing motivation May 22 '12 at 18:58
  • 2
    For "respective jaggedness" you could use Box Counting method used on fractal shapes: this http://en.wikipedia.org/wiki/Box_counting or this http://en.wikipedia.org/wiki/Fractal_dimension_on_networks – Vitaliy Kaurov May 22 '12 at 19:35
  • @J.M. are you sure your routine sorts them anticlockwise (or clockwise)? Also, even if you sort by the polar angle, that does not guarantee that a ray coming out of the centre of the shape will not intersect it twice. – acl May 22 '12 at 19:46
  • @acl: hmm, yes; that is indeed a caveat of sorting by polar angle. I don't have a better approach at the moment. – J. M.'s missing motivation May 22 '12 at 19:49
  • @J.M. I think you also need to use ArcTan[x,y] rather than the form you are currently using – acl May 22 '12 at 19:50
  • @acl: I am indeed implicitly using the two-argument form; Apply[ArcTan, {x, y}] is a valid application, yes? – J. M.'s missing motivation May 22 '12 at 19:51
  • @J.M. let's go to chat – acl May 22 '12 at 20:01
  • 3
    Could you tell us your requirements regarding your "algebraic fit"? – Sjoerd C. de Vries May 22 '12 at 20:20
  • A piecewise and/or parametric polynomial would be fine. I certainly could do it "manually" with some thinking, but I aim at automating the process as much as possible. Could this be done simpler by taking the shape as a polar plot and converting it to cartesian ? – William Briand May 22 '12 at 21:55
  • @WilliamBriand it has overhangs, so it doesn't work; try curvLoc = (Position[ImageData[img], 0]); centredpts = Map[# - N@Mean[curvLoc] &, curvLoc]; topolars = {Sqrt[#[[1]]^2 + #[[2]]^2], ArcTan[#[[1]], #[[2]]]} &; frompolars = {#[[1]]*Cos[#[[2]]], #[[1]]*Sin[#[[2]]]} &; polarpts = Reverse[SortBy[topolars /@ centredpts, Last], 2]; fint = Interpolation[polarpts]; PolarPlot[fint[\[Theta]], {\[Theta], -\[Pi], \[Pi]}] to see what I mean. this happens because ListPlot[polarpts[[All, 2]]] – acl May 22 '12 at 22:39
  • (the first component of polarpts is the polar angle, the second the radius; if you plot the radius, you see there are "overhangs") – acl May 22 '12 at 22:40
  • @William Re: "A piecewise and/or parametric polynomial would be fine." How about just using an interpolating function? Both for this and generally for this task, reducing the number of points using the Ramer-Douglas-Peucker algorithm might be useful. I have a slow but working implementation lying around somewhere, just ping me if you want it. – Szabolcs May 23 '12 at 07:31
  • @Szabolcs I'd be interested in that... – acl May 23 '12 at 11:02
  • If you can order the vertices consistently, that is, so they form a traversal either CW or CCW of the region, then you can separately make interpolations for the x and y values. This will give a parametrization in terms of the two interpolations. – Daniel Lichtblau May 23 '12 at 22:34

1 Answers1

24

You can use Mathematica's image processing functions for questions 1 and 3. Here's how:

1: Area

img = Binarize@Import["https://i.stack.imgur.com/I2gkK.jpg"] ~Erosion~ 1;
(m = MorphologicalComponents[img]) // Colorize

enter image description here

To get the area of the pink part in sq. pixels, use ComponentMeasurements:

2 /. ComponentMeasurements[{m, img}, "Area"]
(* 25168.1 *)

3: Jaggedness

Given two shapes with a similar area, the more jagged one will have a larger perimeter. So a simple way to approximate "jaggedness" would be to use the same function, ComponentMeasurements and compare the perimeters (provided the areas are similar)

2 /. ComponentMeasurements[{m, img}, "PerimeterLength"]
(* 1352 *)

If the areas are not the same, you could find the "EquivalentDiskRadius" (which gives you the radius of the circle with the same area as the shape) for each shape and see by what percentage the corresponding perimeter lengths are off from that of the respective circles.

rm -rf
  • 88,781
  • 21
  • 293
  • 472
  • 3
    +1, though I'm not sure about "Given two shapes with a similar area, the more jagged one will have a larger perimeter." What about two rectangles with the same area but different circumference? Intuitively I'd say both have the same jaggedness. – Sjoerd C. de Vries May 22 '12 at 20:55
  • 1
    @SjoerdC.deVries I see your point. I took the spirit of the question to be along the lines of "Here's this irregular blob and here's another irregular blob. How can I tell which is approximately smoother, for some definition of 'smooth'." – rm -rf May 22 '12 at 21:01
  • 1
    +1, and I think ComponentMeasurements can answer question 2, too, with the Length/Width/SemiAxis/Orientation/Elongation/Eccentricity measurements – Niki Estner May 23 '12 at 09:23
  • How to mark the number of each area if I have multiple shapes in a picture? When I do the similar thing, I can't tell the correspondence between different area and shape. – BNHSX Jul 12 '16 at 18:16
  • +1 I get slightly different results with 10.4 – bobbym Aug 08 '16 at 18:23