Matplotlib/Basemap display (pixel) coordinates are wrong values

Matplotlib/Basemap display (pixel) coordinates are wrong values.

GOAL
Compute display (pixel) coordinates from lat-long input and draw marker on a
map image
using Javascript. Use matplotlib/basemap to compute display coords, return
values to
web browser for drawing on top of map.

PROBLEM

The display (pixel) coordinates computed from geo coordinates by
matplotlib/basemap
are not always correct. It depends on the version of matplotlib and basemap
and it
also depends on the figure size!
Only old versions of matplotlib (0.99.1.1) and basemap (0.99.4) are always
correct!
Newer versions are wrong unless figure size is small!

NOTES: I only found out about pyproj after I started this investigtion, so
my code
example does not use pyproj, but I'm told that shouldn't matter.

My example code (at bottom) prints out geo, projection, display and axes
coordinates.
It also displays a map, using matplotlib, with markers plotted either from
geo,
projection, or axes coordinates;
I have not figured out how to correctly plot markers with display
coordinates
using matplotlib/basemap, although I can plot markers with display
coordinates using javascript.

SCENARIO

I am using Basemap with Miiler projection centered at 180 longitude and a
global map.
The map image is written as PNG and displayed in web browser.
The user can enter lat-long coordinates in HTML FORM and should then see a
marker placed
on the map. This is NOT done by redrawing the map, instead it is done by
computing display coordinates (pixels) and using HTML CANVAS draw
capabilities.

From the HTML FORM, I send lat-long coordinates to server and have

matplotlib/basemap
compute as follows:
geo coordinates (lat-long) to projection coordinates (meters) to display
coordinates (pixels).
Then I return display coordinates to web browser and have javascript place a
symbol
on top of the map image by drawing on the canvas.
This works fine in one environment, but fails on others. The failure is
really strange ...

Correct with Python 2.6.6, matpllotlib 0.99.1.1; mpl_toolkits.basemap 0.99.4
Wrong with two other environments:
  Python 2.6.6; matplotlib 1.4.2; mpl_toolkits.basemap 1.0.7
  Python 2.7.11; matplotlib 1.5.0l; mpl_toolskits.basemap 1.0.7
  
The display (pixel) x-coordinate is correct in all cases.
The display y-coordinate is wrong and appears to be a scaling issue:
it is larger than it should be - the y-coordinate plots at a larger
latitude,
except 0 degrees latitude is correct.

Here's where it gets really strange: the discrepancy depends on the figure
size
and it disappears when the figure is small enough!

Stranger still is that I can take the display coordinates (pixels) and
transform
to axes coordinates (0 to 1), and even though two different matplotlib
environments
result in unequal display coordinates, the axes coordinates are identical!

In all cases I'm using DPI=100.0

For example, here are coordinates for longitude=180, latitude=60,30, with
figsize=(10,10)
      Good Bad Diff (in pixels)
60 647.6 707.8 60.2
30 566.6 593.7 27.1
Diff 81.0 114.1

Step down to (10,8):
      Good Bad Diff
60 547.6 566.3 18.7
30 466.6 475.0 8.4
Diff 81.0 91.3

Shrink to figsize=(10,7) and all is well!
And it will be OK at any height smaller than 7.
      Good Bad Diff
30 415.6 415.6 0.0
60 495.5 495.5 0.0
Diff 79.9 79.9

I also experimented with figsize=(8,y) and found that there was a
discrepancy
when y >= 6 and no discrepancy when y <= 5.

EXAMPLE

This example is incomplete in that it does not include javascript and AJAX
code
that I am using. However, if this program is run in a matplotlib 0.99
environment
and then in a more recent matplotlib (and Basemap), you can see the
discrepancies
in the printed output. The star markers are plotted with display coords and
should
fall on the 180 longitude line (except 180,0 has a small offset).

# coding: utf-8
import matplotlib
#matplotlib.use('GTKAgg')
matplotlib.use('TKAgg')
#matplotlib.use('Agg')
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np
grid_incr = 1.0
lats0 = np.arange(-90.0,90.1,grid_incr)
lons0 = np.arange(0.0,359.9,grid_incr)
lons,lats = np.meshgrid(lons0, lats0)

plt.figure(figsize=(10,8))
map = Basemap(projection='mill', lon_0=180)
#map = Basemap(projection='moll', lon_0=180)
map.drawcoastlines()
map.drawparallels(np.arange(-90,90,30),labels=[1,0,0,0])
map.drawmeridians(np.arange(lons0[0],lons0[-1]+30,60),labels=[0,0,0,1])
map.drawmapboundary(fill_color='grey')

fig = plt.gcf() # get current figure
fig.set_frameon(True)
fig.set_dpi(100.0)
width,height = fig.canvas.get_width_height()
print 'canvas size = {0:.1f},{1:.1f}'.format(width,height)
print 'DPI = %f' %(fig.get_dpi())
print 'fig size pixels = %s' %str((fig.transFigure.transform([1,1])))
fig.patch.set_visible(True)
fig.patch.set_facecolor('green') # just want to know what part of display
is "figure"

