Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

R_Extracting coordinates from SpatialPolygonsDataFrame

Tags:

r

slot

geospatial

Is it only me who have the problem with extracting coordinates of a polygon from SpatialPolygonsDataFrame object? I am able to extract other slots of the object (ID,plotOrder) but not coordinates (coords). I don't know what I am doing wrong. Please find below my R session where bdryData being the SpatialPolygonsDataFrame object with two polygons.

> bdryData
An object of class "SpatialPolygonsDataFrame"
Slot "data":
  ID GRIDCODE
0  1        0
1  2        0

Slot "polygons":
[[1]]
An object of class "Polygons"
Slot "Polygons":
[[1]]
An object of class "Polygon"
Slot "labpt":
[1] 415499.1 432781.7

Slot "area":
[1] 0.6846572

Slot "hole":
[1] FALSE

Slot "ringDir":
[1] 1

Slot "coords":
         [,1]     [,2]
[1,] 415499.6 432781.2
[2,] 415498.4 432781.5
[3,] 415499.3 432782.4
[4,] 415499.6 432781.2



Slot "plotOrder":
[1] 1

Slot "labpt":
[1] 415499.1 432781.7

Slot "ID":
[1] "0"

Slot "area":
[1] 0.6846572


[[2]]
An object of class "Polygons"
Slot "Polygons":
[[1]]
An object of class "Polygon"
Slot "labpt":
[1] 415587.3 432779.4

Slot "area":
[1] 20712.98

Slot "hole":
[1] FALSE

Slot "ringDir":
[1] 1

Slot "coords":
           [,1]     [,2]
  [1,] 415499.6 432781.2
  [2,] 415505.0 432781.8
  [3,] 415506.5 432792.6
  [4,] 415508.9 432792.8
  [5,] 415515.0 432791.5
  [6,] 415517.7 432795.6
  [7,] 415528.6 432797.7
  [8,] 415538.8 432804.2
  [9,] 415543.2 432805.8
 [10,] 415545.1 432803.6
 [11,] 415547.1 432804.7
 [12,] 415551.7 432805.8
 [13,] 415557.5 432812.3
 [14,] 415564.2 432817.1
 [15,] 415568.5 432823.9
 [16,] 415571.0 432826.8
 [17,] 415573.2 432828.7
 [18,] 415574.1 432829.7
 [19,] 415576.2 432830.7
 [20,] 415580.2 432833.8
 [21,] 415589.6 432836.0
 [22,] 415593.1 432841.0
 [23,] 415592.2 432843.7
 [24,] 415590.6 432846.6
 [25,] 415589.0 432853.3
 [26,] 415584.8 432855.3
 [27,] 415579.7 432859.8
 [28,] 415577.7 432866.2
 [29,] 415575.6 432868.1
 [30,] 415566.7 432880.7
 [31,] 415562.7 432887.5
 [32,] 415559.2 432889.1
 [33,] 415561.5 432890.7
 [34,] 415586.2 432889.7
 [35,] 415587.1 432888.6
 [36,] 415588.5 432890.2
 [37,] 415598.2 432888.7
 [38,] 415599.1 432887.7
 [39,] 415601.2 432886.7
 [40,] 415603.1 432885.7
 [41,] 415605.2 432884.7
 [42,] 415606.1 432882.7
 [43,] 415607.2 432880.7
 [44,] 415608.3 432878.3
 [45,] 415612.2 432874.8
 [46,] 415614.7 432871.9
 [47,] 415617.1 432870.7
 [48,] 415622.4 432868.2
 [49,] 415622.0 432862.4
 [50,] 415624.2 432855.4
 [51,] 415633.2 432845.3
 [52,] 415639.0 432841.1
 [53,] 415642.8 432832.9
 [54,] 415647.5 432828.7
 [55,] 415654.3 432820.3
 [56,] 415654.1 432816.5
 [57,] 415658.2 432812.8
 [58,] 415661.9 432808.6
 [59,] 415663.5 432808.7
 [60,] 415668.1 432803.5
 [61,] 415676.5 432801.3
 [62,] 415679.1 432802.7
 [63,] 415680.1 432802.7
 [64,] 415681.1 432802.7
 [65,] 415682.2 432802.7
 [66,] 415685.8 432804.7
 [67,] 415691.8 432802.2
 [68,] 415693.6 432798.9
 [69,] 415696.2 432777.0
 [70,] 415689.8 432773.5
 [71,] 415683.7 432771.6
 [72,] 415680.2 432766.7
 [73,] 415679.0 432765.6
 [74,] 415676.8 432753.7
 [75,] 415671.4 432747.7
 [76,] 415662.7 432747.2
 [77,] 415658.7 432750.0
 [78,] 415657.0 432746.3
 [79,] 415654.1 432743.7
 [80,] 415652.3 432739.8
 [81,] 415649.6 432739.6
 [82,] 415648.0 432739.7
 [83,] 415641.9 432736.4
 [84,] 415633.4 432736.9
 [85,] 415630.2 432734.7
 [86,] 415622.3 432733.6
 [87,] 415614.4 432726.5
 [88,] 415617.1 432719.1
 [89,] 415612.5 432718.1
 [90,] 415610.0 432720.9
 [91,] 415606.2 432716.6
 [92,] 415603.2 432713.9
 [93,] 415601.4 432710.0
 [94,] 415580.3 432708.7
 [95,] 415545.1 432709.7
 [96,] 415543.5 432711.5
 [97,] 415534.0 432715.7
 [98,] 415527.1 432713.7
 [99,] 415521.1 432711.6
