Basemap toolkit: problem with polarplot example

David Huard wrote:

Hi Jeff,

Here is the output. Thanks for looking into this.

In [1]: from matplotlib.toolkits.basemap import pyproj

In [2]: pyproj.test()
Trying:
    from pyproj import Geod
Expecting nothing
ok
Trying:
    g = Geod(ellps='clrk66') # Use Clarke 1966 ellipsoid.
Expecting nothing
ok
Trying:
    boston_lat = 42.+(15./60.); boston_lon = -71.-(7./60.)
Expecting nothing
ok
Trying:
    portland_lat = 45.+(31./60.); portland_lon = -123.-(41./60.)
Expecting nothing
ok
Trying:
    newyork_lat = 40.+(47./60.); newyork_lon = -73.-(58./60.)
Expecting nothing
ok
Trying:
    london_lat = 51.+(32./60.); london_lon = -(5./60.)
Expecting nothing
ok
Trying:
    az12,az21,dist = g.inv(boston_lon,boston_lat,portland_lon,portland_lat)
Expecting nothing
**********************************************************************
File "/usr/local/lib/python2.4/site-packages/matplotlib/toolkits/basemap/pyproj.py", line 401, in matplotlib.toolkits.basemap.pyproj.Geod.__new__
Failed example:
    az12,az21,dist = g.inv(boston_lon,boston_lat,portland_lon,portland_lat)
Exception raised:
    Traceback (most recent call last):
      File " doctest.py", line 1248, in __run
        compileflags, 1) in test.globs
      File "<doctest matplotlib.toolkits.basemap.pyproj.Geod.__new__[6]>", line 1, in ?
        az12,az21,dist = g.inv(boston_lon,boston_lat,portland_lon,portland_lat)
      File "/usr/local/lib/python2.4/site-packages/matplotlib/toolkits/basemap/pyproj.py", line 478, in inv
        _Geod._inv(self, inx, iny, inz, ind, radians=radians)
      File "_geod.pyx", line 123, in _geod.Geod._inv
    ValueError: undefined inverse geodesic (may be an antipodal point)
Trying:
    print "%7.3f %6.3f %12.3f" % (az12,az21,dist)
Expecting:
    -66.531 75.654 4164192.708
**********************************************************************
File "/usr/local/lib/python2.4/site-packages/matplotlib/toolkits/basemap/pyproj.py", line 402, in matplotlib.toolkits.basemap.pyproj.Geod.__new__
Failed example:
    print "%7.3f %6.3f %12.3f" % (az12,az21,dist)
Exception raised:
    Traceback (most recent call last):
      File "doctest.py", line 1248, in __run
        compileflags, 1) in test.globs
      File "<doctest matplotlib.toolkits.basemap.pyproj.Geod.__new_ _[7]>", line 1, in ?
        print "%7.3f %6.3f %12.3f" % (az12,az21,dist)
    NameError: name 'az12' is not defined
Trying:
    endlon, endlat, backaz = g.fwd(boston_lon, boston_lat, az12, dist)
Expecting nothing
**********************************************************************
File "/usr/local/lib/python2.4/site-packages/matplotlib/toolkits/basemap/pyproj.py", line 406, in matplotlib.toolkits.basemap.pyproj.Geod.__new_ _
Failed example:
    endlon, endlat, backaz = g.fwd(boston_lon, boston_lat, az12, dist)
Exception raised:
    Traceback (most recent call last):
      File "doctest.py", line 1248, in __run
        compileflags, 1) in test.globs
      File "<doctest matplotlib.toolkits.basemap.pyproj.Geod.__new__[8]>", line 1, in ?
        endlon, endlat, backaz = g.fwd(boston_lon, boston_lat, az12, dist)
    NameError: name 'az12' is not defined
Trying:
    print "%6.3f %6.3f %13.3f" % (endlat,endlon,backaz)
Expecting:
    45.517 -123.683 75.654
