4.3 Extract information from raster data object
4.3.1 Get CRS
You often need to extract the CRS of a raster object before you interact it with vector data (e.g., extracting values from a raster layer to vector data, or cropping a raster layer to the spatial extent of vector data), which can be done using terra::crs()
:
::crs(IA_cdl_2015_sr) terra
[1] "BOUNDCRS[\n SOURCECRS[\n PROJCRS[\"unnamed\",\n BASEGEOGCRS[\"GRS 1980(IUGG, 1980)\",\n DATUM[\"unknown\",\n ELLIPSOID[\"GRS80\",6378137,298.257222101,\n LENGTHUNIT[\"metre\",1,\n ID[\"EPSG\",9001]]]],\n PRIMEM[\"Greenwich\",0,\n ANGLEUNIT[\"degree\",0.0174532925199433,\n ID[\"EPSG\",9122]]]],\n CONVERSION[\"Albers Equal Area\",\n METHOD[\"Albers Equal Area\",\n ID[\"EPSG\",9822]],\n PARAMETER[\"Latitude of false origin\",23,\n ANGLEUNIT[\"degree\",0.0174532925199433],\n ID[\"EPSG\",8821]],\n PARAMETER[\"Longitude of false origin\",-96,\n ANGLEUNIT[\"degree\",0.0174532925199433],\n ID[\"EPSG\",8822]],\n PARAMETER[\"Latitude of 1st standard parallel\",29.5,\n ANGLEUNIT[\"degree\",0.0174532925199433],\n ID[\"EPSG\",8823]],\n PARAMETER[\"Latitude of 2nd standard parallel\",45.5,\n ANGLEUNIT[\"degree\",0.0174532925199433],\n ID[\"EPSG\",8824]],\n PARAMETER[\"Easting at false origin\",0,\n LENGTHUNIT[\"metre\",1],\n ID[\"EPSG\",8826]],\n PARAMETER[\"Northing at false origin\",0,\n LENGTHUNIT[\"metre\",1],\n ID[\"EPSG\",8827]]],\n CS[Cartesian,2],\n AXIS[\"easting\",east,\n ORDER[1],\n LENGTHUNIT[\"metre\",1,\n ID[\"EPSG\",9001]]],\n AXIS[\"northing\",north,\n ORDER[2],\n LENGTHUNIT[\"metre\",1,\n ID[\"EPSG\",9001]]]]],\n TARGETCRS[\n GEOGCRS[\"WGS 84\",\n DATUM[\"World Geodetic System 1984\",\n ELLIPSOID[\"WGS 84\",6378137,298.257223563,\n LENGTHUNIT[\"metre\",1]]],\n PRIMEM[\"Greenwich\",0,\n ANGLEUNIT[\"degree\",0.0174532925199433]],\n CS[ellipsoidal,2],\n AXIS[\"geodetic latitude (Lat)\",north,\n ORDER[1],\n ANGLEUNIT[\"degree\",0.0174532925199433]],\n AXIS[\"geodetic longitude (Lon)\",east,\n ORDER[2],\n ANGLEUNIT[\"degree\",0.0174532925199433]],\n USAGE[\n SCOPE[\"Horizontal component of 3D system.\"],\n AREA[\"World.\"],\n BBOX[-90,-180,90,180]],\n ID[\"EPSG\",4326]]],\n ABRIDGEDTRANSFORMATION[\"Transformation to WGS84\",\n METHOD[\"Position Vector transformation (geog2D domain)\",\n ID[\"EPSG\",9606]],\n PARAMETER[\"X-axis translation\",0,\n ID[\"EPSG\",8605]],\n PARAMETER[\"Y-axis translation\",0,\n ID[\"EPSG\",8606]],\n PARAMETER[\"Z-axis translation\",0,\n ID[\"EPSG\",8607]],\n PARAMETER[\"X-axis rotation\",0,\n ID[\"EPSG\",8608]],\n PARAMETER[\"Y-axis rotation\",0,\n ID[\"EPSG\",8609]],\n PARAMETER[\"Z-axis rotation\",0,\n ID[\"EPSG\",8610]],\n PARAMETER[\"Scale difference\",1,\n ID[\"EPSG\",8611]]]]"
4.3.2 Subset
You can access specific layers in a multi-layer raster object by indexing:
#--- index ---#
2]] # (originally IA_cdl_2016.tif) IA_cdl_stack_sr[[
class : SpatRaster
dimensions : 11671, 17795, 1 (nrow, ncol, nlyr)
resolution : 30, 30 (x, y)
extent : -52095, 481755, 1938165, 2288295 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=aea +lat_0=23 +lon_0=-96 +lat_1=29.5 +lat_2=45.5 +x_0=0 +y_0=0 +ellps=GRS80 +units=m +no_defs
source : IA_cdl_2016.tif
name : IA_cdl_2016
min value : 0
max value : 241
4.3.3 Get cell values
You can access the values stored in a SpatRaster
object using the terra::values()
function:
#--- terra::values ---#
<- terra::values(IA_cdl_stack_sr)
values_from_rs
#--- take a look ---#
head(values_from_rs)
Layer_1 Layer_1
[1,] 0 0
[2,] 0 0
[3,] 0 0
[4,] 0 0
[5,] 0 0
[6,] 0 0
The returned values come in a matrix form of two columns because we are getting values from a two-layer SpatRaster
object (one column for each layer). In general, terra::values()
returns a \(X\) by \(n\) matrix, where \(X\) is the number of cells and \(n\) is the number of layers.