User Tools

Site Tools


en:learning:schools:s02:code-examples:sm-ce-05-01

C03-1: Reprojection

Introduction

Projections are a very important topic when working with spatial data (either in R or in any other GIS software). R, in contrast to proper GIS software, is not able to do on the fly projections. Therefore, you must be careful in each case, that all of your data feature the same reference system

library(raster)
## Loading required package: sp
library(rgdal)
## rgdal: version: 1.0-4, (SVN revision 548)
##  Geospatial Data Abstraction Library extensions to R successfully loaded
##  Loaded GDAL runtime: GDAL 1.11.2, released 2015/02/10
##  Path to GDAL shared files: /usr/share/gdal/1.11
##  Loaded PROJ.4 runtime: Rel. 4.8.0, 6 March 2012, [PJ_VERSION: 480]
##  Path to PROJ.4 shared files: (autodetected)
##  Linking to sp version: 1.1-1

Where to find information about the reference system

For all spatial objects, the reference system is shown when the object is printed:

## class       : RasterLayer 
## dimensions  : 354, 360, 127440  (nrow, ncol, ncell)
## resolution  : 89.7, 92.2  (x, y)
## extent      : 102329.4, 134621.4, 12550.29, 45189.09  (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=lcc +lat_1=15 +lat_2=16.66666666666667 +lat_0=15.83333333333333 +lon_0=-24 +x_0=161587.83 +y_0=128511.202 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0 
## data source : /home/hanna/Documents/Lehre/kapVerde/dem_fogo.tif 
## names       : dem_fogo 
## values      : -0.8345093, 2763.871  (min, max)

Define Projections

If, for any reason, the reference system is known, but is not specified in the Spatial object, you need to assign it to the object. You can find the information for your reference system on www.spatialreference.org. Each reference system is defined by a unique EPSG code. You can either assign the correct projection using the proj4string which decribes the reference system or simply use the epsg code. Be careful: Assigning a projection does not change the geometry! The coordinates remain the same, only the information about the coordinates is changing.

  ###use the proj4string
  proj4string(dem) <- CRS("+proj=lcc +lat_1=15 +lat_2=16.66666666666667 +lat_0=15.83333333333333 +lon_0=-24 +x_0=161587.83 +y_0=128511.202 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0")
  ###or use the epsg code:
  proj4string(dem) <- CRS("+init=epsg:4826")

Reproject Spatial Data

If you like to change the reference system of your data use the “projectRaster” command for raster data or the “spTransform” command for vector data.

  #e.g with nearest neighbor method:
  newproj="+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
  dem_proj <- projectRaster(dem, crs=newproj,method="ngb")
  print(dem_proj)
## class       : RasterLayer 
## dimensions  : 365, 371, 135415  (nrow, ncol, ncell)
## resolution  : 0.000834, 0.000833  (x, y)
## extent      : -24.55537, -24.24595, 14.78025, 15.08429  (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs +towgs84=0,0,0 
## data source : in memory
## names       : dem_fogo 
## values      : -0.8345093, 2763.871  (min, max)

Note that the coordinates changed, as well as the information of the reference system.

en/learning/schools/s02/code-examples/sm-ce-05-01.txt · Last modified: 2015/09/22 19:32 (external edit)