Revision: 7571
http://matplotlib.svn.sourceforge.net/matplotlib/?rev=7571&view=rev
Author: jswhit
Date: 2009-08-26 01:12:32 +0000 (Wed, 26 Aug 2009)
Log Message:
-----------
added nightshade method
Modified Paths:
--------------
trunk/toolkits/basemap/Changelog
trunk/toolkits/basemap/examples/daynight.py
trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py
trunk/toolkits/basemap/setup.py
Added Paths:
-----------
trunk/toolkits/basemap/lib/mpl_toolkits/basemap/solar.py
Modified: trunk/toolkits/basemap/Changelog
===================================================================
--- trunk/toolkits/basemap/Changelog 2009-08-25 20:06:13 UTC (rev 7570)
+++ trunk/toolkits/basemap/Changelog 2009-08-26 01:12:32 UTC (rev 7571)
@@ -1,3 +1,6 @@
+version 0.99.5 (not yet released)
+ * added 'nightshade' method to shade night regions on a map.
+ * added lonmin, lonmax instance variables.
version 0.99.4 (svn revision 7332)
* ax.frame replaced with ax.spines to maintain compatibility
with matplotlib spines support.
Modified: trunk/toolkits/basemap/examples/daynight.py
===================================================================
--- trunk/toolkits/basemap/examples/daynight.py 2009-08-25 20:06:13 UTC (rev
7570)
+++ trunk/toolkits/basemap/examples/daynight.py 2009-08-26 01:12:32 UTC (rev
7571)
@@ -1,93 +1,16 @@
import numpy as np
-from mpl_toolkits.basemap import Basemap, netcdftime
+from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
from datetime import datetime
-from numpy import ma
# example showing how to compute the day/night terminator and shade nightime
# areas on a map.
-
-def epem(date):
- """
- input: date - datetime object (assumed UTC)
- ouput: gha - Greenwich hour angle, the angle between the Greenwich
- meridian and the meridian containing the subsolar point.
- dec - solar declination.
- """
- dg2rad = np.pi/180.
- rad2dg = 1./dg2rad
- # compute julian day from UTC datetime object.
- jday = netcdftime.JulianDayFromDate(date)
- jd = np.floor(jday) # truncate to integer.
- # utc hour.
- ut = date.hour + date.minute/60. + date.second/3600.
- # calculate number of centuries from J2000
- t = (jd + (ut/24.) - 2451545.0) / 36525.
- # mean longitude corrected for aberration
- l = (280.460 + 36000.770 * t) % 360
- # mean anomaly
- g = 357.528 + 35999.050 * t
- # ecliptic longitude
- lm = l + 1.915 * np.sin(g*dg2rad) + 0.020 * np.sin(2*g*dg2rad)
- # obliquity of the ecliptic
- ep = 23.4393 - 0.01300 * t
- # equation of time
- eqtime = -1.915*np.sin(g*dg2rad) - 0.020*np.sin(2*g*dg2rad) \
- + 2.466*np.sin(2*lm*dg2rad) - 0.053*np.sin(4*lm*dg2rad)
- # Greenwich hour angle
- gha = 15*ut - 180 + eqtime
- # declination of sun
- dec = np.arcsin(np.sin(ep*dg2rad) * np.sin(lm*dg2rad)) * rad2dg
- return gha, dec
-
-def daynightgrid(date, nlons):
- """
- date is datetime object (assumed UTC).
- nlons is # of longitudes used to compute terminator."""
- dg2rad = np.pi/180.
- lons = np.linspace(-180,180,nlons).astype(np.float32)
- # compute greenwich hour angle and solar declination
- # from datetime object (assumed UTC).
- tau, dec = epem(date)
- # compute day/night terminator from hour angle, declination.
- longitude = lons + tau
- lats = np.arctan(-np.cos(longitude*dg2rad)/np.tan(dec*dg2rad))/dg2rad
- # create day/night grid (1 for night, 0 for day)
- nlats = ((nlons-1)/2)+1
- lons2 = np.linspace(-180,180,nlons).astype(np.float32)
- lats2 = np.linspace(-90,90,nlats).astype(np.float32)
- lons2, lats2 = np.meshgrid(lons2,lats2)
- lats = lats[np.newaxis,:]*np.ones((nlats,nlons),dtype=np.float32)
- daynight = np.ones(lons2.shape, np.int8)
- daynight = np.where(lats2>lats,0,daynight)
- daynight = ma.array(daynight,mask=1-daynight) # mask day areas.
- return lons2,lats2,daynight
-
-class Basemap2(Basemap):
- def nightshade(self,date,color="k",nlons=1441,alpha=0.5,zorder=None):
- # create grid of day=0, night=1
- # 1441 means use a 0.25 degree grid.
- lons,lats,daynight = daynightgrid(date,nlons)
- x,y = self(lons,lats)
- # contour the day-night grid, coloring the night area
- # with the specified color and transparency.
- CS = map.contourf(x,y,daynight,1,colors=[color],alpha=alpha)
- # set zorder on ContourSet collections show night shading
- # is on top.
- for c in CS.collections:
- if zorder is None:
- c.set_zorder(c.get_zorder()+1)
- else:
- c.set_zorder(zorder)
- return CS
-
-
# miller projection
-map = Basemap2(projection='mill',lon_0=0)
+map = Basemap(projection='mill',lon_0=180)
# plot coastlines, draw label meridians and parallels.
map.drawcoastlines()
map.drawparallels(np.arange(-90,90,30),labels=[1,0,0,0])
-map.drawmeridians(np.arange(-180,180,60),labels=[0,0,0,1])
+map.drawmeridians(np.arange(map.lonmin,map.lonmax+30,60),labels=[0,0,0,1])
# fill continents 'coral' (with zorder=0), color wet areas 'aqua'
map.drawmapboundary(fill_color='aqua')
map.fillcontinents(color='coral',lake_color='aqua')
Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py
===================================================================
--- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2009-08-25
20:06:13 UTC (rev 7570)
+++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2009-08-26
01:12:32 UTC (rev 7571)
@@ -53,7 +53,7 @@
else:
basemap_datadir = os.sep.join([os.path.dirname(__file__), 'data'])
-__version__ = '0.99.4'
+__version__ = '0.99.5'
# supported map projections.
_projnames = {'cyl' : 'Cylindrical Equidistant',
@@ -3421,6 +3421,42 @@
raise KeyError("barstyle must be 'simple' or 'fancy'")
return rets
+ def nightshade(self,date,color="k",delta=0.25,alpha=0.5,ax=None,zorder=2):
+ """
+ Shade the regions of the map that are in darkness at the time
+ specifed by ``date``. ``date`` is a datetime instance,
+ assumed to be UTC.
+
+ .. tabularcolumns:: |l|L|
+
+ ============== ====================================================
+ Keywords Description
+ ============== ====================================================
+ color color to shade night regions (default black).
+ delta day/night terminator is computed with a
+ a resolution of ``delta`` degrees (default 0.25).
+ alpha alpha transparency for shading (default 0.5, so
+ map background shows through).
+ zorder zorder for shading (default 2).
+ ============== ====================================================
+
+ Extra keyword ``ax`` can be used to override the default axis instance.
+
+ returns a matplotlib.contour.ContourSet instance.
+ """
+ from solar import daynight_grid
+ # create grid of day=0, night=1
+ lons,lats,daynight = daynight_grid(date,delta,self.lonmin,self.lonmax)
+ x,y = self(lons,lats)
+ # contour the day-night grid, coloring the night area
+ # with the specified color and transparency.
+ CS = self.contourf(x,y,daynight,1,colors=[color],alpha=alpha,ax=ax)
+ # set zorder on ContourSet collections show night shading
+ # is on top.
+ for c in CS.collections:
+ c.set_zorder(zorder)
+ return CS
+
def _check_ax(self):
"""
Returns the axis on which to draw.
Added: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/solar.py
===================================================================
--- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/solar.py
(rev 0)
+++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/solar.py 2009-08-26
01:12:32 UTC (rev 7571)
@@ -0,0 +1,66 @@
+# some simple functions to calculate solar position, day-night terminator
+import numpy as np
+from numpy import ma
+import netcdftime
+
+def epem(date):
+ """
+ input: date - datetime object (assumed UTC)
+ ouput: gha - Greenwich hour angle, the angle between the Greenwich
+ meridian and the meridian containing the subsolar point.
+ dec - solar declination.
+ """
+ dg2rad = np.pi/180.
+ rad2dg = 1./dg2rad
+ # compute julian day from UTC datetime object.
+ jday = netcdftime.JulianDayFromDate(date)
+ jd = np.floor(jday) # truncate to integer.
+ # utc hour.
+ ut = date.hour + date.minute/60. + date.second/3600.
+ # calculate number of centuries from J2000
+ t = (jd + (ut/24.) - 2451545.0) / 36525.
+ # mean longitude corrected for aberration
+ l = (280.460 + 36000.770 * t) % 360
+ # mean anomaly
+ g = 357.528 + 35999.050 * t
+ # ecliptic longitude
+ lm = l + 1.915 * np.sin(g*dg2rad) + 0.020 * np.sin(2*g*dg2rad)
+ # obliquity of the ecliptic
+ ep = 23.4393 - 0.01300 * t
+ # equation of time
+ eqtime = -1.915*np.sin(g*dg2rad) - 0.020*np.sin(2*g*dg2rad) \
+ + 2.466*np.sin(2*lm*dg2rad) - 0.053*np.sin(4*lm*dg2rad)
+ # Greenwich hour angle
+ gha = 15*ut - 180 + eqtime
+ # declination of sun
+ dec = np.arcsin(np.sin(ep*dg2rad) * np.sin(lm*dg2rad)) * rad2dg
+ return gha, dec
+
+def daynight_terminator(date, delta, lonmin, lonmax):
+ """
+ date is datetime object (assumed UTC).
+ nlons is # of longitudes used to compute terminator."""
+ dg2rad = np.pi/180.
+ lons = np.arange(lonmin,lonmax+0.5*delta,delta,dtype=np.float32)
+ # compute greenwich hour angle and solar declination
+ # from datetime object (assumed UTC).
+ tau, dec = epem(date)
+ # compute day/night terminator from hour angle, declination.
+ longitude = lons + tau
+ lats = np.arctan(-np.cos(longitude*dg2rad)/np.tan(dec*dg2rad))/dg2rad
+ return lons, lats
+
+def daynight_grid(date, delta, lonmin, lonmax):
+ """
+ date is datetime object (assumed UTC).
+ delta is the grid interval (in degrees) used to compute terminator."""
+ lons, lats = daynight_terminator(date, delta, lonmin, lonmax)
+ # create day/night grid (1 for night, 0 for day)
+ lats2 = np.arange(-90,90+0.5*delta,delta,dtype=np.float32)
+ nlons = len(lons); nlats = len(lats2)
+ lons2, lats2 = np.meshgrid(lons,lats2)
+ lats = lats[np.newaxis,:]*np.ones((nlats,nlons),dtype=np.float32)
+ daynight = np.ones(lons2.shape, np.int8)
+ daynight = np.where(lats2>lats,0,daynight)
+ daynight = ma.array(daynight,mask=1-daynight) # mask day areas.
+ return lons2,lats2,daynight
Modified: trunk/toolkits/basemap/setup.py
===================================================================
--- trunk/toolkits/basemap/setup.py 2009-08-25 20:06:13 UTC (rev 7570)
+++ trunk/toolkits/basemap/setup.py 2009-08-26 01:12:32 UTC (rev 7571)
@@ -225,7 +225,7 @@
package_data = {'mpl_toolkits.basemap':pyproj_datafiles+basemap_datafiles}
setup(
name = "basemap",
- version = "0.99.4",
+ version = "0.99.5",
description = "Plot data on map projections with matplotlib",
long_description = """
An add-on toolkit for matplotlib that lets you plot data
This was sent by the SourceForge.net collaborative development platform, the
world's largest Open Source development site.
------------------------------------------------------------------------------
Let Crystal Reports handle the reporting - Free Crystal Reports 2008 30-Day
trial. Simplify your report design, integration and deployment - and focus on
what you do best, core application coding. Discover what's new with
Crystal Reports now. http://p.sf.net/sfu/bobj-july
_______________________________________________
Matplotlib-checkins mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/matplotlib-checkins
