Sure it helps! Thanks a lot! It works like a charm.

Mathieu

2011/6/9 Roger Bivand <Roger.Bivand@nhh.no>

> On Thu, 9 Jun 2011, Mathieu Rajerison wrote:
>
>  Hi,
>>
>>
>> Thanks a lot for taking my problem into account.
>>
>> I've put a data subset as attached file. It contains two files: roads.shp
>> and cities.shp.
>>
>> routes is the result of an intersection between initial complete roads
>> file
>> and cities. Features in actual roads have been split up, therefore
>> gContains
>> test can directly be executed.
>>
>
> This works, but could be tidied up in a loop:
>
> library(rgeos)
> library(rgdal)
> cities <- readOGR(".", "cities")
> routes <- readOGR(".", "routes")
> candidates <- gIntersects(routes, cities, byid=TRUE)
> i1 <- gIntersection(routes[candidates[1,],], cities[1,])
> i2 <- gIntersection(routes[candidates[2,],], cities[2,])
> i3 <- gIntersection(routes[candidates[3,],], cities[3,])
> ii <- rbind(i1, i2, i3, makeUniqueIDs=TRUE)
> row.names(ii) <- row.names(cities)
> SpatialLinesLengths(ii)
>
> With a loop:
>
> library(rgdal)
> cities <- readOGR(".", "cities")
> routes <- readOGR(".", "routes")
> candidates <- gIntersects(routes, cities, byid=TRUE)
> nc <- nrow(candidates)
> out <- vector(mode="list", length=nc)
> for (i in seq(along=out)) {
>  out[[i]] <- gIntersection(routes[candidates[i,],], cities[i,])
>  row.names(out[[i]]) <- as.character(i)
> }
> ii <- do.call("rbind", out)
> row.names(ii) <- row.names(cities)
> SpatialLinesLengths(ii)
>
> This gives the total lengths by city that you need. Using a candidate
> matrix avoids trying to do intersections on objects that do not intersect.
>
> Hope this helps,
>
> Roger
>
>
>
>> In roads, a column named "contained" indicates whether the road is
>> considered contained or not after my gContains test
>>
>>
>> Mathieu
>>
>> 2011/6/8 Roger Bivand <Roger.Bivand@nhh.no>
>>
>>  On Wed, 8 Jun 2011, Mathieu Rajerison wrote:
>>>
>>>  Hi,
>>>
>>>>
>>>>
>>>> First, thanks to all. It's the first time I use rgeos.
>>>>
>>>> I finally managed to get a layer of roads (routes.bdc object) that are
>>>> the
>>>> result of the intersection between roads and cities (bdc object).
>>>>
>>>> But the problem now is that some intersected roads are not considered
>>>> being
>>>> contained in any of the cities, although they should be, as the result
>>>> of
>>>> an
>>>> intersection. I put two images of them as attached files (city contours
>>>> in
>>>> light orange, roads in blue, and some selected roads in yellow)
>>>>
>>>> http://hpics.li/750b6bf
>>>> http://hpics.li/578ab3d
>>>>
>>>>
>>> What is needed is a reproducible example, so that the classes of the
>>> output
>>> objects become evident. Can you provide such an example with publically
>>> available data, or put a relevant data subset on a link?
>>>
>>> Roger
>>>
>>>
>>>
>>>  Here is the code, for the moment:
>>>> ####
>>>> library(rgeos)
>>>>
>>>> pls <- slot(bdc, "polygons")
>>>> pls1 <- lapply(pls, checkPolygonsHoles)
>>>> slot(bdc, "polygons") <- pls1
>>>>
>>>> routes.bdc = gIntersection(routes, bdc)
>>>>
>>>> routes.bdc = slot(routes.bdc, "lineobj")
>>>>
>>>> bdc.contains = gContains(bdc, routes.bdc, byid=TRUE)
>>>>
>>>> nbdc = apply(bdc.contains, 1, function(x) as.numeric(which(x)))
>>>>
>>>> apply(bdc.contains, 1, function(x) any(x))
>>>> ###
>>>>
>>>> apply(bdc.contains, 1, function(x) any(x))
>>>> mentions many FALSE values which indicate that some roads are not
>>>> contained
>>>> in any city...(?)
>>>>
>>>>
>>>> Any help would be greatly appreciated!
>>>>
>>>> Thanks again
>>>>
>>>> 2011/6/7 Roger Bivand <Roger.Bivand@nhh.no>
>>>>
>>>>  On Tue, 7 Jun 2011, Mathieu Rajerison wrote:
>>>>
>>>>>
>>>>>  Thanks for the answer!
>>>>>
>>>>>
>>>>>> I think I surely have to use gIntersection then gLength for this
>>>>>> purpose
>>>>>>
>>>>>> I've tried gIntersection like this:
>>>>>> routes2 = gIntersection(routes, bdc)
>>>>>>
>>>>>> But I get the following message:
>>>>>> Erreur dans createPolygonsComment(p) :
>>>>>>  rgeos_PolyCreateComment: orphaned hole, cannot find containing
>>>>>> polygon
>>>>>> for
>>>>>> hole at index 1
>>>>>>
>>>>>> Do you have an idea why I get this message and how to correct it?
>>>>>>
>>>>>>
>>>>>>  The polygon object is malformed, and at least one polygon claims only
>>>>> to
>>>>> be
>>>>> an interior ring, which is impossible. You need to look at how you are
>>>>> creating your SpatialPolygons object, and possibly use
>>>>> checkPolygonsHoles()
>>>>> in maptools to repair damaged objects. Typically running:
>>>>>
>>>>> library(maptools)
>>>>> pls <- slot(bdc, "polygons")
>>>>> pls1 <- lapply(pls, checkPolygonsHoles)
>>>>> slot(bdc, "polygons") <- pls1
>>>>>
>>>>> may help. As Colin says, posting a link to the polygons object, and to
>>>>> code
>>>>> used to create it, may be helpful.
>>>>>
>>>>> Roger
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>> 2011/6/7 Roman Lu??trik <roman.lustrik@gmail.com>
>>>>>>
>>>>>>  Package rgeos has all sorts of functions that might come handy.
>>>>>>
>>>>>>
>>>>>>> Cheers,
>>>>>>> Roman
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> On Tue, Jun 7, 2011 at 6:13 PM, Mathieu Rajerison <
>>>>>>> mathieu.rajerison@gmail.com> wrote:
>>>>>>>
>>>>>>>  Hi,
>>>>>>>
>>>>>>>
>>>>>>>>
>>>>>>>> I have a SpatialLines object representing roads and a
>>>>>>>> SpatialPolygons
>>>>>>>> object
>>>>>>>> containing cities.
>>>>>>>>
>>>>>>>> I'd like to know how to overlay a SpatialLines object by
>>>>>>>> aSpatialPolygons
>>>>>>>> object.
>>>>>>>>
>>>>>>>> I'd like to know the total distance of roads that cover each of my
>>>>>>>> cities.
>>>>>>>>
>>>>>>>> Is it possible?
>>>>>>>>
>>>>>>>>
>>>>>>>> Thanks,
>>>>>>>>
>>>>>>>> Mathieu
>>>>>>>>
>>>>>>>>     [[alternative HTML version deleted]]
>>>>>>>>
>>>>>>>> _______________________________________________
>>>>>>>> R-sig-Geo mailing list
>>>>>>>> R-sig-Geo@r-project.org
>>>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>> --
>>>>>>> In God we trust, all others bring data.
>>>>>>>
>>>>>>>
>>>>>>>       [[alternative HTML version deleted]]
>>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>>  --
>>>>>>
>>>>> Roger Bivand
>>>>> Economic Geography Section, Department of Economics, Norwegian School
>>>>> of
>>>>> Economics and Business Administration, Helleveien 30, N-5045 Bergen,
>>>>> Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
>>>>> e-mail: Roger.Bivand@nhh.no
>>>>>
>>>>>
>>>>>
>>>>>        [[alternative HTML version deleted]]
>>>>
>>>> _______________________________________________
>>>> R-sig-Geo mailing list
>>>> R-sig-Geo@r-project.org
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>
>>>>
>>>>  --
>>> Roger Bivand
>>> Economic Geography Section, Department of Economics, Norwegian School of
>>> Economics and Business Administration, Helleveien 30, N-5045 Bergen,
>>> Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
>>> e-mail: Roger.Bivand@nhh.no
>>>
>>>
>>>
>>
> --
> Roger Bivand
> Economic Geography Section, Department of Economics, Norwegian School of
> Economics and Business Administration, Helleveien 30, N-5045 Bergen,
> Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
> e-mail: Roger.Bivand@nhh.no
>
>

	[[alternative HTML version deleted]]

