我有一个Geotiff栅格图像,其坐标系的经度从0到360。
图像的水平中心是180度。参见下图:



我想将其转换为EPSG:4326 SRS,经度范围为-180 180。
我希望图像的中心位于格林威治子午线(0)。我想这个SRS用途非常广泛。我希望结果看起来像这样:



所以我使用gdalwarp命令重新投影:

gdalwarp -s_srs '+proj=latlong +datum=WGS84 +pm=180dW' -t_srs EPSG:4326 test_col.tif test_4326.tif


但是我只会看到尺寸更大(像素更多)和EPSG:4326元数据的tiff。图像本身看起来与初始图像相同。但是,我希望它能交换半球。

我如何对图像进行严格扭曲,使其经度严格定为-180 180 EPSG:4326,并且中心位于经度为0?

这是gdalinfo我的初始文件的内容:

Origin = (-0.102272598067084,89.946211604095552)
Pixel Size = (0.204545196134167,-0.204423208191126)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (  -0.1022726,  89.9462116) (  0d 6' 8.18"W, 89d56'46.36"N)
Lower Left  (  -0.1022726, -89.9462116) (  0d 6' 8.18"W, 89d56'46.36"S)
Upper Right (     359.897,      89.946) (359d53'50.18"E, 89d56'46.36"N)
Lower Right (     359.897,     -89.946) (359d53'50.18"E, 89d56'46.36"S)
Center      ( 179.8975000,  -0.0000000) (179d53'51.00"E,  0d 0' 0.00"S)


这是gdalwarp之后的gdalinfo

Origin = (-180.102727401932952,89.946211604095552)
Pixel Size = (0.091397622896436,-0.091420837939082)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (-180.1027274,  89.9462116) (180d 6' 9.82"W, 89d56'46.36"N)
Lower Left  (-180.1027274, -89.9699975) (180d 6' 9.82"W, 89d58'11.99"S)
Upper Right ( 179.8211116,  89.9462116) (179d49'16.00"E, 89d56'46.36"N)
Lower Right ( 179.8211116, -89.9699975) (179d49'16.00"E, 89d58'11.99"S)
Center      (  -0.1408079,  -0.0118929) (  0d 8'26.91"W,  0d 0'42.81"S)


#1 楼

您可以使用目标范围选项将输出坐标范围显式设置为gdalwarp(即“ -te -180 -90 180 90”),但也可以使用CENTER_LONG配置选项来强制重新包装新的中心经度。这样的事情:
  gdalwarp -s_srs "+proj=longlat +ellps=WGS84" -t_srs WGS84 ~/0_360.tif 180.tif  -wo SOURCE_EXTRA=1000 \
           --config CENTER_LONG 0

还要注意“ SOURCE_EXTRA = 1000”变形选项。重新包装时,源矩形的计算会混淆经度中断,并丢失一些图像。此选项表示需要额外投入。没有它,您将在本初子午线附近看到一个数据间隙。
我认为像您一样设置本初子午线180dW并不是一个好主意。

评论


嗯,--config CENTER_LONG 0不执行任何操作,结果是相同的栅格。我在这里想念什么吗?在GDAL 2.2.3版上运行。

– jurajb
19年7月8日在9:56

#2 楼

基本上,您需要将栅格切成两部分,并用新的偏移/比例将它们重新拼合。

这里有一个示例,说明如何使用gdal_translate和VRT驱动程序从[-180,180]到[0,360]进行操作:
http://trac.osgeo.org/gdal/wiki / UserDocs / RasterProcTutorial

向下扫描到“ 5分钟教程”,详细信息位于“虚拟文件”下。修改示例以使其适合应该足够简单。

#3 楼

可以使用带有rotate软件包的raster函数在R中使用一行代码来完成此操作。

首先加载栅格数据包。

library(raster)
your_raster <- raster("path/to/raster.tif")
rotated_raster <- rotate(your_raster)


#4 楼

如果只想在QGIS中查看栅格,则可以使用参数+ lon_wrap = 180设置自定义投影。

我对此的理解是,默认情况下,proj4会将纬度从0-> 360包裹到-180->180。+ lon_wrap = 180将有效地取消这种包裹,并显示180至西半球360。

选项+ over应该完全禁用环绕,但是-至少就我而言-使用该选项时栅格无法正确显示。

有关更多信息,请参见http://proj4.org/parameters.html#lon-wrap-over-longitude-wrapping。

#5 楼

这是我构建的函数,用于使用JavaScript从0-360到-180-180重新投影网格值的单个暗淡数组。