Take the 2-minute tour ×

I have a set of line features inside a particular polygonal boundary. For each line, I'd like to generate a polygon inside which every possible point is closer to the given line than to any other line in the layer. I've done this in the past for point input features using Delaunay triangulation, but if there's a similar process for doing it with line features I haven't been able to find it.

ETA: Geogeek's solution had occurred to me, but in straighter sections where the input lines have fewer vertices, the resulting polygons get way too close (even overlapping) a line that they shouldn't. Here, the red lines are my inputs, you can see the vertices and the Thiessen polygons generated from them.

enter image description here

Perhaps a quick and (very) dirty solution might be to convert each line to a plentiful set of evenly-spaced points (rather than the line's vertices only), generate Thiessen polygons from those, then dissolve them based on the originating line ID.

share|improve this question
Voronoi diagrams that include line segments along with points are not composed of "polygons"; rather, their cells have boundaries that can include portions of parabolas. For this reason, one of the most efficient and accurate ways to create Voronoi tessellations is to use a raster representation. ESRI calls this procedure Euclidean Allocation. –  whuber Mar 4 '13 at 17:06

2 Answers 2

up vote 6 down vote accepted

To illustrate a raster/image processing solution, I began with the posted image. It is of much lower quality than the original data, due to the superposition of blue dots, gray lines, colored regions, and text; and the thickening of the original red lines. As such it presents a challenge: nevertheless, we can still obtain Voronoi cells with high accuracy.

I extracted the visible parts of the red linear features by subtracting the green from the red channel and then dilating and eroding the brightest parts by three pixels. This was used as the base for a Euclidean distance calculation:

i = Import["http://i.stack.imgur.com/y8xlS.png"];
{r, g, b} = ColorSeparate[i];
string = With[{n = 3}, Erosion[Dilation[Binarize[ImageSubtract[r, g]], n], n]];

Relief plot

(All code shown here is Mathematica 8.)

Identifying the evident "ridges"--which must include all points that separate two adjacent Voronoi cells--and re-combining them with the line layer provides most of what we need to proceed:

ridges = Binarize[ColorNegate[
   LaplacianGaussianFilter[DistanceTransform[ColorNegate[string]], 2] // ImageAdjust], .65];
ColorCombine[{ridges, string}]

Combined images

The red band represents what I could save of the line and the cyan band shows the ridges in the distance transform. (There is still a lot of junk due to the breaks in the original line itself.) These ridges need to be cleaned and closed up through a further dilation--two pixels will do--and then we can identify the connected regions determined by the original lines and the ridges between them (some of which need explicitly to be recombined):

  ColorNegate[ImageAdd[ridges, Dilation[string, 2]]]] /. {2 -> 5, 8 -> 0, 4 -> 3} // Colorize, 2]

What this has accomplished, in effect, is to identify five oriented linear features. We can see three separate linear features emanating from a point of confluence. Each has two sides. I have considered the right side of the two rightmost features as being the same, but have otherwise distinguished everything else, giving the five features. The colored areas show the Voronoi diagram from these five features.


A Euclidean Allocation command based on a layer that distinguishes the three linear features (which I did not have available for this illustration) would not distinguish the different sides of each linear feature, and so it would combine the green and orange regions flanking the leftmost line; it would split the rightmost teal feature into two; and it would combine those split pieces with the corresponding beige and magenta features on their other sides.

Evidently, this raster approach has the power to construct Voronoi tessellations of arbitrary features--points, linear pieces, and even polygons, regardless of their shapes--and it can distinguish sides of linear features.

share|improve this answer
A similar solution is illustrated at mathematica.stackexchange.com/questions/20696/…. –  whuber Mar 5 '13 at 21:14

i think you can:

  • Convert line vertices to points(line_points).
  • make voronoi polygons using the points(line_points).
  • dissolve the resulted polygons using either an id attribute which has been saved from line layer , or by a spatial join with line layer.

i wish i have really understood your question , if not you can you provide a drawing, to explain more your needs.

share|improve this answer
I think you understood it, and that solution did occur to me, but you run into problems where the lines have fewer vertices. I'll update my question with a screenshot. –  Dan C Mar 4 '13 at 17:01
This would work alright if you made the points denser along the line. Although a raster based approach (as whuber mentions in the comments on the question) I would suspect would be much more efficient than this. –  Andy W Mar 4 '13 at 18:08

Your Answer