**********************************************************************
File "/usr/local/lib/python2.4/site-packages/matplotlib/toolkits/basemap/pyproj.py", line 407, in matplotlib.toolkits.basemap.pyproj.Geod.__new__
Failed example:
    print "%6.3f %6.3f %13.3f" % (endlat,endlon,backaz)
Exception raised:
    Traceback (most recent call last):
      File "doctest.py ", line 1248, in __run
        compileflags, 1) in test.globs
      File "<doctest matplotlib.toolkits.basemap.pyproj.Geod.__new__[9]>", line 1, in ?
        print "%6.3f %6.3f %13.3f" % (endlat,endlon,backaz)
    NameError: name 'endlat' is not defined
Trying:
    lons1 = 3*[newyork_lon]; lats1 = 3*[newyork_lat]
Expecting nothing
ok
Trying:
    lons2 = [boston_lon, portland_lon, london_lon]
Expecting nothing
ok
Trying:
    lats2 = [boston_lat, portland_lat, london_lat]
Expecting nothing
ok
Trying:
    az12,az21,dist = g.inv(lons1,lats1,lons2,lats2)
Expecting nothing
**********************************************************************
File "/usr/local/lib/python2.4/site-packages/matplotlib/toolkits/basemap/pyproj.py", line 414, in matplotlib.toolkits.basemap.pyproj.Geod.__new__
Failed example:
    az12,az21,dist = g.inv(lons1,lats1,lons2,lats2)
Exception raised:
    Traceback (most recent call last):
      File "doctest.py", line 1248, in __run
        compileflags, 1) in test.globs
      File "<doctest matplotlib.toolkits.basemap.pyproj.Geod.__new_ _[13]>", line 1, in ?
        az12,az21,dist = g.inv(lons1,lats1,lons2,lats2)
      File "/usr/local/lib/python2.4/site-packages/matplotlib/toolkits/basemap/pyproj.py", line 478, in inv
        _Geod._inv(self, inx, iny, inz, ind, radians=radians)
      File "_geod.pyx", line 123, in _geod.Geod._inv
    ValueError: undefined inverse geodesic (may be an antipodal point)
Trying:
    for faz,baz,d in zip(az12,az21,dist): print "%7.3f %7.3f %9.3f" % (faz,baz,d)
Expecting:
     54.663 -123.448 288303.720
    -65.463 79.342 4013037.318
     51.254 -71.576 5579916.649
**********************************************************************
File "/usr/local/lib/python2.4/site-packages/matplotlib/toolkits/basemap/pyproj.py", line 415, in matplotlib.toolkits.basemap.pyproj.Geod.__new__
Failed example:
    for faz,baz,d in zip(az12,az21,dist): print "%7.3f %7.3f %9.3f" % (faz,baz,d)
Exception raised:
    Traceback (most recent call last):
      File "doctest.py", line 1248, in __run
        compileflags, 1) in test.globs
      File "<doctest matplotlib.toolkits.basemap.pyproj.Geod.__new__[14]>", line 1, in ?
        for faz,baz,d in zip(az12,az21,dist): print "%7.3f %7.3f %9.3f" % (faz,baz,d)
    NameError: name 'az12' is not defined
Trying:
    from pyproj import Geod
Expecting nothing
ok
Trying:
    g = Geod(ellps='clrk66') # Use Clarke 1966 ellipsoid.
Expecting nothing
ok
Trying:
    boston_lat = 42.+(15./60.); boston_lon = -71.-(7./60.)
Expecting nothing
ok
Trying:
    portland_lat = 45.+(31./60.); portland_lon = -123.-(41./60.)
Expecting nothing
ok
Trying:
    lonlats = g.npts(boston_lon,boston_lat,portland_lon,portland_lat,10)
Expecting nothing
ok
Trying:
    for lon,lat in lonlats: print '%6.3f %7.3f' % (lat, lon)
