如何用R语言raster包从GeoTIFF提取像素RGB分量？

Got it, let's work through how to extract RGB components from the Globcover 2009 categorical raster using the raster package. Since this is a land cover classification dataset, each pixel's value corresponds to a class ID—we can't pull RGB bands directly like a regular image. Instead, we'll map each class ID to its matching RGB color from the legend, then generate separate red, green, and blue raster layers.

Step 1: Load Required Packages & Read the Raster #

First, make sure you have the raster package installed, then load it and read your downloaded Globcover TIFF:

# Install package if needed
# install.packages("raster")

library(raster)

# Read the unzipped Globcover raster file
globcover <- raster("Globcover2009_V2.3_Global_.tif")

Step 2: Define the Class-to-RGB Mapping #

You mentioned you have the RGB legend—below is the standard Globcover 2009 legend mapping (adjust this to match your exact legend if needed):

# Create a data frame linking class IDs to their RGB values (0-255 range)
rgb_mapping <- data.frame(
  class_id = c(11, 14, 20, 30, 40, 50, 60, 70, 90, 100, 110, 120, 130, 140, 150, 160, 170, 180, 190, 200, 210, 220),
  red = c(0, 0, 0, 100, 170, 255, 230, 210, 160, 100, 255, 255, 200, 170, 80, 240, 245, 250, 0, 150, 255, 255),
  green = c(0, 100, 150, 170, 220, 255, 240, 230, 160, 100, 255, 190, 110, 90, 130, 240, 250, 250, 0, 150, 255, 200),
  blue = c(255, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 150, 250, 250, 255, 150, 0, 0)
)

Step 3: Generate RGB Component Rasters #

Use the reclassify() function to replace each class ID with its corresponding R, G, or B value:

# Extract red component
red_raster <- reclassify(globcover, rcl = cbind(rgb_mapping$class_id, rgb_mapping$red))

# Extract green component
green_raster <- reclassify(globcover, rcl = cbind(rgb_mapping$class_id, rgb_mapping$green))

# Extract blue component
blue_raster <- reclassify(globcover, rcl = cbind(rgb_mapping$class_id, rgb_mapping$blue))

Step 4: (Optional) Combine & Visualize the RGB Stack #

If you want to verify the result, combine the three components into an RGB stack and plot it:

# Create an RGB raster stack
rgb_stack <- stack(red_raster, green_raster, blue_raster)

# Plot the RGB land cover map
plotRGB(rgb_stack, r=1, g=2, b=3, scale=255)

Step 5: Save the RGB Components #

Save each individual RGB raster as a GeoTIFF for later use:

# Save red component
writeRaster(red_raster, "globcover_red.tif", format = "GTiff", overwrite = TRUE)

# Save green component
writeRaster(green_raster, "globcover_green.tif", format = "GTiff", overwrite = TRUE)

# Save blue component
writeRaster(blue_raster, "globcover_blue.tif", format = "GTiff", overwrite = TRUE)

Important Notes #

• Match Class IDs: Double-check that every class ID in your raster exists in the rgb_mapping data frame—unmatched IDs will be set to NA in the output rasters.
• Memory Management: The full Globcover raster is large. If you run into memory issues, crop it to your area of interest first using crop(globcover, your_extent) before processing.
• Color Scales: If your legend uses 0-1 normalized color values instead of 0-255, adjust the rgb_mapping values accordingly.

内容的提问来源于stack exchange，提问作者NM_