Re: Concentric contours about the centroid, having the same length, and interior to an initial contour.
- To: mathgroup at smc.vnet.net
- Subject: [mg44115] Re: [mg43954] Concentric contours about the centroid, having the same length, and interior to an initial contour.
- From: Daniel Lichtblau <danl at wolfram.com>
- Date: Thu, 23 Oct 2003 07:15:12 -0400 (EDT)
- References: <200310150859.EAA26388@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
Gilmar Rodr?guez Pierluissi wrote:
>
> I have a basic contour defined by:
>
> In[1]:
>
> contour = {{10545, 28012}, {10232, 28166}, {10027, 28218},
> {9714, 28115}, {9354, 28166}, {9042, 28115}, {8785, 28064},
> {8575, 27858}, {8160, 27858}, {7954, 27546}, {8108, 27238},
> {8421, 27238}, {8836, 27079}, {9042, 26874}, {9247, 26510},
> {9406, 26197}, {9611, 25992}, {9662, 25627}, {9560, 25268},
> {9196, 25007}, {8939, 24750}, {8678, 24750}, {8421, 24386},
> {8216, 24125}, {8160, 23919}, {8006, 23607}, {8006, 23401},
> {7852, 23037}, {7749, 22832}, {7539, 22678}, {7179, 22678},
> {6918, 22570}, {6713, 22468}, {6451, 22365}, {6036, 22365},
> {5728, 22468}, {5677, 22104}, {5625, 21693}, {5625, 21483},
> {5364, 21483}, {5159, 21950}, {5000, 22211}, {4846, 22519},
> {4538, 22570}, {4641, 22160}, {4949, 21950}, {4949, 21642},
> {5159, 21278}, {5313, 21068}, {5569, 20862}, {5569, 20708},
> {5518, 20447}, {5415, 20190}, {5415, 19929}, {5313, 19565},
> {5000, 19308}, {4743, 19308}, {4589, 19411}, {4431, 19565},
> {4225, 19672}, {4020, 19621}, {3758, 19621}, {3758, 19359},
> {4020, 19154}, {4328, 18944}, {4482, 18585}, {4897, 18739},
> {5210, 18841}, {5364, 18426}, {5467, 18118}, {5625, 17754},
> {5831, 17339}, {5625, 17026}, {5364, 16821}, {5056, 16667},
> {4897, 16200}, {5000, 15887}, {5159, 15579}, {5107, 15318},
> {4897, 15005}, {4538, 15005}, {4174, 15215}, {4020, 15420},
> {4071, 15112}, {4277, 14851}, {4589, 14697}, {4897, 14594},
> {4897, 14179}, {4897, 13866}, {4743, 13507}, {4795, 13143},
> {4795, 12830}, {4692, 12363}, {4589, 12107}, {4589, 11794},
> {4328, 11486}, {3964, 11327}, {3861, 10968}, {4020, 10604},
> {3861, 10086}, {3810, 9670}, {3707, 9255}, {3446, 8737},
> {3292, 8219}, {3343, 7804}, {3502, 7337}, {3189, 6977},
> {2928, 6926}, {2620, 6870}, {2466, 6613}, {2256, 6147},
> {2050, 6352}, {1948, 6613}, {1686, 6767}, {1481, 6716},
> {1117, 6716}, {804, 6926}, {650, 6819}, {547, 6562},
> {445, 6352}, {393, 5993}, {235, 5941}, {81, 5629},
> {286, 5475}, {701, 5423}, {911, 5372}, {963, 5111},
> {860, 4905}, {650, 4695}, {650, 4280}, {1014, 4126},
> {1271, 3869}, {1686, 3659}, {2050, 3505}, {2466, 3402},
> {2774, 3351}, {3086, 3454}, {3446, 3659}, {3861, 3710},
> {4174, 3659}, {4589, 3710}, {5056, 3608}, {5313, 3244},
> {5467, 2726}, {5677, 2469}, {5933, 2418}, {6195, 2208},
> {6451, 1900}, {6503, 1484}, {6713, 1171}, {6815, 812},
> {7179, 653}, {7488, 551}, {7800, 345}, {8057, 551},
> {8267, 397}, {8575, 187}, {8836, 33}, {8836, 238},
> {8626, 499}, {8524, 705}, {8472, 966}, {8626, 1171},
> {8267, 1274}, {8216, 1069}, {8057, 966}, {7800, 1069},
> {7539, 1171}, {7385, 1330}, {7282, 1690}, {7590, 1741},
> {7642, 1900}, {7642, 2208}, {7800, 2469}, {7590, 2674},
> {7282, 2572}, {7021, 2418}, {6815, 2418}, {6815, 2833},
> {6918, 3038}, {7123, 3090}, {6918, 3244}, {6713, 3454},
> {6605, 3556}, {6451, 3869}, {6713, 3972}, {7021, 4126},
> {7385, 4228}, {6918, 4331}, {6713, 4438}, {6451, 4592},
> {6195, 4228}, {5882, 4387}, {5728, 4695}, {5728, 5008},
> {5625, 5475}, {5728, 5890}, {6143, 5890}, {6661, 5890},
> {7072, 5941}, {6815, 6198}, {6400, 6198}, {6087, 6147},
> {6246, 6511}, {6451, 6819}, {6713, 7080}, {6969, 7495},
> {7231, 7495}, {7800, 7598}, {7436, 7752}, {7179, 7855},
> {7072, 8168}, {6764, 8480}, {6713, 9050}, {6918, 9465},
> {7123, 9824}, {7282, 10240}, {7590, 10968}, {8057, 11379},
> {8318, 11794}, {8785, 12312}, {8990, 12727}, {9406, 13091},
> {9457, 13661}, {9714, 14333}, {10027, 14748}, {10339, 15420},
> {10909, 15579}, {11632, 15733}, {11945, 15630}, {12355, 15318},
> {12981, 15215}, {13443, 14800}, {13653, 14333}, {13756, 13918},
> {13910, 13451}, {14222, 13091}, {14586, 12830}, {14689, 12415},
> {14638, 11948}, {14535, 11430}, {14274, 11071}, {13961, 10707},
> {13807, 10655}, {13602, 10394}, {13961, 10394}, {14274, 10553},
> {14428, 10342}, {14376, 9983}, {14535, 9773}, {14843, 9824},
> {14946, 9568}, {15104, 9358}, {15258, 9101}, {15361, 8788},
> {15566, 8840}, {15982, 8947}, {16243, 8891}, {16603, 8788},
> {16864, 8532}, {17279, 8373}, {17587, 8270}, {17849, 7911},
> {18208, 7649}, {18315, 7080}, {18469, 6767}, {18782, 6665},
> {18936, 6408}, {19193, 6044}, {19557, 5834}, {19972, 6044},
> {20178, 6352}, {20644, 6511}, {20747, 6870}, {20798, 7234},
> {20542, 7547}, {20388, 8014}, {20229, 8322}, {20075, 8634},
> {19921, 9152}, {20229, 9465}, {20490, 9722}, {20747, 9876},
> {20490, 9876}, {20178, 9824}, {19972, 9619}, {19660, 9619},
> {19608, 9876}, {19300, 10086}, {18885, 10137}, {18675, 10342},
> {18623, 10604}, {18782, 10809}, {18936, 11019}, {18572, 11173},
> {18469, 10912}, {18315, 10655}, {18208, 10291}, {17900, 10394},
> {17587, 10394}, {17177, 10342}, {16761, 10394}, {16551, 10450},
> {16295, 10604}, {15879, 10707}, {15566, 11019}, {15622, 11379},
> {15828, 11794}, {16033, 12107}, {15725, 12522}, {15622, 12779},
> {15515, 13143}, {15776, 13451}, {15776, 13764}, {15464, 14025},
> {15258, 14333}, {15048, 14543}, {14894, 14748}, {14843, 14902},
> {14689, 15164}, {14689, 15420}, {14586, 15836}, {14376, 16200},
> {14222, 16615}, {14120, 16872}, {13910, 16923}, {13756, 17026},
> {13499, 17236}, {13186, 17236}, {12981, 17390}, {12668, 17493},
> {12304, 17651}, {12047, 17703}, {11786, 17600}, {11427, 17544},
> {11011, 17600}, {10493, 17600}, {9924, 17651}, {9611, 17908},
> {9457, 18169}, {9714, 18585}, {10181, 18944}, {10391, 19462},
> {10339, 19878}, {10232, 20293}, {10181, 20657}, {10699, 20862},
> {11063, 21278}, {11375, 21847}, {11735, 22160}, {12201, 22468},
> {12514, 22314}, {12719, 21898}, {13032, 21586}, {13084, 22052},
> {13238, 22314}, {13238, 22519}, {12873, 22570}, {12719, 22883},
> {12719, 23350}, {13032, 23765}, {12925, 24125}, {12719, 24437},
> {12463, 24801}, {12355, 25058}, {12047, 25525}, {11837, 25940},
> {11632, 26253}, {11319, 26253}, {11011, 26253}, {10750, 26458},
> {10596, 26822}, {10596, 27182}, {10596, 27494}, {10545, 28012}};
>
> The length of the contour is given by:
>
> In[2]:
>
> Length[contour]
>
> Out[2]:
>
> 375
>
> A plot of this contour is given by:
>
> In[3]:
>
> plt1 = ListPlot[contour, PlotStyle -> {Hue[0], PointSize[.01]}]
>
> Out[3]: graphics
>
> The centroid of the contour is given by:
>
> In[4]:
>
> centroid = N[Plus @@ contour/Length[contour]]
>
> Out[4]: {9426.41, 12877.6}
>
> My question is:
>
> How can I generate (say) five contours, so that:
>
> (1.) Each contour (set) have the same length as the original contour
>
> (in this case, each contour should contain 375 {x,y} points).
>
> (2.) Each contour is equidistant, and concentric (with respect to
>
> the centroid) to the previous contour.
>
> (3.) The five contours are inside the original contour defined above?
>
> Thank you for your help!
In general it is not possible to achieve simultaneously all aspects of
what you request above. I can see at least three alternatives.
(1) You might simply rescale the contour around the centroid. This is
easy to do, as below.
scaledPoints[pts_,frac_] := Module[
{centroid=Apply[Plus,N[pts]]/Length[pts]},
Map[(#+centroid)&, Map[(#-centroid)&,pts]*frac]
]
plots = Table[ListPlot[scaledPoints[contour,j/5],
PlotStyle->{Hue[(j-1)/5],PointSize[.01]},
PlotJoined->True, DisplayFunction->Identity], {j,1,5}];
Show[plots, DisplayFunction->$DisplayFunctionnedit ]
If you execute the above you will see that the "contours" overlap and in
general do not appear (to me, at any rate) that they would make for a
good animation of a "shrinking" contour. This is probably not what you
want.
(2) As Jens-Peer Kuska suggests (see URL below), you might apply a
curvature-based gradient flow to shrink the curve to a point.
http://forums.wolfram.com/mathgroup/archive/2003/Oct/msg00271.html
This too will not in general prevent time snapshot "contours" from
intersecting. it will, however, retain one closed curve (this result due
to M. Grayson in 1987; it does not hold in higher dimensions). This may
or may not be important for your purposes. If it is, then as Kuska
suggested you will want to look into level set methods. See various
works by James Sethian and coauthors.
(3) Since the contour is closed and connected with no self-intersection,
you can regard it as enclosing a region. You can then create new
contours that move equidistantly (rather than at a curvature-dependent
rate) into that region. The drawback is that now the region will tend to
disconnect. If this is acceptable, then the way to proceed is via a
special case of level set methods known as "fast marching" (also due to
Sethian and collaborators).
As it happens, a notebook with Mathematica code for this already exists
and may be found under Downloads at the URL below.
http://library.wolfram.com/infocenter/Demos/4928/
To get started one must preprocess the data. I did this as follows.
ppp = ListPlot[contour, PlotJoined -> True, Axes -> False];
(* Fill in the closed curve. *)
qqq = ppp /. Line -> Polygon;
(* Rasterize the image. *)
Show[contourRasterized =
ImportString[ExportString[qqq, "PNG", ImageResolution -> 125],
"PNG"]];
(* Pull out a simple matrix of values. *)
pixels = contourRasterized[[1, 1]][[All, All, 1]];
(* Remove the boundary pixels as they are not necessarily correct for
our purposes. *)
pixels = Take[pixels, {2, -2}, {2, -2}];
(* Reverse so that the inside is light and outside is dark. This is to
accomodate the convention used in the fast marching code, in order to
proceed only inward with contours.*)
pixels2 = pixels /. {(n_ /; n < 128) :> 255., (n_ /; n >= 128) :> 0.};
(* You can get a good picture of the raster array using
ListContourPlot[pixels2, Contours -> {128}]
Had we not done all this work, we would get equidistant contours in both
directions, i.e. shrinking and expanding. *)
(* Now, using code in the notebook, you can do *)
F[x_, y_] = 1.;
Timing[distances = findDistances[pixels2];]
(* This takes about a minute to run to completion on my Linux machine
(1.5 GHz). You can get a look at the contours using ListContourPlot: *)
ListContourPlot[distances, ContourShading -> False, PlotRange -> All,
ColorFunctionScaling -> False, Contours -> {0, 5, 10, 15, 20, 25}];
(* Alternatively, find the hidden cell (near the end of the notebook)
with the version of mathematica code that uses Compile. Then do *)
Timing[distancesC = findDistancesC[pixels2];]
ListContourPlot[distancesC, ContourShading -> False, PlotRange -> All,
ColorFunctionScaling -> False, Contours -> {0, 5, 10, 15, 20, 25}];
(* This takes about 6 seconds. *)
( *For greater resolution, redo using: *)
Show[contourRasterized =
ImportString[ExportString[qqq, "PNG", ImageResolution -> 180],
"PNG"]];
(* Repeat above preprocessing and finding of distances (takes about 12
seconds using the Compile version), and plot using: *)
ListContourPlot[distancesC, ContourShading -> False, PlotRange -> All,
ColorFunctionScaling -> False, Contours -> {0, 8, 16, 24, 32, 40}];
As for conveniently locating a number of points on the various contours
equal to the number in the original contour, I have no immediate
suggestions.
Daniel Lichtblau
Wolfram Research
- References:
- Concentric contours about the centroid, having the same length, and interior to an initial contour.
- From: gilmar.rodriguez@nwfwmd.state.fl.us (Gilmar Rodr?guez Pierluissi)
- Concentric contours about the centroid, having the same length, and interior to an initial contour.