[R-sig-Geo] Find nearest downstream value of a river network

Jon Olav Skoien jon.skoien at jrc.ec.europa.eu
Wed Apr 13 09:59:16 CEST 2011


Hi,

I have a SpatialLinesDataFrame with predictions on different locations 
of a river network that I would like to plot. However, I have many more 
line segments in the network than I have predictions, so most of the 
data.frame has NA-values. Does anyone know a simple way of finding the 
nearest downstream value for a line segment with an NA-value, so that I 
can get a continuous river network with predictions?

Below is a simple example of what I want to do, based on a shapefile 
("nres") that can be downloaded from:
http://intamap.geo.uu.nl/~jon/sarp/tempaacf/

nres = readOGR(".","nres")
nin = nres
nin$pred[!(nin$OBJECTID %in% c(1015, 1020, 1366, 1369, 4981))] = NA

"nres" is the SpatialLinesDataFrame that I would like to have as result 
after associating all NA-values to the nearest downstream value. "nin" 
is a simplified version of the SpatialLinesDataFrame that I have after 
predicting, and that I would like to use for creating "nres".  "pred" is 
the prediction column in these SLDFs. The following plot shows the 
results (with colors) and thick black lines for the segments where I 
have predictions.

spplot(nres, "pred", col.regions = bpy.colors(),
    panel = function(x,y, ...) {
     panel.polygonsplot(x,y, ...)
     sp.lines(nin[!is.na(nin$pred),], col = "black", lwd = 2)
    })

So far I have got the result in this plot from an iterative procedure 
using the FROMJCT and TOJCT (from and to junction) IDs of the data.frame:

ichange = 1
while (ichange > 0) {
   ichange = 0
   for (i in 1:dim(nin)[1]) {
     if (!is.na(nin$pred[i])) {
       tt = which(nin$TOJCT == nin$FROMJCT[i])
       if (length(tt) > 0 && is.na(nin$pred[tt])) {
         nin$pred[tt] = nin$pred[i]
         ichange = ichange + 1
       }
     }
   }
}

But I think it should be either possible to simplify this (quite slow 
for a large river network), or preferably take advantage of topology of 
SpatialLines objects? Any clues?

Thanks,
Jon



More information about the R-sig-Geo mailing list