2

I'm trying to reproject the Natural Earth I 1:50m raster into an Azimuthal Equidistant projection centered on my hometown. As you may know, the source raster covers from -180° to 180° latitude, and -90° to 90° Here's the command I'm using:

gdalwarp -multi -dstnodata 0 -ot Byte -t_srs "+proj=aeqd +lat_0=48.037 +lon_0=-122.408 +ellps=sphere +units=m +no_defs" -r cubic -of GTiff -ts 10800 10800 NE1_50M_SR_W.tif natural_earth_10800.tif

This works pretty OK, except that it leaves a big, ugly gap west of 180°. Also, it cuts off some of each side of the output raster.

Azimuthal equidistant map of the world, physical theme, with gap along International Dateline

I've experimented with setting a radius instead of the +ellps parameter, with no change. Using WGS84 gives only a small slice of the Earth as output.

Here's what I got when I tried +R_A instead of +ellps=sphere:

Physical map of less than a quarter of the Earth, showing North America with parts of South America and Asia

The result with WGS84 was much the same. Also tried reprojecting to Spherical Mercator as an intermediate projection before reprojecting to my custom AEqD, same artifacts appear.

So, how do I fill in the gap? The pixels are definitely there in the input raster, they're just not making it into the output.

  • You might get these artefacts while reprojecting from the WGS84 ellipsoid to the sphere. – AndreJ Dec 1 '15 at 7:20
  • I did try using Spherical Mercator as an intermediate, thinking that would eliminate any WGS84=>sphere problems, but the same artifacts appeared. – dericke Dec 1 '15 at 7:54
4

I got rid of the date line gap with these options:

gdal_translate -a_srs "+proj=longlat +ellps=sphere" NE1_50M_SR_W.tif ne.tif
gdalwarp -dstnodata 0 -ot Byte -wo SOURCE_EXTRA=1000 -t_srs "+proj=aeqd +lat_0=48.037 +lon_0=-122.408 +ellps=sphere +units=m +no_defs" -of GTiff  ne.tif neaeqd.tif

The SOURCE_EXTRA write option ensures an overlap:

SOURCE_EXTRA: This is a number of extra pixels added around the source window for a given request, and by default it is 1 to take care of rounding error. Setting this larger will increase the amount of data that needs to be read, but can avoid missing source data.

I still get a tolerance condition error, which might come from the antipodal point, reprojected to infinity. The cuts on the sides remain:

enter image description here

  • That did it, cheers! – dericke Dec 2 '15 at 0:16

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.