ax = fig.gca() # get current axes

# We have map projection coords (meters) and geographic coords (degrees).
# You can transform between these two by using the Basemap object.
# map_xy = my_basemap(lon,lat)
# lon,lat = my_basemap(x,y,inverse=True)
# But how to translate into axes?
# Maybe this:
https://scipy.github.io/old-wiki/pages/Cookbook/Matplotlib/Transformations.html
# And there is some gold here, talking about chaining transforms:
https://sourceforge.net/p/matplotlib/mailman/message/19700132/
from matplotlib.transforms import IdentityTransform
# Function performs various coordinate transforms and prints values
def plot_ll(lon,lat,ctype='proj'):
    '''
    Start with lon-lat coordinate and convert to map projection (meters),
    then to display (pixel) coordinates. Finally convert display coordinates
    back to Axes coordinates (0.0 to 1.0). Print all of these.
    Then use ctype to determine which set of coordinates to use to plot
    marker on map.
    '''
    global map
    global fig
    global ax
    global disp_org
    xy_format = '({0:.1f}, {1:.1f})'
    geo_xy = (lon,lat)
    geo_xy_s = '({0:5.1f}, {1:5.1f})'.format(lon,lat)
    proj_xy = map(lon,lat)
    proj_xy_s = '({0:11.1f}, {1:11.1f})'.format(proj_xy[0],proj_xy[1])
    disp_xy = ax.transData.transform(proj_xy)
    disp_xy_s = '({0:6.1f}, {1:6.1f})'.format(disp_xy[0],disp_xy[1])
    inv = ax.transAxes.inverted()
    axes_xy = inv.transform(disp_xy)
    axes_xy_s = '({0:.3f}, {1:.3f})'.format(axes_xy[0],axes_xy[1])
    #print 'plot_ll: geo=%s, proj=%s, disp=%s, axes=%s'
%(str(geo_xy),proj_xy_s,disp_xy_s,axes_xy_s)
    print 'plot_ll: geo=%s, proj=%s, disp=%s, axes=%s'
%(geo_xy_s,proj_xy_s,disp_xy_s,axes_xy_s)
    if ctype == 'proj':
        mcolor = 'r' # red for projection coordinates
        marker = 'o'
        msize = 8
        # Use map object to plot proj_xy with no transform
        map.plot(proj_xy[0],proj_xy[1],color=mcolor,marker=marker,ms=msize)
    elif ctype == 'disp':
        mcolor = 'b' # blue for display coords
        marker = 'o'
        msize = 8
        # Next does not work - I thuoght it would plot '*' marker using
display coords,
        # but the position is not correct.
       
ax.plot(disp_xy[0]-disp_org[0],disp_xy[1]-disp_org[1],color=mcolor,marker='*',ms=msize,
transform=IdentityTransform())
        # Instead invert back to proj_xy and plot
        inv1 = ax.transData.inverted()
        pxy = inv1.transform(disp_xy)
        # Use ax object to plot proj_xy with transform=ax.transData
        ax.plot(pxy[0],pxy[1],color=mcolor,marker=marker,ms=msize,
transform=ax.transData)
    elif ctype == 'axes':
        mcolor = 'g' # green for axes coords
        marker = 'o'
        msize = 8
       
ax.plot(axes_xy[0],axes_xy[1],color=mcolor,marker=marker,ms=msize,transform=ax.transAxes)

# plot points on map, use different coordinates when plotting
proj_org = map(0.0,-90.0) # map origin in projection coords
disp_org = ax.transData.transform(proj_org) # map origin in display coords
print 'display origin: ({0:6.1f}, {1:6.1f})'.format(disp_org[0],disp_org[1])

plot_ll(180.0,0.0,ctype='proj')
plot_ll(180.0,89.0,ctype='proj')
plot_ll(180.0,-89.0,ctype='proj')

plot_ll(175.0,0.0,ctype='axes')
plot_ll(175.0,89.0,ctype='axes')
plot_ll(175.0,-89.0,ctype='axes')
plot_ll(5.0,0.0,ctype='axes')
plot_ll(5.0,89.0,ctype='axes')
plot_ll(5.0,-89.0,ctype='axes')

plot_ll(185.0,0.0,ctype='disp')
plot_ll(185.0,89.0,ctype='disp')
plot_ll(185.0,-89.0,ctype='disp')
plot_ll(355.0,0.0,ctype='disp')
plot_ll(355.0,89.0,ctype='disp')
plot_ll(355.0,-89.0,ctype='disp')
plot_ll(180.0,30.0,ctype='disp')
plot_ll(180.0,60.0,ctype='disp')
plot_ll(180.0,-30.0,ctype='disp')
plot_ll(180.0,-60.0,ctype='disp')

plt.show()

···

--
View this message in context: http://matplotlib.1069221.n5.nabble.com/Matplotlib-Basemap-display-pixel-coordinates-are-wrong-values-tp46900.html
Sent from the matplotlib - users mailing list archive at Nabble.com.