[Matplotlib-users] basemap and fill_between()

2012-05-13 Thread David Craig
Hi, I'm having a problem usinf fill_between() with basemap. I plot two
great circles and want to shade the region between them. My code is below,
it doesnt give any error just creates the plot without filling the area.
Does anyone know if it's possible to do this or should I try a different
method?
Thanks,
David

from mpl_toolkits.basemap import Basemap
from pylab import *

### PARAMETERS FOR MATPLOTLIB :
import matplotlib as mpl
rcParams['font.size'] = 10.
rcParams['font.family'] = 'Comic Sans MS'
rcParams['axes.labelsize'] = 8.
rcParams['xtick.labelsize'] = 6.
rcParams['ytick.labelsize'] = 6.

def shoot(lon, lat, azimuth, maxdist=None):
Shooter Function
Original javascript on http://williams.best.vwh.net/gccalc.htm
Translated to python by Thomas Lecocq

glat1 = lat * pi / 180.
glon1 = lon * pi / 180.
s = maxdist / 1.852
faz = azimuth * pi / 180.

EPS= 0.005
if ((abs(cos(glat1))EPS) and not (abs(sin(faz))EPS)):
alert(Only N-S courses are meaningful, starting at a pole!)

a=6378.13/1.852
f=1/298.257223563
r = 1 - f
tu = r * tan(glat1)
sf = sin(faz)
cf = cos(faz)
if (cf==0):
b=0.
else:
b=2. * arctan2 (tu, cf)

cu = 1. / sqrt(1 + tu * tu)
su = tu * cu
sa = cu * sf
c2a = 1 - sa * sa
x = 1. + sqrt(1. + c2a * (1. / (r * r) - 1.))
x = (x - 2.) / x
c = 1. - x
c = (x * x / 4. + 1.) / c
d = (0.375 * x * x - 1.) * x
tu = s / (r * a * c)
y = tu
c = y + 1
while (abs (y - c)  EPS):

sy = sin(y)
cy = cos(y)
cz = cos(b + y)
e = 2. * cz * cz - 1.
c = y
x = e * cy
y = e + e - 1.
y = (((sy * sy * 4. - 3.) * y * cz * d / 6. + x) *
  d / 4. - cz) * sy * d + tu

b = cu * cy * cf - su * sy
c = r * sqrt(sa * sa + b * b)
d = su * cy + cu * sy * cf
glat2 = (arctan2(d, c) + pi) % (2*pi) - pi
c = cu * cy - su * sy * cf
x = arctan2(sy * sf, c)
c = ((-3. * c2a + 4.) * f + 4.) * c2a * f / 16.
d = ((e * cy * c + cz) * sy * c + y) * sa
glon2 = ((glon1 + x - (1. - c) * d * f + pi) % (2*pi)) - pi

baz = (arctan2(sa, b) + pi) % (2 * pi)

glon2 *= 180./pi
glat2 *= 180./pi
baz *= 180./pi

return (glon2, glat2, baz)

#Create a basemap around N. Atlantic
m = Basemap(llcrnrlon=-45.0,llcrnrlat=30.0,urcrnrlon=15.0,urcrnrlat=75.0,
resolution='i',projection='merc',lon_0=-17.5,lat_0=60.0)


m.drawcountries(linewidth=0.5)
m.drawcoastlines(linewidth=0.5)
m.bluemarble()
m.drawparallels(arange(40.,75.,10.),labels=[1,0,0,0],color='black',dashes=[1,0],labelstyle='+/-',linewidth=0.2)
# draw parallels
m.drawmeridians(arange(-45.,15.,10.),labels=[0,0,0,1],color='black',dashes=[1,0],labelstyle='+/-',linewidth=0.2)
# draw meridians

# Shade region defined by great circles.
x1, y1 = -9.1676613, 51.602
az1 = 270.
az2 = 290.
maxdist = 2000
x2, y2, baz = shoot(x1, y1, az1, maxdist)
x3, y3, baz = shoot(x1, y1, az2, maxdist)

m.drawgreatcircle(x1, y1, x2, y2, del_s=10, color='gray', lw=1.)
m.drawgreatcircle(x1, y1, x3, y3, del_s=10, color='gray', lw=1.)
a=linspace(x3,x1)
b=linspace(y2,y1)
c=linspace(y3,y1)
fill_between(a, b, c, where=None, alpha=0.2)
show()
--
Live Security Virtual Conference
Exclusive live event will cover all the ways today's security and 
threat landscape has changed and how IT managers can respond. Discussions 
will include endpoint security, mobile security and the latest in malware 
threats. http://www.accelacomm.com/jaw/sfrnl04242012/114/50122263/___
Matplotlib-users mailing list
Matplotlib-users@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/matplotlib-users


Re: [Matplotlib-users] basemap and fill_between()

2012-05-13 Thread Jerzy Karczmarczuk
Le 13/05/2012 11:34, David Craig a écrit :
 Hi, I'm having a problem usinf fill_between() with basemap. I plot two 
 great circles and want to shade the region between them. My code is 
 below, it doesnt give any error just creates the plot without filling 
 the area. Does anyone know if it's possible to do this or should I try 
 a different method?
Hello, Boss.
Perhaps I am telling rubbish, but fill_between demands the coordinates 
within your axes, not the basemap angles. Your a,b,c contain some small 
numers (negative x, say -38) won't do any good, x should be, say 3001129 
etc. Use m(...) to convert.

Jerzy K


--
Live Security Virtual Conference
Exclusive live event will cover all the ways today's security and 
threat landscape has changed and how IT managers can respond. Discussions 
will include endpoint security, mobile security and the latest in malware 
threats. http://www.accelacomm.com/jaw/sfrnl04242012/114/50122263/
___
Matplotlib-users mailing list
Matplotlib-users@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/matplotlib-users


Re: [Matplotlib-users] basemap and fill_between()

2012-05-13 Thread Jerzy Karczmarczuk
Appendix.
(and excuses for my approximate syntax in the first message; I was doing 
three things simultaneously).

Maestro David Craig, your filling will be lousy anyway.
If your x and y arrays come from linespace, you will surely get a 
triangle, and not the curvilinear area between GreatCircles. You must 
transform the coordinates in a more clever (but equally easy...) way.

Jerzy K.


--
Live Security Virtual Conference
Exclusive live event will cover all the ways today's security and 
threat landscape has changed and how IT managers can respond. Discussions 
will include endpoint security, mobile security and the latest in malware 
threats. http://www.accelacomm.com/jaw/sfrnl04242012/114/50122263/
___
Matplotlib-users mailing list
Matplotlib-users@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/matplotlib-users


Re: [Matplotlib-users] basemap and fill_between()

2012-05-13 Thread Jeff Whitaker

On 5/13/12 3:34 AM, David Craig wrote:
Hi, I'm having a problem usinf fill_between() with basemap. I plot two 
great circles and want to shade the region between them. My code is 
below, it doesnt give any error just creates the plot without filling 
the area. Does anyone know if it's possible to do this or should I try 
a different method?

Thanks,
David



David Try replacing

m.drawgreatcircle(x1, y1, x2, y2, del_s=10, color='gray', lw=1.)
m.drawgreatcircle(x1, y1, x3, y3, del_s=10, color='gray', lw=1.)
a=linspace(x3,x1)
b=linspace(y2,y1)
c=linspace(y3,y1)
fill_between(a, b, c, where=None, alpha=0.2)

with

xx1,yy1 = m.gcpoints(x1,y1,x2,y2,100)
xx2,yy2 = m.gcpoints(x1,y1,x3,y3,100)
fill_between(xx1, yy1, yy2)

-Jeff



from mpl_toolkits.basemap import Basemap
from pylab import *

### PARAMETERS FOR MATPLOTLIB :
import matplotlib as mpl
rcParams['font.size'] = 10.
rcParams['font.family'] = 'Comic Sans MS'
rcParams['axes.labelsize'] = 8.
rcParams['xtick.labelsize'] = 6.
rcParams['ytick.labelsize'] = 6.

def shoot(lon, lat, azimuth, maxdist=None):
Shooter Function
Original javascript on http://williams.best.vwh.net/gccalc.htm
Translated to python by Thomas Lecocq

glat1 = lat * pi / 180.
glon1 = lon * pi / 180.
s = maxdist / 1.852
faz = azimuth * pi / 180.

EPS= 0.005
if ((abs(cos(glat1))EPS) and not (abs(sin(faz))EPS)):
alert(Only N-S courses are meaningful, starting at a pole!)

a=6378.13/1.852
f=1/298.257223563
r = 1 - f
tu = r * tan(glat1)
sf = sin(faz)
cf = cos(faz)
if (cf==0):
b=0.
else:
b=2. * arctan2 (tu, cf)

cu = 1. / sqrt(1 + tu * tu)
su = tu * cu
sa = cu * sf
c2a = 1 - sa * sa
x = 1. + sqrt(1. + c2a * (1. / (r * r) - 1.))
x = (x - 2.) / x
c = 1. - x
c = (x * x / 4. + 1.) / c
d = (0.375 * x * x - 1.) * x
tu = s / (r * a * c)
y = tu
c = y + 1
while (abs (y - c)  EPS):

sy = sin(y)
cy = cos(y)
cz = cos(b + y)
e = 2. * cz * cz - 1.
c = y
x = e * cy
y = e + e - 1.
y = (((sy * sy * 4. - 3.) * y * cz * d / 6. + x) *
  d / 4. - cz) * sy * d + tu

b = cu * cy * cf - su * sy
c = r * sqrt(sa * sa + b * b)
d = su * cy + cu * sy * cf
glat2 = (arctan2(d, c) + pi) % (2*pi) - pi
c = cu * cy - su * sy * cf
x = arctan2(sy * sf, c)
c = ((-3. * c2a + 4.) * f + 4.) * c2a * f / 16.
d = ((e * cy * c + cz) * sy * c + y) * sa
glon2 = ((glon1 + x - (1. - c) * d * f + pi) % (2*pi)) - pi

baz = (arctan2(sa, b) + pi) % (2 * pi)

glon2 *= 180./pi
glat2 *= 180./pi
baz *= 180./pi

return (glon2, glat2, baz)

#Create a basemap around N. Atlantic
m = Basemap(llcrnrlon=-45.0,llcrnrlat=30.0,urcrnrlon=15.0,urcrnrlat=75.0,
resolution='i',projection='merc',lon_0=-17.5,lat_0=60.0)


m.drawcountries(linewidth=0.5)
m.drawcoastlines(linewidth=0.5)
m.bluemarble()
m.drawparallels(arange(40.,75.,10.),labels=[1,0,0,0],color='black',dashes=[1,0],labelstyle='+/-',linewidth=0.2) 
# draw parallels
m.drawmeridians(arange(-45.,15.,10.),labels=[0,0,0,1],color='black',dashes=[1,0],labelstyle='+/-',linewidth=0.2) 
# draw meridians


# Shade region defined by great circles.
x1, y1 = -9.1676613, 51.602
az1 = 270.
az2 = 290.
maxdist = 2000
x2, y2, baz = shoot(x1, y1, az1, maxdist)
x3, y3, baz = shoot(x1, y1, az2, maxdist)

m.drawgreatcircle(x1, y1, x2, y2, del_s=10, color='gray', lw=1.)
m.drawgreatcircle(x1, y1, x3, y3, del_s=10, color='gray', lw=1.)
a=linspace(x3,x1)
b=linspace(y2,y1)
c=linspace(y3,y1)
fill_between(a, b, c, where=None, alpha=0.2)
show()


--
Live Security Virtual Conference
Exclusive live event will cover all the ways today's security and
threat landscape has changed and how IT managers can respond. Discussions
will include endpoint security, mobile security and the latest in malware
threats. http://www.accelacomm.com/jaw/sfrnl04242012/114/50122263/


___
Matplotlib-users mailing list
Matplotlib-users@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/matplotlib-users


--
Live Security Virtual Conference
Exclusive live event will cover all the ways today's security and 
threat landscape has changed and how IT managers can respond. Discussions 
will include endpoint security, mobile security and the latest in malware 
threats. http://www.accelacomm.com/jaw/sfrnl04242012/114/50122263/___
Matplotlib-users mailing list
Matplotlib-users@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/matplotlib-users