R光栅`projectExtent`“不能做这个转换”答案

【问题标题】：R raster `projectExtent` "cannot do this transformation"R光栅`projectExtent`“不能做这个转换”
【发布时间】：2021-06-11 12:18:42
【问题描述】：

我遇到了这个问题。从下面的一个大栅格中，我想切出 100 个框，存储为 extent 对象。范围来自不同的 CRS，因此它们的单位错误。这是裁剪到法国的 CLC 土地利用栅格：

> CLCF
class      : RasterLayer 
dimensions : 11000, 12000, 1.32e+08  (nrow, ncol, ncell)
resolution : 100, 100  (x, y)
extent     : 3100000, 4300000, 2100000, 3200000  (xmin, xmax, ymin, ymax)
crs        : +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs 
source     : /private/var/folders/ry/d__dq69j6w31gmbnhz7ybj013sn99p/T/Rtmpd4xShx/raster/r_tmp_2021-06-11_133707_3696_56233.grd 
names      : U2018_CLC2018_V2020_20u1 
values     : 1, 128  (min, max)

我的边界框的 CRS 是 WGS84。

> raster::extent(ms[[ix]]$inverse_mask)
class      : Extent 
xmin       : 371250 
xmax       : 381250 
ymin       : 5838250 
ymax       : 5845250 
> raster::crs(ms[[ix]]$inverse_mask)
CRS arguments:
 +proj=moll +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs

遵循常见的建议，这就是我所做的：

empty = raster::raster()
raster::extent(empty) <- raster::extent(ms[[ix]]$inverse_mask)

> empty
class      : RasterLayer 
dimensions : 180, 360, 64800  (nrow, ncol, ncell)
resolution : 27.77778, 38.88889  (x, y)
extent     : 371250, 381250, 5838250, 5845250  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs 

> raster::projectExtent(empty,crs = raster::projection(CLCF))
Error in raster::projectExtent(empty, crs = raster::projection(CLCF)) : 
  cannot do this transformation
In addition: Warning messages:
1: In showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj = prefer_proj) :
  Discarded datum Unknown based on GRS80 ellipsoid in Proj4 definition
2: In rgdal::rawTransform(projfrom, projto, nrow(xy), xy[, 1], xy[,  :
  663 projected point(s) not finite

我在这里做错了什么？谢谢

【问题讨论】：

标签： r r-raster

【解决方案1】：

感谢Transform result of `st_bbox()` to other CRS找到答案

bb = sf::st_bbox(ms[[1]]$inverse_mask)  # sf to the rescue
            bb_ll = sf::st_bbox(
                sf::st_transform( 
                    sf::st_as_sfc(bb), crs = "+proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs")
            )
            e = raster::extent(c(bb_ll["xmin"],bb_ll["xmax"],bb_ll["ymin"],bb_ll["ymax"]))
            OL[[ix]] <- raster::crop(CLCF, raster::extent(ms[[ix]]$inverse_mask))

【讨论】：