如何在没有投影经纬度坐标的情况下从栅格中提取 xy 坐标的数据点

Question

您好，我正在尝试从样本栅格堆栈 ras_dt 中提取 xy 点的值。 ras_dt 是 EQUATES data，网格坐标在 -115.00,38.00,-110.05,45.00 的范围内。如何更改此栅格堆栈的投影和 lat lon 坐标，以便我可以提取点 xy 的数据，如下面的代码所示。

library(raster)
dturl<- "https://www.dropbox.com/s/ztxqpszjfjhpavz/EQUATES_ACONC_O3_SAM.nc?dl=1"
download.file(dturl, "EQUATES_ACONC_O3_SAM.nc")
ras_dt <- raster::stack("EQUATES_ACONC_O3_SAM.nc",varname = "O3") 
ras_dt

# the data domain is -115.00,38.00,-110.05,45.00
plot(ras_dt)

xy <- data.frame(lon=-113.0,lat=40.0)
coordinates(xy) <-  ~lon + lat
extr_dt <- raster::extract(ras_dt, xy) # how to get O3 values for xy here?
extr_dt

ras_dt 在“lambertConformalProjection”中，信息如下：

char LambertConformalProjection;
  :grid_mapping_name = "lambert_conformal_conic";
  :latitude_of_projection_origin = 40.0; // double
  :longitude_of_central_meridian = -97.0; // double
  :standard_parallel = 33.0, 45.0; // double
  :earth_radius = 6370000.0; // double
  :_CoordinateTransformType = "Projection";
  :_CoordinateAxes = "x y";

Answer 1

您需要设置范围和坐标参考系统(CRS)。然后将您的 lon/lat 点转换为该 CRS 并使用提取。根据您编辑的问题，我们现在有了 CRS。

你指定了“域”，但你需要在实际的CRS中进行坐标。我不知道这些是什么，所以我会猜测它，但它会不正确。

数据

library(terra)
dturl<- "https://www.dropbox.com/s/ztxqpszjfjhpavz/EQUATES_ACONC_O3_SAM.nc?dl=1"
download.file(dturl, "EQUATES_ACONC_O3_SAM.nc", mode="wb")
ras_dt <- rast("EQUATES_ACONC_O3_SAM.nc", "O3") 
pcrs <- "+proj=lcc +lon_0=-97 +lat_0=40 +lat_1=33 +lat_2=45 +r =6370000.0"

粗略估计范围：

v <- vect(rbind(c(-115, 38), c(-115, 45), c(-110.05, 38), c(-110.05, 45)), crs="+proj=longlat")
p <- project(v, pcrs)
e <- crds(p) |> apply(2, range) |> as.vector()
e
#[1] -1562368.5 -1025418.3  -139104.3   693662.9

设置范围和crs

ext(ras_dt) <- e
crs(ras_dt) <- pcrs

将点转换为栅格的 crs

xy <- vect(cbind(lon=-113.0,lat=40.0), crs="+proj=longlat")
pxy <- project(xy, pcrs)

并提取

extract(ras_dt, pxy)
#  ID O3_LAY=1_TSTEP=1 O3_LAY=1_TSTEP=2
#1  1         41.88482         40.99662

所以我们的栅格现在看起来像这样。空间分辨率应该是一个整数（12,000？）

ras_dt
#class       : SpatRaster 
#dimensions  : 70, 45, 2  (nrow, ncol, nlyr)
#resolution  : 11932.23, 11896.67  (x, y)
#extent      : -1562369, -1025418, -139104.3, 693662.9  (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=lcc +lat_0=40 +lon_0=-97 +lat_1=33 +lat_2=45 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs 
#source      : EQUATES_ACONC_O3_SAM.nc:O3 
#varname     : O3 (O3              ) 
#names       : O3_LAY=1_TSTEP=1, O3_LAY=1_TSTEP=2 
#unit        : ppbV            , ppbV

您能否在文档的某处找到并提供正确的范围？也许要求数据提供者遵循 NetCDF 约定，以便使用数据所需的所有元数据都存储在文件中。

如何在没有投影经纬度坐标的情况下从栅格中提取 xy 坐标的数据点

How to extract data point for a xy coordiante from raster without projected lat lopn coordinates

r

raster

netcdf

r-raster