[100,] 415505.6 432710.6
[101,] 415501.3 432710.9
[102,] 415499.3 432708.7
[103,] 415495.6 432711.6
[104,] 415482.6 432726.2
[105,] 415477.2 432734.0
[106,] 415478.1 432737.7
[107,] 415479.2 432739.7
[108,] 415480.9 432743.4
[109,] 415486.5 432751.2
[110,] 415493.2 432760.7
[111,] 415494.1 432762.7
[112,] 415498.1 432767.9
[113,] 415497.2 432770.7
[114,] 415490.6 432773.2
[115,] 415493.2 432775.6
[116,] 415496.0 432778.7
[117,] 415499.2 432779.7
[118,] 415499.6 432781.2



Slot "plotOrder":
[1] 1

Slot "labpt":
[1] 415587.3 432779.4

Slot "ID":
[1] "1"

Slot "area":
[1] 20712.98



Slot "plotOrder":
[1] 2 1

Slot "bbox":
       min      max
x 415477.2 415696.2
y 432708.7 432890.7

Slot "proj4string":
CRS arguments:
 +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000
+datum=OSGB36 +units=m +no_defs +ellps=airy
+towgs84=446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894 

Subsetting second polygon from bdryData

> bdryData@polygons[[2]]
An object of class "Polygons"
Slot "Polygons":
[[1]]
An object of class "Polygon"
Slot "labpt":
[1] 415587.3 432779.4

Slot "area":
[1] 20712.98

Slot "hole":
[1] FALSE

Slot "ringDir":
[1] 1

