From: <js...@us...> - 2008-07-31 18:39:00
|
Revision: 5938 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=5938&view=rev Author: jswhit Date: 2008-07-31 18:38:56 +0000 (Thu, 31 Jul 2008) Log Message: ----------- use 'import numpy as np' Modified Paths: -------------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/proj.py Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py 2008-07-31 15:37:56 UTC (rev 5937) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py 2008-07-31 18:38:56 UTC (rev 5938) @@ -1,7 +1,8 @@ """ Performs conversions of netCDF time coordinate data to/from datetime objects. """ -import math, numpy, re, time +import math, re, time +import numpy as np from datetime import datetime as real_datetime _units = ['days','hours','minutes','seconds','day','hour','minute','second'] @@ -629,7 +630,7 @@ except: isscalar = True if not isscalar: - date = numpy.array(date) + date = np.array(date) shape = date.shape if self.calendar in ['julian','standard','gregorian','proleptic_gregorian']: if isscalar: @@ -656,7 +657,7 @@ else: jdelta = [_360DayFromDate(d)-self._jd0 for d in date.flat] if not isscalar: - jdelta = numpy.array(jdelta) + jdelta = np.array(jdelta) # convert to desired units, add time zone offset. if self.units in ['second','seconds']: jdelta = jdelta*86400. + self.tzoffset*60. @@ -669,7 +670,7 @@ if isscalar: return jdelta else: - return numpy.reshape(jdelta,shape) + return np.reshape(jdelta,shape) def num2date(self,time_value): """ @@ -681,8 +682,8 @@ Resolution is 1 second. -Works for scalars, sequences and numpy arrays. -Returns a scalar if input is a scalar, else returns a numpy array. +Works for scalars, sequences and np arrays. +Returns a scalar if input is a scalar, else returns a np array. The datetime instances returned by C{num2date} are 'real' python datetime objects if the date falls in the Gregorian calendar (i.e. @@ -699,7 +700,7 @@ except: isscalar = True if not isscalar: - time_value = numpy.array(time_value) + time_value = np.array(time_value) shape = time_value.shape # convert to desired units, remove time zone offset. if self.units in ['second','seconds']: @@ -734,7 +735,7 @@ if isscalar: return date else: - return numpy.reshape(numpy.array(date),shape) + return np.reshape(np.array(date),shape) def _parse_date(origin): """Parses a date string and returns a tuple Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/proj.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/proj.py 2008-07-31 15:37:56 UTC (rev 5937) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/proj.py 2008-07-31 18:38:56 UTC (rev 5938) @@ -1,4 +1,4 @@ -import numpy as npy +import numpy as np import pyproj import math @@ -92,15 +92,15 @@ self._proj4 = pyproj.Proj(projparams) # find major and minor axes of ellipse defining map proj region. delta = 0.01 - lats = npy.arange(0,90,delta) + lats = np.arange(0,90,delta) lon_0 = projparams['lon_0'] - lons = lon_0*npy.ones(len(lats),'d') + lons = lon_0*np.ones(len(lats),'d') x, y = self._proj4(lons, lats) yi = (y > 1.e20).tolist() ny = yi.index(1)-1 height = y[ny] - lons = npy.arange(lon_0,lon_0+90,delta) - lats = npy.zeros(len(lons),'d') + lons = np.arange(lon_0,lon_0+90,delta) + lats = np.zeros(len(lons),'d') x, y = self(lons, lats) xi = (x > 1.e20).tolist() nx = xi.index(1)-1 @@ -264,8 +264,8 @@ """ dx = (self.urcrnrx-self.llcrnrx)/(nx-1) dy = (self.urcrnry-self.llcrnry)/(ny-1) - x = self.llcrnrx+dx*npy.indices((ny,nx),npy.float32)[1,:,:] - y = self.llcrnry+dy*npy.indices((ny,nx),npy.float32)[0,:,:] + x = self.llcrnrx+dx*np.indices((ny,nx),np.float32)[1,:,:] + y = self.llcrnry+dy*np.indices((ny,nx),np.float32)[0,:,:] lons, lats = self(x, y, inverse=True) if returnxy: return lons, lats, x, y @@ -280,9 +280,9 @@ """ dx = (self.urcrnrx-self.llcrnrx)/(nx-1) dy = (self.urcrnry-self.llcrnry)/(ny-1) - xy = npy.empty((ny,nx,2), npy.float64) - xy[...,0] = self.llcrnrx+dx*npy.indices((ny,nx),npy.float32)[1,:,:] - xy[...,1] = self.llcrnry+dy*npy.indices((ny,nx),npy.float32)[0,:,:] + xy = np.empty((ny,nx,2), np.float64) + xy[...,0] = self.llcrnrx+dx*np.indices((ny,nx),np.float32)[1,:,:] + xy[...,1] = self.llcrnry+dy*np.indices((ny,nx),np.float32)[0,:,:] lonlat = self(xy, inverse=True) if returnxy: return lonlat, xy @@ -329,9 +329,9 @@ t2 = time.clock() print 'compute lats/lons for all points on AWIPS 221 grid (%sx%s)' %(nx,ny) print 'max/min lons' - print min(npy.ravel(lons)),max(npy.ravel(lons)) + print min(np.ravel(lons)),max(np.ravel(lons)) print 'max/min lats' - print min(npy.ravel(lats)),max(npy.ravel(lats)) + print min(np.ravel(lats)),max(np.ravel(lats)) print 'took',t2-t1,'secs' print 'Same thing but with a single 3-D array' t1 = time.clock() This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2008-08-08 20:22:58
|
Revision: 6007 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=6007&view=rev Author: jswhit Date: 2008-08-08 20:22:50 +0000 (Fri, 08 Aug 2008) Log Message: ----------- update pupynere to 1.0.1 Modified Paths: -------------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py Added Paths: ----------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2008-08-08 12:19:47 UTC (rev 6006) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2008-08-08 20:22:50 UTC (rev 6007) @@ -36,7 +36,7 @@ import numpy as np import numpy.ma as ma from shapelib import ShapeFile -import _geoslib, pupynere, netcdftime +import _geoslib, netcdf, netcdftime # basemap data files now installed in lib/matplotlib/toolkits/basemap/data basemap_datadir = os.sep.join([os.path.dirname(__file__), 'data']) @@ -3640,7 +3640,7 @@ else: return corners -def NetCDFFile(file, maskandscale=True): +def NetCDFFile(file, mode='r', maskandscale=True): """NetCDF File reader. API is the same as Scientific.IO.NetCDF. If ``file`` is a URL that starts with `http`, it is assumed to be a remote OPenDAP dataset, and the python dap client is used @@ -3656,9 +3656,9 @@ ``maskandscale`` keyword to False. """ if file.startswith('http'): - return pupynere._RemoteFile(file,maskandscale) + return netcdf._RemoteFile(file,maskandscale=maskandscale) else: - return pupynere._LocalFile(file,maskandscale) + return netcdf.netcdf_file(file,mode=mode,maskandscale=maskandscale) def num2date(times,units='days since 0001-01-01 00:00:00',calendar='proleptic_gregorian'): """ Added: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py (rev 0) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py 2008-08-08 20:22:50 UTC (rev 6007) @@ -0,0 +1,76 @@ +from numpy import ma, squeeze +from pupynere import netcdf_file, _maskandscale +from dap.client import open as open_remote +from dap.dtypes import ArrayType, GridType, typemap + +_typecodes = dict([[_v,_k] for _k,_v in typemap.items()]) + +class _RemoteFile(object): + """A NetCDF file reader. API is the same as Scientific.IO.NetCDF.""" + + def __init__(self, file, maskandscale=False): + self._buffer = open_remote(file) + self._maskandscale = maskandscale + self._parse() + + def read(self, size=-1): + """Alias for reading the file buffer.""" + return self._buffer.read(size) + + def _parse(self): + """Initial parsing of the header.""" + # Read header info. + self._dim_array() + self._gatt_array() + self._var_array() + + def _dim_array(self): + """Read a dict with dimensions names and sizes.""" + self.dimensions = {} + self._dims = [] + for k,d in self._buffer.iteritems(): + if (isinstance(d, ArrayType) or isinstance(d, GridType)) and len(d.shape) == 1 and k == d.dimensions[0]: + name = k + length = len(d) + self.dimensions[name] = length + self._dims.append(name) # preserve dim order + + def _gatt_array(self): + """Read global attributes.""" + self.__dict__.update(self._buffer.attributes) + + def _var_array(self): + """Read all variables.""" + # Read variables. + self.variables = {} + for k,d in self._buffer.iteritems(): + if isinstance(d, GridType) or isinstance(d, ArrayType): + name = k + self.variables[name] = _RemoteVariable(d,self._maskandscale) + + def close(self): + # this is a no-op provided for compatibility + pass + +class _RemoteVariable(object): + def __init__(self, var, maskandscale): + self._var = var + self._maskandscale = maskandscale + self.dtype = var.type + self.shape = var.shape + self.dimensions = var.dimensions + self.__dict__.update(var.attributes) + + def __getitem__(self, index): + datout = squeeze(self._var.__getitem__(index)) + # automatically + # - remove singleton dimensions + # - create a masked array using missing_value or _FillValue attribute + # - apply scale_factor and add_offset to packed integer data + if self._maskandscale: + return _maskandscale(self,datout) + else: + return datout + + def typecode(self): + return _typecodes[self.dtype] Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py 2008-08-08 12:19:47 UTC (rev 6006) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py 2008-08-08 20:22:50 UTC (rev 6007) @@ -1,45 +1,89 @@ -"""NetCDF reader. +""" +NetCDF reader/writer module. -Pupynere implements a PUre PYthon NEtcdf REader. +This module implements the Scientific.IO.NetCDF API to read and create +NetCDF files. The same API is also used in the PyNIO and pynetcdf +modules, allowing these modules to be used interchangebly when working +with NetCDF files. The major advantage of ``scipy.io.netcdf`` over other +modules is that it doesn't require the code to be linked to the NetCDF +libraries as the other modules do. -Copyright (c) 2003-2006 Roberto De Almeida <ro...@py...> +The code is based on the `NetCDF file format specification +<http://www.unidata.ucar.edu/software/netcdf/guide_15.html>`_. A NetCDF +file is a self-describing binary format, with a header followed by +data. The header contains metadata describing dimensions, variables +and the position of the data in the file, so access can be done in an +efficient manner without loading unnecessary data into memory. We use +the ``mmap`` module to create Numpy arrays mapped to the data on disk, +for the same purpose. -Permission is hereby granted, free of charge, to any person obtaining -a copy of this software and associated documentation files (the -"Software"), to deal in the Software without restriction, including -without limitation the rights to use, copy, modify, merge, publish, -distribute, sublicense, and/or sell copies of the Software, and to -permit persons to whom the Software is furnished to do so, subject -to the following conditions: +The structure of a NetCDF file is as follows: -The above copyright notice and this permission notice shall be -included in all copies or substantial portions of the Software. + C D F <VERSION BYTE> <NUMBER OF RECORDS> + <DIMENSIONS> <GLOBAL ATTRIBUTES> <VARIABLES METADATA> + <NON-RECORD DATA> <RECORD DATA> -THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, -EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF -MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. -IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR -ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF -CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION -WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. +Record data refers to data where the first axis can be expanded at +will. All record variables share a same dimension at the first axis, +and they are stored at the end of the file per record, ie + + A[0], B[0], ..., A[1], B[1], ..., etc, + +so that new data can be appended to the file without changing its original +structure. Non-record data are padded to a 4n bytes boundary. Record data +are also padded, unless there is exactly one record variable in the file, +in which case the padding is dropped. All data is stored in big endian +byte order. + +The Scientific.IO.NetCDF API allows attributes to be added directly to +instances of ``netcdf_file`` and ``netcdf_variable``. To differentiate +between user-set attributes and instance attributes, user-set attributes +are automatically stored in the ``_attributes`` attribute by overloading +``__setattr__``. This is the reason why the code sometimes uses +``obj.__dict__['key'] = value``, instead of simply ``obj.key = value``; +otherwise the key would be inserted into userspace attributes. + +To create a NetCDF file:: + + >>> import time + >>> f = netcdf_file('simple.nc', 'w') + >>> f.history = 'Created for a test' + >>> f.createDimension('time', 10) + >>> time = f.createVariable('time', 'i', ('time',)) + >>> time[:] = range(10) + >>> time.units = 'days since 2008-01-01' + >>> f.close() + +To read the NetCDF file we just created:: + + >>> f = netcdf_file('simple.nc', 'r') + >>> print f.history + Created for a test + >>> time = f.variables['time'] + >>> print time.units + days since 2008-01-01 + >>> print time.shape + (10,) + >>> print time[-1] + 9 + >>> f.close() + +TODO: properly implement ``_FillValue``. """ -__author__ = "Roberto De Almeida <ro...@py...>" +__all__ = ['netcdf_file', 'netcdf_variable'] -import struct -import itertools -import mmap +from operator import mul +from mmap import mmap, ACCESS_READ -from numpy import ndarray, empty, array, ma, squeeze, zeros -import numpy +from numpy import fromstring, ndarray, dtype, empty, array, asarray, squeeze, zeros, ma +from numpy import little_endian as LITTLE_ENDIAN -from dap.client import open as open_remote -from dap.dtypes import ArrayType, GridType, typemap -ABSENT = '\x00' * 8 -ZERO = '\x00' * 4 -NC_BYTE = '\x00\x00\x00\x01' +ABSENT = '\x00\x00\x00\x00\x00\x00\x00\x00' +ZERO = '\x00\x00\x00\x00' +NC_BYTE = '\x00\x00\x00\x01' NC_CHAR = '\x00\x00\x00\x02' NC_SHORT = '\x00\x00\x00\x03' NC_INT = '\x00\x00\x00\x04' @@ -49,360 +93,531 @@ NC_VARIABLE = '\x00\x00\x00\x0b' NC_ATTRIBUTE = '\x00\x00\x00\x0c' -_typecodes = dict([[_v,_k] for _k,_v in typemap.items()]) -# default _FillValue for netcdf types (apply also to corresponding -# DAP types). -_default_fillvals = {'c':'\0', - 'S':"", - 'b':-127, - 'B':-127, - 'h':-32767, - 'H':65535, - 'i':-2147483647L, - 'L':4294967295L, - 'q':-2147483647L, - 'f':9.9692099683868690e+36, - 'd':9.9692099683868690e+36} -def NetCDFFile(file, maskandscale=True): - """NetCDF File reader. API is the same as Scientific.IO.NetCDF. - If 'file' is a URL that starts with 'http', it is assumed - to be a remote OPenDAP dataset, and the python dap client is used - to retrieve the data. Only the OPenDAP Array and Grid data - types are recognized. If file does not start with 'http', it - is assumed to be a local netCDF file. Data will - automatically be converted to masked arrays if the variable has either - a 'missing_value' or '_FillValue' attribute, and some data points - are equal to the value specified by that attribute. In addition, - variables stored as integers that have the 'scale_factor' and - 'add_offset' attribute will automatically be rescaled to floats when - read. To suppress these automatic conversions, set the - maskandscale keyword to False. +TYPEMAP = { NC_BYTE: ('b', 1), + NC_CHAR: ('c', 1), + NC_SHORT: ('h', 2), + NC_INT: ('i', 4), + NC_FLOAT: ('f', 4), + NC_DOUBLE: ('d', 8) } + +REVERSE = { 'b': NC_BYTE, + 'c': NC_CHAR, + 'h': NC_SHORT, + 'i': NC_INT, + 'f': NC_FLOAT, + 'd': NC_DOUBLE, + + # these come from asarray(1).dtype.char and asarray('foo').dtype.char, + # used when getting the types from generic attributes. + 'l': NC_INT, + 'S': NC_CHAR } + + +class netcdf_file(object): """ - if file.startswith('http'): - return _RemoteFile(file,maskandscale) - else: - return _LocalFile(file,maskandscale) - -def _maskandscale(var,datout): - totalmask = zeros(datout.shape,numpy.bool) - fillval = None - if hasattr(var, 'missing_value') and (datout == var.missing_value).any(): - fillval = var.missing_value - totalmask += datout==fillval - if hasattr(var, '_FillValue') and (datout == var._FillValue).any(): - if fillval is None: - fillval = var._FillValue - totalmask += datout==var._FillValue - elif (datout == _default_fillvals[var.typecode()]).any(): - if fillval is None: - fillval = _default_fillvals[var.typecode()] - totalmask += datout==_default_fillvals[var.dtype] - if fillval is not None: - datout = ma.masked_array(datout,mask=totalmask,fill_value=fillval) - try: - datout = var.scale_factor*datout + var.add_offset - except: - pass - return datout + A ``netcdf_file`` object has two standard attributes: ``dimensions`` and + ``variables``. The values of both are dictionaries, mapping dimension + names to their associated lengths and variable names to variables, + respectively. Application programs should never modify these + dictionaries. -class _RemoteFile(object): - """A NetCDF file reader. API is the same as Scientific.IO.NetCDF.""" + All other attributes correspond to global attributes defined in the + NetCDF file. Global file attributes are created by assigning to an + attribute of the ``netcdf_file`` object. - def __init__(self, file, maskandscale): - self._buffer = open_remote(file) + """ + def __init__(self, filename, mode='r', maskandscale=False): + self._maskandscale = maskandscale - self._parse() - def read(self, size=-1): - """Alias for reading the file buffer.""" - return self._buffer.read(size) + self.filename = filename - def _parse(self): - """Initial parsing of the header.""" - # Read header info. - self._dim_array() - self._gatt_array() - self._var_array() + assert mode in 'rw', "Mode must be either 'r' or 'w'." + self.mode = mode - def _dim_array(self): - """Read a dict with dimensions names and sizes.""" self.dimensions = {} + self.variables = {} + self._dims = [] - for k,d in self._buffer.iteritems(): - if (isinstance(d, ArrayType) or isinstance(d, GridType)) and len(d.shape) == 1 and k == d.dimensions[0]: - name = k - length = len(d) - self.dimensions[name] = length - self._dims.append(name) # preserve dim order + self._recs = 0 + self._recsize = 0 - def _gatt_array(self): - """Read global attributes.""" - self.__dict__.update(self._buffer.attributes) + self.fp = open(self.filename, '%sb' % mode) - def _var_array(self): - """Read all variables.""" - # Read variables. - self.variables = {} - for k,d in self._buffer.iteritems(): - if isinstance(d, GridType) or isinstance(d, ArrayType): - name = k - self.variables[name] = _RemoteVariable(d,self._maskandscale) + self._attributes = {} + if mode is 'r': + self._read() + + def __setattr__(self, attr, value): + # Store user defined attributes in a separate dict, + # so we can save them to file later. + try: + self._attributes[attr] = value + except AttributeError: + pass + self.__dict__[attr] = value + def close(self): - # this is a no-op provided for compatibility - pass + if not self.fp.closed: + try: + self.flush() + finally: + self.fp.close() + __del__ = close + def createDimension(self, name, length): + self.dimensions[name] = length + self._dims.append(name) -class _RemoteVariable(object): - def __init__(self, var, maskandscale): - self._var = var - self._maskandscale = maskandscale - self.dtype = var.type - self.shape = var.shape - self.dimensions = var.dimensions - self.__dict__.update(var.attributes) + def createVariable(self, name, type, dimensions): + shape = tuple([self.dimensions[dim] for dim in dimensions]) + shape_ = tuple([dim or 0 for dim in shape]) # replace None with 0 for numpy - def __getitem__(self, index): - datout = squeeze(self._var.__getitem__(index)) - # automatically - # - remove singleton dimensions - # - create a masked array using missing_value or _FillValue attribute - # - apply scale_factor and add_offset to packed integer data - if self._maskandscale: - return _maskandscale(self,datout) + if isinstance(type, basestring): type = dtype(type) + typecode, size = type.char, type.itemsize + dtype_ = '>%s' % typecode + if size > 1: dtype_ += str(size) + + data = empty(shape_, dtype=dtype_) + self.variables[name] = netcdf_variable(data, typecode, shape, dimensions) + return self.variables[name] + + def flush(self): + if self.mode is 'w': + self._write() + sync = flush + + def _write(self): + self.fp.write('CDF') + + self.__dict__['version_byte'] = 1 + self.fp.write(array(1, '>b').tostring()) + + # Write headers and data. + self._write_numrecs() + self._write_dim_array() + self._write_gatt_array() + self._write_var_array() + + def _write_numrecs(self): + # Get highest record count from all record variables. + for var in self.variables.values(): + if not var.shape[0] and len(var.data) > self._recs: + self.__dict__['_recs'] = len(var.data) + self._pack_int(self._recs) + + def _write_dim_array(self): + if self.dimensions: + self.fp.write(NC_DIMENSION) + self._pack_int(len(self.dimensions)) + for name, length in self.dimensions.items(): + self._pack_string(name) + self._pack_int(length or 0) # replace None with 0 for record dimension else: - return datout + self.fp.write(ABSENT) - def typecode(self): - return _typecodes[self.dtype] + def _write_gatt_array(self): + self._write_att_array(self._attributes) + def _write_att_array(self, attributes): + if attributes: + self.fp.write(NC_ATTRIBUTE) + self._pack_int(len(attributes)) + for name, values in attributes.items(): + self._pack_string(name) + self._write_values(values) + else: + self.fp.write(ABSENT) -class _LocalFile(object): - """A NetCDF file reader. API is the same as Scientific.IO.NetCDF.""" + def _write_var_array(self): + if self.variables: + self.fp.write(NC_VARIABLE) + self._pack_int(len(self.variables)) - def __init__(self, file, maskandscale): - self._buffer = open(file, 'rb') - self._maskandscale = maskandscale - self._parse() + # Sort variables non-recs first, then recs. + variables = self.variables.items() + variables.sort(key=lambda (k, v): v.shape and v.shape[0] is not None) + variables.reverse() + variables = [k for (k, v) in variables] - def read(self, size=-1): - """Alias for reading the file buffer.""" - return self._buffer.read(size) + # Set the metadata for all variables. + for name in variables: + self._write_var_metadata(name) + # Now that we have the metadata, we know the vsize of + # each record variable, so we can calculate recsize. + self.__dict__['_recsize'] = sum([ + var._vsize for var in self.variables.values() + if var.shape[0] is None]) + # Set the data for all variables. + for name in variables: + self._write_var_data(name) + else: + self.fp.write(ABSENT) - def _parse(self): - """Initial parsing of the header.""" - # Check magic bytes. - assert self.read(3) == 'CDF' + def _write_var_metadata(self, name): + var = self.variables[name] - # Read version byte. - byte = self.read(1) - self.version_byte = struct.unpack('>b', byte)[0] + self._pack_string(name) + self._pack_int(len(var.dimensions)) + for dimname in var.dimensions: + dimid = self._dims.index(dimname) + self._pack_int(dimid) - # Read header info. - self._numrecs() - self._dim_array() - self._gatt_array() - self._var_array() + self._write_att_array(var._attributes) - def _numrecs(self): - """Read number of records.""" - self._nrecs = self._unpack_int() + nc_type = REVERSE[var.typecode()] + self.fp.write(nc_type) - def _dim_array(self): - """Read a dict with dimensions names and sizes.""" - assert self.read(4) in [ZERO, NC_DIMENSION] + if var.shape[0]: + vsize = var.data.size * var.data.itemsize + vsize += -vsize % 4 + else: # record variable + vsize = var.data[0].size * var.data.itemsize + rec_vars = len([var for var in self.variables.values() + if var.shape[0] is None]) + if rec_vars > 1: + vsize += -vsize % 4 + self.variables[name].__dict__['_vsize'] = vsize + self._pack_int(vsize) + + # Pack a bogus begin, and set the real value later. + self.variables[name].__dict__['_begin'] = self.fp.tell() + self._pack_begin(0) + + def _write_var_data(self, name): + var = self.variables[name] + + # Set begin in file header. + the_beguine = self.fp.tell() + self.fp.seek(var._begin) + self._pack_begin(the_beguine) + self.fp.seek(the_beguine) + + # Write data. + if var.shape[0]: + self.fp.write(var.data.tostring()) + count = var.data.size * var.data.itemsize + self.fp.write('0' * (var._vsize - count)) + else: # record variable + # Handle rec vars with shape[0] < nrecs. + if self._recs > len(var.data): + shape = (self._recs,) + var.data.shape[1:] + var.data.resize(shape) + + pos0 = pos = self.fp.tell() + for rec in var.data: + # Apparently scalars cannot be converted to big endian. If we + # try to convert a ``=i4`` scalar to, say, '>i4' the dtype + # will remain as ``=i4``. + if not rec.shape and (rec.dtype.byteorder == '<' or + (rec.dtype.byteorder == '=' and LITTLE_ENDIAN)): + rec = rec.byteswap() + self.fp.write(rec.tostring()) + # Padding + count = rec.size * rec.itemsize + self.fp.write('0' * (var._vsize - count)) + pos += self._recsize + self.fp.seek(pos) + self.fp.seek(pos0 + var._vsize) + + def _write_values(self, values): + values = asarray(values) + values = values.astype(values.dtype.newbyteorder('>')) + + nc_type = REVERSE[values.dtype.char] + self.fp.write(nc_type) + + if values.dtype.char == 'S': + nelems = values.itemsize + else: + nelems = values.size + self._pack_int(nelems) + + if not values.shape and (values.dtype.byteorder == '<' or + (values.dtype.byteorder == '=' and LITTLE_ENDIAN)): + values = values.byteswap() + self.fp.write(values.tostring()) + count = values.size * values.itemsize + self.fp.write('0' * (-count % 4)) # pad + + def _read(self): + # Check magic bytes and version + assert self.fp.read(3) == 'CDF', "Error: %s is not a valid NetCDF 3 file" % self.filename + self.__dict__['version_byte'] = fromstring(self.fp.read(1), '>b')[0] + + # Read file headers and set data. + self._read_numrecs() + self._read_dim_array() + self._read_gatt_array() + self._read_var_array() + + def _read_numrecs(self): + self.__dict__['_recs'] = self._unpack_int() + + def _read_dim_array(self): + assert self.fp.read(4) in [ZERO, NC_DIMENSION] count = self._unpack_int() - self.dimensions = {} - self._dims = [] for dim in range(count): - name = self._read_string() - length = self._unpack_int() - if length == 0: length = None # record dimension + name = self._unpack_string() + length = self._unpack_int() or None # None for record dimension self.dimensions[name] = length - self._dims.append(name) # preserve dim order + self._dims.append(name) # preserve order - def _gatt_array(self): - """Read global attributes.""" - self.attributes = self._att_array() + def _read_gatt_array(self): + for k, v in self._read_att_array().items(): + self.__setattr__(k, v) - # Update __dict__ for compatibility with S.IO.N - self.__dict__.update(self.attributes) - - def _att_array(self): - """Read a dict with attributes.""" - assert self.read(4) in [ZERO, NC_ATTRIBUTE] + def _read_att_array(self): + assert self.fp.read(4) in [ZERO, NC_ATTRIBUTE] count = self._unpack_int() - # Read attributes. attributes = {} - for attribute in range(count): - name = self._read_string() - nc_type = self._unpack_int() - n = self._unpack_int() + for attr in range(count): + name = self._unpack_string() + attributes[name] = self._read_values() + return attributes - # Read value for attributes. - attributes[name] = self._read_values(n, nc_type) + def _read_var_array(self): + assert self.fp.read(4) in [ZERO, NC_VARIABLE] - return attributes + begin = 0 + dtypes = {'names': [], 'formats': []} + rec_vars = [] + count = self._unpack_int() + for var in range(count): + name, dimensions, shape, attributes, typecode, size, dtype_, begin_, vsize = self._read_var() + if shape and shape[0] is None: + rec_vars.append(name) + self.__dict__['_recsize'] += vsize + if begin == 0: begin = begin_ + dtypes['names'].append(name) + dtypes['formats'].append(str(shape[1:]) + dtype_) - def _var_array(self): - """Read all variables.""" - assert self.read(4) in [ZERO, NC_VARIABLE] + # Handle padding with a virtual variable. + if typecode in 'bch': + actual_size = reduce(mul, (1,) + shape[1:]) * size + padding = -actual_size % 4 + if padding: + dtypes['names'].append('_padding_%d' % var) + dtypes['formats'].append('(%d,)>b' % padding) - # Read size of each record, in bytes. - self._read_recsize() + # Data will be set later. + data = None + else: + mm = mmap(self.fp.fileno(), begin_+vsize, access=ACCESS_READ) + data = ndarray.__new__(ndarray, shape, dtype=dtype_, + buffer=mm, offset=begin_, order=0) - # Read variables. - self.variables = {} - count = self._unpack_int() - for variable in range(count): - name = self._read_string() - self.variables[name] = self._read_var() + # Add variable. + self.variables[name] = netcdf_variable( + data, typecode, shape, dimensions, attributes) - def _read_recsize(self): - """Read all variables and compute record bytes.""" - pos = self._buffer.tell() - - recsize = 0 - count = self._unpack_int() - for variable in range(count): - name = self._read_string() - n = self._unpack_int() - isrec = False - for i in range(n): - dimid = self._unpack_int() - name = self._dims[dimid] - dim = self.dimensions[name] - if dim is None and i == 0: - isrec = True - attributes = self._att_array() - nc_type = self._unpack_int() - vsize = self._unpack_int() - begin = [self._unpack_int, self._unpack_int64][self.version_byte-1]() + if rec_vars: + # Remove padding when only one record variable. + if len(rec_vars) == 1: + dtypes['names'] = dtypes['names'][:1] + dtypes['formats'] = dtypes['formats'][:1] - if isrec: recsize += vsize + # Build rec array. + mm = mmap(self.fp.fileno(), begin+self._recs*self._recsize, access=ACCESS_READ) + rec_array = ndarray.__new__(ndarray, (self._recs,), dtype=dtypes, + buffer=mm, offset=begin, order=0) - self._recsize = recsize - self._buffer.seek(pos) + for var in rec_vars: + self.variables[var].__dict__['data'] = rec_array[var] def _read_var(self): + name = self._unpack_string() dimensions = [] shape = [] - n = self._unpack_int() - isrec = False - for i in range(n): + dims = self._unpack_int() + + for i in range(dims): dimid = self._unpack_int() - name = self._dims[dimid] - dimensions.append(name) - dim = self.dimensions[name] - if dim is None and i == 0: - dim = self._nrecs - isrec = True + dimname = self._dims[dimid] + dimensions.append(dimname) + dim = self.dimensions[dimname] shape.append(dim) dimensions = tuple(dimensions) shape = tuple(shape) - attributes = self._att_array() - nc_type = self._unpack_int() + attributes = self._read_att_array() + nc_type = self.fp.read(4) vsize = self._unpack_int() - - # Read offset. begin = [self._unpack_int, self._unpack_int64][self.version_byte-1]() - return _LocalVariable(self._buffer.fileno(), nc_type, vsize, begin, shape, dimensions, attributes, isrec, self._recsize, maskandscale=self._maskandscale) + typecode, size = TYPEMAP[nc_type] + if typecode is 'c': + dtype_ = '>c' + else: + dtype_ = '>%s' % typecode + if size > 1: dtype_ += str(size) - def _read_values(self, n, nc_type): - bytes = [1, 1, 2, 4, 4, 8] - typecodes = ['b', 'c', 'h', 'i', 'f', 'd'] - - count = n * bytes[nc_type-1] - values = self.read(count) - padding = self.read((4 - (count % 4)) % 4) - - typecode = typecodes[nc_type-1] - if nc_type != 2: # not char - values = struct.unpack('>%s' % (typecode * n), values) - values = array(values, dtype=typecode) + return name, dimensions, shape, attributes, typecode, size, dtype_, begin, vsize + + def _read_values(self): + nc_type = self.fp.read(4) + n = self._unpack_int() + + typecode, size = TYPEMAP[nc_type] + + count = n*size + values = self.fp.read(count) + self.fp.read(-count % 4) # read padding + + if typecode is not 'c': + values = fromstring(values, dtype='>%s%d' % (typecode, size)) + if values.shape == (1,): values = values[0] else: - # Remove EOL terminator. - if values.endswith('\x00'): values = values[:-1] - + values = values.rstrip('\x00') return values + def _pack_begin(self, begin): + if self.version_byte == 1: + self._pack_int(begin) + elif self.version_byte == 2: + self._pack_int64(begin) + + def _pack_int(self, value): + self.fp.write(array(value, '>i').tostring()) + _pack_int32 = _pack_int + def _unpack_int(self): - return struct.unpack('>i', self.read(4))[0] + return fromstring(self.fp.read(4), '>i')[0] _unpack_int32 = _unpack_int + def _pack_int64(self, value): + self.fp.write(array(value, '>q').tostring()) + def _unpack_int64(self): - return struct.unpack('>q', self.read(8))[0] + return fromstring(self.fp.read(8), '>q')[0] - def _read_string(self): - count = struct.unpack('>i', self.read(4))[0] - s = self.read(count) - # Remove EOL terminator. - if s.endswith('\x00'): s = s[:-1] - padding = self.read((4 - (count % 4)) % 4) + def _pack_string(self, s): + count = len(s) + self._pack_int(count) + self.fp.write(s) + self.fp.write('0' * (-count % 4)) # pad + + def _unpack_string(self): + count = self._unpack_int() + s = self.fp.read(count).rstrip('\x00') + self.fp.read(-count % 4) # read padding return s - def close(self): - self._buffer.close() +class netcdf_variable(object): + """ + ``netcdf_variable`` objects are constructed by calling the method + ``createVariable`` on the netcdf_file object. -class _LocalVariable(object): - def __init__(self, fileno, nc_type, vsize, begin, shape, dimensions, attributes, isrec=False, recsize=0, maskandscale=True): - self._nc_type = nc_type - self._vsize = vsize - self._begin = begin - self.shape = shape + ``netcdf_variable`` objects behave much like array objects defined in + Numpy, except that their data resides in a file. Data is read by + indexing and written by assigning to an indexed subset; the entire + array can be accessed by the index ``[:]`` or using the methods + ``getValue`` and ``assignValue``. ``netcdf_variable`` objects also + have attribute ``shape`` with the same meaning as for arrays, but + the shape cannot be modified. There is another read-only attribute + ``dimensions``, whose value is the tuple of dimension names. + + All other attributes correspond to variable attributes defined in + the NetCDF file. Variable attributes are created by assigning to an + attribute of the ``netcdf_variable`` object. + + """ + def __init__(self, data, typecode, shape, dimensions, attributes=None, maskandscale=True): + self.data = data + self._typecode = typecode + self._shape = shape self.dimensions = dimensions - self.attributes = attributes # for ``dap.plugins.netcdf`` - self.__dict__.update(attributes) - self._is_record = isrec self._maskandscale = maskandscale - # Number of bytes and type. - self._bytes = [1, 1, 2, 4, 4, 8][self._nc_type-1] - type_ = ['i', 'S', 'i', 'i', 'f', 'f'][self._nc_type-1] - dtype = '>%s%d' % (type_, self._bytes) - bytes = self._begin + self._vsize + self._attributes = attributes or {} + for k, v in self._attributes.items(): + self.__dict__[k] = v - if isrec: - # Record variables are not stored contiguosly on disk, so we - # need to create a separate array for each record. - self.__array_data__ = empty(shape, dtype) - bytes += (shape[0] - 1) * recsize - for n in range(shape[0]): - offset = self._begin + (n * recsize) - mm = mmap.mmap(fileno, bytes, access=mmap.ACCESS_READ) - self.__array_data__[n] = ndarray.__new__(ndarray, shape[1:], dtype=dtype, buffer=mm, offset=offset, order=0) - else: - # Create buffer and data. - mm = mmap.mmap(fileno, bytes, access=mmap.ACCESS_READ) - self.__array_data__ = ndarray.__new__(ndarray, shape, dtype=dtype, buffer=mm, offset=self._begin, order=0) + def __setattr__(self, attr, value): + # Store user defined attributes in a separate dict, + # so we can save them to file later. + try: + self._attributes[attr] = value + except AttributeError: + pass + self.__dict__[attr] = value - # N-D array interface - self.__array_interface__ = {'shape' : shape, - 'typestr': dtype, - 'data' : self.__array_data__, - 'version': 3, - } + @property + def shape(self): + return self.data.shape + def getValue(self): + return self.data.item() + + def assignValue(self, value): + self.data.itemset(value) + + def typecode(self): + return self._typecode + def __getitem__(self, index): - datout = squeeze(self.__array_data__.__getitem__(index)) - # automatically - # - remove singleton dimensions - # - create a masked array using missing_value or _FillValue attribute - # - apply scale_factor and add_offset to packed integer data + datout = squeeze(self.data[index]) if self._maskandscale: return _maskandscale(self,datout) else: return datout - def getValue(self): - """For scalars.""" - return self.__array_data__.item() + def __setitem__(self, index, data): + # Expand data for record vars? + if not self._shape[0]: + if isinstance(index, tuple): + rec_index = index[0] + else: + rec_index = index + if isinstance(rec_index, slice): + recs = (rec_index.start or 0) + len(data) + else: + recs = rec_index + 1 + if recs > len(self.data): + shape = (recs,) + self._shape[1:] + self.data.resize(shape) + self.data[index] = data - def typecode(self): - return ['b', 'c', 'h', 'i', 'f', 'd'][self._nc_type-1] + +NetCDFFile = netcdf_file +NetCDFVariable = netcdf_variable + +# default _FillValue for netcdf types (apply also to corresponding +# DAP types). +_default_fillvals = {'c':'\0', + 'S':"", + 'b':-127, + 'B':-127, + 'h':-32767, + 'H':65535, + 'i':-2147483647L, + 'L':4294967295L, + 'q':-2147483647L, + 'f':9.9692099683868690e+36, + 'd':9.9692099683868690e+36} +def _maskandscale(var,datout): + totalmask = zeros(datout.shape,bool) + fillval = None + if hasattr(var, 'missing_value') and (datout == var.missing_value).any(): + fillval = var.missing_value + totalmask += datout==fillval + if hasattr(var, '_FillValue') and (datout == var._FillValue).any(): + if fillval is None: + fillval = var._FillValue + totalmask += datout==var._FillValue + elif (datout == _default_fillvals[var.typecode()]).any(): + if fillval is None: + fillval = _default_fillvals[var.typecode()] + totalmask += datout==_default_fillvals[var.dtype] + if fillval is not None: + datout = ma.masked_array(datout,mask=totalmask,fill_value=fillval) + try: + datout = var.scale_factor*datout + var.add_offset + except: + pass + return datout This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2008-08-09 12:57:57
|
Revision: 6010 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=6010&view=rev Author: jswhit Date: 2008-08-09 12:57:54 +0000 (Sat, 09 Aug 2008) Log Message: ----------- update to pupynere 1.0.2 Modified Paths: -------------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2008-08-08 20:52:12 UTC (rev 6009) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2008-08-09 12:57:54 UTC (rev 6010) @@ -3641,19 +3641,25 @@ return corners def NetCDFFile(file, mode='r', maskandscale=True): - """NetCDF File reader. API is the same as Scientific.IO.NetCDF. + """NetCDF File reader/writer. API is the same as Scientific.IO.NetCDF. + If ``file`` is a URL that starts with `http`, it is assumed - to be a remote OPenDAP dataset, and the python dap client is used + to be a remote OPenDAP dataset, and pydap is used to retrieve the data. Only the OPenDAP Array and Grid data types are recognized. If file does not start with `http`, it - is assumed to be a local netCDF file. Data will - automatically be converted to masked arrays if the variable has either - a ``missing_value`` or ``_FillValue`` attribute, and some data points - are equal to the value specified by that attribute. In addition, - variables stored as integers that have the ``scale_factor`` and - ``add_offset`` attribute will automatically be rescaled to floats when - read. To suppress these automatic conversions, set the - ``maskandscale`` keyword to False. + is assumed to be a local netCDF file and is read + with scipy.io.netcdf. Both pydap and scipy.io.netcdf are written + by Roberto De Almeida. + + Data will + automatically be converted to and from masked arrays if the variable + has either a ``missing_value`` or ``_FillValue`` attribute, and + some data points are equal to the value specified by that attribute. + In addition, variables that have the ``scale_factor`` and ``add_offset`` + attribute will automatically be converted to and from short integers. + To suppress these automatic conversions, set the ``maskandscale`` + keyword to False. + """ if file.startswith('http'): return netcdf._RemoteFile(file,maskandscale=maskandscale) Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py 2008-08-08 20:52:12 UTC (rev 6009) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py 2008-08-09 12:57:54 UTC (rev 6010) @@ -1,5 +1,5 @@ from numpy import ma, squeeze -from pupynere import netcdf_file, _maskandscale +from pupynere import netcdf_file, _unmaskandscale from dap.client import open as open_remote from dap.dtypes import ArrayType, GridType, typemap @@ -68,7 +68,7 @@ # - create a masked array using missing_value or _FillValue attribute # - apply scale_factor and add_offset to packed integer data if self._maskandscale: - return _maskandscale(self,datout) + return _unmaskandscale(self,datout) else: return datout Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py 2008-08-08 20:52:12 UTC (rev 6009) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py 2008-08-09 12:57:54 UTC (rev 6010) @@ -128,9 +128,7 @@ """ def __init__(self, filename, mode='r', maskandscale=False): - self._maskandscale = maskandscale - self.filename = filename assert mode in 'rw', "Mode must be either 'r' or 'w'." @@ -181,7 +179,7 @@ if size > 1: dtype_ += str(size) data = empty(shape_, dtype=dtype_) - self.variables[name] = netcdf_variable(data, typecode, shape, dimensions) + self.variables[name] = netcdf_variable(data, typecode, shape, dimensions, maskandscale=self._maskandscale) return self.variables[name] def flush(self): @@ -204,7 +202,7 @@ def _write_numrecs(self): # Get highest record count from all record variables. for var in self.variables.values(): - if not var.shape[0] and len(var.data) > self._recs: + if not var._shape[0] and len(var.data) > self._recs: self.__dict__['_recs'] = len(var.data) self._pack_int(self._recs) @@ -238,7 +236,7 @@ # Sort variables non-recs first, then recs. variables = self.variables.items() - variables.sort(key=lambda (k, v): v.shape and v.shape[0] is not None) + variables.sort(key=lambda (k, v): v._shape and v._shape[0] is not None) variables.reverse() variables = [k for (k, v) in variables] @@ -249,7 +247,7 @@ # each record variable, so we can calculate recsize. self.__dict__['_recsize'] = sum([ var._vsize for var in self.variables.values() - if var.shape[0] is None]) + if var._shape[0] is None]) # Set the data for all variables. for name in variables: self._write_var_data(name) @@ -270,13 +268,13 @@ nc_type = REVERSE[var.typecode()] self.fp.write(nc_type) - if var.shape[0]: + if var._shape[0]: vsize = var.data.size * var.data.itemsize vsize += -vsize % 4 else: # record variable vsize = var.data[0].size * var.data.itemsize rec_vars = len([var for var in self.variables.values() - if var.shape[0] is None]) + if var._shape[0] is None]) if rec_vars > 1: vsize += -vsize % 4 self.variables[name].__dict__['_vsize'] = vsize @@ -296,7 +294,7 @@ self.fp.seek(the_beguine) # Write data. - if var.shape[0]: + if var._shape[0]: self.fp.write(var.data.tostring()) count = var.data.size * var.data.itemsize self.fp.write('0' * (var._vsize - count)) @@ -413,7 +411,7 @@ # Add variable. self.variables[name] = netcdf_variable( - data, typecode, shape, dimensions, attributes) + data, typecode, shape, dimensions, attributes, maskandscale=self._maskandscale) if rec_vars: # Remove padding when only one record variable. @@ -527,7 +525,7 @@ attribute of the ``netcdf_variable`` object. """ - def __init__(self, data, typecode, shape, dimensions, attributes=None, maskandscale=True): + def __init__(self, data, typecode, shape, dimensions, attributes=None, maskandscale=False): self.data = data self._typecode = typecode self._shape = shape @@ -561,13 +559,15 @@ return self._typecode def __getitem__(self, index): - datout = squeeze(self.data[index]) + data = squeeze(self.data[index]) if self._maskandscale: - return _maskandscale(self,datout) + return _unmaskandscale(self,data) else: - return datout + return data def __setitem__(self, index, data): + if self._maskandscale: + data = _maskandscale(self,data) # Expand data for record vars? if not self._shape[0]: if isinstance(index, tuple): @@ -600,24 +600,55 @@ 'q':-2147483647L, 'f':9.9692099683868690e+36, 'd':9.9692099683868690e+36} -def _maskandscale(var,datout): - totalmask = zeros(datout.shape,bool) - fillval = None - if hasattr(var, 'missing_value') and (datout == var.missing_value).any(): - fillval = var.missing_value - totalmask += datout==fillval - if hasattr(var, '_FillValue') and (datout == var._FillValue).any(): - if fillval is None: +def _unmaskandscale(var,data): + # if _maskandscale mode set to True, perform + # automatic unpacking using scale_factor/add_offset + # and automatic conversion to masked array using + # missing_value/_Fill_Value. + totalmask = zeros(data.shape, bool) + fill_value = None + if hasattr(var, 'missing_value') and (data == var.missing_value).any(): + mask=data==var.missing_value + fill_value = var.missing_value + totalmask += mask + if hasattr(var, '_FillValue') and (data == var._FillValue).any(): + mask=data==var._FillValue + if fill_value is None: + fill_value = var._FillValue + totalmask += mask + else: + fillval = _default_fillvals[var.typecode()] + if (data == fillval).any(): + mask=data==fillval + if fill_value is None: + fill_value = fillval + totalmask += mask + # all values where data == missing_value or _FillValue are + # masked. fill_value set to missing_value if it exists, + # otherwise _FillValue. + if fill_value is not None: + data = ma.masked_array(data,mask=totalmask,fill_value=fill_value) + # if variable has scale_factor and add_offset attributes, rescale. + if hasattr(var, 'scale_factor') and hasattr(var, 'add_offset'): + data = var.scale_factor*data + var.add_offset + return data + +def _maskandscale(var,data): + # if _maskandscale mode set to True, perform + # automatic packing using scale_factor/add_offset + # and automatic filling of masked arrays using + # missing_value/_Fill_Value. + # use missing_value as fill value. + # if no missing value set, use _FillValue. + if hasattr(data,'mask'): + if hasattr(var, 'missing_value'): + fillval = var.missing_value + elif hasattr(var, '_FillValue'): fillval = var._FillValue - totalmask += datout==var._FillValue - elif (datout == _default_fillvals[var.typecode()]).any(): - if fillval is None: + else: fillval = _default_fillvals[var.typecode()] - totalmask += datout==_default_fillvals[var.dtype] - if fillval is not None: - datout = ma.masked_array(datout,mask=totalmask,fill_value=fillval) - try: - datout = var.scale_factor*datout + var.add_offset - except: - pass - return datout + data = data.filled(fill_value=fillval) + # pack using scale_factor and add_offset. + if hasattr(var, 'scale_factor') and hasattr(var, 'add_offset'): + data = (data - var.add_offset)/var.scale_factor + return data This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2008-09-05 02:02:52
|
Revision: 6066 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=6066&view=rev Author: jswhit Date: 2008-09-05 02:02:49 +0000 (Fri, 05 Sep 2008) Log Message: ----------- allow for username and password for remote opendap datasets. Modified Paths: -------------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2008-09-04 19:56:37 UTC (rev 6065) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2008-09-05 02:02:49 UTC (rev 6066) @@ -3646,7 +3646,7 @@ else: return corners -def NetCDFFile(file, mode='r', maskandscale=True): +def NetCDFFile(file, mode='r', maskandscale=True, username=None, password=None): """NetCDF File reader/writer. API is the same as Scientific.IO.NetCDF. If ``file`` is a URL that starts with `http`, it is assumed @@ -3668,7 +3668,7 @@ """ if file.startswith('http'): - return netcdf._RemoteFile(file,maskandscale=maskandscale) + return netcdf._RemoteFile(file,maskandscale=maskandscale,username=username,password=password) else: return netcdf.netcdf_file(file,mode=mode,maskandscale=maskandscale) Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py 2008-09-04 19:56:37 UTC (rev 6065) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py 2008-09-05 02:02:49 UTC (rev 6066) @@ -8,8 +8,8 @@ class _RemoteFile(object): """A NetCDF file reader. API is the same as Scientific.IO.NetCDF.""" - def __init__(self, file, maskandscale=False): - self._buffer = open_remote(file) + def __init__(self, file, maskandscale=False, username=None, password=None): + self._buffer = open_remote(file,username=username,password=password) self._maskandscale = maskandscale self._parse() This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2008-09-05 15:37:09
|
Revision: 6070 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=6070&view=rev Author: jswhit Date: 2008-09-05 15:37:01 +0000 (Fri, 05 Sep 2008) Log Message: ----------- allow for 'cache' option when accessing remote datasets. Modified Paths: -------------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2008-09-05 13:23:21 UTC (rev 6069) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2008-09-05 15:37:01 UTC (rev 6070) @@ -3646,7 +3646,8 @@ else: return corners -def NetCDFFile(file, mode='r', maskandscale=True, username=None, password=None): +def NetCDFFile(file, mode='r', maskandscale=True, cache=None,\ + username=None, password=None, verbose=False): """NetCDF File reader/writer. API is the same as Scientific.IO.NetCDF. If ``file`` is a URL that starts with `http`, it is assumed @@ -3666,9 +3667,17 @@ To suppress these automatic conversions, set the ``maskandscale`` keyword to False. + The keywords ``cache``, ``username``, ``password`` and ``verbose`` are only + valid for remote OPenDAP datasets. ``username`` and ``password`` are used + to access OPenDAP datasets that require authentication. ``verbose=True`` + will make the pydap client print out the URLs being accessed. + ``cache`` is a location (a directory) for caching data, so that repeated + accesses to the same URL avoid the network. + """ if file.startswith('http'): - return netcdf._RemoteFile(file,maskandscale=maskandscale,username=username,password=password) + return netcdf._RemoteFile(file,maskandscale=maskandscale,\ + cache=cache,username=username,password=password,verbose=verbose) else: return netcdf.netcdf_file(file,mode=mode,maskandscale=maskandscale) Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py 2008-09-05 13:23:21 UTC (rev 6069) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py 2008-09-05 15:37:01 UTC (rev 6070) @@ -8,8 +8,10 @@ class _RemoteFile(object): """A NetCDF file reader. API is the same as Scientific.IO.NetCDF.""" - def __init__(self, file, maskandscale=False, username=None, password=None): - self._buffer = open_remote(file,username=username,password=password) + def __init__(self, file, maskandscale=False, cache=None,\ + username=None, password=None, verbose=False): + self._buffer = open_remote(file,cache=cache,\ + username=username,password=password,verbose=verbose) self._maskandscale = maskandscale self._parse() This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2008-10-18 01:15:22
|
Revision: 6254 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=6254&view=rev Author: jswhit Date: 2008-10-18 01:15:15 +0000 (Sat, 18 Oct 2008) Log Message: ----------- move tests to separate file. Modified Paths: -------------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py Added Paths: ----------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/test.py Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2008-10-17 20:24:50 UTC (rev 6253) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2008-10-18 01:15:15 UTC (rev 6254) @@ -3871,57 +3871,3 @@ """ cdftime = netcdftime.utime(units,calendar=calendar) return cdftime.date2num(dates) - - - -# beginnings of a test suite. - -from numpy.testing import NumpyTestCase,assert_almost_equal -class TestRotateVector(NumpyTestCase): - def make_array(self): - lat = np.array([0, 45, 75, 90]) - lon = np.array([0,90,180,270]) - u = np.ones((len(lat), len(lon))) - v = np.zeros((len(lat), len(lon))) - return u,v,lat,lon - - def test_cylindrical(self): - # Cylindrical case - B = Basemap() - u,v,lat,lon=self.make_array() - ru, rv = B.rotate_vector(u,v, lon, lat) - - # Check that the vectors are identical. - assert_almost_equal(ru, u) - assert_almost_equal(rv, v) - - def test_nan(self): - B = Basemap() - u,v,lat,lon=self.make_array() - # Set one element to 0, so that the vector magnitude is 0. - u[1,1] = 0. - ru, rv = B.rotate_vector(u,v, lon, lat) - assert not np.isnan(ru).any() - assert_almost_equal(u, ru) - assert_almost_equal(v, rv) - - def test_npstere(self): - # NP Stereographic case - B=Basemap(projection='npstere', boundinglat=50., lon_0=0.) - u,v,lat,lon=self.make_array() - v = np.ones((len(lat), len(lon))) - - ru, rv = B.rotate_vector(u,v, lon, lat) - - assert_almost_equal(ru[2, :],[1,-1,-1,1], 6) - assert_almost_equal(rv[2, :],[1,1,-1,-1], 6) - -def test(): - """ - Run some tests. - """ - import unittest - suite = unittest.makeSuite(TestRotateVector,'test') - runner = unittest.TextTestRunner() - runner.run(suite) - Added: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/test.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/test.py (rev 0) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/test.py 2008-10-18 01:15:15 UTC (rev 6254) @@ -0,0 +1,56 @@ +from mpl_toolkits.basemap import Basemap +import numpy as np + +# beginnings of a test suite. + +from numpy.testing import NumpyTestCase,assert_almost_equal +class TestRotateVector(NumpyTestCase): + def make_array(self): + lat = np.array([0, 45, 75, 90]) + lon = np.array([0,90,180,270]) + u = np.ones((len(lat), len(lon))) + v = np.zeros((len(lat), len(lon))) + return u,v,lat,lon + + def test_cylindrical(self): + # Cylindrical case + B = Basemap() + u,v,lat,lon=self.make_array() + ru, rv = B.rotate_vector(u,v, lon, lat) + + # Check that the vectors are identical. + assert_almost_equal(ru, u) + assert_almost_equal(rv, v) + + def test_nan(self): + B = Basemap() + u,v,lat,lon=self.make_array() + # Set one element to 0, so that the vector magnitude is 0. + u[1,1] = 0. + ru, rv = B.rotate_vector(u,v, lon, lat) + assert not np.isnan(ru).any() + assert_almost_equal(u, ru) + assert_almost_equal(v, rv) + + def test_npstere(self): + # NP Stereographic case + B=Basemap(projection='npstere', boundinglat=50., lon_0=0.) + u,v,lat,lon=self.make_array() + v = np.ones((len(lat), len(lon))) + + ru, rv = B.rotate_vector(u,v, lon, lat) + + assert_almost_equal(ru[2, :],[1,-1,-1,1], 6) + assert_almost_equal(rv[2, :],[1,1,-1,-1], 6) + +def test(): + """ + Run some tests. + """ + import unittest + suite = unittest.makeSuite(TestRotateVector,'test') + runner = unittest.TextTestRunner() + runner.run(suite) + +if __name__ == '__main__': + test() This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2009-01-26 19:48:09
|
Revision: 6830 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=6830&view=rev Author: jswhit Date: 2009-01-26 19:48:04 +0000 (Mon, 26 Jan 2009) Log Message: ----------- update netcdftime to version 0.7, add date2index function. Modified Paths: -------------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2009-01-26 16:39:14 UTC (rev 6829) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2009-01-26 19:48:04 UTC (rev 6830) @@ -17,6 +17,8 @@ :func:`num2date`: convert from a numeric time value to a datetime object. :func:`date2num`: convert from a datetime object to a numeric time value. + +:func:`date2index`: compute a time variable index corresponding to a date. """ from matplotlib import __version__ as _matplotlib_version from matplotlib.cbook import is_scalar, dedent @@ -3941,6 +3943,42 @@ cdftime = netcdftime.utime(units,calendar=calendar) return cdftime.date2num(dates) +def date2index(dates, nctime, calendar=None): + """ + Return indices of a netCDF time variable corresponding to the given dates. + + .. tabularcolumns:: |l|L| + + ============== ==================================================== + Arguments Description + ============== ==================================================== + dates A datetime object or a sequence of datetime objects. + The datetime objects should not include a + time-zone offset. + nctime A netCDF time variable object. The nctime object + must have a ``units`` attribute. + ============== ==================================================== + + .. tabularcolumns:: |l|L| + + ============== ==================================================== + Keywords Description + ============== ==================================================== + calendar describes the calendar used in the time + calculations. All the values currently defined in + the CF metadata convention + (http://cf-pcmdi.llnl.gov/documents/cf-conventions/) + are supported. + Valid calendars ``standard``, ``gregorian``, + ``proleptic_gregorian``, ``noleap``, ``365_day``, + ``julian``, ``all_leap``, ``366_day``. + Default is ``proleptic_gregorian``. + ============== ==================================================== + + Returns an index or a sequence of indices. + """ + return netcdftime.date2index(dates, nctime, calendar=None) + def maskoceans(lonsin,latsin,datain,inlands=False): """ mask data (``datain``), defined on a grid with latitudes ``latsin`` Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py 2009-01-26 16:39:14 UTC (rev 6829) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py 2009-01-26 19:48:04 UTC (rev 6830) @@ -1,14 +1,13 @@ """ Performs conversions of netCDF time coordinate data to/from datetime objects. """ -import math, re, time -import numpy as np +import math, numpy, re, time from datetime import datetime as real_datetime _units = ['days','hours','minutes','seconds','day','hour','minute','second'] _calendars = ['standard','gregorian','proleptic_gregorian','noleap','julian','all_leap','365_day','366_day','360_day'] -__version__ = '0.6' +__version__ = '0.7' class datetime: """ @@ -467,8 +466,8 @@ The B{C{calendar}} keyword describes the calendar used in the time calculations. All the values currently defined in the U{CF metadata convention -<http://www.cgd.ucar.edu/cms/eaton/cf-metadata/CF-1.0.html#time>} are -accepted. The default is C{'standard'}, which corresponds to the mixed +<http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.1/cf-conventions.html#time-coordinate>} +are accepted. The default is C{'standard'}, which corresponds to the mixed Gregorian/Julian calendar used by the C{udunits library}. Valid calendars are: @@ -534,8 +533,7 @@ C{'standard'} or C{'gregorian'} calendars. An exception will be raised if you pass a 'datetime-like' object in that range to the C{L{date2num}} class method. -Words of Wisdom from the British MetOffice concerning reference dates -U{http://www.metoffice.com/research/hadleycentre/models/GDT/ch26.html}: +Words of Wisdom from the British MetOffice concerning reference dates: "udunits implements the mixed Gregorian/Julian calendar system, as followed in England, in which dates prior to 1582-10-15 are assumed to use @@ -560,8 +558,8 @@ @keyword calendar: describes the calendar used in the time calculations. All the values currently defined in the U{CF metadata convention -<http://www.cgd.ucar.edu/cms/eaton/cf-metadata/CF-1.0.html#time>} are -accepted. The default is C{'standard'}, which corresponds to the mixed +<http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.1/cf-conventions.html#time-coordinate>} +are accepted. The default is C{'standard'}, which corresponds to the mixed Gregorian/Julian calendar used by the C{udunits library}. Valid calendars are: - C{'gregorian'} or C{'standard'} (default): @@ -630,7 +628,7 @@ except: isscalar = True if not isscalar: - date = np.array(date) + date = numpy.array(date) shape = date.shape if self.calendar in ['julian','standard','gregorian','proleptic_gregorian']: if isscalar: @@ -657,7 +655,7 @@ else: jdelta = [_360DayFromDate(d)-self._jd0 for d in date.flat] if not isscalar: - jdelta = np.array(jdelta) + jdelta = numpy.array(jdelta) # convert to desired units, add time zone offset. if self.units in ['second','seconds']: jdelta = jdelta*86400. + self.tzoffset*60. @@ -670,7 +668,7 @@ if isscalar: return jdelta else: - return np.reshape(jdelta,shape) + return numpy.reshape(jdelta,shape) def num2date(self,time_value): """ @@ -682,8 +680,8 @@ Resolution is 1 second. -Works for scalars, sequences and np arrays. -Returns a scalar if input is a scalar, else returns a np array. +Works for scalars, sequences and numpy arrays. +Returns a scalar if input is a scalar, else returns a numpy array. The datetime instances returned by C{num2date} are 'real' python datetime objects if the date falls in the Gregorian calendar (i.e. @@ -700,7 +698,7 @@ except: isscalar = True if not isscalar: - time_value = np.array(time_value) + time_value = numpy.array(time_value, dtype='d') shape = time_value.shape # convert to desired units, remove time zone offset. if self.units in ['second','seconds']: @@ -735,7 +733,7 @@ if isscalar: return date else: - return np.reshape(np.array(date),shape) + return numpy.reshape(numpy.array(date),shape) def _parse_date(origin): """Parses a date string and returns a tuple @@ -852,3 +850,147 @@ for site in sites: s = s[:site] + syear + s[site+4:] return s + +def date2num(dates,units,calendar='standard'): + """ +date2num(dates,units,calendar='standard') + +Return numeric time values given datetime objects. The units +of the numeric time values are described by the L{units} argument +and the L{calendar} keyword. The datetime objects must +be in UTC with no time-zone offset. If there is a +time-zone offset in C{units}, it will be applied to the +returned numeric values. + +Like the matplotlib C{date2num} function, except that it allows +for different units and calendars. Behaves the same if +C{units = 'days since 0001-01-01 00:00:00'} and +C{calendar = 'proleptic_gregorian'}. + +@param dates: A datetime object or a sequence of datetime objects. + The datetime objects should not include a time-zone offset. + +@param units: a string of the form C{'B{time units} since B{reference time}}' + describing the time units. B{C{time units}} can be days, hours, minutes + or seconds. B{C{reference time}} is the time origin. A valid choice + would be units=C{'hours since 1800-01-01 00:00:00 -6:00'}. + +@param calendar: describes the calendar used in the time calculations. + All the values currently defined in the U{CF metadata convention + <http://cf-pcmdi.llnl.gov/documents/cf-conventions/>} are supported. + Valid calendars C{'standard', 'gregorian', 'proleptic_gregorian' + 'noleap', '365_day', '360_day', 'julian', 'all_leap', '366_day'}. + Default is C{'standard'}, which is a mixed Julian/Gregorian calendar. + +@return: a numeric time value, or an array of numeric time values. + +The maximum resolution of the numeric time values is 1 second. + """ + cdftime = utime(units,calendar=calendar) + return cdftime.date2num(dates) + +def num2date(times,units,calendar='standard'): + """ +num2date(times,units,calendar='standard') + +Return datetime objects given numeric time values. The units +of the numeric time values are described by the C{units} argument +and the C{calendar} keyword. The returned datetime objects represent +UTC with no time-zone offset, even if the specified +C{units} contain a time-zone offset. + +Like the matplotlib C{num2date} function, except that it allows +for different units and calendars. Behaves the same if +C{units = 'days since 001-01-01 00:00:00'} and +C{calendar = 'proleptic_gregorian'}. + +@param times: numeric time values. Maximum resolution is 1 second. + +@param units: a string of the form C{'B{time units} since B{reference time}}' +describing the time units. B{C{time units}} can be days, hours, minutes +or seconds. B{C{reference time}} is the time origin. A valid choice +would be units=C{'hours since 1800-01-01 00:00:00 -6:00'}. + +@param calendar: describes the calendar used in the time calculations. +All the values currently defined in the U{CF metadata convention +<http://cf-pcmdi.llnl.gov/documents/cf-conventions/>} are supported. +Valid calendars C{'standard', 'gregorian', 'proleptic_gregorian' +'noleap', '365_day', '360_day', 'julian', 'all_leap', '366_day'}. +Default is C{'standard'}, which is a mixed Julian/Gregorian calendar. + +@return: a datetime instance, or an array of datetime instances. + +The datetime instances returned are 'real' python datetime +objects if the date falls in the Gregorian calendar (i.e. +C{calendar='proleptic_gregorian'}, or C{calendar = 'standard'} or C{'gregorian'} +and the date is after 1582-10-15). Otherwise, they are 'phony' datetime +objects which support some but not all the methods of 'real' python +datetime objects. This is because the python datetime module cannot +the uses the C{'proleptic_gregorian'} calendar, even before the switch +occured from the Julian calendar in 1582. The datetime instances +do not contain a time-zone offset, even if the specified C{units} +contains one. + """ + cdftime = utime(units,calendar=calendar) + return cdftime.num2date(times) + + +def _check_index(indices, dates, nctime, calendar): + """Assert that the time indices given correspond to the given dates.""" + t = nctime[indices] + assert numpy.all( num2date(t, nctime.units, calendar) == dates) + + +def date2index(dates, nctime, calendar=None): + """ + date2index(dates, nctime, calendar=None) + + Return indices of a netCDF time variable corresponding to the given dates. + + @param dates: A datetime object or a sequence of datetime objects. + The datetime objects should not include a time-zone offset. + + @param nctime: A netCDF time variable object. The nctime object must have a + C{units} attribute. + + @param calendar: Describes the calendar used in the time calculation. + Valid calendars C{'standard', 'gregorian', 'proleptic_gregorian' + 'noleap', '365_day', '360_day', 'julian', 'all_leap', '366_day'}. + Default is C{'standard'}, which is a mixed Julian/Gregorian calendar + If C{calendar} is None, its value is given by C{nctime.calendar} or + C{standard} if no such attribute exists. + """ + # Setting the calendar. + if calendar is None: + calendar = getattr(nctime, 'calendar', 'standard') + + num = numpy.atleast_1d(date2num(dates, nctime.units, calendar)) + + index = numpy.empty(numpy.alen(dates), int) + + # Trying to infer the correct index from the starting time and the stride. + try: + t0, t1 = nctime[:2] + dt = t1 - t0 + index[:] = (num-t0)/dt + + # convert numpy scalars or single element arrays to python ints. + if not len(index.shape) or index.shape == (1,): + index = index.item() + + # Checking that the index really corresponds to the given date. + _check_index(index, dates, nctime, calendar) + + except AssertionError: + # If check fails, use brute force method. + index[:] = numpy.digitize(num, nctime[:]) - 1 + + # convert numpy scalars or single element arrays to python ints. + if not len(index.shape) or index.shape == (1,): + index = index.item() + + # Perform check again. + _check_index(index, dates, nctime, calendar) + + return index + This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2009-09-09 16:39:34
|
Revision: 7726 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=7726&view=rev Author: jswhit Date: 2009-09-09 16:39:26 +0000 (Wed, 09 Sep 2009) Log Message: ----------- added __len__ methods (required for date2index) Modified Paths: -------------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py 2009-09-09 16:28:43 UTC (rev 7725) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py 2009-09-09 16:39:26 UTC (rev 7726) @@ -74,5 +74,8 @@ else: return datout + def __len__(self): + return self.shape[0] + def typecode(self): return _typecodes[self.dtype] Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py 2009-09-09 16:28:43 UTC (rev 7725) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py 2009-09-09 16:39:26 UTC (rev 7726) @@ -580,6 +580,9 @@ def shape(self): return self.data.shape + def __len__(self): + return self.data.shape[0] + def getValue(self): return self.data.item() This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <ef...@us...> - 2009-12-08 02:02:59
|
Revision: 8012 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=8012&view=rev Author: efiring Date: 2009-12-08 02:02:49 +0000 (Tue, 08 Dec 2009) Log Message: ----------- Remove spurious whitespace Modified Paths: -------------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/proj.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pupynere.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/pyproj.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/solar.py Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2009-12-07 01:17:27 UTC (rev 8011) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2009-12-08 02:02:49 UTC (rev 8012) @@ -1,7 +1,7 @@ """ Module for plotting data on maps with matplotlib. -Contains the :class:`Basemap` class (which does most of the +Contains the :class:`Basemap` class (which does most of the heavy lifting), and the following functions: :func:`NetCDFFile`: Read local and remote NetCDF datasets. @@ -26,8 +26,8 @@ _mpl_required_version = '0.98' if _matplotlib_version < _mpl_required_version: msg = dedent(""" - your matplotlib is too old - basemap requires version %s or - higher, you have version %s""" % + your matplotlib is too old - basemap requires version %s or + higher, you have version %s""" % (_mpl_required_version,_matplotlib_version)) raise ImportError(msg) from matplotlib import rcParams, is_interactive, _pylab_helpers @@ -180,12 +180,12 @@ lat_0 center of desired map domain (in degrees). ============== ==================================================== - For ``sinu``, ``moll``, ``npstere``, ``spstere``, ``nplaea``, ``splaea``, + For ``sinu``, ``moll``, ``npstere``, ``spstere``, ``nplaea``, ``splaea``, ``npaeqd``, ``spaeqd``, ``robin`` or ``mbtfpq``, the values of llcrnrlon, llcrnrlat, urcrnrlon, urcrnrlat, width and height are ignored (because either they are computed internally, or entire globe is - always plotted). - + always plotted). + For the cylindrical projections (``cyl``, ``merc``, ``mill`` and ``gall``), the default is to use llcrnrlon=-180,llcrnrlat=-90, urcrnrlon=180 and urcrnrlat=90). For all other @@ -193,7 +193,7 @@ corners or width and height must be specified by the user. For ``ortho`` and ``geos``, the lat/lon values of the corners may be specified, - or the x/y values of the corners (llcrnrx,llcrnry,urcrnrx,urcrnry) in the + or the x/y values of the corners (llcrnrx,llcrnry,urcrnrx,urcrnry) in the coordinate system of the global projection (with x=0,y=0 at the center of the global projection). If the corners are not specified, the entire globe is plotted. @@ -206,7 +206,7 @@ Keyword Description ============== ==================================================== resolution resolution of boundary database to use. Can be ``c`` - (crude), ``l`` (low), ``i`` (intermediate), ``h`` + (crude), ``l`` (low), ``i`` (intermediate), ``h`` (high), ``f`` (full) or None. If None, no boundary data will be read in (and class methods such as drawcoastlines will raise an @@ -217,40 +217,40 @@ (http://www.soest.hawaii.edu/wessel/gshhs/gshhs.html). State, country and river datasets from the Generic Mapping Tools (http://gmt.soest.hawaii.edu). - area_thresh coastline or lake with an area smaller than - area_thresh in km^2 will not be plotted. + area_thresh coastline or lake with an area smaller than + area_thresh in km^2 will not be plotted. Default 10000,1000,100,10,1 for resolution ``c``, ``l``, ``i``, ``h``, ``f``. rsphere radius of the sphere used to define map projection (default 6370997 meters, close to the arithmetic mean - radius of the earth). If given as a sequence, the + radius of the earth). If given as a sequence, the first two elements are interpreted as the radii - of the major and minor axes of an ellipsoid. - Note: sometimes an ellipsoid is specified by the + of the major and minor axes of an ellipsoid. + Note: sometimes an ellipsoid is specified by the major axis and an inverse flattening parameter (if). The minor axis (b) can be computed from the major - axis (a) and the inverse flattening parameter using + axis (a) and the inverse flattening parameter using the formula if = a/(a-b). suppress_ticks suppress automatic drawing of axis ticks and labels - in map projection coordinates. Default False, + in map projection coordinates. Default False, so parallels and meridians can be labelled instead. If parallel or meridian labelling is requested (using drawparallels and drawmeridians methods), automatic tick labelling will be supressed even if suppress_ticks=False. suppress_ticks=False is useful if you want to use your own custom tick - formatter, or if you want to let matplotlib label + formatter, or if you want to let matplotlib label the axes in meters using map projection coordinates. fix_aspect fix aspect ratio of plot to match aspect ratio of map projection region (default True). anchor determines how map is placed in axes rectangle - (passed to axes.set_aspect). Default is ``C``, + (passed to axes.set_aspect). Default is ``C``, which means map is centered. - Allowed values are - ``C``, ``SW``, ``S``, ``SE``, ``E``, ``NE``, + Allowed values are + ``C``, ``SW``, ``S``, ``SE``, ``E``, ``NE``, ``N``, ``NW``, and ``W``. - ax set default axes instance + ax set default axes instance (default None - matplotlib.pyplot.gca() may be used to get the current axes instance). If you don``t want matplotlib.pyplot to be imported, @@ -258,10 +258,10 @@ instance, or use the ``ax`` keyword in each Basemap method call that does drawing. In the first case, all Basemap method calls will draw to the same axes - instance. In the second case, you can draw to + instance. In the second case, you can draw to different axes with the same Basemap instance. You can also use the ``ax`` keyword in individual - method calls to selectively override the default + method calls to selectively override the default axes instance. ============== ==================================================== @@ -279,36 +279,36 @@ and mercator projections. default is lat_0 for stereographic projection. default is 0 for mercator projection. - lat_1 first standard parallel for lambert conformal, + lat_1 first standard parallel for lambert conformal, albers equal area and equidistant conic. - Latitude of one of the two points on the projection - centerline for oblique mercator. If lat_1 is not given, but - lat_0 is, lat_1 is set to lat_0 for lambert + Latitude of one of the two points on the projection + centerline for oblique mercator. If lat_1 is not given, but + lat_0 is, lat_1 is set to lat_0 for lambert conformal, albers equal area and equidistant conic. - lat_2 second standard parallel for lambert conformal, + lat_2 second standard parallel for lambert conformal, albers equal area and equidistant conic. Latitude of one of the two points on the projection - centerline for oblique mercator. If lat_2 is not - given it is set to lat_1 for lambert conformal, + centerline for oblique mercator. If lat_2 is not + given it is set to lat_1 for lambert conformal, albers equal area and equidistant conic. lon_1 Longitude of one of the two points on the projection centerline for oblique mercator. lon_2 Longitude of one of the two points on the projection centerline for oblique mercator. no_rot only used by oblique mercator. - If set to True, the map projection coordinates will + If set to True, the map projection coordinates will not be rotated to true North. Default is False (projection coordinates are automatically rotated). - lat_0 central latitude (y-axis origin) - used by all - projections. + lat_0 central latitude (y-axis origin) - used by all + projections. lon_0 central meridian (x-axis origin) - used by all projections. boundinglat bounding latitude for pole-centered projections - (npstere,spstere,nplaea,splaea,npaeqd,spaeqd). + (npstere,spstere,nplaea,splaea,npaeqd,spaeqd). These projections are square regions centered on the north or south pole. The longitude lon_0 is at 6-o'clock, and the - latitude circle boundinglat is tangent to the edge + latitude circle boundinglat is tangent to the edge of the map at lon_0. satellite_height height of satellite (in m) above equator - only relevant for geostationary projections @@ -325,11 +325,11 @@ projection map projection. Print the module variable ``supported_projections`` to see a list of allowed values. - aspect map aspect ratio + aspect map aspect ratio (size of y dimension / size of x dimension). llcrnrlon longitude of lower left hand corner of the selected map domain. - llcrnrlat latitude of lower left hand corner of the + llcrnrlat latitude of lower left hand corner of the selected map domain. urcrnrlon longitude of upper right hand corner of the selected map domain. @@ -337,7 +337,7 @@ selected map domain. llcrnrx x value of lower left hand corner of the selected map domain in map projection coordinates. - llcrnry y value of lower left hand corner of the + llcrnry y value of lower left hand corner of the selected map domain in map projection coordinates. urcrnrx x value of upper right hand corner of the selected map domain in map projection coordinates. @@ -345,8 +345,8 @@ selected map domain in map projection coordinates. rmajor equatorial radius of ellipsoid used (in meters). rminor polar radius of ellipsoid used (in meters). - resolution resolution of boundary dataset being used (``c`` - for crude, ``l`` for low, etc.). + resolution resolution of boundary dataset being used (``c`` + for crude, ``l`` for low, etc.). If None, no boundary dataset is associated with the Basemap instance. proj4string the string describing the map projection that is @@ -372,7 +372,7 @@ self.llcrnrlon and self.llcrnrlat. Input arguments lon, lat can be either scalar floats, sequences - or numpy arrays. + or numpy arrays. **Example Usage:** @@ -650,7 +650,7 @@ self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat elif projection in _cylproj: if projection == 'merc': - if lat_ts is None: + if lat_ts is None: lat_ts = 0. projparams['lat_ts']=lat_ts if not using_corners: @@ -795,7 +795,7 @@ ind.append(len(xd)) for i in ind: # don't add empty lists. - if len(range(iprev,i)): + if len(range(iprev,i)): coastsegs.append(zip(x[iprev:i],y[iprev:i])) iprev = i else: @@ -887,7 +887,7 @@ raise ValueError('%s projection cannot cross pole'%(self.projection)) # make sure orthographic or gnomonic projection has containsPole=True # we will compute the intersections in stereographic - # coordinates, then transform to orthographic. This is + # coordinates, then transform to orthographic. This is # because these projections are only defined on a hemisphere, and # some boundary features (like Eurasia) would be undefined otherwise. if self.projection in ['ortho','gnom'] and name == 'gshhs': @@ -895,7 +895,7 @@ lon_0=self.projparams['lon_0'] lat_0=self.projparams['lat_0'] re = self.projparams['R'] - # center of stereographic projection restricted to be + # center of stereographic projection restricted to be # nearest one of 6 points on the sphere (every 90 deg lat/lon). lon0 = 90.*(np.around(lon_0/90.)) lat0 = 90.*(np.around(lat_0/90.)) @@ -958,7 +958,7 @@ else: poly = Shape(b) # this is a workaround to avoid - # "GEOS_ERROR: TopologyException: + # "GEOS_ERROR: TopologyException: # found non-noded intersection between ..." # with geos 3.0.0 if _geoslib.__geos_major_version__ > 2: @@ -1035,7 +1035,7 @@ # create a GEOS geometry object. poly = Shape(b) # this is a workaround to avoid - # "GEOS_ERROR: TopologyException: + # "GEOS_ERROR: TopologyException: # found non-noded intersection between ..." # with geos 3.0.0 if _geoslib.__geos_major_version__ > 2: @@ -1212,7 +1212,7 @@ linewidth line width for boundary (default 1.) color color of boundary line (default black) fill_color fill the map region background with this - color (default is no fill or fill with axis + color (default is no fill or fill with axis background color). zorder sets the zorder for filling map background (default 0). @@ -1411,7 +1411,7 @@ antialiased antialiasing switch for coastlines (default True). ax axes instance (overrides default axes instance) zorder sets the zorder for the coastlines (if not specified, - uses default zorder for + uses default zorder for matplotlib.patches.LineCollections). ============== ==================================================== @@ -1443,11 +1443,11 @@ ============== ==================================================== linewidth country boundary line width (default 0.5) color country boundary line color (default black) - antialiased antialiasing switch for country boundaries (default + antialiased antialiasing switch for country boundaries (default True). ax axes instance (overrides default axes instance) zorder sets the zorder for the country boundaries (if not - specified uses default zorder for + specified uses default zorder for matplotlib.patches.LineCollections). ============== ==================================================== @@ -1486,7 +1486,7 @@ antialiased antialiasing switch for state boundaries (default True). ax axes instance (overrides default axes instance) - zorder sets the zorder for the state boundaries (if not + zorder sets the zorder for the state boundaries (if not specified, uses default zorder for matplotlib.patches.LineCollections). ============== ==================================================== @@ -1523,10 +1523,10 @@ ============== ==================================================== linewidth river boundary line width (default 0.5) color river boundary line color (default black) - antialiased antialiasing switch for river boundaries (default + antialiased antialiasing switch for river boundaries (default True). ax axes instance (overrides default axes instance) - zorder sets the zorder for the rivers (if not + zorder sets the zorder for the rivers (if not specified uses default zorder for matplotlib.patches.LineCollections). ============== ==================================================== @@ -1588,25 +1588,25 @@ Argument Description ============== ==================================================== shapefile path to shapefile components. Example: - shapefile='/home/jeff/esri/world_borders' assumes - that world_borders.shp, world_borders.shx and + shapefile='/home/jeff/esri/world_borders' assumes + that world_borders.shp, world_borders.shx and world_borders.dbf live in /home/jeff/esri. name name for Basemap attribute to hold the shapefile - vertices or points in map projection - coordinates. Class attribute name+'_info' is a list - of dictionaries, one for each shape, containing + vertices or points in map projection + coordinates. Class attribute name+'_info' is a list + of dictionaries, one for each shape, containing attributes of each shape from dbf file, For example, if name='counties', self.counties - will be a list of x,y vertices for each shape in + will be a list of x,y vertices for each shape in map projection coordinates and self.counties_info will be a list of dictionaries with shape - attributes. Rings in individual Polygon + attributes. Rings in individual Polygon shapes are split out into separate polygons, and additional keys 'RINGNUM' and 'SHAPENUM' are added to the shape attribute dictionary. ============== ==================================================== - The following optional keyword arguments are only relevant for Polyline + The following optional keyword arguments are only relevant for Polyline and Polygon shape types, for Point and MultiPoint shapes they are ignored. @@ -1615,13 +1615,13 @@ ============== ==================================================== Keyword Description ============== ==================================================== - drawbounds draw boundaries of shapes (default True). - zorder shape boundary zorder (if not specified, - default for mathplotlib.lines.LineCollection + drawbounds draw boundaries of shapes (default True). + zorder shape boundary zorder (if not specified, + default for mathplotlib.lines.LineCollection is used). linewidth shape boundary line width (default 0.5) color shape boundary line color (default black) - antialiased antialiasing switch for shape boundaries + antialiased antialiasing switch for shape boundaries (default True). ax axes instance (overrides default axes instance) ============== ==================================================== @@ -1632,7 +1632,7 @@ the SHPT* constants defined in the shapelib module, see http://shapelib.maptools.org/shp_api.html) and min and max are 4-element lists with the minimum and maximum values of the - vertices. If ``drawbounds=True`` a + vertices. If ``drawbounds=True`` a matplotlib.patches.LineCollection object is appended to the tuple. """ if not os.path.exists('%s.shp'%shapefile): @@ -1723,7 +1723,7 @@ fmt='%g',xoffset=None,yoffset=None,ax=None,latmax=None, **kwargs): """ - Draw and label parallels (latitude lines) for values (in degrees) + Draw and label parallels (latitude lines) for values (in degrees) given in the sequence ``circles``. .. tabularcolumns:: |l|L| @@ -1749,7 +1749,7 @@ labelled with "N" and "S". fmt a format string to format the parallel labels (default '%g') **or** a function that takes a - latitude value in degrees as it's only argument + latitude value in degrees as it's only argument and returns a formatted string. xoffset label offset from edge of map in x-direction (default is 0.01 times width of map in map @@ -1761,13 +1761,13 @@ latmax absolute value of latitude to which meridians are drawn (default is 80). \**kwargs additional keyword arguments controlling text - for labels that are passed on to + for labels that are passed on to the text method of the axes instance (see matplotlib.pyplot.text documentation). ============== ==================================================== returns a dictionary whose keys are the parallel values, and - whose values are tuples containing lists of the + whose values are tuples containing lists of the matplotlib.lines.Line2D and matplotlib.text.Text instances associated with each parallel. """ @@ -1811,7 +1811,7 @@ x,y = self(lons,lats) # remove points outside domain. # leave a little slop around edges (3*xdelta) - # don't really know why, but this appears to be needed to + # don't really know why, but this appears to be needed to # or lines sometimes don't reach edge of plot. testx = np.logical_and(x>=self.xmin-3*xdelta,x<=self.xmax+3*xdelta) x = np.compress(testx, x) @@ -2001,7 +2001,7 @@ labelled with "E" and "W". fmt a format string to format the meridian labels (default '%g') **or** a function that takes a - longitude value in degrees as it's only argument + longitude value in degrees as it's only argument and returns a formatted string. xoffset label offset from edge of map in x-direction (default is 0.01 times width of map in map @@ -2013,13 +2013,13 @@ latmax absolute value of latitude to which meridians are drawn (default is 80). \**kwargs additional keyword arguments controlling text - for labels that are passed on to + for labels that are passed on to the text method of the axes instance (see matplotlib.pyplot.text documentation). ============== ==================================================== returns a dictionary whose keys are the meridian values, and - whose values are tuples containing lists of the + whose values are tuples containing lists of the matplotlib.lines.Line2D and matplotlib.text.Text instances associated with each meridian. """ @@ -2049,7 +2049,7 @@ x,y = self(lons,lats) # remove points outside domain. # leave a little slop around edges (3*xdelta) - # don't really know why, but this appears to be needed to + # don't really know why, but this appears to be needed to # or lines sometimes don't reach edge of plot. testx = np.logical_and(x>=self.xmin-3*xdelta,x<=self.xmax+3*xdelta) x = np.compress(testx, x) @@ -2218,7 +2218,7 @@ Draw a polygon centered at ``lon_0,lat_0``. The polygon approximates a circle on the surface of the earth with radius ``radius_deg`` degrees latitude along longitude ``lon_0``, - made up of ``npts`` vertices. + made up of ``npts`` vertices. The polygon represents a Tissot's indicatrix (http://en.wikipedia.org/wiki/Tissot's_Indicatrix), which when drawn on a map shows the distortion @@ -2308,7 +2308,7 @@ Interpolate a scalar field (``datin``) from a lat/lon grid with longitudes = ``lons`` and latitudes = ``lats`` to a ``ny`` by ``nx`` map projection grid. Typically used to transform data to - map projection coordinates for plotting on a map with + map projection coordinates for plotting on a map with the :meth:`imshow`. .. tabularcolumns:: |l|L| @@ -2320,7 +2320,7 @@ lons, lats rank-1 arrays containing longitudes and latitudes (in degrees) of input data in increasing order. For non-cylindrical projections (those other than - ``cyl``, ``merc``, ``gall`` and ``mill``) lons must + ``cyl``, ``merc``, ``gall`` and ``mill``) lons must fit within range -180 to 180. nx, ny The size of the output regular grid in map projection coordinates @@ -2391,12 +2391,12 @@ lons, lats rank-1 arrays containing longitudes and latitudes (in degrees) of input data in increasing order. For non-cylindrical projections (those other than - ``cyl``, ``merc``, ``gall`` and ``mill``) lons must + ``cyl``, ``merc``, ``gall`` and ``mill``) lons must fit within range -180 to 180. nx, ny The size of the output regular grid in map projection coordinates ============== ==================================================== - + .. tabularcolumns:: |l|L| ============== ==================================================== @@ -2442,7 +2442,7 @@ def rotate_vector(self,uin,vin,lons,lats,returnxy=False): """ Rotate a vector field (``uin,vin``) on a rectilinear grid - with longitudes = ``lons`` and latitudes = ``lats`` from + with longitudes = ``lons`` and latitudes = ``lats`` from geographical (lat/lon) into map projection (x/y) coordinates. Differs from transform_vector in that no interpolation is done. @@ -2463,13 +2463,13 @@ lons, lats Arrays containing longitudes and latitudes (in degrees) of input data in increasing order. For non-cylindrical projections (those other than - ``cyl``, ``merc``, ``gall`` and ``mill``) lons must + ``cyl``, ``merc``, ``gall`` and ``mill``) lons must fit within range -180 to 180. ============== ==================================================== Returns ``uout, vout`` (rotated vector field). - If the optional keyword argument - ``returnxy`` is True (default is False), + If the optional keyword argument + ``returnxy`` is True (default is False), returns ``uout,vout,x,y`` (where ``x,y`` are the map projection coordinates of the grid defined by ``lons,lats``). """ @@ -2480,7 +2480,7 @@ if lons.ndim == lats.ndim == 1 and uin.ndim == vin.ndim == 2 and\ uin.shape[1] == vin.shape[1] == lons.shape[0] and\ uin.shape[0] == vin.shape[0] == lats.shape[0]: - lons, lats = np.meshgrid(lons, lats) + lons, lats = np.meshgrid(lons, lats) else: if not lons.shape == lats.shape == uin.shape == vin.shape: raise TypeError("shapes of lons,lats and uin,vin don't match") @@ -2493,41 +2493,41 @@ vin = vin.filled(1) else: masked = False - + # Map the (lon, lat) vector in the complex plane. uvc = uin + 1j*vin uvmag = np.abs(uvc) theta = np.angle(uvc) - - # Define a displacement (dlon, dlat) that moves all - # positions (lons, lats) a small distance in the - # direction of the original vector. + + # Define a displacement (dlon, dlat) that moves all + # positions (lons, lats) a small distance in the + # direction of the original vector. dc = 1E-5 * np.exp(theta*1j) dlat = dc.imag * np.cos(np.radians(lats)) - dlon = dc.real - + dlon = dc.real + # Deal with displacements that overshoot the North or South Pole. farnorth = np.abs(lats+dlat) >= 90.0 somenorth = farnorth.any() if somenorth: dlon[farnorth] *= -1.0 dlat[farnorth] *= -1.0 - + # Add displacement to original location and find the native coordinates. lon1 = lons + dlon lat1 = lats + dlat xn, yn = self(lon1, lat1) - - # Determine the angle of the displacement in the native coordinates. + + # Determine the angle of the displacement in the native coordinates. vecangle = np.arctan2(yn-y, xn-x) if somenorth: vecangle[farnorth] += np.pi - + # Compute the x-y components of the original vector. uvcout = uvmag * np.exp(1j*vecangle) uout = uvcout.real vout = uvcout.imag - + if masked: uout = ma.array(uout, mask=mask) vout = ma.array(vout, mask=mask) @@ -2581,7 +2581,7 @@ """ ax = kwargs.pop('ax', None) or self._check_ax() # if ax kwarg not supplied, and ax attribute not set, import pyplot. - if self.ax is None and kwargs.pop('ax', None) is None: + if self.ax is None and kwargs.pop('ax', None) is None: import matplotlib.pyplot as plt # allow callers to override the hold state by passing hold=True|False b = ax.ishold() @@ -2612,7 +2612,7 @@ def plot(self, *args, **kwargs): """ - Draw lines and/or markers on the map + Draw lines and/or markers on the map (see matplotlib.pyplot.plot documentation). Extra keyword ``ax`` can be used to override the default axis instance. @@ -2641,7 +2641,7 @@ def imshow(self, *args, **kwargs): """ - Display an image over the map + Display an image over the map (see matplotlib.pyplot.imshow documentation). ``extent`` and ``origin`` keywords set automatically so image @@ -2655,7 +2655,7 @@ """ ax = kwargs.pop('ax', None) or self._check_ax() # if ax kwarg not supplied, and ax attribute not set, import pyplot. - if self.ax is None and kwargs.pop('ax', None) is None: + if self.ax is None and kwargs.pop('ax', None) is None: import matplotlib.pyplot as plt kwargs['extent']=(self.llcrnrx,self.urcrnrx,self.llcrnry,self.urcrnry) # use origin='lower', unless overridden. @@ -2703,7 +2703,7 @@ """ ax = kwargs.pop('ax', None) or self._check_ax() # if ax kwarg not supplied, and ax attribute not set, import pyplot. - if self.ax is None and kwargs.pop('ax', None) is None: + if self.ax is None and kwargs.pop('ax', None) is None: import matplotlib.pyplot as plt # make x,y masked arrays # (masked where data is outside of projection limb) @@ -2747,7 +2747,7 @@ """ ax = kwargs.pop('ax', None) or self._check_ax() # if ax kwarg not supplied, and ax attribute not set, import pyplot. - if self.ax is None and kwargs.pop('ax', None) is None: + if self.ax is None and kwargs.pop('ax', None) is None: import matplotlib.pyplot as plt # allow callers to override the hold state by passing hold=True|False b = ax.ishold() @@ -2778,7 +2778,7 @@ def contour(self,x,y,data,*args,**kwargs): """ - Make a contour plot over the map + Make a contour plot over the map (see matplotlib.pyplot.contour documentation). Extra keyword ``ax`` can be used to override the default axis instance. @@ -2787,7 +2787,7 @@ """ ax = kwargs.pop('ax', None) or self._check_ax() # if ax kwarg not supplied, and ax attribute not set, import pyplot. - if self.ax is None and kwargs.pop('ax', None) is None: + if self.ax is None and kwargs.pop('ax', None) is None: import matplotlib.pyplot as plt # make sure x is monotonically increasing - if not, # print warning suggesting that the data be shifted in longitude @@ -2862,7 +2862,7 @@ """ ax = kwargs.pop('ax', None) or self._check_ax() # if ax kwarg not supplied, and ax attribute not set, import pyplot. - if self.ax is None and kwargs.pop('ax', None) is None: + if self.ax is None and kwargs.pop('ax', None) is None: import matplotlib.pyplot as plt # make sure x is monotonically increasing - if not, # print warning suggesting that the data be shifted in longitude @@ -2942,7 +2942,7 @@ """ ax = kwargs.pop('ax', None) or self._check_ax() # if ax kwarg not supplied, and ax attribute not set, import pyplot. - if self.ax is None and kwargs.pop('ax', None) is None: + if self.ax is None and kwargs.pop('ax', None) is None: import matplotlib.pyplot as plt # allow callers to override the hold state by passing hold=True|False b = ax.ishold() @@ -2972,8 +2972,8 @@ Other \*args and \**kwargs passed on to matplotlib.pyplot.barbs - Returns two matplotlib.axes.Barbs instances, one for the Northern - Hemisphere and one for the Southern Hemisphere. + Returns two matplotlib.axes.Barbs instances, one for the Northern + Hemisphere and one for the Southern Hemisphere. """ if _matplotlib_version < '0.98.3': msg = dedent(""" @@ -2982,7 +2982,7 @@ raise NotImplementedError(msg) ax = kwargs.pop('ax', None) or self._check_ax() # if ax kwarg not supplied, and ax attribute not set, import pyplot. - if self.ax is None and kwargs.pop('ax', None) is None: + if self.ax is None and kwargs.pop('ax', None) is None: import matplotlib.pyplot as plt # allow callers to override the hold state by passing hold=True|False b = ax.ishold() @@ -3025,7 +3025,7 @@ ============== ==================================================== Keywords Description ============== ==================================================== - land_color desired land color (color name or rgba tuple). + land_color desired land color (color name or rgba tuple). Default gray ("0.8"). ocean_color desired ocean color (color name or rgba tuple). Default white. @@ -3042,7 +3042,7 @@ lsmask_lats 1d array of latitudes for lsmask (ignored if lsmask is None). Latitudes must be ordered from -90 S northward. - \**kwargs extra keyword arguments passed on to + \**kwargs extra keyword arguments passed on to :meth:`imshow` ============== ==================================================== @@ -3089,7 +3089,7 @@ # redo the interpolation (unless a new land-sea mask is passed # in via the lsmask, lsmask_lons, lsmask_lats keywords). - # is it a cylindrical projection whose limits lie + # is it a cylindrical projection whose limits lie # outside the limits of the image? cylproj = self.projection in _cylproj and \ (self.urcrnrlon > lsmask_lons[-1] or \ @@ -3172,12 +3172,12 @@ If image is a URL (starts with 'http'), it is downloaded to a temp file using urllib.urlretrieve. - Default (if ``image`` not specified) is to display + Default (if ``image`` not specified) is to display 'blue marble next generation' image from http://visibleearth.nasa.gov/. Specified image must have pixels covering the whole globe in a regular lat/lon grid, starting and -180W and the South Pole. - Works with the global images from + Works with the global images from http://earthobservatory.nasa.gov/Features/BlueMarble/BlueMarble_monthlies.php. The ``scale`` keyword can be used to downsample (rescale) the image. @@ -3232,7 +3232,7 @@ delta = 360./float(nlons) self._bm_lons = np.arange(-180.+0.5*delta,180.,delta) self._bm_lats = np.arange(-90.+0.5*delta,90.,delta) - # is it a cylindrical projection whose limits lie + # is it a cylindrical projection whose limits lie # outside the limits of the image? cylproj = self.projection in _cylproj and \ (self.urcrnrlon > self._bm_lons[-1] or \ @@ -3257,13 +3257,13 @@ if newfile or not hasattr(self,'_bm_rgba_warped'): # transform to nx x ny regularly spaced native # projection grid. - # nx and ny chosen to have roughly the + # nx and ny chosen to have roughly the # same horizontal res as original image. if self.projection != 'cyl': dx = 2.*np.pi*self.rmajor/float(nlons) nx = int((self.xmax-self.xmin)/dx)+1 ny = int((self.ymax-self.ymin)/dx)+1 - else: + else: dx = 360./float(nlons) nx = int((self.urcrnrlon-self.llcrnrlon)/dx)+1 ny = int((self.urcrnrlat-self.llcrnrlat)/dx)+1 @@ -3322,9 +3322,9 @@ fontcolor='k',fillcolor1='w',fillcolor2='k',ax=None,\ format='%d'): """ - Draw a map scale at ``lon,lat`` of length ``length`` + Draw a map scale at ``lon,lat`` of length ``length`` representing distance in the map - projection coordinates at ``lon0,lat0``. + projection coordinates at ``lon0,lat0``. .. tabularcolumns:: |l|L| @@ -3334,24 +3334,24 @@ units the units of the length argument (Default km). barstyle ``simple`` or ``fancy`` (roughly corresponding to the styles provided by Generic Mapping Tools). - Default ``simple``. + Default ``simple``. fontsize for map scale annotations, default 9. color for map scale annotations, default black. labelstype ``simple`` (default) or ``fancy``. For ``fancy`` the map scale factor (ratio betwee the actual distance and map projection distance - at lon0,lat0) and the value of lon0,lat0 are also - displayed on the top of the scale bar. For + at lon0,lat0) and the value of lon0,lat0 are also + displayed on the top of the scale bar. For ``simple``, just the units are display on top and the distance below the scale bar. If equal to False, plot an empty label. format a string formatter to format numeric values yoffset yoffset controls how tall the scale bar is, and how far the annotations are offset from the - scale bar. Default is 0.02 times the height of + scale bar. Default is 0.02 times the height of the map (0.02*(self.ymax-self.ymin)). fillcolor1(2) colors of the alternating filled regions - (default white and black). Only relevant for + (default white and black). Only relevant for 'fancy' barstyle. ============== ==================================================== @@ -3365,7 +3365,7 @@ # convert length to meters lenlab = length if units == 'km': - length = length*1000 + length = length*1000 elif units == 'mi': length = length*1609.344 elif units == 'nmi': @@ -3592,7 +3592,7 @@ ============== ==================================================== Arguments Description ============== ==================================================== - datain a rank-2 array with 1st dimension corresponding to + datain a rank-2 array with 1st dimension corresponding to y, 2nd dimension x. xin, yin rank-1 arrays containing x and y of datain grid in increasing order. @@ -3606,15 +3606,15 @@ ============== ==================================================== checkbounds If True, values of xout and yout are checked to see that they lie within the range specified by xin - and xin. - If False, and xout,yout are outside xin,yin, + and xin. + If False, and xout,yout are outside xin,yin, interpolated values will be clipped to values on boundary of input grid (xin,yin) Default is False. masked If True, points outside the range of xin and yin - are masked (in a masked array). + are masked (in a masked array). If masked is set to a number, then - points outside the range of xin and yin will be + points outside the range of xin and yin will be set to that number. Default False. order 0 for nearest-neighbor interpolation, 1 for bilinear interpolation (default 1). @@ -3827,16 +3827,16 @@ In addition, variables that have the ``scale_factor`` and ``add_offset`` attribute will automatically be converted to and from short integers. To suppress these automatic conversions, set the ``maskandscale`` - keyword to False. + keyword to False. The keywords ``cache``, ``username``, ``password`` and ``verbose`` are only - valid for remote OPenDAP datasets. ``username`` and ``password`` are used + valid for remote OPenDAP datasets. ``username`` and ``password`` are used to access OPenDAP datasets that require authentication. ``verbose=True`` will make the pydap client print out the URLs being accessed. ``cache`` is a location (a directory) for caching data, so that repeated - accesses to the same URL avoid the network. + accesses to the same URL avoid the network. - The keyword ``mmap`` is only valid for local netCDF files. When + The keyword ``mmap`` is only valid for local netCDF files. When ``mmap=True`` (default), the mmap module is used to access the data. This may be slow for very large netCDF variables. """ @@ -3851,8 +3851,8 @@ """ Return datetime objects given numeric time values. The units of the numeric time values are described by the ``units`` argument - and the ``calendar`` keyword. The returned datetime objects represent - UTC with no time-zone offset, even if the specified + and the ``calendar`` keyword. The returned datetime objects represent + UTC with no time-zone offset, even if the specified units contain a time-zone offset. Default behavior is the same as the matplotlib.dates.num2date function @@ -3872,14 +3872,14 @@ ============== ==================================================== Keywords Description ============== ==================================================== - units a string of the form '<time units> since + units a string of the form '<time units> since <reference time>' describing the units and origin of the time coordinate. <time units> can be days, hours, minutes - or seconds. <reference time> is the time origin. + or seconds. <reference time> is the time origin. Default is 'days since 0001-01-01 00:00:00'. calendar describes the calendar used in the time - calculations. All the values currently defined in + calculations. All the values currently defined in the CF metadata convention (http://cf-pcmdi.llnl.gov/documents/cf-conventions/) are supported. @@ -3891,11 +3891,11 @@ Returns a datetime instance, or an array of datetime instances. - The datetime instances returned are 'real' python datetime - objects if the date falls in the Gregorian calendar (i.e. + The datetime instances returned are 'real' python datetime + objects if the date falls in the Gregorian calendar (i.e. calendar=``proleptic_gregorian``, or calendar = ``standard`` or ``gregorian`` and the date is after 1582-10-15). - Otherwise, they are 'phony' datetime + Otherwise, they are 'phony' datetime objects which support some but not all the methods of 'real' python datetime objects. The datetime instances do not contain a time-zone offset, even if the specified units contains one. @@ -3908,7 +3908,7 @@ Return numeric time values given datetime objects. The units of the numeric time values are described by the ``units`` argument and the ``calendar`` keyword. The datetime objects must - be in UTC with no time-zone offset. If there is a + be in UTC with no time-zone offset. If there is a time-zone offset in units, it will be applied to the returned numeric values. @@ -3931,14 +3931,14 @@ ============== ==================================================== Keywords Description ============== ==================================================== - units a string of the form '<time units> since + units a string of the form '<time units> since <reference time>' describing the units and origin of the time coordinate. <time units> can be days, hours, minutes - or seconds. <reference time> is the time origin. + or seconds. <reference time> is the time origin. Default is 'days since 0001-01-01 00:00:00'. calendar describes the calendar used in the time - calculations. All the values currently defined in + calculations. All the values currently defined in the CF metadata convention (http://cf-pcmdi.llnl.gov/documents/cf-conventions/) are supported. @@ -3977,7 +3977,7 @@ Keywords Description ============== ==================================================== calendar describes the calendar used in the time - calculations. All the values currently defined in + calculations. All the values currently defined in the CF metadata convention (http://cf-pcmdi.llnl.gov/documents/cf-conventions/) are supported. @@ -3985,16 +3985,16 @@ ``proleptic_gregorian``, ``noleap``, ``365_day``, ``julian``, ``all_leap``, ``366_day``. If ``calendar=None``, will use ``calendar`` attribute - of ``nctime`` object, and if that attribute does - not exist calendar is set to ``standard``. + of ``nctime`` object, and if that attribute does + not exist calendar is set to ``standard``. Default is ``None``. select The index selection method. ``exact`` will return the - indices perfectly matching the dates given. + indices perfectly matching the dates given. ``before`` and ``after`` will return the indices corresponding to the dates just before or after the given dates if an exact match cannot be found. - ``nearest`` will return the indices that - correspond to the closest dates. Default ``exact``. + ``nearest`` will return the indices that + correspond to the closest dates. Default ``exact``. ============== ==================================================== Returns an index or a sequence of indices. Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py 2009-12-07 01:17:27 UTC (rev 8011) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py 2009-12-08 02:02:49 UTC (rev 8012) @@ -49,7 +49,7 @@ creates a Julian Day from a 'datetime-like' object. Returns the fractional Julian Day (resolution 1 second). -if calendar='standard' or 'gregorian' (default), Julian day follows Julian +if calendar='standard' or 'gregorian' (default), Julian day follows Julian Calendar on and before 1582-10-5, Gregorian calendar after 1582-10-15. if calendar='proleptic_gregorian', Julian Day follows gregorian calendar. @@ -62,7 +62,7 @@ Virginia. p. 63 """ - + # based on redate.py by David Finlayson. year=date.year; month=date.month; day=date.day @@ -74,13 +74,13 @@ if (month < 3): month = month + 12 year = year - 1 - + A = int(year/100) jd = int(365.25 * (year + 4716)) + int(30.6001 * (month + 1)) + \ day - 1524.5 - # optionally adjust the jd for the switch from + # optionally adjust the jd for the switch from # the Julian to Gregorian Calendar # here assumed to have occurred the day after 1582 October 4 if calendar in ['standard','gregorian']: @@ -98,21 +98,21 @@ B = 0 else: raise ValueError, 'unknown calendar, must be one of julian,standard,gregorian,proleptic_gregorian, got %s' % calendar - + # adjust for Julian calendar if necessary jd = jd + B - - return jd + return jd + def _NoLeapDayFromDate(date): """ -creates a Julian Day for a calendar with no leap years from a datetime +creates a Julian Day for a calendar with no leap years from a datetime instance. Returns the fractional Julian Day (resolution 1 second). """ - + year=date.year; month=date.month; day=date.day hour=date.hour; minute=date.minute; second=date.second # Convert time to fractions of a day @@ -122,12 +122,12 @@ if (month < 3): month = month + 12 year = year - 1 - + jd = int(365. * (year + 4716)) + int(30.6001 * (month + 1)) + \ day - 1524.5 - - return jd + return jd + def _AllLeapFromDate(date): """ @@ -137,7 +137,7 @@ Returns the fractional Julian Day (resolution 1 second). """ - + year=date.year; month=date.month; day=date.day hour=date.hour; minute=date.minute; second=date.second # Convert time to fractions of a day @@ -147,12 +147,12 @@ if (month < 3): month = month + 12 year = year - 1 - + jd = int(366. * (year + 4716)) + int(30.6001 * (month + 1)) + \ day - 1524.5 - - return jd + return jd + def _360DayFromDate(date): """ @@ -162,23 +162,23 @@ Returns the fractional Julian Day (resolution 1 second). """ - + year=date.year; month=date.month; day=date.day hour=date.hour; minute=date.minute; second=date.second # Convert time to fractions of a day day = day + hour/24.0 + minute/1440.0 + second/86400.0 jd = int(360. * (year + 4716)) + int(30. * (month - 1)) + day - - return jd + return jd + def DateFromJulianDay(JD,calendar='standard'): """ -returns a 'datetime-like' object given Julian Day. Julian Day is a +returns a 'datetime-like' object given Julian Day. Julian Day is a fractional day with a resolution of 1 second. -if calendar='standard' or 'gregorian' (default), Julian day follows Julian +if calendar='standard' or 'gregorian' (default), Julian day follows Julian Calendar on and before 1582-10-5, Gregorian calendar after 1582-10-15. if calendar='proleptic_gregorian', Julian Day follows gregorian calendar. @@ -200,7 +200,7 @@ """ # based on redate.py by David Finlayson. - + if JD < 0: raise ValueError, 'Julian Day must be positive' @@ -250,18 +250,18 @@ leap = 1 if calendar == 'proleptic_gregorian' or \ (calendar in ['standard','gregorian'] and JD >= 2299160.5): - if year % 100 == 0 and year % 400 != 0: + if year % 100 == 0 and year % 400 != 0: print year % 100, year % 400 leap = 0 if leap and month > 2: dayofyr = dayofyr + leap - - # Convert fractions of a day to time + + # Convert fractions of a day to time (dfrac, days) = math.modf(day/1.0) (hfrac, hours) = math.modf(dfrac * 24.0) (mfrac, minutes) = math.modf(hfrac * 60.0) seconds = round(mfrac * 60.0) # seconds are rounded - + if seconds > 59: seconds = 0 minutes = minutes + 1 @@ -271,7 +271,7 @@ if hours > 23: hours = 0 days = days + 1 - + # return a 'real' datetime instance if calendar is gregorian. if calendar == 'proleptic_gregorian' or \ (calendar in ['standard','gregorian'] and JD >= 2299160.5): @@ -283,13 +283,13 @@ def _DateFromNoLeapDay(JD): """ -returns a 'datetime-like' object given Julian Day for a calendar with no leap +returns a 'datetime-like' object given Julian Day for a calendar with no leap days. Julian Day is a fractional day with a resolution of 1 second. """ # based on redate.py by David Finlayson. - + if JD < 0: raise ValueError, 'Julian Day must be positive' @@ -318,13 +318,13 @@ year = C - 4716 else: year = C - 4715 - - # Convert fractions of a day to time + + # Convert fractions of a day to time (dfrac, days) = math.modf(day/1.0) (hfrac, hours) = math.modf(dfrac * 24.0) (mfrac, minutes) = math.modf(hfrac * 60.0) seconds = round(mfrac * 60.0) # seconds are rounded - + if seconds > 59: seconds = 0 minutes = minutes + 1 @@ -334,7 +334,7 @@ if hours > 23: hours = 0 days = days + 1 - + return datetime(year,month,int(days),int(hours),int(minutes),int(seconds), dayofwk, dayofyr) def _DateFromAllLeap(JD): @@ -347,7 +347,7 @@ """ # based on redate.py by David Finlayson. - + if JD < 0: raise ValueError, 'Julian Day must be positive' @@ -378,13 +378,13 @@ year = C - 4716 else: year = C - 4715 - - # Convert fractions of a day to time + + # Convert fractions of a day to time (dfrac, days) = math.modf(day/1.0) (hfrac, hours) = math.modf(dfrac * 24.0) (mfrac, minutes) = math.modf(hfrac * 60.0) seconds = round(mfrac * 60.0) # seconds are rounded - + if seconds > 59: seconds = 0 minutes = minutes + 1 @@ -394,7 +394,7 @@ if hours > 23: hours = 0 days = days + 1 - + return datetime(year,month,int(days),int(hours),int(minutes),int(seconds), dayofwk, dayofyr) def _DateFrom360Day(JD): @@ -412,16 +412,16 @@ #jd = int(360. * (year + 4716)) + int(30. * (month - 1)) + day (F, Z) = math.modf(JD) year = int((Z-0.5)/360.) - 4716 - dayofyr = JD - (year+4716)*360 + dayofyr = JD - (year+4716)*360 month = int((dayofyr-0.5)/30)+1 - day = dayofyr - (month-1)*30 + F - - # Convert fractions of a day to time + day = dayofyr - (month-1)*30 + F + + # Convert fractions of a day to time (dfrac, days) = math.modf(day/1.0) (hfrac, hours) = math.modf(dfrac * 24.0) (mfrac, minutes) = math.modf(hfrac * 60.0) seconds = round(mfrac * 60.0) # seconds are rounded - + if seconds > 59: seconds = 0 minutes = minutes + 1 @@ -431,7 +431,7 @@ if hours > 23: hours = 0 days = days + 1 - + return datetime(year,month,int(days),int(hours),int(minutes),int(seconds),-1, int(dayofyr)) def _dateparse(timestr): @@ -455,20 +455,20 @@ To initialize: C{t = utime(unit_string,calendar='standard')} -where +where B{C{unit_string}} is a string of the form C{'time-units since <time-origin>'} defining the time units. -Valid time-units are days, hours, minutes and seconds (the singular forms -are also accepted). An example unit_string would be C{'hours +Valid time-units are days, hours, minutes and seconds (the singular forms +are also accepted). An example unit_string would be C{'hours since 0001... [truncated message content] |
From: <js...@us...> - 2011-02-12 13:56:51
|
Revision: 8977 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=8977&view=rev Author: jswhit Date: 2011-02-12 13:56:44 +0000 (Sat, 12 Feb 2011) Log Message: ----------- remove more vestiges of netcdf support Modified Paths: -------------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py Removed Paths: ------------- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2011-02-12 13:54:34 UTC (rev 8976) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/__init__.py 2011-02-12 13:56:44 UTC (rev 8977) @@ -11,12 +11,6 @@ :func:`shiftgrid`: shifts global lat/lon grids east or west. :func:`addcyclic`: Add cyclic (wraparound) point in longitude. - -:func:`num2date`: convert from a numeric time value to a datetime object. - -:func:`date2num`: convert from a datetime object to a numeric time value. - -:func:`date2index`: compute a time variable index corresponding to a date. """ from matplotlib import __version__ as _matplotlib_version from matplotlib.cbook import is_scalar, dedent @@ -38,7 +32,7 @@ import numpy as np import numpy.ma as ma from shapelib import ShapeFile -import _geoslib, netcdftime +import _geoslib # basemap data files now installed in lib/matplotlib/toolkits/basemap/data # check to see if environment variable BASEMAPDATA set to a directory, @@ -3974,160 +3968,6 @@ else: return corners -def num2date(times,units='days since 0001-01-01 00:00:00',calendar='proleptic_gregorian'): - """ - Return datetime objects given numeric time values. The units - of the numeric time values are described by the ``units`` argument - and the ``calendar`` keyword. The returned datetime objects represent - UTC with no time-zone offset, even if the specified - units contain a time-zone offset. - - Default behavior is the same as the matplotlib.dates.num2date function - but the reference time and calendar can be changed via the - ``units`` and ``calendar`` keywords. - - .. tabularcolumns:: |l|L| - - ============== ==================================================== - Arguments Description - ============== ==================================================== - times numeric time values. Maximum resolution is 1 second. - ============== ==================================================== - - .. tabularcolumns:: |l|L| - - ============== ==================================================== - Keywords Description - ============== ==================================================== - units a string of the form '<time units> since - <reference time>' describing the units and - origin of the time coordinate. - <time units> can be days, hours, minutes - or seconds. <reference time> is the time origin. - Default is 'days since 0001-01-01 00:00:00'. - calendar describes the calendar used in the time - calculations. All the values currently defined in - the CF metadata convention - (http://cf-pcmdi.llnl.gov/documents/cf-conventions/) - are supported. - Valid calendars ``standard``, ``gregorian``, - ``proleptic_gregorian``, ``noleap``, ``365_day``, - ``julian``, ``all_leap``, ``366_day``. - Default is ``proleptic_gregorian``. - ============== ==================================================== - - Returns a datetime instance, or an array of datetime instances. - - The datetime instances returned are 'real' python datetime - objects if the date falls in the Gregorian calendar (i.e. - calendar=``proleptic_gregorian``, or calendar = ``standard`` - or ``gregorian`` and the date is after 1582-10-15). - Otherwise, they are 'phony' datetime - objects which support some but not all the methods of 'real' python - datetime objects. The datetime instances do not contain - a time-zone offset, even if the specified units contains one. - """ - cdftime = netcdftime.utime(units,calendar=calendar) - return cdftime.num2date(times) - -def date2num(dates,units='days since 0001-01-01 00:00:00',calendar='proleptic_gregorian'): - """ - Return numeric time values given datetime objects. The units - of the numeric time values are described by the ``units`` argument - and the ``calendar`` keyword. The datetime objects must - be in UTC with no time-zone offset. If there is a - time-zone offset in units, it will be applied to the - returned numeric values. - - Default behavior is the same as the matplotlib.dates.date2num function - but the reference time and calendar can be changed via the - ``units`` and ``calendar`` keywords. - - .. tabularcolumns:: |l|L| - - ============== ==================================================== - Arguments Description - ============== ==================================================== - dates A datetime object or a sequence of datetime objects. - The datetime objects should not include a - time-zone offset. - ============== ==================================================== - - .. tabularcolumns:: |l|L| - - ============== ==================================================== - Keywords Description - ============== ==================================================== - units a string of the form '<time units> since - <reference time>' describing the units and - origin of the time coordinate. - <time units> can be days, hours, minutes - or seconds. <reference time> is the time origin. - Default is 'days since 0001-01-01 00:00:00'. - calendar describes the calendar used in the time - calculations. All the values currently defined in - the CF metadata convention - (http://cf-pcmdi.llnl.gov/documents/cf-conventions/) - are supported. - Valid calendars ``standard``, ``gregorian``, - ``proleptic_gregorian``, ``noleap``, ``365_day``, - ``julian``, ``all_leap``, ``366_day``. - Default is ``proleptic_gregorian``. - ============== ==================================================== - - Returns a numeric time value, or an array of numeric time values. - - The maximum resolution of the numeric time values is 1 second. - """ - cdftime = netcdftime.utime(units,calendar=calendar) - return cdftime.date2num(dates) - -def date2index(dates, nctime, calendar=None,select='exact'): - """ - Return indices of a netCDF time variable corresponding to the given dates. - - .. tabularcolumns:: |l|L| - - ============== ==================================================== - Arguments Description - ============== ==================================================== - dates A datetime object or a sequence of datetime objects. - The datetime objects should not include a - time-zone offset. - nctime A netCDF time variable object. The nctime object - must have a ``units`` attribute. - ============== ==================================================== - - .. tabularcolumns:: |l|L| - - ============== ==================================================== - Keywords Description - ============== ==================================================== - calendar describes the calendar used in the time - calculations. All the values currently defined in - the CF metadata convention - (http://cf-pcmdi.llnl.gov/documents/cf-conventions/) - are supported. - Valid calendars ``standard``, ``gregorian``, - ``proleptic_gregorian``, ``noleap``, ``365_day``, - ``julian``, ``all_leap``, ``366_day``. - If ``calendar=None``, will use ``calendar`` attribute - of ``nctime`` object, and if that attribute does - not exist calendar is set to ``standard``. - Default is ``None``. - select The index selection method. ``exact`` will return the - indices perfectly matching the dates given. - ``before`` and ``after`` will return the indices - corresponding to the dates just before or after - the given dates if an exact match cannot be found. - ``nearest`` will return the indices that - correspond to the closest dates. Default ``exact``. - ============== ==================================================== - - Returns an index or a sequence of indices. - """ - return netcdftime.date2index(dates, nctime, calendar=calendar, select=select) - def maskoceans(lonsin,latsin,datain,inlands=False): """ mask data (``datain``), defined on a grid with latitudes ``latsin`` Deleted: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py 2011-02-12 13:54:34 UTC (rev 8976) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdf.py 2011-02-12 13:56:44 UTC (rev 8977) @@ -1,81 +0,0 @@ -from numpy import ma, squeeze -from pupynere import netcdf_file, _unmaskandscale -from dap.client import open as open_remote -from dap.dtypes import ArrayType, GridType, typemap - -_typecodes = dict([[_v,_k] for _k,_v in typemap.items()]) - -class _RemoteFile(object): - """A NetCDF file reader. API is the same as Scientific.IO.NetCDF.""" - - def __init__(self, file, maskandscale=False, cache=None,\ - username=None, password=None, verbose=False): - self._buffer = open_remote(file,cache=cache,\ - username=username,password=password,verbose=verbose) - self._maskandscale = maskandscale - self._parse() - - def read(self, size=-1): - """Alias for reading the file buffer.""" - return self._buffer.read(size) - - def _parse(self): - """Initial parsing of the header.""" - # Read header info. - self._dim_array() - self._gatt_array() - self._var_array() - - def _dim_array(self): - """Read a dict with dimensions names and sizes.""" - self.dimensions = {} - self._dims = [] - for k,d in self._buffer.iteritems(): - if (isinstance(d, ArrayType) or isinstance(d, GridType)) and len(d.shape) == 1 and k == d.dimensions[0]: - name = k - length = len(d) - self.dimensions[name] = length - self._dims.append(name) # preserve dim order - - def _gatt_array(self): - """Read global attributes.""" - self.__dict__.update(self._buffer.attributes) - - def _var_array(self): - """Read all variables.""" - # Read variables. - self.variables = {} - for k,d in self._buffer.iteritems(): - if isinstance(d, GridType) or isinstance(d, ArrayType): - name = k - self.variables[name] = _RemoteVariable(d,self._maskandscale) - - def close(self): - # this is a no-op provided for compatibility - pass - -class _RemoteVariable(object): - def __init__(self, var, maskandscale): - self._var = var - self._maskandscale = maskandscale - self.dtype = var.type - self.shape = var.shape - self.dimensions = var.dimensions - self.__dict__.update(var.attributes) - - def __getitem__(self, index): - datout = squeeze(self._var.__getitem__(index)) - # automatically - # - remove singleton dimensions - # - create a masked array using missing_value or _FillValue attribute - # - apply scale_factor and add_offset to packed integer data - if self._maskandscale: - return _unmaskandscale(self,datout) - else: - return datout - - def __len__(self): - return self.shape[0] - - def typecode(self): - return _typecodes[self.dtype] Deleted: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py =================================================================== --- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py 2011-02-12 13:54:34 UTC (rev 8976) +++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/netcdftime.py 2011-02-12 13:56:44 UTC (rev 8977) @@ -1,1050 +0,0 @@ -""" -Performs conversions of netCDF time coordinate data to/from datetime objects. -""" -import math, numpy, re, time -from datetime import datetime as real_datetime - -_units = ['days','hours','minutes','seconds','day','hour','minute','second'] -_calendars = ['standard','gregorian','proleptic_gregorian','noleap','julian','all_leap','365_day','366_day','360_day'] - -__version__ = '0.7.1' - -class datetime: - """ -Phony datetime object which mimics the python datetime object, -but allows for dates that don't exist in the proleptic gregorian calendar. -Doesn't do timedelta operations, doesn't overload + and -. - -Has strftime, timetuple and __repr__ methods. The format -of the string produced by __repr__ is controlled by self.format -(default %Y-%m-%d %H:%M:%S). - -Instance variables are year,month,day,hour,minute,second,dayofwk,dayofyr -and format. - """ - def __init__(self,year,month,day,hour=0,minute=0,second=0,dayofwk=-1,dayofyr=1): - """dayofyr set to 1 by default - otherwise time.strftime will complain""" - self.year=year - self.month=month - self.day=day - self.hour=hour - self.minute=minute - self.dayofwk=dayofwk - self.dayofyr=dayofyr - self.second=second - self.format='%Y-%m-%d %H:%M:%S' - def strftime(self,format=None): - if format is None: - format = self.format - return _strftime(self,format) - def timetuple(self): - return (self.year,self.month,self.day,self.hour,self.minute,self.second,self.dayofwk,self.dayofyr,-1) - def __repr__(self): - return self.strftime(self.format) - -def JulianDayFromDate(date,calendar='standard'): - - """ - -creates a Julian Day from a 'datetime-like' object. Returns the fractional -Julian Day (resolution 1 second). - -if calendar='standard' or 'gregorian' (default), Julian day follows Julian -Calendar on and before 1582-10-5, Gregorian calendar after 1582-10-15. - -if calendar='proleptic_gregorian', Julian Day follows gregorian calendar. - -if calendar='julian', Julian Day follows julian calendar. - -Algorithm: - -Meeus, Jean (1998) Astronomical Algorithms (2nd Edition). Willmann-Bell, -Virginia. p. 63 - - """ - - # based on redate.py by David Finlayson. - - year=date.year; month=date.month; day=date.day - hour=date.hour; minute=date.minute; second=date.second - # Convert time to fractions of a day - day = day + hour/24.0 + minute/1440.0 + second/86400.0 - - # Start Meeus algorithm (variables are in his notation) - if (month < 3): - month = month + 12 - year = year - 1 - - A = int(year/100) - - jd = int(365.25 * (year + 4716)) + int(30.6001 * (month + 1)) + \ - day - 1524.5 - - # optionally adjust the jd for the switch from - # the Julian to Gregorian Calendar - # here assumed to have occurred the day after 1582 October 4 - if calendar in ['standard','gregorian']: - if jd >= 2299170.5: - # 1582 October 15 (Gregorian Calendar) - B = 2 - A + int(A/4) - elif jd < 2299160.5: - # 1582 October 5 (Julian Calendar) - B = 0 - else: - raise ValueError, 'impossible date (falls in gap between end of Julian calendar and beginning of Gregorian calendar' - elif calendar == 'proleptic_gregorian': - B = 2 - A + int(A/4) - elif calendar == 'julian': - B = 0 - else: - raise ValueError, 'unknown calendar, must be one of julian,standard,gregorian,proleptic_gregorian, got %s' % calendar - - # adjust for Julian calendar if necessary - jd = jd + B - - return jd - -def _NoLeapDayFromDate(date): - - """ - -creates a Julian Day for a calendar with no leap years from a datetime -instance. Returns the fractional Julian Day (resolution 1 second). - - """ - - year=date.year; month=date.month; day=date.day - hour=date.hour; minute=date.minute; second=date.second - # Convert time to fractions of a day - day = day + hour/24.0 + minute/1440.0 + second/86400.0 - - # Start Meeus algorithm (variables are in his notation) - if (month < 3): - month = month + 12 - year = year - 1 - - jd = int(365. * (year + 4716)) + int(30.6001 * (month + 1)) + \ - day - 1524.5 - - return jd - -def _AllLeapFromDate(date): - - """ - -creates a Julian Day for a calendar where all years have 366 days from -a 'datetime-like' object. -Returns the fractional Julian Day (resolution 1 second). - - """ - - year=date.year; month=date.month; day=date.day - hour=date.hour; minute=date.minute; second=date.second - # Convert time to fractions of a day - day = day + hour/24.0 + minute/1440.0 + second/86400.0 - - # Start Meeus algorithm (variables are in his notation) - if (month < 3): - month = month + 12 - year = year - 1 - - jd = int(366. * (year + 4716)) + int(30.6001 * (month + 1)) + \ - day - 1524.5 - - return jd - -def _360DayFromDate(date): - - """ - -creates a Julian Day for a calendar where all months have 30 daysfrom -a 'datetime-like' object. -Returns the fractional Julian Day (resolution 1 second). - - """ - - year=date.year; month=date.month; day=date.day - hour=date.hour; minute=date.minute; second=date.second - # Convert time to fractions of a day - day = day + hour/24.0 + minute/1440.0 + second/86400.0 - - jd = int(360. * (year + 4716)) + int(30. * (month - 1)) + day - - return jd - -def DateFromJulianDay(JD,calendar='standard'): - """ - -returns a 'datetime-like' object given Julian Day. Julian Day is a -fractional day with a resolution of 1 second. - -if calendar='standard' or 'gregorian' (default), Julian day follows Julian -Calendar on and before 1582-10-5, Gregorian calendar after 1582-10-15. - -if calendar='proleptic_gregorian', Julian Day follows gregorian calendar. - -if calendar='julian', Julian Day follows julian calendar. - -The datetime object is a 'real' datetime object if the date falls in -the Gregorian calendar (i.e. calendar='proleptic_gregorian', or -calendar = 'standard'/'gregorian' and the date is after 1582-10-15). -Otherwise, it's a 'phony' datetime object which is actually an instance -of netcdftime.datetime. - - -Algorithm: - -Meeus, Jean (1998) Astronomical Algorithms (2nd Edition). Willmann-Bell, -Virginia. p. 63 - - """ - - # based on redate.py by David Finlayson. - - if JD < 0: - raise ValueError, 'Julian Day must be positive' - - dayofwk = int(math.fmod(int(JD + 1.5),7)) - (F, Z) = math.modf(JD + 0.5) - Z = int(Z) - if calendar in ['standard','gregorian']: - if JD < 2299160.5: - A = Z - else: - alpha = int((Z - 1867216.25)/36524.25) - A = Z + 1 + alpha - int(alpha/4) - - elif calendar == 'proleptic_gregorian': - alpha = int((Z - 1867216.25)/36524.25) - A = Z + 1 + alpha - int(alpha/4) - elif calendar == 'julian': - A = Z - else: - raise ValueError, 'unknown calendar, must be one of julian,standard,gregorian,proleptic_gregorian, got %s' % calendar - - B = A + 1524 - C = int((B - 122.1)/365.25) - D = int(365.25 * C) - E = int((B - D)/30.6001) - - # Convert to date - day = B - D - int(30.6001 * E) + F - nday = B-D-123 - if nday <= 305: - dayofyr = nday+60 - else: - dayofyr = nday-305 - if E < 14: - month = E - 1 - else: - month = E - 13 - - if month > 2: - year = C - 4716 - else: - year = C - 4715 - - # a leap year? - leap = 0 - if year % 4 == 0: - leap = 1 - if calendar == 'proleptic_gregorian' or \ - (calendar in ['standard','gregorian'] and JD >= 2299160.5): - if year % 100 == 0 and year % 400 != 0: - print year % 100, year % 400 - leap = 0 - if leap and month > 2: - dayofyr = dayofyr + leap - - # Convert fractions of a day to time - (dfrac, days) = math.modf(day/1.0) - (hfrac, hours) = math.modf(dfrac * 24.0) - (mfrac, minutes) = math.modf(hfrac * 60.0) - seconds = round(mfrac * 60.0) # seconds are rounded - - if seconds > 59: - seconds = 0 - minutes = minutes + 1 - if minutes > 59: - minutes = 0 - hours = hours + 1 - if hours > 23: - hours = 0 - days = days + 1 - - # return a 'real' datetime instance if calendar is gregorian. - if calendar == 'proleptic_gregorian' or \ - (calendar in ['standard','gregorian'] and JD >= 2299160.5): - return real_datetime(year,month,int(days),int(hours),int(minutes),int(seconds)) - else: - # or else, return a 'datetime-like' instance. - return datetime(year,month,int(days),int(hours),int(minutes),int(seconds),dayofwk,dayofyr) - -def _DateFromNoLeapDay(JD): - """ - -returns a 'datetime-like' object given Julian Day for a calendar with no leap -days. Julian Day is a fractional day with a resolution of 1 second. - - """ - - # based on redate.py by David Finlayson. - - if JD < 0: - raise ValueError, 'Julian Day must be positive' - - dayofwk = int(math.fmod(int(JD + 1.5),7)) - (F, Z) = math.modf(JD + 0.5) - Z = int(Z) - A = Z - B = A + 1524 - C = int((B - 122.1)/365.) - D = int(365. * C) - E = int((B - D)/30.6001) - - # Convert to date - day = B - D - int(30.6001 * E) + F - nday = B-D-123 - if nday <= 305: - dayofyr = nday+60 - else: - dayofyr = nday-305 - if E < 14: - month = E - 1 - else: - month = E - 13 - - if month > 2: - year = C - 4716 - else: - year = C - 4715 - - # Convert fractions of a day to time - (dfrac, days) = math.modf(day/1.0) - (hfrac, hours) = math.modf(dfrac * 24.0) - (mfrac, minutes) = math.modf(hfrac * 60.0) - seconds = round(mfrac * 60.0) # seconds are rounded - - if seconds > 59: - seconds = 0 - minutes = minutes + 1 - if minutes > 59: - minutes = 0 - hours = hours + 1 - if hours > 23: - hours = 0 - days = days + 1 - - return datetime(year,month,int(days),int(hours),int(minutes),int(seconds), dayofwk, dayofyr) - -def _DateFromAllLeap(JD): - """ - -returns a 'datetime-like' object given Julian Day for a calendar where all -years have 366 days. -Julian Day is a fractional day with a resolution of 1 second. - - """ - - # based on redate.py by David Finlayson. - - if JD < 0: - raise ValueError, 'Julian Day must be positive' - - dayofwk = int(math.fmod(int(JD + 1.5),7)) - (F, Z) = math.modf(JD + 0.5) - Z = int(Z) - A = Z - B = A + 1524 - C = int((B - 122.1)/366.) - D = int(366. * C) - E = int((B - D)/30.6001) - - # Convert to date - day = B - D - int(30.6001 * E) + F - nday = B-D-123 - if nday <= 305: - dayofyr = nday+60 - else: - dayofyr = nday-305 - if E < 14: - month = E - 1 - else: - month = E - 13 - if month > 2: - dayofyr = dayofyr+1 - - if month > 2: - year = C - 4716 - else: - year = C - 4715 - - # Convert fractions of a day to time - (dfrac, days) = math.modf(day/1.0) - (hfrac, hours) = math.modf(dfrac * 24.0) - (mfrac, minutes) = math.modf(hfrac * 60.0) - seconds = round(mfrac * 60.0) # seconds are rounded - - if seconds > 59: - seconds = 0 - minutes = minutes + 1 - if minutes > 59: - minutes = 0 - hours = hours + 1 - if hours > 23: - hours = 0 - days = days + 1 - - return datetime(year,month,int(days),int(hours),int(minutes),int(seconds), dayofwk, dayofyr) - -def _DateFrom360Day(JD): - """ - -returns a 'datetime-like' object given Julian Day for a calendar where all -months have 30 days. -Julian Day is a fractional day with a resolution of 1 second. - - """ - - if JD < 0: - raise ValueError, 'Julian Day must be positive' - - #jd = int(360. * (year + 4716)) + int(30. * (month - 1)) + day - (F, Z) = math.modf(JD) - year = int((Z-0.5)/360.) - 4716 - dayofyr = JD - (year+4716)*360 - month = int((dayofyr-0.5)/30)+1 - day = dayofyr - (month-1)*30 + F - - # Convert fractions of a day to time - (dfrac, days) = math.modf(day/1.0) - (hfrac, hours) = math.modf(dfrac * 24.0) - (mfrac, minutes) = math.modf(hfrac * 60.0) - seconds = round(mfrac * 60.0) # seconds are rounded - - if seconds > 59: - seconds = 0 - minutes = minutes + 1 - if minutes > 59: - minutes = 0 - hours = hours + 1 - if hours > 23: - hours = 0 - days = days + 1 - - return datetime(year,month,int(days),int(hours),int(minutes),int(seconds),-1, int(dayofyr)) - -def _dateparse(timestr): - """parse a string of the form time-units since yyyy-mm-dd hh:mm:ss - return a tuple (units, datetimeinstance)""" - timestr_split = timestr.split() - units = timestr_split[0].lower() - if units not in _units: - raise ValueError,"units must be one of 'seconds', 'minutes', 'hours' or 'days' (or singular version of these), got '%s'" % units - if timestr_split[1].lower() != 'since': - raise ValueError,"no 'since' in unit_string" - # parse the date string. - n = timestr.find('since')+6 - year,month,day,hour,minute,second,utc_offset = _parse_date(timestr[n:]) - return units, utc_offset, datetime(year, month, day, hour, minute, second) - -class utime: - """ -Performs conversions of netCDF time coordinate -data to/from datetime objects. - -To initialize: C{t = utime(unit_string,calendar='standard')} - -where - -B{C{unit_string}} is a string of the form -C{'time-units since <time-origin>'} defining the time units. - -Valid time-units are days, hours, minutes and seconds (the singular forms -are also accepted). An example unit_string would be C{'hours -since 0001-01-01 00:00:00'}. - -The B{C{calendar}} keyword describes the calendar used in the time calculations. -All the values currently defined in the U{CF metadata convention -<http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.1/cf-conventions.html#time-coordinate>} -are accepted. The default is C{'standard'}, which corresponds to the mixed -Gregorian/Julian calendar used by the C{udunits library}. Valid calendars -are: - -C{'gregorian'} or C{'standard'} (default): - -Mixed Gregorian/Julian calendar as defined by udunits. - -C{'proleptic_gregorian'}: - -A Gregorian calendar extended to dates before 1582-10-15. That is, a year -is a leap year if either (i) it is divisible by 4 but not by 100 or (ii) -it is divisible by 400. - -C{'noleap'} or C{'365_day'}: - -Gregorian calendar without leap years, i.e., all years are 365 days long. -all_leap or 366_day Gregorian calendar with every year being a leap year, -i.e., all years are 366 days long. - -C{'360_day'}: - -All years are 360 days divided into 30 day months. - -C{'julian'}: - -Proleptic Julian calendar, extended to dates after 1582-10-5. A year is a -leap year if it is divisible by 4. - -The C{L{num2date}} and C{L{date2num}} class methods can used to convert datetime -instances to/from the specified time units using the specified calendar. - -The datetime instances returned by C{num2date} are 'real' python datetime -objects if the date falls in the Gregorian calendar (i.e. -C{calendar='proleptic_gregorian', 'standard'} or C{'gregorian'} and -the date is after 1582-10-15). Otherwise, they are 'phony' datetime -objects which are actually instances of C{L{netcdftime.datetime}}. This is -because the python datetime module cannot handle the weird dates in some -calendars (such as C{'360_day'} and C{'all_leap'}) which don't exist in any real -world calendar. - - -Example usage: - ->>> from netcdftime import utime ->>> from datetime import datetime ->>> cdftime = utime('hours since 0001-01-01 00:00:00') ->>> date = datetime.now() ->>> print date -2006-03-17 16:04:02.561678 ->>> ->>> t = cdftime.date2num(date) ->>> print t -17577328.0672 ->>> ->>> date = cdftime.num2date(t) ->>> print date -2006-03-17 16:04:02 ->>> - -The resolution of the transformation operation is 1 second. - -Warning: Dates between 1582-10-5 and 1582-10-15 do not exist in the -C{'standard'} or C{'gregorian'} calendars. An exception will be raised if you pass -a 'datetime-like' object in that range to the C{L{date2num}} class method. - -Words of Wisdom from the British MetOffice concerning reference dates: - -"udunits implements the mixed Gregorian/Julian calendar system, as -followed in England, in which dates prior to 1582-10-15 are assumed to use -the Julian calendar. Other software cannot be relied upon to handle the -change of calendar in the same way, so for robustness it is recommended -that the reference date be later than 1582. If earlier dates must be used, -it should be noted that udunits treats 0 AD as identical to 1 AD." - -@ivar origin: datetime instance defining the origin of the netCDF time variable. -@ivar calendar: the calendar used (as specified by the C{calendar} keyword). -@ivar unit_string: a string defining the the netCDF time variable. -@ivar units: the units part of C{unit_string} (i.e. 'days', 'hours', 'seconds'). - """ - def __init__(self,unit_string,calendar='standard'): - """ -@param unit_string: a string of the form -C{'time-units since <time-origin>'} defining the time units. - -Valid time-units are days, hours, minutes and seconds (the singular forms -are also accepted). An example unit_string would be C{'hours -since 0001-01-01 00:00:00'}. - -@keyword calendar: describes the calendar used in the time calculations. -All the values currently defined in the U{CF metadata convention -<http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.1/cf-conventions.html#time-coordinate>} -are accepted. The default is C{'standard'}, which corresponds to the mixed -Gregorian/Julian calendar used by the C{udunits library}. Valid calendars -are: - - C{'gregorian'} or C{'standard'} (default): - Mixed Gregorian/Julian calendar as defined by udunits. - - C{'proleptic_gregorian'}: - A Gregorian calendar extended to dates before 1582-10-15. That is, a year - is a leap year if either (i) it is divisible by 4 but not by 100 or (ii) - it is divisible by 400. - - C{'noleap'} or C{'365_day'}: - Gregorian calendar without leap years, i.e., all years are 365 days long. - - C{'all_leap'} or C{'366_day'}: - Gregorian calendar with every year being a leap year, i.e., - all years are 366 days long. - -C{'360_day'}: - All years are 360 days divided into 30 day months. - -C{'julian'}: - Proleptic Julian calendar, extended to dates after 1582-10-5. A year is a - leap year if it is divisible by 4. - -@returns: A class instance which may be used for converting times from netCDF -units to datetime objects. - """ - if calendar in _calendars: - self.calendar = calendar - else: - raise ValueError, "calendar must be one of %s, got '%s'" % (str(_calendars),calendar) - units, tzoffset, self.origin = _dateparse(unit_string) - self.tzoffset = tzoffset # time zone offset in minutes - self.units = units - self.unit_string = unit_string - if self.calendar in ['noleap','365_day'] and self.origin.month == 2 and self.origin.day == 29: - raise ValueError, 'cannot specify a leap day as the reference time with the noleap calendar' - if self.calendar == '360_day' and self.origin.day > 30: - raise ValueError, 'there are only 30 days in every month with the 360_day calendar' - if self.calendar in ['noleap','365_day']: - self._jd0 = _NoLeapDayFromDate(self.origin) - elif self.calendar in ['all_leap','366_day']: - self._jd0 = _AllLeapFromDate(self.origin) - elif self.calendar == '360_day': - self._jd0 = _360DayFromDate(self.origin) - else: - self._jd0 = JulianDayFromDate(self.origin,calendar=self.calendar) - - def date2num(self,date): - """ -Returns C{time_value} in units described by L{unit_string}, using -the specified L{calendar}, given a 'datetime-like' object. - -The datetime object must represent UTC with no time-zone offset. -If there is a time-zone offset implied by L{unit_string}, it will -be applied to the returned numeric values. - -Resolution is 1 second. - -If C{calendar = 'standard'} or C{'gregorian'} (indicating -that the mixed Julian/Gregorian calendar is to be used), an -exception will be raised if the 'datetime-like' object describes -a date between 1582-10-5 and 1582-10-15. - -Works for scalars, sequences and numpy arrays. -Returns a scalar if input is a scalar, else returns a numpy array. - """ - isscalar = False - try: - date[0] - except: - isscalar = True - if not isscalar: - date = numpy.array(date) - shape = date.shape - if self.calendar in ['julian','standard','gregorian','proleptic_gregorian']: - if isscalar: - jdelta = JulianDayFromDate(date,self.calendar)-self._jd0 - else: - jdelta = [JulianDayFromDate(d,self.calendar)-self._jd0 for d in date.flat] - elif self.calendar in ['noleap','365_day']: - if date.month == 2 and date.day == 29: - raise ValueError, 'there is no leap day in the noleap calendar' - if isscalar: - jdelta = _NoLeapDayFromDate(date) - self._jd0 - else: - jdelta = [_NoLeapDayFromDate(d)-self._jd0 for d in date.flat] - elif self.calendar in ['all_leap','366_day']: - if isscalar: - jdelta = _AllLeapFromDate(date) - self._jd0 - else: - jdelta = [_AllLeapFromDate(d)-self._jd0 for d in date.flat] - elif self.calendar == '360_day': - if self.calendar == '360_day' and date.day > 30: - raise ValueError, 'there are only 30 days in every month with the 360_day calendar' - if isscalar: - jdelta = _360DayFromDate(date) - self._jd0 - else: - jdelta = [_360DayFromDate(d)-self._jd0 for d in date.flat] - if not isscalar: - jdelta = numpy.array(jdelta) - # convert to desired units, add time zone offset. - if self.units in ['second','seconds']: - jdelta = jdelta*86400. + self.tzoffset*60. - elif self.units in ['minute','minutes']: - jdelta = jdelta*1440. + self.tzoffset - elif self.units in ['hour','hours']: - jdelta = jdelta*24. + self.tzoffset/60. - elif self.units in ['day','days']: - jdelta = jdelta + self.tzoffset/1440. - if isscalar: - return jdelta - else: - return numpy.reshape(jdelta,shape) - - def num2date(self,time_value): - """ -Return a 'datetime-like' object given a C{time_value} in units -described by L{unit_string}, using L{calendar}. - -dates are in UTC with no offset, even if L{unit_string} contains -a time zone offset from UTC. - -Resolution is 1 second. - -Works for scalars, sequences and numpy arrays. -Returns a scalar if input is a scalar, else returns a numpy array. - -The datetime instances returned by C{num2date} are 'real' python datetime -objects if the date falls in the Gregorian calendar (i.e. -C{calendar='proleptic_gregorian'}, or C{calendar = 'standard'/'gregorian'} and -the date is after 1582-10-15). Otherwise, they are 'phony' datetime -objects which are actually instances of netcdftime.datetime. This is -because the python datetime module cannot handle the weird dates in some -calendars (such as C{'360_day'} and C{'all_leap'}) which -do not exist in any real world calendar. - """ - isscalar = False - try: - time_value[0] - except: - isscalar = True - if not isscalar: - time_value = numpy.array(time_value, dtype='d') - shape = time_value.shape - # convert to desired units, remove time zone offset. - if self.units in ['second','seconds']: - jdelta = time_value/86400. - self.tzoffset/1440. - elif self.units in ['minute','minutes']: - jdelta = time_value/1440. - self.tzoffset/1440. - elif self.units in ['hour','hours']: - jdelta = time_value/24. - self.tzoffset/1440. - elif self.units in ['day','days']: - jdelta = time_value - self.tzoffset/1440. - jd = self._jd0 + jdelta - if self.calendar in ['julian','standard','gregorian','proleptic_gregorian']: - if not isscalar: - date = [DateFromJulianDay(j,self.calendar) for j in jd.flat] - else: - date = DateFromJulianDay(jd,self.calendar) - elif self.calendar in ['noleap','365_day']: - if not isscalar: - date = [_DateFromNoLeapDay(j) for j in jd.flat] - else: - date = _DateFromNoLeapDay(jd) - elif self.calendar in ['all_leap','366_day']: - if not isscalar: - date = [_DateFromAllLeap(j) for j in jd.flat] - else: - date = _DateFromAllLeap(jd) - elif self.calendar == '360_day': - if not isscalar: - date = [_DateFrom360Day(j) for j in jd.flat] - else: - date = _DateFrom360Day(jd) - if isscalar: - return date - else: - return numpy.reshape(numpy.array(date),shape) - -def _parse_date(origin): - """Parses a date string and returns a tuple - (year,month,day,hour,minute,second,utc_offset). - utc_offset is in minutes. - - This function parses the 'origin' part of the time unit. It should be - something like:: - - 2004-11-03 14:42:27.0 +2:00 - - Lots of things are optional; just the date is mandatory. - - by Roberto De Almeida - - excerpted from coards.py - http://cheeseshop.python.org/pypi/coards/ - """ - # yyyy-mm-dd [hh:mm:ss[.s][ [+-]hh[:][mm]]] - p = re.compile( r'''(?P<year>\d{1,4}) # yyyy - - # - (?P<month>\d{1,2}) # mm or m - - # - (?P<day>\d{1,2}) # dd or d - # - (?: # [optional time and timezone] - \s # - (?P<hour>\d{1,2}) # hh or h - : # - (?P<min>\d{1,2}) # mm or m - (?: - \: - (?P<sec>\d{1,2}) # ss or s (optional) - )? - # - (?: # [optional decisecond] - \. # . - (?P<dsec>\d) # s - )? # - (?: # [optional timezone] - \s # - (?P<ho>[+-]? # [+ or -] - \d{1,2}) # hh or h - :? # [:] - (?P<mo>\d{2})? # [mm] - )? # - )? # - $ # EOL - ''', re.VERBOSE) - - m = p.match(origin.strip()) - if m: - c = m.groupdict(0) - # UTC offset. - offset = int(c['ho'])*60 + int(c['mo']) - return int(c['year']),int(c['month']),int(c['day']),int(c['hour']),int(c['min']),int(c['sec']),offset - - raise Exception('Invalid date origin: %s' % origin) - -# remove the unsupposed "%s" command. But don't -# do it if there's an even number of %s before the s -# because those are all escaped. Can't simply -# remove the s because the result of -# %sY -# should be %Y if %s isn't supported, not the -# 4 digit year. -_illegal_s = re.compile(r"((^|[^%])(%%)*%s)") - -def _findall(text, substr): - # Also finds overlaps - sites = [] - i = 0 - while 1: - j = text.find(substr, i) - if j == -1: - break - sites.append(j) - i=j+1 - return sites - -# Every 28 years the calendar repeats, except through century leap -# years where it's 6 years. But only if you're using the Gregorian -# calendar. ;) - -def _strftime(dt, fmt): - if _illegal_s.search(fmt): - raise TypeError("This strftime implementation does not handle %s") - # don't use strftime method at all. - #if dt.year > 1900: - # return dt.strftime(fmt) - - year = dt.year - # For every non-leap year century, advance by - # 6 years to get into the 28-year repeat cycle - delta = 2000 - year - off = 6*(delta // 100 + delta // 400) - year = year + off - - # Move to around the year 2000 - year = year + ((2000 - year)//28)*28 - timetuple = dt.timetuple() - s1 = time.strftime(fmt, (year,) + timetuple[1:]) - sites1 = _findall(s1, str(year)) - - s2 = time.strftime(fmt, (year+28,) + timetuple[1:]) - sites2 = _findall(s2, str(year+28)) - - sites = [] - for site in sites1: - if site in sites2: - sites.append(site) - - s = s1 - syear = "%4d" % (dt.year,) - for site in sites: - s = s[:site] + syear + s[site+4:] - return s - -def date2num(dates,units,calendar='standard'): - """ -date2num(dates,units,calendar='standard') - -Return numeric time values given datetime objects. The units -of the numeric time values are described by the L{units} argument -and the L{calendar} keyword. The datetime objects must -be in UTC with no time-zone offset. If there is a -time-zone offset in C{units}, it will be applied to the -returned numeric values. - -Like the matplotlib C{date2num} function, except that it allows -for different units and calendars. Behaves the same if -C{units = 'days since 0001-01-01 00:00:00'} and -C{calendar = 'proleptic_gregorian'}. - -@param dates: A datetime object or a sequence of datetime objects. - The datetime objects should not include a time-zone offset. - -@param units: a string of the form C{'B{time units} since B{reference time}}' - describing the time units. B{C{time units}} can be days, hours, minutes - or seconds. B{C{reference time}} is the time origin. A valid choice - would be units=C{'hours since 1800-01-01 00:00:00 -6:00'}. - -@param calendar: describes the calendar used in the time calculations. - All the values currently defined in the U{CF metadata convention - <http://cf-pcmdi.llnl.gov/documents/cf-conventions/>} are supported. - Valid calendars C{'standard', 'gregorian', 'proleptic_gregorian' - 'noleap', '365_day', '360_day', 'julian', 'all_leap', '366_day'}. - Default is C{'standard'}, which is a mixed Julian/Gregorian calendar. - -@return: a numeric time value, or an array of numeric time values. - -The maximum resolution of the numeric time values is 1 second. - """ - cdftime = utime(units,calendar=calendar) - return cdftime.date2num(dates) - -def num2date(times,units,calendar='standard'): - """ -num2date(times,units,calendar='standard') - -Return datetime objects given numeric time values. The units -of the numeric time values are described by the C{units} argument -and the C{calendar} keyword. The returned datetime objects represent -UTC with no time-zone offset, even if the specified -C{units} contain a time-zone offset. - -Like the matplotlib C{num2date} function, except that it allows -for different units and calendars. Behaves the same if -C{units = 'days since 001-01-01 00:00:00'} and -C{calendar = 'proleptic_gregorian'}. - -@param times: numeric time values. Maximum resolution is 1 second. - -@param units: a string of the form C{'B{time units} since B{reference time}}' -describing the time units. B{C{time units}} can be days, hours, minutes -or seconds. B{C{reference time}} is the time origin. A valid choice -would be units=C{'hours since 1800-01-01 00:00:00 -6:00'}. - -@param calendar: describes the calendar used in the time calculations. -All the values currently defined in the U{CF metadata convention -<http://cf-pcmdi.llnl.gov/documents/cf-conventions/>} are supported. -Valid calendars C{'standard', 'gregorian', 'proleptic_gregorian' -'noleap', '365_day', '360_day', 'julian', 'all_leap', '366_day'}. -Default is C{'standard'}, which is a mixed Julian/Gregorian calendar. - -@return: a datetime instance, or an array of datetime instances. - -The datetime instances returned are 'real' python datetime -objects if the date falls in the Gregorian calendar (i.e. -C{calendar='proleptic_gregorian'}, or C{calendar = 'standard'} or C{'gregorian'} -and the date is after 1582-10-15). Otherwise, they are 'phony' datetime -objects which support some but not all the methods of 'real' python -datetime objects. This is because the python datetime module cannot -the uses the C{'proleptic_gregorian'} calendar, even before the switch -occured from the Julian calendar in 1582. The datetime instances -do not contain a time-zone offset, even if the specified C{units} -contains one. - """ - cdftime = utime(units,calendar=calendar) - return cdftime.num2date(times) - -def _check_index(indices, dates, nctime, calendar): - """Return True if the time indices given correspond to the given dates, - False otherwise. - - Parameters: - - indices : sequence of integers - Positive integers indexing the time variable. - - dates : sequence of datetime objects - Reference dates. - - nctime : netCDF Variable object - NetCDF time object. - - calendar : string - Calendar of nctime. - """ - if (indices <0).any(): - return False - - if (indices >= nctime.shape[0]).any(): - return False - -# t = nctime[indices] -# fancy indexing not available, fall back on this. - t=[] - for ind in indices: - t.append(nctime[ind]) - return numpy.all( num2date(t, nctime.units, calendar) == dates) - - -def date2index(dates, nctime, calendar=None, select='exact'): - """ - date2index(dates, nctime, calendar=None, select='exact') - - Return indices of a netCDF time variable corresponding to the given dates. - - @param dates: A datetime object or a sequence of datetime objects. - The datetime objects should not include a time-zone offset. - - @param nctime: A netCDF time variable object. The nctime object must have a - C{units} attribute. The entries are assumed to be stored in increasing - order. - - @param calendar: Describes the calendar used in the time calculation. - Valid calendars C{'standard', 'gregorian', 'proleptic_gregorian' - 'noleap', '365_day', '360_day', 'julian', 'all_leap', '366_day'}. - Default is C{'standard'}, which is a mixed Julian/Gregorian calendar - If C{calendar} is None, its value is given by C{nctime.calendar} or - C{standard} if no such attribute exists. - - @param select: C{'exact', 'before', 'after', 'nearest'} - The index selection method. C{exact} will return the indices perfectly - matching the dates given. C{before} and C{after} will return the indices - corresponding to the dates just before or just after the given dates if - an exact match cannot be found. C{nearest} will return the indices that - correpond to the closest dates. - """ - # Setting the calendar. - - if calendar == None: - calendar = getattr(nctime, 'calendar', 'standard') - - num = numpy.atleast_1d(date2num(dates, nctime.units, calendar)) - - # Trying to infer the correct index from the starting time and the stride. - # This assumes that the times are increasing uniformly. - t0, t1 = nctime[:2] - dt = t1 - t0 - index = numpy.array((num-t0)/dt, int) - - # Checking that the index really corresponds to the given date. - # If the times do not correspond, then it means that the times - # are not increasing uniformly and we try the bisection method. - if not _check_index(index, dates, nctime, calendar): - - # Use the bisection method. Assumes the dates are ordered. - import bisect - - index = numpy.array([bisect.bisect_left(nctime, n) for n in num], int) - - # Find the dates for which the match is not perfect. - # Use list comprehension instead of the simpler `nctime[index]` since - # not all time objects support numpy integer indexing (eg dap). - ncnum = numpy.squeeze([nctime[i] for i in index]) - mismatch = numpy.nonzero(ncnum != num)[0] - - if select == 'exact': - if len(mismatch) > 0: - raise ValueError, 'Some dates not found.' - - elif select == 'before': - index[mismatch] -= 1 - - elif select == 'after': - pass - - elif select == 'nearest': - nearest_to_left = num[mismatch] < numpy.array( [nctime[i-1] + nctime[i] for i in index[mismatch]]) / 2. - index[mismatch] = index[mismatch] - 1 * nearest_to_left - - else: - raise ValueError("%s is not an option for the `select` argument."%select) - - # convert numpy scalars or single element arrays to python ints. - return _toscalar(index) - - -def _toscalar(a): - if a.shape in [(),(1,)]: - return a.item() - else: - return a This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |