对齐 Natural Earth Geojson 和 Raster 以在 D3 中呈现

Aligning Natural Earth Geojson and Raster to render in D3

我正在尝试使用 D3 渲染带有高程数据的世界地图。

为此,我使用 Natural Earth 50m land geojson: https://github.com/martynafford/natural-earth-geojson/tree/master/50m/physical

和自然地球高程栅格数据: https://www.naturalearthdata.com/downloads/50m-raster-data/50m-shaded-relief/

我正在使用本教程:https://datawanderings.com/2020/08/08/raster-backgrounds/

所以我首先找到了 geojson 的边界:

ogrinfo ne_50m_land.json -so -al | grep Extent  
Extent: (-180.000000, -89.998926) - (180.000000, 83.599609)

然后我剪切了 TIF 文件:

gdal_translate -projwin -180.000000 83.599609 180.000000 -89.998926 SR_50M.tif SR_50M-box.tif

并投影到墨卡托投影中:

gdalwarp -overwrite -s_srs "+proj=longlat +datum=WGS84 +no_defs" -t_srs EPSG:3395 SR_50M-box.tif SR_50M-proj.tif

最后我导出到 PNG :

gdal_translate -of PNG SR_50M-proj.tif SR_50M.png

然后我使用 D3 渲染所有内容

   this.projection = d3.geoMercator()
        .translate([0, 0])
        .scale(1)


   this.rasterImage = new Image();
   this.rasterImage.src = raster;

   this.path = d3.geoPath(this.projection, this.ctx);
   this.bb = this.path.bounds(this.land);

   const s = 1 / Math.max((this.bb[1][0] - this.bb[0][0]) / this.getWidth(), (this.bb[1][1] - this.bb[0][1]) / this.getHeight());
    // transform
   const t = [(this.getWidth() - s * (this.bb[1][0] + this.bb[0][0])) / 2, (this.getHeight() - s * (this.bb[1][1] + this.bb[0][1])) / 2];

    // update projection
   this.projection
       .scale(s)
       .translate(t);

   this.raster_width = (this.bb[1][0] - this.bb[0][0]) * s;
   this.raster_height = (this.bb[1][1] - this.bb[0][1]) * s;

   this.rtranslate_x = (this.getWidth() - this.raster_width) / 2;
   this.rtranslate_y = (this.getHeight() - this.raster_height) / 2;

   this.ctx.beginPath();
   this.path(this.land);
   this.ctx.fill();

   this.ctx.save();
   this.ctx.globalAlpha = 0.8;
   this.ctx.translate(this.rtranslate_x, this.rtranslate_y);
   this.ctx.drawImage(this.rasterImage, 0, 0, this.raster_width, this.raster_height);
   this.ctx.restore();

然而最后,geojson 土地和栅格没有对齐:

我试过先投影再剪切,或者投影时使用伪墨卡托或墨卡托,但都没有用..任何人有想法吗?

墨卡托通常被剪裁成大约 85 度 North/South (~85.05113 N/S) - 比这更远你会得到一张比宽高的地图,并且得到更多范围内每增加一个学位 north/south 就会高很多..

使用此限制的 D3 剪辑功能:

The spherical Mercator projection. Defines a default projection.clipExtent such that the world is projected to a square, clipped to approximately ±85° latitude.

北部边界很好,但是 geojson 的南部边界是 -89.998926 度,您可以使用它来切割图像。但是当 D3 剪辑 geojson 时,与 geojson 相比,您将图像拉伸不同的量,因此您会看到这个问题。

解决方案应该是将图像裁剪到一个边界,该边界代表 D3 将为墨卡托(85.05113 度以南)呈现的限制,而不是数据本身的限制。

我没有查过 gdal 是如何忠实地实现 EPSG:3395 的,因为定义提供了南纬 80 度和北纬 84 度的投影边界 - 虽然看图片,这并不看起来不是问题。

您还可以对 D3 投影 (d3v4+) 使用更清晰的 fitSize 方法:

 projection.fitSize([width,height],geojsonObject)

这将根据提供的 width/height.

为您设置比例和翻译