Slot "coords":
           [,1]     [,2]
  [1,] 415499.6 432781.2
  [2,] 415505.0 432781.8
  [3,] 415506.5 432792.6
  [4,] 415508.9 432792.8
  [5,] 415515.0 432791.5
  [6,] 415517.7 432795.6
  [7,] 415528.6 432797.7
  [8,] 415538.8 432804.2
  [9,] 415543.2 432805.8
 [10,] 415545.1 432803.6
 [11,] 415547.1 432804.7
 [12,] 415551.7 432805.8
 [13,] 415557.5 432812.3
 [14,] 415564.2 432817.1
 [15,] 415568.5 432823.9
 [16,] 415571.0 432826.8
 [17,] 415573.2 432828.7
 [18,] 415574.1 432829.7
 [19,] 415576.2 432830.7
 [20,] 415580.2 432833.8
 [21,] 415589.6 432836.0
 [22,] 415593.1 432841.0
 [23,] 415592.2 432843.7
 [24,] 415590.6 432846.6
 [25,] 415589.0 432853.3
 [26,] 415584.8 432855.3
 [27,] 415579.7 432859.8
 [28,] 415577.7 432866.2
 [29,] 415575.6 432868.1
 [30,] 415566.7 432880.7
 [31,] 415562.7 432887.5
 [32,] 415559.2 432889.1
 [33,] 415561.5 432890.7
 [34,] 415586.2 432889.7
 [35,] 415587.1 432888.6
 [36,] 415588.5 432890.2
 [37,] 415598.2 432888.7
 [38,] 415599.1 432887.7
 [39,] 415601.2 432886.7
 [40,] 415603.1 432885.7
 [41,] 415605.2 432884.7
 [42,] 415606.1 432882.7
 [43,] 415607.2 432880.7
 [44,] 415608.3 432878.3
 [45,] 415612.2 432874.8
 [46,] 415614.7 432871.9
 [47,] 415617.1 432870.7
 [48,] 415622.4 432868.2
 [49,] 415622.0 432862.4
 [50,] 415624.2 432855.4
 [51,] 415633.2 432845.3
 [52,] 415639.0 432841.1
 [53,] 415642.8 432832.9
 [54,] 415647.5 432828.7
 [55,] 415654.3 432820.3
 [56,] 415654.1 432816.5
 [57,] 415658.2 432812.8
 [58,] 415661.9 432808.6
 [59,] 415663.5 432808.7
 [60,] 415668.1 432803.5
 [61,] 415676.5 432801.3
 [62,] 415679.1 432802.7
 [63,] 415680.1 432802.7
 [64,] 415681.1 432802.7
 [65,] 415682.2 432802.7
 [66,] 415685.8 432804.7
 [67,] 415691.8 432802.2
 [68,] 415693.6 432798.9
 [69,] 415696.2 432777.0
 [70,] 415689.8 432773.5
 [71,] 415683.7 432771.6
 [72,] 415680.2 432766.7
 [73,] 415679.0 432765.6
 [74,] 415676.8 432753.7
 [75,] 415671.4 432747.7
 [76,] 415662.7 432747.2
 [77,] 415658.7 432750.0
 [78,] 415657.0 432746.3
 [79,] 415654.1 432743.7
 [80,] 415652.3 432739.8
 [81,] 415649.6 432739.6
 [82,] 415648.0 432739.7
 [83,] 415641.9 432736.4
 [84,] 415633.4 432736.9
 [85,] 415630.2 432734.7
 [86,] 415622.3 432733.6
 [87,] 415614.4 432726.5
 [88,] 415617.1 432719.1
 [89,] 415612.5 432718.1
 [90,] 415610.0 432720.9
 [91,] 415606.2 432716.6
 [92,] 415603.2 432713.9
 [93,] 415601.4 432710.0
 [94,] 415580.3 432708.7
 [95,] 415545.1 432709.7
 [96,] 415543.5 432711.5
 [97,] 415534.0 432715.7
 [98,] 415527.1 432713.7
 [99,] 415521.1 432711.6
[100,] 415505.6 432710.6
[101,] 415501.3 432710.9
[102,] 415499.3 432708.7
[103,] 415495.6 432711.6
[104,] 415482.6 432726.2
[105,] 415477.2 432734.0
[106,] 415478.1 432737.7
[107,] 415479.2 432739.7
[108,] 415480.9 432743.4
[109,] 415486.5 432751.2
[110,] 415493.2 432760.7
[111,] 415494.1 432762.7
[112,] 415498.1 432767.9
[113,] 415497.2 432770.7
[114,] 415490.6 432773.2
[115,] 415493.2 432775.6
[116,] 415496.0 432778.7
[117,] 415499.2 432779.7
[118,] 415499.6 432781.2



Slot "plotOrder":
[1] 1

Slot "labpt":
[1] 415587.3 432779.4

Slot "ID":
[1] "1"

Slot "area":
[1] 20712.98

Extracting slots

> bdryData@polygons[[2]]@ID 
[1] "1"

> bdryData@polygons[[2]]@plotOrder
[1] 1

But problem with coordinates

> bdryData@polygons[[2]]@coords
Error: no slot of name "coords" for this object of class "Polygons"

Any help is really appreciated. Thanks.

like image 220
rm167 Avatar asked Apr 22 '15 16:04

rm167


4 Answers

This question was also addressed on gis.stackexchange, here. I made an example below testing all the options mentioned here by @mdsumner. Also have a look here

library(sp)
library(sf)
#> Warning: package 'sf' was built under R version 3.5.3
#> Linking to GEOS 3.6.1, GDAL 2.2.3, PROJ 4.9.3
library(raster)
library(spbabel)
#> Warning: package 'spbabel' was built under R version 3.5.3
library(tmap)

