Albert Heim
Albert Heim

Reputation: 13

Failing to assign projection to sp object SpatialPointsDataFrame

I have a SpatialPointsDataFrame called johnny, created from a vanilla dataframe by assigning coordinates. These coordinates are in coordinate system EPSG 4326 (the standard GPS geographic coordinate system), but johnny does not know that. So ,I am trying to assign EPSG 4326 to johnny, essentially as in this earlier question data projection in R using package SP . I, too, am using sp. My ultimate goal is to project johnny to projected_johnny. However, I can't seem to assign the existing projection correctly first. Who sees my mistake?

library(sp)
x <- seq(80,90,by=1)
y <- seq(40,50,by=1)
value <- seq(10,20,by=1)
johnny <- data.frame(cbind(x,y,value))
coordinates(johnny) <- ~x+y
class(johnny)

[1] "SpatialPointsDataFrame"
attr(,"package")
[1] "sp"

proj4string(johnny) <- CRS("+init=epsg:4326")

Error in if (is.na(get("has_proj_def.dat", envir = .RGDAL_CACHE))) { : argument is of length zero

I have considered and rejected the following possible solutions after trying them out:

  1. Adding library rdgal directly
  2. using CRS("+proj=longlat +datum=WGS84") instead of CRS("+init=epsg:4326")

I am using R 3.6.0 and sp 1.3-1. The rgdal version loaded via sp is 1.5-15. Any ideas are welcome. This should be such a simple action...

Upvotes: 0

Views: 1116

Answers (1)

Gray
Gray

Reputation: 1388

I looked over your code and guessed what you are probably trying to accomplish. But the way you are going about things is more different than it needs to be. There is a simple way to accomplished this. By far, the easiest way to accomplish this is by using those tools found in the R, sf package. Know that the sf package is a newer package than the sp package. And the sf package provides easy to use tools for accomplishing these tasks.

The code below is somewhat different from your code. A two column matrix was used instead of your three column data frame.

The simple feature geometry points were created from the matrix. Then the simple feature column object was created from the geometry points. Then the plot was created.

Code:

# Create matrix

x <- seq(80,90,by=1)
y <- seq(40,50,by=1)
# value <- seq(10,20,by=1)
#johnny <- data.frame(cbind(x,y))
jm <- matrix(data = c(x,y), nrow = 11, ncol = 2)

# coordinates(johnny) <- ~x+y
# class(johnny)
# johnny

Create sf multipoint geometry:

jm.sfg <- st_multipoint(jm)
jm.sfg

Create sf column object:

jm.sfc <- st_sfc(jm.sfg, crs = 4326)
jm.sfc

Plot

plot(jm.sfc, axes = TRUE)

The plot can be viewed from link below.

enter image description here

Upvotes: 1

Related Questions