3

I am unable to get centroids from polygons in R. The confusing bit is, it used to work in the past and has suddenly stopped working! This is the piece of code:

polygon.centroids <- SpatialPointsDataFrame(gCentroid(polygon, byid = TRUE), 
                                      polygon@data, match.ID = FALSE)

I can also confirm that the proj4string for the polygon data already exists and is as follows:

proj4string:    +proj=utm +zone=55 +south +ellps=GRS80 +units=m +no_defs

The error message that I am getting upon executing the centroids code is as follows:

   Error in TopologyFunc(spgeom, id, byid, "rgeos_getcentroid") : 
  no slot of name "proj4string" for this object of class "sf" 

Any idea what's going on? As mentioned earlier, this code HAS WORKED OKAY in the past.

Thanks

SymbolixAU
  • 25,502
  • 4
  • 67
  • 139
VGu
  • 386
  • 5
  • 23
  • 3
    You appear to be mixing `sf` objects with `sp` objects. what is the output of `class(polygon)` ? Perhaps you're after the `st_centroid(polygon)` function from `library(sf)` ? – SymbolixAU Mar 08 '18 at 03:27
  • @SymbolixAU- Thanks for your response. I checked and it is an sf object and st_centroid seems to work. However, I have further processing to do in which I have to perform an intersect with an sp object. How do I go about converting sp to sf or vice versa whichever is easier? Thanks once again. – VGu Mar 08 '18 at 03:36
  • 1
    use `sf::st_as_sf` to go from `sp` to `sf` – SymbolixAU Mar 08 '18 at 03:46
  • Hi @SymbolixAU- sp to sf does not work for me unfortunately. It returns an error: invalid crs: +init=epsg:28355 +proj=utm +zone=55 +south +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs, reason: no system list, errno: 22. This occurs even if I assign the proj4string. – VGu Mar 08 '18 at 04:23
  • Could be something to do with a broken `rgdal` installation. I'm just guessing though – SymbolixAU Mar 08 '18 at 04:57

1 Answers1

5

Use the as() function from the methods package to convert your sf object to a SpatialPolygonsDataFrame object - as mentioned in Reading, Writing and Converting Simple Features.

Here, I use polygons from the City of Chicago's 77 community areas to reproduce this suggestion.

Update

Rather than transforming your entire sf object for the sole purpose of calculating centroids, you can do the conversion inside of gCentroid().

# load neccessary packages
library( sf )
library( rgeos )

# load data
comarea606 <-
    read_sf( "https://data.cityofchicago.org/api/geospatial/cauq-8yn6?method=export&format=GeoJSON" )

# find the centroid of each polygon
gCentroid( spgeom = comarea606 )
# Error in TopologyFunc(spgeom, id, byid, "rgeos_getcentroid") : 
#     no slot of name "proj4string" for this object of class "sf"

# convert from sf
# to spatial polygon data frame
# to find the centroid of each polygon
centroids <- 
    gCentroid( 
        spgeom = methods::as( object = comarea606, Class = "Spatial" )
        , byid = TRUE 
        )

# view results
centroids
# SpatialPoints:
#     x        y
# 1  -87.61868 41.83512
# 2  -87.60322 41.82375
# 3  -87.63242 41.80909
# 4  -87.61786 41.81295
# 5  -87.59618 41.80892
# 6  -87.68752 41.97517
# 7  -87.61793 41.79236
# 8  -87.59231 41.79409
# 9  -87.59493 41.77888
# 10 -87.67017 42.00962
# 11 -87.77117 41.97883
# 12 -87.75836 41.99393
# 13 -87.72358 41.98364
# 14 -87.72156 41.96807
# 15 -87.76340 41.95403
# 16 -87.72345 41.95358
# 17 -87.80602 41.94651
# 18 -87.79803 41.92930
# 19 -87.76550 41.92726
# 20 -87.69501 42.00157
# 21 -87.73474 41.92435
# 22 -87.71121 41.93867
# 23 -87.69916 41.92276
# 24 -87.72092 41.90007
# 25 -87.67636 41.90121
# 26 -87.76311 41.89410
# 27 -87.73023 41.87859
# 28 -87.70590 41.87891
# 29 -87.66352 41.87401
# 30 -87.71722 41.86019
# 31 -87.65588 41.96581
# 32 -87.71400 41.83909
# 33 -87.66757 41.85027
# 34 -87.62033 41.85718
# 35 -87.63397 41.84208
# 36 -87.80345 41.98524
# 37 -87.63331 41.89960
# 38 -87.62519 41.87887
# 39 -87.57278 41.76158
# 40 -87.61597 41.74021
# 41 -87.58635 41.74420
# 42 -87.55143 41.74124
# 43 -87.59648 41.72818
# 44 -87.57272 41.72968
# 45 -87.62337 41.70659
# 46 -87.68384 41.94779
# 47 -87.59825 41.70613
# 48 -87.57006 41.69064
# 49 -87.53490 41.70731
# 50 -87.63574 41.67382
# 51 -87.60285 41.66014
# 52 -87.54658 41.66053
# 53 -87.76424 41.79619
# 54 -87.72636 41.81088
# 55 -87.69886 41.81737
# 56 -87.67250 41.82992
# 57 -87.65600 41.94423
# 58 -87.64879 41.83615
# 59 -87.65917 41.80902
# 60 -87.72421 41.79298
# 61 -87.69644 41.79543
# 62 -87.76851 41.77958
# 63 -87.72693 41.76978
# 64 -87.69567 41.77185
# 65 -87.66660 41.77593
# 66 -87.64250 41.77720
# 67 -87.61613 41.76325
# 68 -87.64949 41.92269
# 69 -87.70837 41.74576
# 70 -87.65631 41.74421
# 71 -87.67508 41.71315
# 72 -87.64890 41.71749
# 73 -87.71319 41.69488
# 74 -87.66905 41.68973
# 75 -87.89370 41.97568
# 76 -87.66342 41.98671
# 77 -87.81378 42.00761
# Coordinate Reference System (CRS) arguments: +proj=longlat
# +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0 

# end of script #
Cristian E. Nuno
  • 2,822
  • 2
  • 19
  • 33
  • Thanks for your response. Convert from sf to sp code does not seem to work. Returns the following error: Error in methods::as(object = union.centroids, class = "Spatial") : unused argument (class = "Spatial"). – VGu Mar 08 '18 at 03:42
  • Try changing it to `methods::as(object = union.centroids, Class = "Spatial") `and see how it works. – Cristian E. Nuno Mar 08 '18 at 03:46
  • Your code is of course correct, so +1. However, I would recommend converting from `sp` to `sf` as the latter has now superseded the former. – SymbolixAU Mar 08 '18 at 03:51
  • 1
    @SymbolixAU thank you for your advice. I will take your word for this as good practice. I'm in the land of `sp` and importing shape files using `rgdal` but am trying to force myself to learn and use `sf` more. – Cristian E. Nuno Mar 08 '18 at 03:57