将兰伯特方位角等面积的 Hydro1k 更改为 WGS 84

Change Hydro1k with Lambert Azimuthal Equal Area to WGS 84

我正在使用 wgs84 投影分析栅格数据,我还需要从 USGS (https://lta.cr.usgs.gov/HYDRO1K) 下载的 Hydro1k 数据。然而,Hydro1k 使用具有特定原点和其他参数的 Lambert Azimuthal Equal Area 投影。所以我决定将 Hydro1k 的投影转换为 WGS 84 以供日后分析。试了网上找的一种方法,改投影好像没问题。

下面是我使用的代码和栅格的属性:

#Import the selected raster after knowing name
As.fa.1k=raster(Path.all.1k[4])

#Show the attribute of the input raster from README file:
Projection used:  Lambert Azimuthal Equal Area

                  Units = meters

                  Pixel Size = 1000 meters

                  Radius of Sphere of Influence = 6,370,997 meters

                  Longitude of Origin = 20 00 00E

                  Latitude of Origin = 55 00 00N

                  False Easting = 0.0

                  False Northing = 0.0

#Custom the projection based on the origin of raster 
As.proj="+proj=laea +lon_0=100 +lat_0=45 +ellps=sphere"

#Assign CRS to layer
crs(As.fa.1k)=As.proj

#Get wgs4 projection
wgs = "+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0"

#Assign wgs84 to raster with lambert projection
As.wgs=projectRaster(As.fa.1k, crs=wgs)

然而,经过投影转换的栅格中的值看起来很奇怪:

#values in raster before being transformed
class       : RasterLayer 

dimensions  : 8384, 9102, 76311168  (nrow, ncol, ncell)

resolution  : 1000, 1000  (x, y)

extent      : -4462500, 4639500, -3999500, 4384500  (xmin, xmax, ymin, ymax)

coord. ref. : +proj=laea +lon_0=-100 +lat_0=45 +ellps=sphere 

names       : na_fd 

values      : -9999, 255  (min, max)

-9999应该是海洋中单元格的值。

#values in raster after after transformed
class       : RasterLayer 

dimensions  : 1910, 5497, 10499270  (nrow, ncol, ncell)

resolution  : 0.0655, 0.04495  (x, y)

extent      : -180.0035, 180.05, -0.5929785, 85.26152  (xmin, xmax, ymin, ymax)

coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 

names       : na_fd 

values      : 1, 55537  (min, max)

而且我提取了栅格变换前的值,最大值也是55537,看来栅格变换前后的值相差很大

如何在转换之前保持栅格的值?

谢谢。

projectRaster 将在重新投影期间扭曲栅格像元,尤其是从投影坐标系(Lambert Azimuthal)到地理坐标系(wgs84),并且需要对值进行一些插值。

您还应该注意到单元格大小的巨大差异,也许您可​​以在 projectRaster 中使用 'res=' 强制分辨率?

此外,尝试 'method = "ngb"' 作为 projectRaster 中的插值方法以保持值相同