Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

I am trying to find out the points of inflection , i.e. points where the curves in a line start and end. If you look at the image, the green line may be a road or a stream, and the black points are the points where the curves start and end. enter image description here

What would be the high level steps to automate the generation of these points? I have ArcGIS desktop and am quite handy with ArcObjects.

share|improve this question
    
The source data is a polyline made of line segments and you want to aproximate it with curves, or does it already has arc segments? –  U2ros Oct 18 '12 at 12:47
    
Currently, it is made of line segments. –  Devdatta Tengshe Oct 18 '12 at 12:53
1  
The illustration in this question looks remarkably like one published at esri.com/news/arcuser/0110/turning.html. –  whuber Sep 23 '13 at 18:05
    
@whuber: Very astute observation. That was exactly the data source I had used to create the image. –  Devdatta Tengshe Sep 24 '13 at 3:33
add comment

3 Answers

up vote 9 down vote accepted

When the curve is comprised of line segments, then all interior points of those segments are inflection points, which is not interesting. Instead, the curve should be thought of as being approximated by the vertices of those segments. By splining a piecewise twice-differentiable curve through those segments, we can then compute the curvature. An inflection point, strictly speaking, is then a place where the curvature is zero.

In the example there are longish stretches where the curvature is nearly zero. This suggests that the indicated points ought to approximate the ends of such stretches of low-curvature regions.

An effective algorithm will therefore spline the vertices, compute the curvature along a dense set of intermediate points, identify ranges of near-zero curvature (using some reasonable estimate of what it means to be "near"), and mark the endpoints of those ranges.

Here is working R code to illustrate these ideas. Let's begin with a line string expressed as a sequence of coordinates:

xy <- matrix(c(5,20, 3,18, 2,19, 1.5,16, 5.5,9, 4.5,8, 3.5,12, 2.5,11, 3.5,3, 
               2,3, 2,6, 0,6, 2.5,-4, 4,-5, 6.5,-2, 7.5,-2.5, 7.7,-3.5, 6.5,-8), ncol=2, byrow=TRUE)

Spline the x and y coordinates separately to achieve a parametrization of the curve. (The parameter will be called time.)

n <- dim(xy)[1]
fx <- splinefun(1:n, xy[,1], method="natural")
fy <- splinefun(1:n, xy[,2], method="natural")

Interpolate the splines for plotting and computation:

time <- seq(1,n,length.out=511)
uv <- sapply(time, function(t) c(fx(t), fy(t)))

We need a function to compute the curvature of a parameterized curve. It needs to estimate the first and second derivatives of the spline. With many splines (such as cubic splines) this is an easy algebraic calculation. R provides the first three derivatives automatically. (In other environments, one might want to compute the derivatives numerically.)

curvature <- function(t, fx, fy) {
  # t is an argument to spline functions fx and fy.
  xp <- fx(t,1); yp <- fy(t,1)            # First derivatives
  xpp <- fx(t,2); ypp <- fy(t,2)          # Second derivatives
  v <- sqrt(xp^2 + yp^2)                  # Speed
  (xp*ypp - yp*xpp) / v^3                 # (Signed) curvature
  # (Left turns have positive curvature; right turns, negative.)
}

kappa <- abs(curvature(time, fx, fy))     # Absolute curvature of the data

I propose to estimate a threshold for zero curvature in terms of the extent of the curve. This at least is a good starting point; it ought to be adjusted according to the tortuosity of the curve (that is, increased for longer curves). This will later be used for coloring the plots according to curvature.

curvature.zero <- 2*pi / max(range(xy[,1]), range(xy[,2])) # A small threshold
i.col <- 1 + floor(127 * curvature.zero/(curvature.zero + kappa)) 
palette(terrain.colors(max(i.col)))                        # Colors

Now that the vertices have been splined and the curvature computed, it remains only to find the inflection points. To show them we may plot the vertices, plot the spline, and mark the inflection points on it.