Expecting:
    43.528 -75.414
    44.637 -79.883
    45.565 -84.512
    46.299 -89.279
    46.830 -94.156
    47.149 -99.112
    47.251 -104.106
    47.136 -109.100
    46.805 -114.051
    46.262 -118.924
ok
Trying:
    from pyproj import Proj
Expecting nothing
ok
Trying:
    p = Proj(proj='utm',zone=10,ellps='WGS84')
Expecting nothing
ok
Trying:
    x,y = p(-120.108, 34.36116666)
Expecting nothing
ok
Trying:
    print 'x=%9.3f y=%11.3f' % (x,y)
Expecting:
    x=765975.641 y=3805993.134
ok
Trying:
    print 'lon=%8.3f lat=%5.3f' % p(x,y,inverse=True)
Expecting:
    lon=-120.108 lat=34.361
ok
Trying:
    lons = (-119.72,-118.40,-122.38)
Expecting nothing
ok
Trying:
    lats = (36.77, 33.93, 37.62 )
Expecting nothing
ok
Trying:
    x,y = p(lons, lats)
Expecting nothing
ok
Trying:
    print 'x: %9.3f %9.3f %9.3f' % x
Expecting:
    x: 792763.863 925321.537 554714.301
ok
Trying:
    print 'y: %9.3f %9.3f %9.3f' % y
Expecting:
    y: 4074377.617 3763936.941 4163835.303
ok
Trying:
    lons, lats = p(x, y, inverse=True) # inverse transform
Expecting nothing
ok
Trying:
    print 'lons: %8.3f %8.3f %8.3f' % lons
Expecting:
    lons: -119.720 -118.400 -122.380
ok
Trying:
    print 'lats: %8.3f %8.3f %8.3f' % lats
Expecting:
    lats: 36.770 33.930 37.620
ok
Trying:
    p1 = Proj(init='epsg:26915')
Expecting nothing
ok
Trying:
    p2 = Proj(init='epsg:26715')
Expecting nothing
ok
Trying:
    x1, y1 = p1(-92.199881,38.56694)
Expecting nothing
ok
Trying:
    x2, y2 = transform(p1,p2,x1,y1)
Expecting nothing
ok
Trying:
    print '%9.3f %11.3f' % (x1,y1)
Expecting:
    569704.566 4269024.671
ok
Trying:
    print '%9.3f %11.3f' % (x2,y2)
Expecting:
    569706.333 4268817.680
ok
Trying:
    print '%8.3f %5.3f' % p2(x2,y2,inverse=True)
Expecting:
     -92.200 38.567
ok
Trying:
    lats = (38.83,39.32,38.75) # Columbia, KC and StL Missouri
Expecting nothing
ok
Trying:
    lons = (-92.22,-94.72,-90.37)
Expecting nothing
ok
Trying:
    x1, y1 = p1(lons,lats)
Expecting nothing
ok
Trying:
    x2, y2 = transform(p1,p2,x1,y1)
Expecting nothing
ok
Trying:
    xy = x1+y1
Expecting nothing
ok
Trying:
    print '%9.3f %9.3f %9.3f %11.3f %11.3f %11.3f' % xy
Expecting:
    567703.344 351730.944 728553.093 4298200.739 4353698.725 4292319.005
ok
Trying:
    xy = x2+y2
Expecting nothing
ok
Trying:
    print '%9.3f %9.3f %9.3f %11.3f %11.3f %11.3f' % xy
Expecting:
    567705.072 351727.113 728558.917 4297993.157 4353490.111 4292111.678
ok
Trying:
    lons, lats = p2(x2,y2,inverse=True)
Expecting nothing
ok
Trying:
    xy = lons+lats
Expecting nothing
ok
Trying:
    print '%8.3f %8.3f %8.3f %5.3f %5.3f %5.3f' % xy
Expecting:
     -92.220 -94.720 -90.370 38.830 39.320 38.750
