Problem projecting 2d arrays of coordinates..

[guziy]

I know we have a lot of issues already, but this one is intriguing and maybe gives an opportunity to learn more about cython...

Basemap works fine in the following situations (below lon2d and lat2d are 2D numpy arrays of latitudes and longitudes):

a)

xx, yy = basemap(lons2d.T, lats2d.T)

b)

xx, yy = basemap(lons2d.flatten(), lats2d.flatten())

c) Looping through each point of the arrays and projecting works fine as well...

BUT

when I try to reproject the unchanged arrays

 xx, yy = basemap(lons2d, lats2d)

I get the exception below:

In [4]: xx, yy = bmp(lons, lats)
---------------------------------------------------------------------------
RuntimeError                              Traceback (most recent call last)
<ipython-input-4-cd4fecf35ebb> in <module>()
----> 1 xx, yy = bmp(lons, lats)

/RESCUE/skynet3_rech1/huziy/python_virtualenvs/py3.4/lib/python3.4/site-packages/mpl_toolkits/basemap/__init__.py in __call__(self, x, y, inverse)
   1146             except TypeError:
   1147                 y = [_dg2rad*yy for yy in y]
-> 1148         xout,yout = self.projtran(x,y,inverse=inverse)
   1149         if self.celestial and inverse:
   1150             try:

/RESCUE/skynet3_rech1/huziy/python_virtualenvs/py3.4/lib/python3.4/site-packages/mpl_toolkits/basemap/proj.py in __call__(self, *args, **kw)
    284             outxy = self._proj4(xy, inverse=inverse)
    285         else:
--> 286             outx,outy = self._proj4(x, y, inverse=inverse)
    287         if inverse:
    288             if self.projection in ['merc','mill','gall']:

/RESCUE/skynet3_rech1/huziy/python_virtualenvs/py3.4/lib/python3.4/site-packages/mpl_toolkits/basemap/pyproj.py in __call__(self, *args, **kw)
    386             _proj.Proj._inv(self, inx, iny, radians=radians, errcheck=errcheck)
    387         else:
--> 388             _proj.Proj._fwd(self, inx, iny, radians=radians, errcheck=errcheck)
    389         # if inputs were lists, tuples or floats, convert back.
    390         outx = _convertback(xisfloat,xislist,xistuple,inx)

_proj.pyx in _proj.Proj._fwd (src/_proj.c:1571)()

RuntimeError:

I am using Python 3.4.3, Basemap 1.0.7

I've pickled the list [lons, lats, basemap] and attached the file here if anoyone wants to try to reproduce (you can download it from here: http://scaweb.sca.uqam.ca/~huziy/lon_lat_bmp_problem.bin).

and then use it as

import pickle
lons, lats, bmp = pickle.load(open("lon_lat_bmp_problem.bin","rb"))
xx, yy = bmp(lons, lats)  # <-- should throw an exception

Cheers

[QuLogic]

Not that I don't trust you, but pickles can be dangerous. You should output a NumPy file with np.savez_compressed if you want people to look at it.

[guziy]

Hi:

Here is the npz file: http://scaweb.sca.uqam.ca/~huziy/test.npz
It can be used as follows:

import numpy as np
lons, lats, bmp = np.load("test.npz")["arr_0"]
bmp(lons, lats)  # <-- Should cause the exception 

PS: Would numpy prevent me from saving to npz if calling bmp was doing something evil?

Cheers

[QuLogic]

PS: Would numpy prevent me from saving to npz if calling bmp was doing something evil?

The danger in pickles doesn't have anything to do with Basemap. Pickles allow arbitrary code execution at load time, and are not necessarily portable across versions. Plus, NumPy's file is smaller.

[QuLogic]

PS: Would numpy prevent me from saving to npz if calling bmp was doing something evil?

Ah, sorry, I didn't see you were referring to the Basemap object also saved there. For objects, I think NumPy ends up deferring to a pickle in the end, so no, it wouldn't. Note though that you can save them as separate arrays (np.savez(lons=lons, lats=lats, bmp=bmp)) so that if someone were especially wary, they would only have to load the arrays.

[guziy]

Thanks for the clarification))

Cheers