plot(xy, asp=1, xlab="x",ylab="y", type="n")
tmp <- sapply(2:length(kappa), function(i) lines(rbind(uv[,i-1],uv[,i]), lwd=2, col=i.col[i]))
points(t(sapply(time[diff(kappa < curvature.zero/2) != 0], 
       function(t) c(fx(t), fy(t)))), pch=19, col="Black")
points(xy)

Plot

The open points are the original vertices in xy and the black points are the inflection points automatically identified with this algorithm. Because curvature cannot reliably be computed at the endpoints of the curve, those points are not specially marked.

share|improve this answer
    
Maybe the terminology I used was wrong. What you assumed is exactly what I wanted. Your answer looks promising, and I will have to work with R to process my Shapefile. –  Devdatta Tengshe Oct 18 '12 at 16:45
add comment

You can use the Densify tool. For this case, you choose to densify by angle, Next, choose the the maximum angle accepted in a straight line. Then apply to result line to the tool Splite line at vertices. Finally, delete the lines having shape_length smaller the minimum road length.

enter image description here

In this picture, we see three steps:

1- Densify line using angle. I have used 10 degrees as the parameter, and we used splitline. In the picture, the curved line is in its initial phase.

arcpy.Densify_edit("line" , "ANGLE" , "","",10)
arcpy.SplitLine_management("line" , "line_split")

2- Select the segments where having shape_length not redundant. As we see in the table, I have not selected those redundant lengths. Then, I select them into a new feature class.

arcpy.Select_analysis("line_split" , "line_split_selected")

3- We have extracted the vertices located in the edges of the lines, which are inflection points.

arcpy.FeatureVerticesToPoints_management("line_split_selected" , "line_split_pnt" , "DANGLE")
share|improve this answer
    
I have the same comments and questions concerning your other answer: it's a nice idea, but at the same time it's not clear that it will produce the desired result, nor how one should choose the threshold angle. Could you provide an illustration of the output so readers can evaluate what this proposal really does? Providing worked examples is especially important when recommending ESRI software as part of a solution, because their algorithms usually are not documented, making it impossible to know exactly what they are doing. –  whuber Oct 18 '12 at 17:15
    
to be very sure that is a working solution , i need to test it, but i can't test it, i'm missing the data,so i suppose that the proposed tools by ESRI gonna work like expected, but this answers need to be tested further. –  geogeek Oct 18 '12 at 17:25
    
we could name them ideas and not answers –  geogeek Oct 18 '12 at 17:25
1  
Would you like me to move them into comments, then? BTW, if you want data to test, you could--for a start--use the coordinates I posted in my answer, because they're close to the illustration in the question. But why not just use any geographic data you have handy? –  whuber Oct 18 '12 at 18:33
2  
yes really this solution is working better at extracting just straight lines. –  geogeek Oct 18 '12 at 20:36
show 2 more comments

You can use the Generalize tool which has the maximum offset from the original line as a parameter, so you can choose the offset which fits your case.

enter image description here

If we name the original line "line_cur", and the generalized one "line_gen", then we could clip "line_cur" by "line_gen". The result will be the straight segment of "line_cur". Then we could clean some very short segment by deleting them with an sql query which selects the Shape_length greater that minimum road length.

share|improve this answer
    
This is a nice idea. It's not clear how well it would work in practice, though. Could you perhaps show an example showing the inflection points that are found? –  whuber Oct 18 '12 at 16:46
    
i have made an edit to include a picture, the picture explains how this tool could make a line sticking with straight segments, so we have to make a clip to the old line, to extract just the old straight lines segments –  geogeek Oct 18 '12 at 16:50
    
is anything not clear i'm available to respond your questions ? –  geogeek Oct 18 '12 at 16:51
    
I don't see any inflection points identified in the illustration. Where would they be, exactly? And how should one choose the tolerance for generalization? –  whuber Oct 18 '12 at 16:53
    
i need to some data to perform the test, but i think we should choose the tolerance by experimentation –  geogeek Oct 18 '12 at 16:55
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.