如何更改范围的坐标系以匹配栅格?

问题描述 投票:0回答:1

这是我当前的 R 代码:

library(ebirdst)
library(raster)
library(sf)
library(terra)

# specify parameters
species = "American Robin"
region_extent <- ext(-76.353957, -75.246633, 44.965633, 45.536983)

# get occurrence data
ebirdst_download_status(species, download_abundance = FALSE, download_occurrence = TRUE, force = FALSE)
occurrence_raster <- load_raster(species, product=c("occurrence"))

# filter to county
occurrence_raster <- project(occurrence_raster, crs(region_extent))
cropped_raster <- crop(occurrence_raster, region_extent)

这可行,但是在最后第二步中将栅格投影到范围的坐标系上非常慢。有没有办法可以改变范围的坐标系?当我打印occurrence_raster时,我得到的是:

class       : SpatRaster 
dimensions  : 5630, 13511, 52  (nrow, ncol, nlyr)
resolution  : 2962.807, 2962.809  (x, y)
extent      : -20015109, 20015371, -6673060, 10007555  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=sinu +lon_0=0 +x_0=0 +y_0=0 +R=6371007.181 +units=m +no_defs 
r raster
1个回答
0
投票

您可以:

  1. 使用定义的CRS创建裁剪对象(此处假设为WGS84/EPSG:4326)
  2. 获取SpatRaster的CRS,并用于投影裁剪对象
  3. 使用投影裁剪对象到
    crop()
    SpatRaster
library(ebirdst)
library(raster)
library(sf)
library(terra)

# Specify parameters
species = "American Robin"

# Set crop extent as SpatVector with defined CRS
region_extent <- vect(ext(-76.353957, -75.246633, 44.965633, 45.536983), crs = "EPSG:4326")

# Get occurrence data
ebirdst_download_status(species,
                        download_abundance = FALSE,
                        download_occurrence = TRUE,
                        force = FALSE)

occurrence_raster <- load_raster(species, product = c("occurrence"))

# Get CRS from occurrence_raster
crs <- crs(occurrence_raster)

# Project region to match CRS of occurrence_raster
region_prj <- project(region_extent, crs)

region_prj
# class       : SpatVector 
# geometry    : polygons 
# dimensions  : 1, 0  (geometries, attributes)
# extent      : -6007065, -5860692, 4999956, 5063487  (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=sinu +lon_0=0 +x_0=0 +y_0=0 +R=6371007.181 +units=m +no_defs 

# Crop occurrence_raster to projected region extent
cropped_raster <- crop(occurrence_raster, ext(region_prj))

cropped_raster
# class       : SpatRaster 
# dimensions  : 21, 49, 52  (nrow, ncol, nlyr)
# resolution  : 2962.807, 2962.809  (x, y)
# extent      : -6006959, -5861782, 5000408, 5062626  (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=sinu +lon_0=0 +x_0=0 +y_0=0 +R=6371007.181 +units=m +no_defs 
# source(s)   : memory
# varname     : amerob_occurrence_median_3km_2022 
# names       : 2022-01-04, 2022-01-11, 2022-01-18, 2022-01-25, 2022-02-01, 2022-02-08, ... 
# min values  :  0.0000000,  0.0000000,  0.0000000,  0.0000000,   0.000000,  0.0000000, ... 
# max values  :  0.7491944,  0.7453837,  0.7386807,  0.7669981,   0.759219,  0.7370626, ...
© www.soinside.com 2019 - 2024. All rights reserved.