Create simple Basemap that corresponds to a specific region on the Earth

Question

I'm am trying to familarize myself with matplotlib and Basemap. As a start, I'm trying to generate an image of greenland that matches a specific grid for which I have data.

The gruesome details below describe my problem: I can't create an image that is the proper size matching the projection/area desired.

The projection and grid that I would like to match: Projection as Proj4 string: "+proj=stere +lat_0=90 +lon_0=-45 +lat_ts=70 +ellps=WGS84 +datum=WGS84 +units=m"

The area defined by a grid is an 800x1400 2000m resolution grid where: Outer edge of LowerLeft Corner(m): -700,000., -3,400,000. Outer edge of UpperRight Corner (m): 900,000., -600,000. => (-700,000 + 2000 * 800, -3,400,000 + 2000 * 1400)

Basemap won't let me specify the corners in meters-xy for a stereographic projection so I have to convert these to lat/lon.

> gdaltransform -s_srs "+proj=stere +lat_0=90 +lon_0=-45 +lat_ts=70 +ellps=WGS84 +datum=WGS84 +units=m" -t_srs "+proj=latlong"`
-700000 -3400000
-56.6336339989404 58.7244253840871 0
900000 -600000
11.3099324740202 80.0389929796586 0

Now I should have all of the information to create an 800x1400 image.

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt

def create_map():
    fig = plt.figure(1, figsize=(8, 14), frameon=False, dpi=100)
    fig.add_axes([0, 0, 1, 1])

    m = Basemap(resolution="i",
                projection='stere', lat_ts=70, lat_0=90., lon_0=-45.,
                llcrnrlon=-56.6336339989404, llcrnrlat=58.7244253840871,
                urcrnrlon=11.3099324740202, urcrnrlat=80.0389929796586,
                rsphere=(6378137.0, 6356752.3142))

    m.drawcoastlines()
    m.fillcontinents(color='#c1c1c1')
    m.drawmapboundary(fill_color='#6587ad', linewidth=0.0)
    plt.savefig('greenland.png', pad_inches=0.0, bbox_inches='tight')

if __name__ == '__main__':
    create_map()

The problem I face is that when I do this, I get an 800x1399 image. If I don't include the bbox_inches='tight' in the plt.savefig command, I get an 800x1400 image, with a single strip of invisible pixels along the (edit) bottom edge (edit).

Can anyone help me so that I can be sure that I am setting up my Basemap properly? I feel like I'm probably just missing a simple trick, but not getting an image the size I expect is strange.

As always, thanks in advance.

Answer 1

It appears that this might be a result of a bug in matplotlib. Jeff Whitaker took a look and said it looked right and I tried to reproduce this behavior without using Basemap and I was able to.

It seems that the aspect of the data values might cause the output image to be the wrong size.

Here's some code that shows the problem. Sorry for the false alarm.

# rectangle.py   --
import matplotlib.pyplot as plt


def create_image():
    fig = plt.figure(1, figsize=(8, 14), frameon=False, dpi=100)
    fig.add_axes([0, 0, 1, 1])
    ax = plt.gca()

    # This isn't necessary to create the issue unless you want to see the
    # transparent pixels at bottom.

    # for spine in ax.spines.values():
    #     spine.set_linewidth(0.0)

    limb = ax.axesPatch
    limb.set_facecolor('#6587ad')

    x1 = 0.0
    y1 = 0.0
    x2 = 16.

    # Use this line and get what I was expecting:
    #    y2 = 27.999999999999994671   # produces 800 x 1400 image

    # Use this line and get the wrong size
    y2 = 27.999999999999994670   # produces (wrong?) 800 x 1399 image

    corners = ((x1, y1), (x2, y2))
    ax.update_datalim(corners)
    ax.set_xlim((x1, x2))
    ax.set_ylim((y1, y2))

    ax.set_aspect('equal', anchor='C')
    ax.set_xticks([])
    ax.set_yticks([])

    plt.savefig('rectangle.png', pad_inches=0.0, bbox_inches='tight')

    # If you use this below, the file size is correct, but there is a single
    # line transparent pixels along the bottom of the image if you set the
    # linewidth to zero...

    #  plt.savefig('rectangle.png', pad_inches=0.0)


if __name__ == '__main__':
    create_image()