library(microbenchmark)
library(ggplot2)

# Prepare data
data(World)
# Convert from sf to sp objects
atf_sf <- World[World$iso_a3 == "ATF", ]
atf_sp <- as(atf_sf, "Spatial")
atf_sp
#> class       : SpatialPolygonsDataFrame 
#> features    : 1 
#> extent      : 5490427, 5660887, -6048972, -5932855  (xmin, xmax, ymin, ymax)
#> coord. ref. : +proj=eck4 +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0 
#> variables   : 15
#> # A tibble: 1 x 15
#>   iso_a3 name  sovereignt continent area  pop_est pop_est_dens economy
#>   <fct>  <fct> <fct>      <fct>     <S3:>   <dbl>        <dbl> <fct>  
#> 1 ATF    Fr. ~ France     Seven se~ 7257~     140       0.0193 6. Dev~
#> # ... with 7 more variables: income_grp <fct>, gdp_cap_est <dbl>,
#> #   life_exp <dbl>, well_being <dbl>, footprint <dbl>, inequality <dbl>,
#> #   HPI <dbl>

# Try various functions:

raster::geom(atf_sp)
#>       object part cump hole       x        y
#>  [1,]      1    1    1    0 5550200 -5932855
#>  [2,]      1    1    1    0 5589907 -5964836
#>  [3,]      1    1    1    0 5660887 -5977490
#>  [4,]      1    1    1    0 5656160 -5996685
#>  [5,]      1    1    1    0 5615621 -6042456
#>  [6,]      1    1    1    0 5490427 -6048972
#>  [7,]      1    1    1    0 5509148 -5995424
#>  [8,]      1    1    1    0 5536900 -5953683
#>  [9,]      1    1    1    0 5550200 -5932855

ggplot2::fortify(atf_sp)
#> Regions defined for each Polygons
#>      long      lat order  hole piece id group
#> 1 5550200 -5932855     1 FALSE     1  8   8.1
#> 2 5589907 -5964836     2 FALSE     1  8   8.1
#> 3 5660887 -5977490     3 FALSE     1  8   8.1
#> 4 5656160 -5996685     4 FALSE     1  8   8.1
#> 5 5615621 -6042456     5 FALSE     1  8   8.1
#> 6 5490427 -6048972     6 FALSE     1  8   8.1
#> 7 5509148 -5995424     7 FALSE     1  8   8.1
#> 8 5536900 -5953683     8 FALSE     1  8   8.1
#> 9 5550200 -5932855     9 FALSE     1  8   8.1

spbabel::sptable(atf_sp)
#> # A tibble: 9 x 6
#>   object_ branch_ island_ order_       x_        y_
#>     <int>   <int> <lgl>    <int>    <dbl>     <dbl>
#> 1       1       1 TRUE         1 5550200. -5932855.
#> 2       1       1 TRUE         2 5589907. -5964836.
#> 3       1       1 TRUE         3 5660887. -5977490.
#> 4       1       1 TRUE         4 5656160. -5996685.
#> 5       1       1 TRUE         5 5615621. -6042456.
#> 6       1       1 TRUE         6 5490427. -6048972.
#> 7       1       1 TRUE         7 5509148. -5995424.
#> 8       1       1 TRUE         8 5536900. -5953683.
#> 9       1       1 TRUE         9 5550200. -5932855.

