0

在写入然后读取具有某些文件格式（ESRI .hdr 标签格式、netCDF 和 ENVI，但不是 GeoTIFF 或 Erdas HFA）的 SpatRaster 时，光栅原点、分辨率和坐标略有变化。这些参数的相等性测试足以返回 FALSE 值。这是一个例子，GeoTiff 然后是 flt。

Nrow=45; Ncol=108
r1 <- rast(nrow=Nrow,ncol=Ncol,vals=1:(Nrow*Ncol),ext=c(-180,180,-60,90))
origin(r1)
# [1] 0 0

writeRaster(r1,"test.tif",overwrite=TRUE)

r2 <- rast("test.tif")
origin(r2)  
# origin is OK, same as r1
# [1] 0 0
 
res(r1)==res(r2)
# [1] TRUE TRUE

writeRaster(r1,"test.flt",overwrite=TRUE)
r3 <- rast("test.flt")
origin(r3)  
# origin has been shifted
# [1] 1.705303e-13 5.684342e-14

res(r1)==res(r3)   
# the lag is large enough for r3 to be considered as a different raster from r1
# [1] FALSE FALSE

这里有三个问题：

. 为什么会发生这样的转变？

. 有没有办法阻止它们？

. 如果没有，阅读后纠正它们的最佳方法是什么？

谢谢。

4

1 回答 1

0

发生偏移的原因是用于存储分数的小数位数不足，存储范围或栅格的错误方法加剧了这种情况。该FLT文件使用“原点/分辨率”方法进行地理配准

readLines("test.hdr")
 [1] "BYTEORDER      I"                 "LAYOUT         BIL"               "NROWS          45"               
 [4] "NCOLS          108"               "NBANDS         1"                 "NBITS          32"               
 [7] "BANDROWBYTES   432"               "TOTALROWBYTES  432"               "PIXELTYPE      FLOAT"            
[10] "ULXMAP         -178.333333333333" "ULYMAP         88.3333333333333"  "XDIM           3.33333333333333" 
[13] "YDIM           3.33333333333333"  "NODATA         1.#QNAN"

可以看到，有ULXMAPand ULYMAP（单元格中心的左上x和y坐标），但是没有存储另一个角的坐标。所以必须计算。

在这种情况下，最小（左）x 坐标是

XDIM <- 3.33333333333333
xmin <- -178.333333333333 - 0.5 * XDIM
print(xmin, 16)
# -179.9999999999997

并且最大（右）x坐标是

NCOLS <- 108
xmax <- xmin + NCOLS * XDIM
print(xmax, 16)
#180

yres = (xmax - xmin) / NCOLS
print(yres, 17)
[1] 3.3333333333333304

它应该在哪里

print(360/NCOLS, 17)
[1] 3.3333333333333335

一个非常小的差异，但可以通过存储最大 x 和 y 坐标而不是 x 和 y 分辨率来轻松避免精度损失。不幸的是，这不是在这里完成的，甚至 GDAL 在内部也使用了这种（我认为有缺陷的）方法。

我认为您看到的文件类型之间的差异是由于存储的小数位数（双精度与单精度）。

如果这会导致问题（通常会忽略如此小的差异），您可以通过分配正确的范围来解决此问题

 ext(r3) <- c(-180, 180, -90, 90)
于 2021-04-10T05:11:02.797 回答