Page Comparison

To subset the extent of a raster to a desired area, we use the 'crop' function in the 'raster' package. This is the most commonly used function to modify the extent of a raster object. The raster function takes two main arguments:

...

Below examples on how to crop raster objects is provided. These examples were taken from the “Effects of Cyclone Yasi on Green Cover at Mission Beach” tutorial. It can be beneficial to put Put the code snippets in context, by looking at a broader section of the R script. Code snippets have a grey background, and outputs have a white background, could be benefitial. Boxes with grey background contain code snippets, and boxes with white background containt code (text) outputs.

In the example, a SpatialPolygons object for the area of interest (i.e. the Study Area) had been previously created to prepare a map showing the Region & Study Area. The raster is in EPSG:3577(=‘GDA94’/‘Australian Albers’), an Australian-specific CRS (further details for the EPSG:3577 CRS in this link). The polygon for the Study Area, on the other hand, is in EPSG:4326 CRS, a general world spatial reference (more details for this CRS in this link. Therefore, we need to transform the CRS for the Study Area spatial polygon form EPSG:4326 to EPSG:3577. We could then directly use this re-projected polygon as an argument in the ‘crop’ function, as the ‘crop’ function would internally extract the extent for this object. However, for clarity we first create an Extent object from the re-projected Study Area, and use this Extent object as an argument in the ‘crop’ function.

...

#==================================================================================
# Subset Rasters to Study Area
#==================================================================================

#----------------------------------------------------------------------------------
# Create Cookie Cutter for the Study Area
#----------------------------------------------------------------------------------
# Change Study Area spatial polygon CRS (EPSG:4326) to the Australian-specific CRS used for raster (EPSG:3577).

# Find Raster CRS
SPGC.2010q3.rl

.

## class       : RasterLayer 
## dimensions  : 135159, 141481, 19122430479  (nrow, ncol, ncell)
## resolution  : 30, 30  (x, y)
## extent      : -1944645, 2299785, -4910195, -855425  (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=aea +lat_1=-18 +lat_2=-36 +lat_0=0 +lon_0=132 +x_0=0 +y_0=0 +ellps=GRS80 +units=m +no_defs 
## data source : C:/Users/uqbblanc/Documents/TERN/04b-DSDP_GitHub/Prep/Landscapes_AusCover-RemoteSensing/YasiEffectsonGCatMB/YasiEffectsonGCatMB_Tutorial/SPGC_2010q3.vrt 
## names       : SPGC_2010q3 
## values      : 0, 255  (min, max)

...

# Double check that this is the EPSG:3577 CRS (= to the raster's CRS)
CRS("+init=epsg:3577")

.

## CRS arguments:
##  +init=epsg:3577 +proj=aea +lat_1=-18 +lat_2=-36 +lat_0=0
## +lon_0=132 +x_0=0 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0
## +units=m +no_defs

...

# Transform the Spatial Polygon CRS (from EPSG:4326 to EPSG:3577)
StudyArea.SP.reprj = spTransform(StudyArea.SP, CRS("+init=epsg:3577"))

# Study Area Reprojected Extent
StudyArea.extent.reprj = extent(StudyArea.SP.reprj) 
StudyArea.extent.reprj

...

#----------------------------------------------------------------------------------
# Subset (i.e. crop) raster layer to area of interest
#----------------------------------------------------------------------------------
# Subset both Year-Season rasters
SPGC.StudyArea.2010q3.rl = crop(SPGC.2010q3.rl, StudyArea.extent.reprj)
SPGC.StudyArea.2010q3.rl

...