2

I have been struggling for a few hours trying to transform SpatialDataFrame coordinates to another type. I have two objects from two different sources, I need to plot them together.

> class(Object1)
[1] "RasterLayer"
attr(,"package")
[1] "raster"

> head(coordinates(Object1))
            x        y
[1,] 13.35000 46.99917
[2,] 13.35083 46.99917
[3,] 13.35167 46.99917
[4,] 13.35250 46.99917
[5,] 13.35333 46.99917

> projection(Object1)
[1] "+proj=tmerc +lat_0=0 +lon_0=15 +k=0.9999 +x_0=500000 +y_0=-5000000 +ellps=bessel +units=m +no_defs"
[6,] 13.35417 46.99917

> class(Object2)
[1] "SpatialPolygonsDataFrame"
attr(,"package")
[1] "sp"

> coordinates(Object2)
          [,1]     [,2]
7657  454777.5 119375.8
32166 539358.5 129256.9
49121 617992.2 151211.0

> projection(Object2)
[1] "+proj=tmerc +lat_0=0 +lon_0=15 +k=0.9999 +x_0=500000 +y_0=-5000000 +ellps=bessel +units=m +no_defs"

If I try to plot them together, only the first Object is displayed.

library(tmap)
tm_shape(Object1) +tm_raster() +tm_shape(Object2) + tm_borders()

One might hope that this would work, but only in my dreams:

> coordinates(Object2) <- coordinates(Object2)/10000
Error in `coordinates<-`(`*tmp*`, value = c(45.4777468962329, 53.9358522795272,  : 
  setting coordinates cannot be done on Spatial objects, where they have already been set

I know that I should use spTransform(), but simply I don't know how to correctly reference the transformation.

4

You are working with different Coordinate Reference System (CRS), so you cannot expect to show both layers at the same time without reprojecting one of them. You have a raster layer and a vector layer, in this case, always is recommend to reproject vector layer (reprojecting a raster layer distort the information).

To use spTransform(), you need to set the desired output projection:

projection(Object1)
[1] "+proj=longlat +ellps=WGS84 +no_defs"

projection(Object2)
[1] "+proj=utm +zone=19 +south +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"

# use crs() to extract CRS from raster
Object2reprojected <- spTransform(Object2,crs(Object1))

library(tmap)
tm_shape(Object1) +tm_raster() +tm_shape(Object2reprojected) + tm_borders()

enter image description here

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.