[R-sig-Geo] fraction of a polygon that is intersected by another polygon

satishr satish.regonda at gmail.com
Thu Dec 22 20:12:59 CET 2011


Hi,
   I need to calculate the fraction of a polygon that is intersected by
another polygon. details of the problem: I have data in two shape files; the
first file corresponds to a rectangular grid with uniform intervals whereas
the second file corresponds to US climate divisions, which are polygons. The
polygons representing climate divisions vary in size, some times these
climate divisions/polygons are so small such that they fit in one
rectangular grid and some times require more than 20 rectangular grids. What
I need is: identify rectangular grids that intersect with a polygon and
calculate the fraction of the intersection between the rectangular grid and
polygon.
    The function 'extract' with 'weights=TRUE' was able to calculate
fractions but it has following issues: (a) the function did not consider a
grid(s) although there is a overlap with a climate division, and (b)
fractional amounts are changed when I changed the rectangular grid size.
Here is my code. Any help is greatly appreciated.
Thanks,
Satish

#          Read two shape files, rasterize rectangular grid and use it with
'extract'

#          First shape file, rectangular grid;

            spatGP=readShapeSpatial("Gaussian_Poylgon_ExtendedCONUS.shp")
            spatGP1=as(spatGP,"SpatialPolygons")
            spatGP2=slot(spatGP1,"polygons")

 #          rasterization

             spatras0 <- raster(ncol=180, nrow=180)
             extent(spatras0)<- extent(spatGP);
             spatras <- rasterize(spatGP, spatras0)

 #          second shapefile, polygons, i.e., climate divisions

             cddata=readShapeSpatial("divisions.shp")
             spcddata=as(cddata,"SpatialPolygons")
             regionscddata=slot(spcddata,"polygons")
 
             ipolygon=256   # say, I am interested in this polygon

             sub.tmp = slot(regionscddata[[ipolygon]],'Polygons')
#          create a polygon
            poly.tmp=cbind(sub.tmp[[1]]@coords[,1],sub.tmp[[1]]@coords[,2])
            polys = SpatialPolygons(list(Polygons(list(Polygon(poly.tmp)),
1)))

#           Calculate fractions
            weights.tmp = extract(spatras, polys, weights=TRUE)

‘weights.tmp’ is a matrix – provides information of cells that are covered
by the polygon, i.e., first column cell positions and second column fraction
of cell.

http://r-sig-geo.2731867.n2.nabble.com/file/n7119685/USClimateDivsions_RectangularGrid.jpeg 


--
View this message in context: http://r-sig-geo.2731867.n2.nabble.com/fraction-of-a-polygon-that-is-intersected-by-another-polygon-tp7119685p7119685.html
Sent from the R-sig-geo mailing list archive at Nabble.com.



More information about the R-sig-Geo mailing list