EGB13

Here are the results of some playing with gdalwarp and NASA’s Blue Marble images of the Earth.  My experiences with gdalwarp and NOAA’s national weather radar mosaic, used to create a time-lapse animation of twelve months of weather, prompted a comment on the video’s YouTube page by someone having difficulty using gdalwarp with Blue Marble images.  That prompted me to write up how I used my limited understanding of gdalwarp and related subjects for those radar images.  But of course I had to try my hand on the Blue Marble next.  What follows is how I went about it and the results I got.

I was initially put off by this post by a subject matter expert, which I definitely am not, explaining why “EPSG:4326 is usually the wrong ‘projection’.”  (See my earlier post that very briefly explains EPSG codes and their use in gdalwarp.)  After trying to follow his reasoning and to find some implementation of his suggested alternative, I went back and tried EPSG:4326 in this context and found that it works fine, at least as far as my eye can tell at these scales.  So for now I have put aside Mr. Nielsen’s concerns to investigate further another day.  In all that follows EPSG:4326 rules for the source projection.

The EPSG code is a shorthand way to define the use of a specific spatial reference system (SRS).  EPSG:4326 appears to work for this case where the source is a rectangular array of pixels with fixed a constant scale.  Each pixel, regardless of where it is in the source image, has a constant degree width, and each pixel in the image has a constant degree height.  Therefore the pixel’s [x,y] position in the image is trivially converted to geographic longitude (x) and latitude (y).

An image of the world from NASA's Blue Marble

Here is a “Blue Marble” image from NASA’s Earth Observatory.  It is a “plate carrée” equirectangular projection that is not really a projection in the geometric sense.  But it has that constancy of scale across the whole map to which we can apply EPSG:4326 and with which we can do some simple arithmetic for related purposes.

The initial image I worked with was scaled to a pixel dimensions of 1024 wide by 512 high for convenience.  The image spans the whole Earth — 360° west-to-east and 180° from pole to pole.  So each pixel is 360°/1024 wide and 180°/512 high, both of which equal 0.3515625° per pixel.  So the world file for this image is trivially:

0.3515625
0.00000
0.00000
-0.3515625
-180.0
90.0

Albers equal area conic projection from page 35 of the proj v. 3 user manual.

The world file calculator will do this math for you and produce the world file in the proper format. Remember that the world file is needed to tell the software where the image fits in the spatial reference system (SRS) being used.  We now have both an SRS (EPSG:4326) and a world file, so we’re ready to try some projections.

Trial projections

To the right are two images. The first is an excerpt with associated text taken from the user manual by Gerald I. Evenden for version 3 of the proj program. I manually colored the map to make it clearer.  Many GIS programs, including gdalwarp, use the same +proj notation for SRS notation, so we can take this example of an Albers equal area conic projection from the proj document and try to reproduce it with gdal using the Blue Marble image.  The command is:

gdalwarp -s_srs EPSG:4326 -t_srs '+proj=aea +lon_0=90w +lat_1=20n +lat_2=60n'  world.jpg aea-world.jpg

Blue Marble in Albers equal area conic projection

The result in the second image is almost identical to the example shown in the manual.  Success!  The gdalwarp program does crop the projection differently.  I manually shaded the portion of the projection encompassed by a triangle between the center of the north polar circle and where the outside of the projection meets the top of the frame.  The remaining color part of the projection appears to exactly reproduce the manual’s illustration.

Excerpt from page 36 of the original "proj" (version 3) user manual showing a polyconic projection.

The next projection that I tried is the “Polyconic (American) Projection” on page 37 of the version 3 manual for proj.  The corresponding gdalwarp command used here is:

gdalwarp -s_srs EPSG:4326 -t_srs '+proj=poly +lon_0=90w'  world.jpg poly-world.jpg

On the left is the example image from the manual, and the result from gdalwarp is show below.  Gdalwarp apparently only shows the results that fit in the center circle of the manual’s example.  I manually colored the example to highlight that portion.  The gdalwarp results again very closely mirror the example, giving confidence that this approach is working correctly.  You will notice that whereas the example in the proj documentation appears perfectly circular, the output of gdalwarp seems slightly elongated, and I have no explanation for that yet.  The overall results appear correct, though — compare Africa, Greenland, Alaska, and Antarctica in the two images.  I will call my initial efforts at projecting the Blue Marble image are a success.

Blue Marble in polyconic projection

One minor note to add about the use of gdalwarp in these projections is that if there is an EPSG equivalent to the target projections, then those codes could be substituted for the target SRS.  I didn’t bother to look those up because I wanted to be sure that I exactly follwed the examples.

My next post will document a more complex and practical use of the Blue Marble image.