as.data.frame(as(as(atf_sp, "SpatialLinesDataFrame"),"SpatialPointsDataFrame"))
#>     iso_a3                   name sovereignt               continent
#> 8      ATF Fr. S. Antarctic Lands     France Seven seas (open ocean)
#> 8.1    ATF Fr. S. Antarctic Lands     France Seven seas (open ocean)
#> 8.2    ATF Fr. S. Antarctic Lands     France Seven seas (open ocean)
#> 8.3    ATF Fr. S. Antarctic Lands     France Seven seas (open ocean)
#> 8.4    ATF Fr. S. Antarctic Lands     France Seven seas (open ocean)
#> 8.5    ATF Fr. S. Antarctic Lands     France Seven seas (open ocean)
#> 8.6    ATF Fr. S. Antarctic Lands     France Seven seas (open ocean)
#> 8.7    ATF Fr. S. Antarctic Lands     France Seven seas (open ocean)
#> 8.8    ATF Fr. S. Antarctic Lands     France Seven seas (open ocean)
#>                area pop_est pop_est_dens              economy
#> 8   7257.455 [km^2]     140   0.01929051 6. Developing region
#> 8.1 7257.455 [km^2]     140   0.01929051 6. Developing region
#> 8.2 7257.455 [km^2]     140   0.01929051 6. Developing region
#> 8.3 7257.455 [km^2]     140   0.01929051 6. Developing region
#> 8.4 7257.455 [km^2]     140   0.01929051 6. Developing region
#> 8.5 7257.455 [km^2]     140   0.01929051 6. Developing region
#> 8.6 7257.455 [km^2]     140   0.01929051 6. Developing region
#> 8.7 7257.455 [km^2]     140   0.01929051 6. Developing region
#> 8.8 7257.455 [km^2]     140   0.01929051 6. Developing region
#>                  income_grp gdp_cap_est life_exp well_being footprint
#> 8   2. High income: nonOECD    114285.7       NA         NA        NA
#> 8.1 2. High income: nonOECD    114285.7       NA         NA        NA
#> 8.2 2. High income: nonOECD    114285.7       NA         NA        NA
#> 8.3 2. High income: nonOECD    114285.7       NA         NA        NA
#> 8.4 2. High income: nonOECD    114285.7       NA         NA        NA
#> 8.5 2. High income: nonOECD    114285.7       NA         NA        NA
#> 8.6 2. High income: nonOECD    114285.7       NA         NA        NA
#> 8.7 2. High income: nonOECD    114285.7       NA         NA        NA
#> 8.8 2. High income: nonOECD    114285.7       NA         NA        NA
#>     inequality HPI Lines.NR Lines.ID Line.NR coords.x1 coords.x2
#> 8           NA  NA        1        8       1   5550200  -5932855
#> 8.1         NA  NA        1        8       1   5589907  -5964836
#> 8.2         NA  NA        1        8       1   5660887  -5977490
#> 8.3         NA  NA        1        8       1   5656160  -5996685
#> 8.4         NA  NA        1        8       1   5615621  -6042456
#> 8.5         NA  NA        1        8       1   5490427  -6048972
#> 8.6         NA  NA        1        8       1   5509148  -5995424
#> 8.7         NA  NA        1        8       1   5536900  -5953683
#> 8.8         NA  NA        1        8       1   5550200  -5932855

# What about speed? raster::geom is the fastest
res <- microbenchmark(raster::geom(atf_sp),
                      ggplot2::fortify(atf_sp),
                      spbabel::sptable(atf_sp),
                      as.data.frame(as(as(atf_sp, "SpatialLinesDataFrame"),
                                       "SpatialPointsDataFrame")))
ggplot2::autoplot(res)
#> Coordinate system already present. Adding new coordinate system, which will replace the existing one.

Created on 2019-03-23 by the reprex package (v0.2.1)

like image 110
Valentin Avatar answered Nov 16 '22 14:11

Valentin


Finally, I figured out that I didn't parse the output correctly. The correct way to do is bdryData@polygons[[2]]@Polygons[[1]]@coords. Mind the difference in command polygons(Polygons and polygons) and it took me ages to find out.

like image 45
rm167 Avatar answered Nov 16 '22 13:11

rm167


Use the coordinates() function from the sp package. It should give you the values in a list format.

You can also get the Polygon attribute from the shapefile.

mfile = readOGR(dsn=dsn,layer=layername)
polys = attr(mfile,'polygons')
npolys = length(polys)
for (i in 1:npolys){
  poly = polys[[i]]
  polys2 = attr(poly,'Polygons')
  npolys2 = length(polys2)
  for (j in 1:npolys2){
     #do stuff with these values
     coords = coordinates(polys2[[j]])
  }
}
like image 20
Max Candocia Avatar answered Nov 16 '22 14:11

Max Candocia


This took me a while to figure out too. The following function I wrote worked for me. sp.df should be SpatialPolygonsDataFrame.

extractCoords <- function(sp.df)
{
    results <- list()
    for(i in 1:length(sp.df@polygons[[1]]@Polygons))
    {
        results[[i]] <- sp.df@polygons[[1]]@Polygons[[i]]@coords
    }
    results <- Reduce(rbind, results)
    results
}
like image 5
forlooper Avatar answered Nov 16 '22 14:11

forlooper