Nyamisi Peter
  • Blog
  • Resources
    • R Weekly
    • R Bloggers
    • Semba-blog
    • Nyamisi ShinyApps
    • Wior WebApp
  • Archive

On this page

  • Introduction
  • Data
  • Defining the extent
  • Cropping
  • Plotting the cropped area
    • Consultated references

Cropping raster image using terra and tidyterra packages in R

terra
visualization
code
Author

Nyamisi Peter

Published

April 15, 2024

Introduction

Raster data is a type of spatial data that is stored as a grid of squares or rectangles called pixels. Each pixel contains a value that represents a specific location on the Earth’s surface (Hijmans, 2024). Terra package provides methods to manipulate spatial data in both raster and vector format (Hijmans, 2024). Terra package supports SpatRaster data – that can handle large raster files and SpatVector data – that supports all types of geometric operations such as intersections in vector files. Tidyterra is a package that uses common methods from the tidyverse for SpatRaster and SpatVector objects created with the terra package (Hernangómez, 2023).

First, we will load packages that we are going to use in this post:

require(tidyterra)
require(terra)
require(leaflet)

Data

We will use satellite chlorophyll-a data measured by the MODIS sensor in Network Common Data Form (NetCDF) with file extension .nc for Tanzania coastal waters. The data will be loaded using rast() function of terra package (Hijmans, 2024).

chl = rast("M2012214-2012244.vw.all_products.MYO.median.01aug12-31aug12.1368530332.data.nc")

chl
class       : SpatRaster 
dimensions  : 2224, 1112, 3  (nrow, ncol, nlyr)
resolution  : 0.004496807, 0.004496807  (x, y)
extent      : 37.99775, 42.9982, -12.99775, -2.996853  (xmin, xmax, ymin, ymax)
coord. ref. : lon/lat WGS 84 
sources     : M2012214-2012244.vw.all_products.MYO.median.01aug12-31aug12.1368530332.data.nc:l2_flags  
              M2012214-2012244.vw.all_products.MYO.median.01aug12-31aug12.1368530332.data.nc:chl_oc3  
              M2012214-2012244.vw.all_products.MYO.median.01aug12-31aug12.1368530332.data.nc:land_mask  
varnames    : l2_flags 
              chl_oc3 (Chlorophyll-a concentration in sea water as measured by the MODIS sensor) 
              land_mask 
names       : l2_flags,       chl_oc3, land_mask 
unit        :         , milligram m-3,           

Our chl data is in SpatRaster with 2224 number of rows, 1112 number of columns and 3 number of layers. Layers include l2_flags, chl_oc3 and land_mask. The spatial resolution of our data is 0.004496807 pixels which is equivalent to 0.4946488 km approximately 0.5 km. We will select (using select() of tidyterra package) and plot (using plot() of terra package) the second layer (chl_oc3) which contain chlorophyll-a values.

chl |> 
  select(chl_oc3) |> 
  filter(chl_oc3 < 1) |>
  terra::plot()

Defining the extent

Suppose our analysis is focus only to a portion of Tanzania coastal area and not the whole area. So we need to crop that area on interest and leave the others. The focus of this post is around the Pemba Channel within longitude 38.9 and 40 East, and latitude -5.7, -4.6 South of the Equator. The first step is to specify the extent of our area using the ext() function of terra package (Hijmans, 2024):

site.extent = ext( 38.9, 40, -5.7, -4.6)

Cropping

Then, we will crop our extent area from chl data using the crop() function of terra package (Hijmans, 2024).

site.crop = chl |> 
  crop(site.extent)

Plotting the cropped area

Lastly we will plot our cropped area and visualize them;

site.crop |> 
  select(chl_oc3) |> 
  filter(chl_oc3 < 1) |> 
  terra::plot()

We may also visualize the cropped area using leaflet package (Cheng et al., 2024);

pal = colorNumeric(c("#7f007f", "#0000ff",  "#007fff", "#00ffff", "#00bf00", "#7fdf00", "#ffff00", "#ff7f00", "#ff3f00", "#ff0000", "#bf0000"),  na.color = "transparent", values(site.crop |> select(chl_oc3) |> filter(chl_oc3 < 1)))              

leaflet::leaflet() |> 
  leaflet::addTiles() |> 
  leaflet::addRasterImage(site.crop |> select(chl_oc3) |> filter(chl_oc3 < 1),
                          colors = pal,
                          opacity = 1) |> 
  leaflet::addLegend(pal = pal,
                     values = values(site.crop |> select(chl_oc3) |> filter(chl_oc3 < 1)),
                     title = "Chlorophyll")

Thank you for following and do not miss the next post!!

Consultated references

Cheng, J., Schloerke, B., Karambelkar, B., Xie, Y., 2024. Leaflet: Create interactive web maps with the JavaScript ’leaflet’ library.
Hernangómez, D., 2023. Using the tidyverse with terra objects: The tidyterra package. Journal of Open Source Software 8, 5751. https://doi.org/10.21105/joss.05751
Hijmans, R.J., 2024. Terra: Spatial data analysis.
 
Cookie Preferences