ok
11 items had no tests:
    matplotlib.toolkits.basemap.pyproj
    matplotlib.toolkits.basemap.pyproj.Geod
    matplotlib.toolkits.basemap.pyproj.Geod.fwd
    matplotlib.toolkits.basemap.pyproj.Geod.inv
    matplotlib.toolkits.basemap.pyproj.Proj
    matplotlib.toolkits.basemap.pyproj.Proj.__call__
    matplotlib.toolkits.basemap.pyproj.Proj.is_geocent
    matplotlib.toolkits.basemap.pyproj.Proj.is_latlong
    matplotlib.toolkits.basemap.pyproj._convertback
    matplotlib.toolkits.basemap.pyproj._copytobuffer
    matplotlib.toolkits.basemap.pyproj.test
3 items passed all tests:
   6 tests in matplotlib.toolkits.basemap.pyproj.Geod.npts
  13 tests in matplotlib.toolkits.basemap.pyproj.Proj.__new__
  18 tests in matplotlib.toolkits.basemap.pyproj.transform
**********************************************************************
1 items had failures:
   6 of 15 in matplotlib.toolkits.basemap.pyproj.Geod.__new__
52 tests in 15 items.
46 passed and 6 failed.
***Test Failed*** 6 failures.

2007/9/12, Jeff Whitaker < jswhit@…196… <mailto:jswhit@…196…>>:

    David Huard wrote:
    > Hi, the pyproj package seems to cause a problem in the polarmap
    > example of the basemap toolkit.
    >
    > Thanks,
    >
    > david
    >
    > huardda@...559...:~/svnrepos/toolkits/basemap/examples$ python
    polarmaps.py
    > min/max etopo20 data:
    > -9026.625 6228.8125
    > plotting North Polar Lambert Azimuthal Equal Area basemap ...
    > plotting North Polar Stereographic basemap ...
    > plotting North Polar Azimuthal Equidistant basemap ...
    > Traceback (most recent call last):
    > File "polarmaps.py", line 51, in ?
    > resolution='c',area_thresh=10000.,lat_0=lat_0,lon_0=lon_0_ortho)
    > File
    >
    "/usr/local/lib/python2.4/site-packages/matplotlib/toolkits/basemap/basemap.py",

    > line 1111, in __init__
    >
    az1,alpha21,dist=gc.inv(lon_0,lat_0,math.radians(lonsnew[0]),math.radians(latsnew[0]),radians=True)
    > File
    >
    "/usr/local/lib/python2.4/site-packages/matplotlib/toolkits/basemap/pyproj.py",

    > line 478, in inv
    > _Geod._inv(self, inx, iny, inz, ind, radians=radians)
    > File "_geod.pyx", line 123, in _geod.Geod._inv
    > ValueError: undefined inverse geodesic (may be an antipodal point)
    >
    > This is from a fresh SVN version of both matplotlib and basemap.
    > Linux, Ubuntu edgy, Xeon-64.

    David: Odd - I can't reproduce that on my mac. Can you try this

    >>> from matplotlib.toolkits.basemap import pyproj
    >>> pyproj.test()

    and let me know if any of the tests fail?

    -Jeff

    –
    Jeffrey S. Whitaker Phone : (303)497-6313
    Meteorologist FAX : (303)497-6449
    NOAA/OAR/PSD R/PSD1 Email : Jeffrey.S.Whitaker@…236…
    <mailto:Jeffrey.S.Whitaker@…236…>
    325 Broadway Office : Skaggs Research Cntr 1D-124
    Boulder, CO, USA 80303-3328 Web : http://tinyurl.com/5telg

David: It might be a 64-bit issue. I'll try to find a 64-bit Xeon box to debug it, and let you know what I find.

-Jeff

···

--
Jeffrey S. Whitaker Phone : (303)497-6313
Meteorologist FAX : (303)497-6449
NOAA/OAR/PSD R/PSD1 Email : Jeffrey.S.Whitaker@...236...
325 Broadway Office : Skaggs Research Cntr 1D-124
Boulder, CO, USA 80303-3328 Web : http://tinyurl.com/5telg