Context Navigation

Virtual.py @ 939

Revision 939, 181.3 KB (checked in by jjr8, 14 months ago)
Fixed/implemented: * #533: Add CoRTAD as a data product * #535: When accessing Aviso NRT datasets and Daily temporal resolution is selected, tools fail with "Please select daily temporal resolution or switch to a Delayed Time (DT) product that offers weekly resolution." * Tweaked documentation and tool names.

Line
1	# Datasets/Virtual.py - Classes representing virtual datasets (e.g.
2	# ClippedGrid).
3	#
4	# Copyright (C) 2010 Jason J. Roberts
5	#
6	# This program is free software; you can redistribute it and/or
7	# modify it under the terms of the GNU General Public License
8	# as published by the Free Software Foundation; either version 2
9	# of the License, or (at your option) any later version.
10	#
11	# This program is distributed in the hope that it will be useful,
12	# but WITHOUT ANY WARRANTY; without even the implied warranty of
13	# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14	# GNU General Public License (available in the file LICENSE.TXT)
15	# for more details.
16	#
17	# You should have received a copy of the GNU General Public License
18	# along with this program; if not, write to the Free Software
19	# Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
20
21	import bisect
22	import datetime
23	import math
24	import os
25
26	from GeoEco.Datasets import DatasetCollection, QueryableAttribute, Grid
27	from GeoEco.DynamicDocString import DynamicDocString
28	from GeoEco.Internationalization import _
29
30
31	class TimeSeriesGridStack(Grid):
32	__doc__ = DynamicDocString()
33
34	def _GetReportProgress(self):
35	return self._ReportProgress
36
37	def _SetReportProgress(self, value):
38	# TOOO: Validation
39	self._ReportProgress = value
40
41	ReportProgress = property(_GetReportProgress, _SetReportProgress, doc=DynamicDocString())
42
43	def __init__(self, collection, expression=None, reportProgress=True, **options):
44	# TODO: self.__class__.__doc__.Obj.ValidateMethodInvocation()
45
46	# TODO: What can we validate about the collection without
47	# accessing it?
48
49	# Validate that the collection has a queryable attribute with
50	# data type DateTimeTypeMetadata.
51
52	from GeoEco.Types import DateTimeTypeMetadata
53
54	dateTimeAttrs = collection.GetQueryableAttributesWithDataType(DateTimeTypeMetadata)
55	if len(dateTimeAttrs) <= 0:
56	raise ValueError(_(u'This dataset collection does not have a queryable attribute defined with the data type DateTimeTypeMetadata. In order to build a TimeSeriesGridStack from it, it must have an attribute with that data type.'))
57	if len(dateTimeAttrs) > 1: # Should never happen; CollectibleObject.__init__ prevents it
58	raise ValueError(_(u'This dataset collection has multiple queryable attributes defined with the data type DateTimeTypeMetadata. In order to build a TimeSeriesGridStack from it, only one queryable attribute of that type must be defined.'))
59
60	# Query the collection for the oldest grid within it.
61
62	self._CachedOldestGrid = collection.GetOldestDataset(expression, **options)
63
64	if self._CachedOldestGrid is None:
65	raise CollectionIsEmptyError(collection.DisplayName, expression)
66
67	if not issubclass(self._CachedOldestGrid.__class__, Grid):
68	raise TypeError(_(u'The dataset collection %(dn)s does not contain Grid datasets. %(cls)s can only be used with dataset collections that contain Grid datasets.') % {u'dn': collection.DisplayName, u'cls': self.__class__.__name__})
69
70	# Copy all of the queryable attributes, except the
71	# DateTimeTypeMetadata one, and their values from the oldest
72	# grid. All of the grids returned by the expression should
73	# have the same values for any given queryable attribute
74	# (except the DateTimeTypeMetadata attribute and attributes
75	# derived from it). We do not verify this, however.
76
77	queryableAttributes = []
78	obj = self._CachedOldestGrid
79	while obj is not None:
80	if obj._QueryableAttributes is not None:
81	for attr in obj._QueryableAttributes:
82	if attr.Name != dateTimeAttrs[0].Name:
83	queryableAttributes.append(QueryableAttribute(attr.Name, attr.DisplayName, attr.DataType, None, attr.DerivedFromAttr, attr.DerivedValueMap, attr.DerivedValueFunc)) # Do not copy attr.DerivedLazyDatasetProps
84	obj = obj.ParentCollection
85
86	queryableAttributeValues = {}
87	for attr in queryableAttributes:
88	if attr.DerivedFromAttr != dateTimeAttrs[0].Name:
89	queryableAttributeValues[attr.Name] = self._CachedOldestGrid.GetQueryableAttributeValue(attr.Name)
90
91	# Initialize our properties.
92
93	self._Collection = collection
94	self._Expression = expression
95	self._ReportProgress = reportProgress
96	self._Options = options
97	self._DateTimeAttrName = dateTimeAttrs[0].Name
98	self._CachedQueryExpression = None
99	self._CachedDatasets = None
100
101	# Initialize the base class.
102
103	super(TimeSeriesGridStack, self).__init__(queryableAttributes=tuple(queryableAttributes), queryableAttributeValues=queryableAttributeValues)
104
105	def _Close(self):
106	if hasattr(self, '_CachedDatasets') and self._CachedDatasets is not None:
107	while len(self._CachedDatasets) > 0:
108	self._CachedDatasets[0].Close()
109	del self._CachedDatasets[0]
110	self._CachedQueryExpression = None
111	self._CachedDatasets = None
112
113	if hasattr(self, '_Collection'):
114	self._Collection.Close()
115
116	super(TimeSeriesGridStack, self)._Close()
117
118	def _GetDisplayName(self):
119	return self._Collection.DisplayName
120
121	def _GetLazyPropertyPhysicalValue(self, name):
122
123	# If the oldest grid does not have a t dimension already, we
124	# are stacking a collection of 2D (yx) or 3D (zyx) grids into
125	# a 3D (tyx) or 4D (tzyx) stack.
126
127	if 't' not in self._CachedOldestGrid.Dimensions:
128
129	# Handle properties that are not identical to but are
130	# easily calculated from the values of the oldest grid.
131
132	if name == 'Dimensions':
133	return 't' + self._CachedOldestGrid.Dimensions
134
135	if name == 'PhysicalDimensions':
136	return 't' + self._CachedOldestGrid.Dimensions # Transposing of the underlying time slices is done when we fetch each slice, thus they are all properly ordered by the time we receive them.
137
138	if name == 'PhysicalDimensionsFlipped':
139	return tuple([False] * (len(self._CachedOldestGrid.Dimensions) + 1)) # Flipping of the underlying time slices is done when we fetch each slice, thus they are all properly oriented by the time we receive them.
140
141	if name == 'CoordDependencies':
142	return tuple([None] + list(self._CachedOldestGrid.CoordDependencies))
143
144	if name == 'CoordIncrements':
145	return tuple([self._CachedOldestGrid.GetLazyPropertyValue('TIncrement')] + list(self._CachedOldestGrid.CoordIncrements))
146
147	if name == 'CornerCoords':
148	return tuple([self._CachedOldestGrid.GetQueryableAttributeValue(self._DateTimeAttrName)] + list(self._CachedOldestGrid.GetLazyPropertyValue('CornerCoords')))
149
150	# Handle the shape. This is more complicated because we
151	# have to determine how many time slices there are.
152
153	if name == 'Shape':
154
155	# If the t increment is not None, we do not need to
156	# retrieve the full list of grids to know how many
157	# time slices there are. Instead, we can calculate how
158	# many must appear between the oldest grid and the
159	# newest grid.
160
161	if self.CoordIncrements[0] is not None:
162	newestGrid = self._Collection.GetNewestDataset(self._Expression, **self._Options)
163	if not issubclass(newestGrid.__class__, Grid):
164	raise TypeError(_(u'The dataset collection %(dn)s does not contain Grid datasets. %(cls)s can only be used with dataset collections that contain Grid datasets.') % {u'dn': self._Collection.DisplayName, u'cls': self.__class__.__name__})
165	newestGridDateTime = newestGrid.GetQueryableAttributeValue(self._DateTimeAttrName)
166
167	# Estimate the number of time slices between the
168	# oldest grid and newest grid.
169
170	delta = newestGridDateTime - self._CachedOldestGrid.GetQueryableAttributeValue(self._DateTimeAttrName)
171	delta = delta.days * 86400. + delta.seconds
172
173	if self.TIncrementUnit == 'year':
174	numTimeSlices = int(1.1 * delta / 86400. / 365.)
175	elif self.TIncrementUnit == 'season':
176	numTimeSlices = int(1.1 * delta / 86400. / 365. * 4.)
177	elif self.TIncrementUnit == 'month':
178	numTimeSlices = int(1.1 * delta / 86400. / 365. * 12.)
179	elif self.TIncrementUnit == 'day':
180	numTimeSlices = int(1.1 * delta / 86400.)
181	elif self.TIncrementUnit == 'hour':
182	numTimeSlices = int(1.1 * delta / 3600.)
183	elif self.TIncrementUnit == 'minute':
184	numTimeSlices = int(1.1 * delta / 60.)
185	elif self.TIncrementUnit == 'second':
186	numTimeSlices = int(1.1 * delta)
187	else:
188	raise NotImplementedError(_(u'Programming error in this tool: the t increment unit \'%(unit)s\' is unknown. Please contact the author of this tool for assistance.') % {u'unit': self.TIncrementUnit})
189
190	# Get a list of t coordinates starting with the
191	# first time slice.
192
193	tCornerCoordType = self.GetLazyPropertyValue('TCornerCoordType').lower()
194	if tCornerCoordType == 'min':
195	fixedIncrementOffset = -0.5
196	elif tCornerCoordType == 'center':
197	fixedIncrementOffset = 0.0
198	elif tCornerCoordType == 'max':
199	fixedIncrementOffset = 0.5
200	else:
201	raise NotImplementedError(_(u'Programming error in this tool: the t corner coordinate type \'%(type)s\' is unknown. Please contact the author of this tool for assistance.') % {u'type': self.GetLazyPropertyValue('TCornerCoordType')})
202
203	tCoords = self._GetTCoordsList(fixedIncrementOffset, numTimeSlices)
204
205	# While the time of the newest grid is newer than
206	# the newest t coordinate, double the size of the
207	# list. This should probably never happen.
208
209	while tCoords[-1] < newestGridDateTime:
210	numTimeSlices *= 2
211	tCoords = self._GetTCoordsList(fixedIncrementOffset, numTimeSlices)
212
213	# Search the t coordinates backwards for the time
214	# of the newest grid. This tells us the number of
215	# time slices. If we do not find it, something odd
216	# is going on; the parsed datetime is inconsistent
217	# with the definition of the dataset.
218
219	i = len(tCoords) - 1
220	while i >= 0:
221	if tCoords[i] == newestGridDateTime:
222	break
223	if tCoords[i] < newestGridDateTime:
224	raise ValueError(_(u'Failed to compute a time coordinate that matched the time of the newest grid in this %(cls)s. The datetime of that grid is %(last)s but the two closest time coordinates are %(dt1)s and %(dt2)s.') % {u'cls': self.__class__.__name__, u'last': newestGridDateTime.strftime('%Y-%m-%d %H:%M:%S'), u'dt1': tCoords[i].strftime('%Y-%m-%d %H:%M:%S'), u'dt2': tCoords[i+1].strftime('%Y-%m-%d %H:%M:%S')})
225	i -= 1
226	if i < 0:
227	raise ValueError(_(u'Programming error in this tool: The datetime of the newest grid in this %(cls)s is %(last)s, which comes before the datetime of the first time coordinate, %(dt1)s. Please contact the author of this tool for assistance') % {u'cls': self.__class__.__name__, u'last': newestGridDateTime.strftime('%Y-%m-%d %H:%M:%S'), u'dt1': tCoords[0].strftime('%Y-%m-%d %H:%M:%S')})
228
229	# Set and return the shape.
230
231	shape = tuple([i+1] + list(self._CachedOldestGrid.Shape))
232	self.SetLazyPropertyValue('Shape', shape)
233	return shape
234
235	# The t increment is None. We need to retrieve the
236	# full list of grids to know the number of time
237	# slices. We do not currently support this. TODO:
238	# implement it.
239
240	raise NotImplementedError(_(u'The dataset collection %(dn)s contains grids that do not have a fixed time increment. The current implementation of %(cls)s does not support grids of this kind.') % {u'dn': self._Collection.DisplayName, u'cls': self.__class__.__name__})
241
242	# If the contained grids have a t dimension already, we are
243	# concatenating a collection of 3D (tyx) or 4D (tzyx) grids.
244	# The stack will have the same dimensions as an individual
245	# grid in the collection. We do not currently support this.
246	# TODO: implement it.
247
248	else:
249	raise NotImplementedError(_(u'The dataset collection %(dn)s contains grids that have a time dimension. The current implementation of %(cls)s does not support grids with a time dimension.') % {u'dn': self._Collection.DisplayName, u'cls': self.__class__.__name__})
250
251	# If we got to here, the caller has requested a lazy property
252	# that is assumed to be the same for all grids in the
253	# collection as well as the stack itself. Return the value
254	# from the oldest grid.
255
256	return self._CachedOldestGrid.GetLazyPropertyValue(name)
257
258	def _GetCoords(self, coord, coordNum, slices, sliceDims, fixedIncrementOffset):
259	return self._Grid._GetCoords(coord, coordNum, slices, sliceDims, fixedIncrementOffset)
260
261	def _ReadNumpyArray(self, sliceList):
262
263	# Get a list of t coordinates for the requested time slices.
264
265	tCornerCoordType = self.GetLazyPropertyValue('TCornerCoordType').lower()
266	if tCornerCoordType == 'min':
267	tCoords = self.MinCoords.__getitem__(('t', sliceList[0]))
268	elif tCornerCoordType == 'center':
269	tCoords = self.CenterCoords.__getitem__(('t', sliceList[0]))
270	elif tCornerCoordType == 'max':
271	tCoords = self.MaxCoords.__getitem__(('t', sliceList[0]))
272	else:
273	raise NotImplementedError(_(u'Programming error in this tool: the t corner coordinate type \'%(type)s\' is unknown. Please contact the author of this tool for assistance.') % {u'type': self.GetLazyPropertyValue('TCornerCoordType')})
274
275	# TODO: handle len(tCoords) == 0
276
277	# Query the collection for a list of datasets with t
278	# coordinates that fall within the min and max coordinates of
279	# the requested time slices.
280
281	expression = self._DateTimeAttrName + ' >= ' + tCoords[0].strftime('#%Y-%m-%d %H:%M:%S# AND Year >= %Y') + ' AND ' + self._DateTimeAttrName + ' <= ' + tCoords[-1].strftime('#%Y-%m-%d %H:%M:%S# AND Year <= %Y')
282	if self._Expression is not None:
283	expression = '(' + expression + ') AND (' + self._Expression + ')'
284
285	if self._CachedQueryExpression is None or self._CachedQueryExpression != expression:
286	if self._CachedDatasets is not None:
287	while len(self._CachedDatasets) > 0:
288	self._CachedDatasets[0].Close()
289	del self._CachedDatasets[0]
290	self._CachedDatasets = None
291	self._CachedQueryExpression = None
292
293	self._CachedDatasets = self._Collection.QueryDatasets(expression, self._ReportProgress and len(tCoords) > 1, **self._Options)
294	self._CachedQueryExpression = expression
295
296	# Most likely, the datasets are already sorted in ascending
297	# time order, but this is not required. Sort them, to be sure.
298
299	self._CachedDatasets.sort(key=lambda ds: ds.GetQueryableAttributeValue(self._DateTimeAttrName))
300
301	# Allocate a numpy array to return.
302
303	import numpy
304	data = numpy.zeros(map(lambda s: s.stop - s.start, sliceList), str(self._CachedOldestGrid.UnscaledDataType))
305
306	# Fill each time slice by retrieving the data from the
307	# corresponding dataset. If we encounter a time slice for
308	# which there is no dataset with a matching time coordinate,
309	# allocate a slice of NoData and report a warning. If there is
310	# no NoData value, leave the values at zero.
311
312	t = 0
313	i = 0
314	while t < len(tCoords):
315	while i < len(self._CachedDatasets) and self._CachedDatasets[i].GetQueryableAttributeValue(self._DateTimeAttrName) < tCoords[t]:
316	i += 1
317
318	if i < len(self._CachedDatasets) and self._CachedDatasets[i].GetQueryableAttributeValue(self._DateTimeAttrName) == tCoords[t]:
319	data[t] = self._CachedDatasets[i].UnscaledData.__getitem__(tuple(sliceList[1:]))
320	elif self._CachedOldestGrid.UnscaledNoDataValue is not None:
321	if tCornerCoordType == 'min':
322	self._LogWarning(_(u'There is no data in %(dn)s for the time slice that starts on %(ts)s.') % {u'dn': self._Collection.DisplayName, u'ts': tCoords[t].strftime('%Y-%m-%d %H:%M:%S')})
323	elif tCornerCoordType == 'center':
324	self._LogWarning(_(u'There is no data in %(dn)s for the time slice %(ts)s.') % {u'dn': self._Collection.DisplayName, u'ts': tCoords[t].strftime('%Y-%m-%d %H:%M:%S')})
325	else:
326	self._LogWarning(_(u'There is no data in %(dn)s for the time slice that ends on %(ts)s.') % {u'dn': self._Collection.DisplayName, u'ts': tCoords[t].strftime('%Y-%m-%d %H:%M:%S')})
327	data[t] += self._CachedOldestGrid.UnscaledNoDataValue
328	else:
329	if tCornerCoordType == 'min':
330	self._LogWarning(_(u'There is no data in %(dn)s for the time slice that starts on %(ts)s but the datasets in this collection do not have a NoData value defined. The values of this time slice will be set to 0.') % {u'dn': self._Collection.DisplayName, u'ts': tCoords[t].strftime('%Y-%m-%d %H:%M:%S')})
331	elif tCornerCoordType == 'center':
332	self._LogWarning(_(u'There is no data in %(dn)s for the time slice %(ts)s but the datasets in this collection do not have a NoData value defined. The values of this time slice will be set to 0.') % {u'dn': self._Collection.DisplayName, u'ts': tCoords[t].strftime('%Y-%m-%d %H:%M:%S')})
333	else:
334	self._LogWarning(_(u'There is no data in %(dn)s for the time slice that ends on %(ts)s but the datasets in this collection do not have a NoData value defined. The values of this time slice will be set to 0.') % {u'dn': self._Collection.DisplayName, u'ts': tCoords[t].strftime('%Y-%m-%d %H:%M:%S')})
335
336	t += 1
337
338	# Return the populated numpy array.
339
340	return data, self._CachedOldestGrid.UnscaledNoDataValue
341
342
343	class CollectionIsEmptyError(Exception):
344
345	def __init__(self, collectionDisplayName=None, expression=None):
346	self.CollectionDisplayName = collectionDisplayName
347	self.Expression = expression
348
349	def __str__(self):
350	if self.CollectionDisplayName is None and self.Expression is None:
351	return _(u'This collection contains no datasets.')
352	if self.CollectionDisplayName is not None and self.Expression is None:
353	return _(u'The %(dn)s contains no datasets.') % {u'dn': self.CollectionDisplayName}
354	return _(u'The %(dn)s contains no datasets matching the expression %(expr)s.') % {u'dn': self.CollectionDisplayName, u'expr': self.Expression}
355
356
357	class GridSlice(Grid):
358	__doc__ = DynamicDocString()
359
360	def __init__(self, grid, tIndex=None, zIndex=None, tQAName=u'DateTime', tQADisplayName=_(u'Date'), tQACoordType=u'min', zQAName=u'Depth', zQADisplayName=_(u'Depth'), zQACoordType=u'center', displayName=None, additionalQueryableAttributeValues=None):
361	# TODO: Validation
362
363	# Perform additional validation.
364
365	if tIndex is None and zIndex is None:
366	raise ValueError(_(u'Both tIndex and zIndex are None. A value must be provided for at least one of them.'))
367
368	if tIndex is not None:
369	if 't' not in grid.Dimensions:
370	raise TypeError(_(u'A value was provided for tIndex but %(dn)s does not have a t dimension.') % {u'dn': grid.DisplayName})
371	if tQAName is None:
372	raise TypeError(_(u'If a value is provided for tIndex, a value must also be provided for tQAName.'))
373	if tQADisplayName is None:
374	raise TypeError(_(u'If a value is provided for tIndex, a value must also be provided for tQADisplayName.'))
375	if tQACoordType is None:
376	raise TypeError(_(u'If a value is provided for tIndex, a value must also be provided for tQACoordType.'))
377
378	if zIndex is not None:
379	if 'z' not in grid.Dimensions:
380	raise TypeError(_(u'A value was provided for zIndex but %(dn)s does not have a z dimension.') % {u'dn': grid.DisplayName})
381	if zQAName is None:
382	raise TypeError(_(u'If a value is provided for zIndex, a value must also be provided for zQAName.'))
383	if zQADisplayName is None:
384	raise TypeError(_(u'If a value is provided for zIndex, a value must also be provided for zQADisplayName.'))
385	if zQACoordType is None:
386	raise TypeError(_(u'If a value is provided for zIndex, a value must also be provided for zQACoordType.'))
387
388	# Validate the provided indices and make them positive.
389
390	if tIndex is not None:
391	if tIndex < 0:
392	if tIndex + grid.Shape[0] < 0:
393	raise IndexError(_(u'tIndex is out of range.'))
394	tIndex += grid.Shape[0]
395	elif tIndex >= grid.Shape[0]:
396	raise IndexError(_(u'tIndex is out of range.'))
397
398	if zIndex is not None:
399	if 'z' not in grid.Dimensions:
400	if zIndex + grid.Shape[grid.Dimensions.index('z')] < 0:
401	raise IndexError(_(u'zIndex is out of range.'))
402	zIndex += grid.Shape[grid.Dimensions.index('z')]
403	elif zIndex >= grid.Shape[grid.Dimensions.index('z')]:
404	raise IndexError(_(u'zIndex is out of range.'))
405
406	# Initialize our properties.
407
408	self._Grid = grid
409	self._TIndex = tIndex
410	self._ZIndex = zIndex
411
412	if self._TIndex is not None:
413	self._TQAName = tQAName
414	self._TQADisplayName = tQADisplayName
415	self._TQACoordType = tQACoordType
416	else:
417	self._TQAName = None
418	self._TQADisplayName = None
419	self._TQACoordType = None
420
421	if self._ZIndex is not None:
422	self._ZQAName = zQAName
423	self._ZQADisplayName = zQADisplayName
424	self._ZQACoordType = zQACoordType
425	else:
426	self._ZQAName = None
427	self._ZQADisplayName = None
428	self._ZQACoordType = None
429
430	if displayName is not None:
431	self._DisplayName = displayName
432	elif self._TIndex is not None and self._ZIndex is not None:
433	self._DisplayName = _(u'%(tdn)s, %(zdn)s slice [%(tIndex)i, %(zIndex)i] of %(dn)s') % {u'tdn': self._TQADisplayName.lower(), u'zdn': self._ZQADisplayName.lower(), u'tIndex': self._TIndex, u'zIndex': self._ZIndex, u'dn': self._Grid.DisplayName}
434	elif self._TIndex is not None:
435	self._DisplayName = _(u'%(tdn)s slice %(tIndex)i of %(dn)s') % {u'tdn': self._TQADisplayName.lower(), u'tIndex': self._TIndex, u'dn': self._Grid.DisplayName}
436	else:
437	self._DisplayName = _(u'%(zdn)s slice %(zIndex)i of %(dn)s') % {u'zdn': self._ZQADisplayName.lower(), u'zIndex': self._ZIndex, u'dn': self._Grid.DisplayName}
438
439	# For our queryable attributes, use all of those of the grid
440	# plus the ones for the t and/or z dimensions.
441
442	queryableAttributes = []
443	queryableAttributes.extend(grid.GetAllQueryableAttributes())
444
445	queryableAttributeValues = {}
446	for qa in queryableAttributes:
447	if qa.DerivedFromAttr is None or qa.DerivedFromAttr not in [self._TQAName, self._ZQAName]:
448	queryableAttributeValues[qa.Name] = grid.GetQueryableAttributeValue(qa.Name)
449
450	if self._TIndex is not None:
451	queryableAttributes.append(QueryableAttribute(self._TQAName, self._TQADisplayName, DateTimeTypeMetadata()))
452	if self._TQACoordType == 'min':
453	queryableAttributeValues[self._TQAName] = self._Grid.MinCoords['t', self._TIndex]
454	elif self._TQACoordType == 'center':
455	queryableAttributeValues[self._TQAName] = self._Grid.CenterCoords['t', self._TIndex]
456	else:
457	queryableAttributeValues[self._TQAName] = self._Grid.MaxCoords['t', self._TIndex]
458
459	if self._ZIndex is not None:
460	queryableAttributes.append(QueryableAttribute(self._ZQAName, self._ZQADisplayName, FloatTypeMetadata()))
461	if self._ZQACoordType == 'min':
462	queryableAttributeValues[self._ZQAName] = self._Grid.MinCoords['z', self._ZIndex]
463	elif self._ZQACoordType == 'center':
464	queryableAttributeValues[self._ZQAName] = self._Grid.CenterCoords['z', self._ZIndex]
465	else:
466	queryableAttributeValues[self._ZQAName] = self._Grid.MaxCoords['z', self._ZIndex]
467
468	if additionalQueryableAttributeValues is not None:
469	queryableAttributeValues.update(additionalQueryableAttributeValues)
470
471	# Our goal is to imitate the contained grid except with fewer
472	# dimensions. In order to do this, we have to override the
473	# dimensions and other lazy properties related to it. Obtain
474	# the indices of the remaining dimensions.
475
476	if 'z' in self._Grid.Dimensions and self._ZIndex is None:
477	self._RemainingDimensionIndices = [1,2,3] # It is a t slice of a tzyx grid, so remaining dimensions are zyx
478	elif 't' in self._Grid.Dimensions and self._TIndex is None:
479	self._RemainingDimensionIndices = [0,2,3] # It is a z slice of a tzyx grid, so remaining dimensions are tyx
480	else:
481	self._RemainingDimensionIndices = range(len(self._Grid.Dimensions))[-2:] # It is either a tz slice of a tzyx grid, a t slice of a tyx grid, or a z slice of a zyx grid, so remaining dimensions are yx
482
483	# Initialize the base class.
484
485	super(GridSlice, self).__init__(queryableAttributes=tuple(queryableAttributes), queryableAttributeValues=queryableAttributeValues)
486
487	def _GetDisplayName(self):
488	return self._DisplayName
489
490	def _GetLazyPropertyPhysicalValue(self, name):
491
492	# If the requested property is sequence related to the
493	# dimensions, get the values from the grid we're slicing but
494	# remove the element corresponding to the sliced dimension(s).
495
496	if name == 'Dimensions':
497	return u''.join([self._Grid.Dimensions[i] for i in self._RemainingDimensionIndices])
498
499	if name == 'Shape':
500	return tuple([self._Grid.Shape[i] for i in self._RemainingDimensionIndices])
501
502	if name == 'CoordDependencies':
503	return tuple([self._Grid.CoordDependencies[i] for i in self._RemainingDimensionIndices]) # TODO: to the contained lists need to be modified?
504
505	if name == 'CoordIncrements':
506	return tuple([self._Grid.CoordIncrements[i] for i in self._RemainingDimensionIndices])
507
508	if name == 'CornerCoords':
509	return tuple([self._Grid.GetLazyPropertyValue('CornerCoords')[i] for i in self._RemainingDimensionIndices])
510
511	# If the requested property is PhysicalDimensions or
512	# PhysicalDimensionsFlipped, return values indicating the
513	# dimensions in the ideal order. The contained grid takes care
514	# of reordering, if needed.
515
516	if name == 'PhysicalDimensions':
517	return self.Dimensions
518
519	if name == 'PhysicalDimensionsFlipped':
520	return tuple([False] * len(self.Dimensions))
521
522	# Otherwise just get the unaltered value from the contained
523	# grid.
524
525	return self._Grid.GetLazyPropertyValue(name)
526
527	@classmethod
528	def _TestCapability(cls, capability):
529	return cls._Grid._TestCapability(capability)
530
531	@classmethod
532	def _GetSRTypeForSetting(cls):
533	return cls._Grid._GetSRTypeForSetting()
534
535	def _SetSpatialReference(self, sr):
536	return self._Grid._SetSpatialReference(sr)
537
538	def _GetUnscaledDataAsArray(self, key):
539	return self._Grid._GetUnscaledDataAsArray(self._AddSlicedDimsToKey(key))
540
541	def _SetUnscaledDataWithArray(self, key, value):
542	return self._Grid._SetUnscaledDataWithArray(self._AddSlicedDimsToKey(key), value)
543
544	def _AddSlicedDimsToKey(self, key):
545
546	# Validate the key. Although we are calling the
547	# _ValidateAndFlipKey function, because our
548	# PhysicalDimensionsFlipped contains only False, none of the
549	# key's indices will be flipped.
550
551	key2 = self._ValidateAndFlipKey(key)
552
553	# The key does not include the dimensions the sliced
554	# dimensions. Add these to the key as single indices.
555
556	if self._ZIndex is not None:
557	key2.insert(0, self._ZIndex)
558
559	if self._TIndex is not None:
560	key2.insert(0, self._TIndex)
561
562	# Return a tuple.
563
564	return tuple(key2)
565
566	def _GetCoords(self, coord, coordNum, slices, sliceDims, fixedIncrementOffset):
567	return self._Grid._GetCoords(coord, coordNum, slices, sliceDims, fixedIncrementOffset)
568
569
570	class GridSliceCollection(DatasetCollection):
571	__doc__ = DynamicDocString()
572
573	def __init__(self, grid, tQAName=u'DateTime', tQADisplayName=_(u'Time'), tQACoordType=u'min', zQAName=u'Depth', zQADisplayName=_(u'Depth'), zQACoordType=u'center', displayName=None, additionalQueryableAttributes=None):
574	# TODO: Validation
575
576	# Perform additional validation.
577
578	if tQAName is not None and tQADisplayName is None:
579	raise TypeError(_(u'If a value is provided for tQAName, a value must also be provided for tQADisplayName.'))
580
581	if tQAName is not None and tQACoordType is None:
582	raise TypeError(_(u'If a value is provided for tQAName, a value must also be provided for tQACoordType.'))
583
584	if zQAName is not None and zQADisplayName is None:
585	raise TypeError(_(u'If a value is provided for zQAName, a value must also be provided for zQADisplayName.'))
586
587	if zQAName is not None and zQACoordType is None:
588	raise TypeError(_(u'If a value is provided for zQAName, a value must also be provided for zQACoordType.'))
589
590	if not ('t' in grid.Dimensions or 'z' in grid.Dimensions):
591	raise ValueError(_(u'Cannot construct a GridSliceCollection from %(dn)s because it does not have a t dimension or a z dimension.') % {u'dn': grid.DisplayName})
592
593	if not ('t' in grid.Dimensions and tQAName is not None or 'z' in grid.Dimensions and zQAName is not None):
594	raise ValueError(_(u'Cannot construct a GridSliceCollection from %(dn)s. Although it has a t dimension and/or z dimension, it was not sliced by that (or those) dimensions. It must be sliced by at least one of them.') % {u'dn': grid.DisplayName})
595
596	# Initialize our properties.
597
598	self._Grid = grid
599
600	if 't' in grid.Dimensions and tQAName is not None:
601	self._TQAName = tQAName
602	self._TQADisplayName = tQADisplayName
603	self._TQACoordType = tQACoordType
604	else:
605	self._TQAName = None
606	self._TQADisplayName = None
607	self._TQACoordType = None
608
609	if 'z' in grid.Dimensions and zQAName is not None:
610	self._ZQAName = zQAName
611	self._ZQADisplayName = zQADisplayName
612	self._ZQACoordType = zQACoordType
613	else:
614	self._ZQAName = None
615	self._ZQADisplayName = None
616	self._ZQACoordType = None
617
618	if displayName is not None:
619	self._DisplayName = displayName
620	elif self._TQAName is not None and self._ZQAName is not None:
621	self._DisplayName = _(u'%(tdn)s and %(zdn)s slices of %(dn)s') % {u'tdn': self._TQADisplayName.lower(), u'zdn': self._ZQADisplayName.lower(), u'dn': self._Grid.DisplayName}
622	elif self._TQAName is not None:
623	self._DisplayName = _(u'%(tdn)s slices of %(dn)s') % {u'tdn': self._TQADisplayName.lower(), u'dn': self._Grid.DisplayName}
624	else:
625	self._DisplayName = _(u'%(zdn)s slices of %(dn)s') % {u'zdn': self._ZQADisplayName.lower(), u'dn': self._Grid.DisplayName}
626
627	# For our queryable attributes, use all of those of the grid
628	# plus the ones for the t and/or z dimensions.
629
630	queryableAttributes = []
631	if grid._QueryableAttributes is not None:
632	queryableAttributes.extend(grid._QueryableAttributes)
633
634	if self._TQAName is not None:
635	queryableAttributes.append(QueryableAttribute(self._TQAName, self._TQADisplayName, DateTimeTypeMetadata()))
636
637	if self._ZQAName is not None:
638	queryableAttributes.append(QueryableAttribute(self._ZQAName, self._ZQADisplayName, FloatTypeMetadata()))
639
640	self._AdditionalQueryableAttributes = None
641	if additionalQueryableAttributes is not None:
642	self._AdditionalQueryableAttributes = tuple(additionalQueryableAttributes)
643	queryableAttributes.extend(additionalQueryableAttributes)
644
645	# Initialize the base class.
646
647	super(GridSliceCollection, self).__init__(parentCollection=grid.ParentCollection, queryableAttributes=tuple(queryableAttributes))
648
649	def _Close(self):
650	if hasattr(self, '_Grid') and self._Grid is not None:
651	self._Grid.Close()
652	super(GridSliceCollection, self)._Close()
653
654	def _GetDisplayName(self):
655	return self._DisplayName
656
657	def _QueryDatasets(self, parsedExpression, progressReporter, options, parentAttrValues):
658	return self._QueryGridSlices(parsedExpression, progressReporter)
659
660	def _GetOldestDataset(self, parsedExpression, options, parentAttrValues, dateTimeAttrName):
661	datasets = self._QueryGridSlices(parsedExpression, numResults=1)
662	if len(datasets) > 0:
663	return datasets[0]
664	return None
665
666	def _GetNewestDataset(self, parsedExpression, options, parentAttrValues, dateTimeAttrName):
667	datasets = self._QueryGridSlices(parsedExpression, numResults=1, reverseOrder=True)
668	if len(datasets) > 0:
669	return datasets[0]
670	return None
671
672	def _QueryGridSlices(self, parsedExpression, progressReporter=None, numResults=None, reverseOrder=False):
673
674	attrValues = {}
675	for attr in self._QueryableAttributes:
676	attrValues[attr.Name] = self._Grid.GetQueryableAttributeValue(attr.Name)
677
678	# Iterate through the slices in the appropriate order, z
679	# changing before t.
680
681	results = []
682
683	if self._TQAName is not None:
684	maxT = self._Grid.Shape[0]
685	if reverseOrder:
686	t = maxT - 1
687	else:
688	t = 0
689	else:
690	t = None
691
692	if self._ZQAName is not None:
693	maxZ = self._Grid.Shape[self._Grid.Dimensions.index('z')]
694	if reverseOrder:
695	z = maxZ - 1
696	else:
697	z = 0
698	else:
699	z = None
700
701	while (numResults is None or len(results) < numResults) and \
702	(self._TQAName is None or (t >= 0 and t < maxT)) and \
703	(self._ZQAName is None or (z >= 0 and z < maxZ)):
704
705	# Set the t and/or z queryable attribute values for this
706	# slice.
707
708	if self._TQAName is not None:
709	if self._TQACoordType == 'min':
710	tCoord = self._Grid.MinCoords['t', t]
711	elif self._TQACoordType == 'center':
712	tCoord = self._Grid.CenterCoords['t', t]
713	else:
714	tCoord = self._Grid.MaxCoords['t', t]
715	attrValues[self._TQAName] = tCoord
716	attrValues['Year'] = tCoord.year
717	attrValues['Month'] = tCoord.month
718	attrValues['Day'] = tCoord.day
719	attrValues['Hour'] = tCoord.hour
720	attrValues['Minute'] = tCoord.minute
721	attrValues['Second'] = tCoord.second
722	attrValues['DayOfYear'] = (datetime.datetime(tCoord.year, tCoord.month, tCoord.day) - datetime.datetime(tCoord.year, 1, 1)).days + 1
723	else:
724	tCoord = None
725
726	if self._ZQAName is not None:
727	if self._ZQACoordType == 'min':
728	zCoord = self._Grid.MinCoords['z', z]
729	elif self._ZQACoordType == 'center':
730	zCoord = self._Grid.CenterCoords['z', z]
731	else:
732	zCoord = self._Grid.MaxCoords['z', z]
733	attrValues[self._ZQAName] = zCoord
734	else:
735	zCoord = None
736
737	# Evaluate the expression for this slice. The only thing
738	# that will be different about this slice compared to
739	# others is the t and/or z queryable attribute values.
740
741	results.extend(self._EvaluateExpressionForSlice(parsedExpression, attrValues, t, z, tCoord, zCoord, progressReporter))
742
743	# Go on to the next slice.
744
745	if self._ZQAName is not None:
746	if reverseOrder:
747	z -= 1
748	else:
749	z += 1
750
751	if self._TQAName is not None and (self._ZQAName is None or z < 0 or z >= maxZ):
752	if reverseOrder:
753	if self._ZQAName is not None:
754	z = maxZ - 1
755	t -= 1
756	else:
757	if self._ZQAName is not None:
758	z = 0
759	t += 1
760
761	return results
762
763	def _EvaluateExpressionForSlice(self, parsedExpression, attrValues, t, z, tCoord, zCoord, progressReporter):
764
765	if parsedExpression is not None:
766	try:
767	result = parsedExpression.eval(attrValues)
768	except Exception, e:
769	return [] # TODO: report better message
770	else:
771	result = True
772
773	if result is None or result:
774	if self._TQAName is not None and self._ZQAName is not None:
775	self._LogDebug(_(u'%(class)s 0x%(id)08X: Query result for t=%(t)s (%(tCoord)s), z=%(z)i (%(zCoord)s): %(result)s'), {u'class': self.__class__.__name__, u'id': id(self), u't': t, u'z': z, u'tCoord': str(tCoord), u'zCoord': repr(zCoord), u'result': repr(result)})
776	elif self._TQAName is not None:
777	self._LogDebug(_(u'%(class)s 0x%(id)08X: Query result for t=%(t)s (%(tCoord)s): %(result)s'), {u'class': self.__class__.__name__, u'id': id(self), u't': t, u'tCoord': str(tCoord), u'result': repr(result)}) # TODO: Re-enable this.
778	else:
779	self._LogDebug(_(u'%(class)s 0x%(id)08X: Query result for z=%(z)i (%(zCoord)s): %(result)s'), {u'class': self.__class__.__name__, u'id': id(self), u'z': z, u'zCoord': repr(zCoord), u'result': repr(result)})
780
781	if result:
782	grids = [GridSlice(self._Grid, tIndex=t, zIndex=z, tQAName=self._TQAName, tQADisplayName=self._TQADisplayName, tQACoordType=self._TQACoordType, zQAName=self._ZQAName, zQADisplayName=self._ZQADisplayName, zQACoordType=self._ZQACoordType)]
783	if progressReporter is not None:
784	progressReporter.ReportProgress()
785	else:
786	grids = []
787
788	return grids
789
790
791	class ClippedGrid(Grid):
792	__doc__ = DynamicDocString()
793
794	def __init__(self, grid, clipBy=u'Cell indices', xMin=None, xMax=None, yMin=None, yMax=None, zMin=None, zMax=None, tMin=None, tMax=None):
795	self.__class__.__doc__.Obj.ValidateMethodInvocation()
796
797	# Validate the provided indices and convert them to a list of
798	# slices with positive indices.
799
800	sliceList = [self._GetSlicesForClippedExtent(grid, 'y', clipBy, yMin, yMax), self._GetSlicesForClippedExtent(grid, 'x', clipBy, xMin, xMax)]
801
802	if 'z' in grid.Dimensions:
803	sliceList = [self._GetSlicesForClippedExtent(grid, 'z', clipBy, zMin, zMax)] + sliceList
804	elif zMin is not None or zMax is not None:
805	raise ValueError(_(u'Values were provided for zMin and/or zMax but %(dn)s does not have a z dimension.') % {u'dn': grid.DisplayName})
806
807	if 't' in grid.Dimensions:
808	sliceList = [self._GetSlicesForClippedExtent(grid, 't', clipBy, tMin, tMax)] + sliceList
809	elif tMin is not None or tMax is not None:
810	raise ValueError(_(u'Values were provided for tMin and/or tMax but %(dn)s does not have a t dimension.') % {u'dn': grid.DisplayName})
811
812	# Initialize our properties.
813
814	self._Grid = grid
815	self._SliceList = sliceList
816
817	sliceListForDisplayName = []
818	for i in range(len(grid.Dimensions)):
819	if sliceList[i].start > 0:
820	sliceListForDisplayName.append(_(u'%(dim)sMin = %(index)i') % {u'dim': grid.Dimensions[i], u'index': sliceList[i].start})
821	if sliceList[i].stop < grid.Shape[i]:
822	sliceListForDisplayName.append(_(u'%(dim)sMax = %(index)i') % {u'dim': grid.Dimensions[i], u'index': sliceList[i].stop - 1})
823
824	if len(sliceListForDisplayName) > 0:
825	self._DisplayName = _(u'%(dn)s, clipped to indices %(indices)s') % {u'dn': self._Grid.DisplayName, u'indices': u', '.join(sliceListForDisplayName)}
826	else:
827	self._DisplayName = self._Grid.DisplayName
828
829	# Initialize the base class.
830
831	super(ClippedGrid, self).__init__(self._Grid.ParentCollection, queryableAttributes=self._Grid._QueryableAttributes, queryableAttributeValues=self._Grid._QueryableAttributeValues)
832
833	# Our goal is to imitate the contained grid except with a
834	# smaller extent. In order to do this, we have to override the
835	# lazy properties for the shape and corner coordinates. This
836	# task is complicated because we have to expose the same
837	# queryable attributes and have the same parent collection,
838	# and those could be used to set the shape and corner
839	# coordinates. Thus we can't just wait to be called at
840	# _GetLazyPropertyPhysicalValue to return the shape and corner
841	# coordinates because if a queryable attribute sets them, we
842	# won't ever be called.
843	#
844	# To work around this, set the shape now (we know it already)
845	# and see if the corner coordinates are available from the
846	# contained grid without accessing physical storage (i.e. they
847	# are already cached by the contained grid or are set by a
848	# queryable attribute of it or its parents). If they are, then
849	# set our modified values now. Otherwise, do nothing; we know
850	# that we'll be called at _GetLazyPropertyPhysicalValue when
851	# they are needed, and we can retrieve and modify the values
852	# from the contained grid at that point.
853
854	self.SetLazyPropertyValue('Shape', tuple(map(lambda s: s.stop - s.start, sliceList)))
855
856	oldCornerCoords = grid.GetLazyPropertyValue('CornerCoords', allowPhysicalValue=False)
857	if oldCornerCoords is not None:
858	self.SetLazyPropertyValue('CornerCoords', self._GetNewCornerCoords(oldCornerCoords, grid, sliceList))
859
860	def _Close(self):
861	if hasattr(self, '_Grid') and self._Grid is not None:
862	self._Grid.Close()
863	super(ClippedGrid, self)._Close()
864
865	def _GetDisplayName(self):
866	return self._DisplayName
867
868	def _GetLazyPropertyPhysicalValue(self, name):
869	if name == 'CornerCoords':
870	return self._GetNewCornerCoords(self._Grid.GetLazyPropertyValue('CornerCoords'), self._Grid, self._SliceList)
871	return self._Grid.GetLazyPropertyValue(name)
872
873	@classmethod
874	def _TestCapability(cls, capability):
875	return cls._Grid._TestCapability(capability)
876
877	@classmethod
878	def _GetSRTypeForSetting(cls):
879	return cls._Grid._GetSRTypeForSetting()
880
881	def _SetSpatialReference(self, sr):
882	return self._Grid._SetSpatialReference(sr)
883
884	def _GetCoords(self, coord, coordNum, slices, sliceDims, fixedIncrementOffset):
885	return self._Grid._GetCoords(coord, coordNum, slices, sliceDims, fixedIncrementOffset)
886
887	@classmethod
888	def _GetSlicesForClippedExtent(cls, grid, dim, clipBy, start, stop):
889	dimNum = grid.Dimensions.index(dim)
890
891	if start is not None:
892	if clipBy == u'cell indices':
893	start = int(start)
894	if start < 0:
895	absStart = grid.Shape[dimNum] + start
896	else:
897	absStart = start
898	if absStart < 0 or absStart > grid.Shape[dimNum] - 1:
899	raise IndexError(_(u'%(dim)sMin (%(value)i) is out of range.') % {u'dim': dim, u'value': start})
900	else:
901	absStart = bisect.bisect_left(grid.CenterCoords[dim], start)
902	if absStart > grid.Shape[dimNum] - 1:
903	raise IndexError(_(u'%(dim)sMin (%(value)s) is out of range. It must be less than or equal to %(max)s, the %(dim)s coordinate of the center of the right-most cell.') % {u'dim': dim, u'value': repr(start), u'max': repr(grid.CenterCoords[dim, -1])})
904	else:
905	absStart = 0
906
907	if stop is not None:
908	if clipBy == u'cell indices':
909	stop = int(stop)
910	if stop < 0:
911	absStop = grid.Shape[dimNum] + stop
912	else:
913	absStop = stop
914	if absStop < 0 or absStop > grid.Shape[dimNum]:
915	raise IndexError(_(u'%(dim)sMax (%(value)i) is out of range.') % {u'dim': dim, u'value': stop})
916	else:
917	if start is not None and stop < start:
918	raise IndexError(_(u'%(dim)sMin (%(value1)s) is greater than %(dim)sMax (%(value2)s). %(dim)sMax must be greater than %(dim)sMin.') % {u'dim': dim, u'value1': repr(start), u'value2': repr(stop)})
919	absStop = bisect.bisect_right(grid.CenterCoords[dim], stop)
920	if absStop == 0:
921	raise IndexError(_(u'%(dim)sMax (%(value)s) is out of range. It must be greater than or equal to %(min)s, the %(dim)s coordinate of the center of the left-most cell.') % {u'dim': dim, u'value': repr(stop), u'min': repr(grid.CenterCoords[dim, 0])})
922	else:
923	absStop = grid.Shape[dimNum]
924
925	if absStart == absStop:
926	if clipBy == u'cell indices':
927	raise IndexError(_(u'%(dim)sMin and %(dim)sMax are the same (%(value)i). %(dim)sMax must be greater than %(dim)sMin.') % {u'dim': dim, u'value': start})
928	else:
929	raise IndexError(_(u'%(dim)sMin and %(dim)sMax do not enclose the center of at least one cell. The clipped grid must have at least one cell along the %(dim)s axis. For this to happen, %(dim)sMin and %(dim)sMax must enclose the center of at least one cell.') % {u'dim': dim})
930	elif absStart > absStop:
931	raise IndexError(_(u'%(dim)sMin (%(value1)i) is greater than %(dim)sMax (%(value2)i). %(dim)sMax must be greater than %(dim)sMin.') % {u'dim': dim, u'value1': start, u'value2': stop})
932
933	return slice(absStart, absStop)
934
935	@classmethod
936	def _GetNewCornerCoords(cls, oldCornerCoords, grid, sliceList): # TODO: Does this need to be a classmethod?
937	newCornerCoords = list(oldCornerCoords)
938
939	for i in range(len(newCornerCoords)):
940	if newCornerCoords[i] is not None:
941	if grid.Dimensions[i] != 't':
942	if grid.CoordDependencies[i] is None:
943	newCornerCoords[i] = grid.CenterCoords[grid.Dimensions[i], sliceList[i].start]
944	else:
945	key = [grid.Dimensions[i]]
946	for d in grid.Dimensions:
947	if d == grid.Dimensions[i] or d in grid.CoordDependencies[i]:
948	key.append(sliceList[i].start)
949	newCornerCoords[i] = grid.CenterCoords.__getitem__(tuple(key))
950	else:
951	tCornerCoordType = grid.GetLazyPropertyValue('TCornerCoordType')
952	if tCornerCoordType == 'min':
953	newCornerCoords[i] = grid.MinCoords['t', sliceList[i].start]
954	elif tCornerCoordType == 'center':
955	newCornerCoords[i] = grid.CenterCoords['t', sliceList[i].start]
956	else:
957	newCornerCoords[i] = grid.MaxCoords['t', sliceList[i].start]
958
959	return tuple(newCornerCoords)
960
961	def _GetCoordsForOffset(self, key, fixedIncrementOffset):
962
963	# Validate the key.
964
965	coord, coordNum, slices, sliceDims = self._GetSlicesForCoordsKey(key)
966
967	# Adjust the slices list, which could be None, or a list of
968	# slices and/or integers.
969
970	if slices is not None:
971	for i in range(0, len(sliceDims)):
972	dimNum = self.Dimensions.index(sliceDims[i])
973
974	if isinstance (slices[i], int):
975	slices[i] = self._AdjustCoord(slices[i], dimNum)
976
977	elif slices[i].start is None and slices[i].stop is None:
978	if slices[i].step is not None and slices[i].step < 0:
979	if self._SliceList[coordNum].start > 0:
980	slices[i] = slice(self._SliceList[coordNum].stop - 1, self._SliceList[coordNum].start - 1, slices[i].step)
981	else:
982	slices[i] = slice(self._SliceList[coordNum].stop - 1, None, slices[i].step)
983	else:
984	slices[i] = slice(self._SliceList[coordNum].start, self._SliceList[coordNum].stop, slices[i].step)
985
986	elif slices[i].start is not None and slices[i].stop is None:
987	if slices[i].step is not None and slices[i].step < 0:
988	if self._SliceList[coordNum].start > 0:
989	slices[i] = slice(self._AdjustCoord(slices[i].start, dimNum), self._SliceList[coordNum].start - 1, slices[i].step)
990	else:
991	slices[i] = slice(self._AdjustCoord(slices[i].start, dimNum), None, slices[i].step)
992	else:
993	slices[i] = slice(self._AdjustCoord(slices[i].start, dimNum), self._SliceList[coordNum].stop, slices[i].step)
994
995	elif slices[i].start is None and slices[i].stop is not None:
996	if slices[i].step is not None and slices[i].step < 0:
997	slices[i] = slice(self._SliceList[coordNum].stop - 1, self._AdjustCoord(slices[i].stop, dimNum), slices[i].step)
998	else:
999	slices[i] = slice(self._SliceList[coordNum].start, self._AdjustCoord(slices[i].stop, dimNum), slices[i].step)
1000
1001	else:
1002	slices[i] = slice(self._AdjustCoord(slices[i].start, dimNum), self._AdjustCoord(slices[i].stop, dimNum), slices[i].step)
1003	else:
1004	slices = []
1005	for i in range(0, len(sliceDims)):
1006	dimNum = self.Dimensions.index(sliceDims[i])
1007	slices.append(self._SliceList[dimNum])
1008
1009	# Get the coordinates from the contained grid.
1010
1011	return self._Grid._GetCoordsForOffset(tuple([coord] + slices), fixedIncrementOffset)
1012
1013	def _AdjustCoord(self, c, coordNum):
1014	if c >= 0:
1015	return c + self._SliceList[coordNum].start
1016	return c - (self._Grid.Shape[coordNum] - self._SliceList[coordNum].stop)
1017
1018	def _ReadNumpyArray(self, sliceList):
1019
1020	# Because the we did not override the PhysicalDimensions lazy
1021	# property, the sliceList is ordered according to
1022	# PhysicalDimensions, not Dimensions. But self._SliceList is
1023	# ordered according to Dimensions, so we have to be careful
1024	# about pulling indices from it when computing the slice to
1025	# retreive from the contained grid.
1026
1027	slicesToRead = []
1028	for i, d in enumerate(self.GetLazyPropertyValue('PhysicalDimensions')):
1029	offset = self._SliceList[self.Dimensions.index(d)].start
1030	slicesToRead.append(slice(sliceList[i].start + offset, sliceList[i].stop + offset, sliceList[i].step))
1031
1032	return self._Grid._ReadNumpyArray(slicesToRead)
1033
1034	def _WriteNumpyArray(self, sliceList, data):
1035
1036	# Because the we did not override the PhysicalDimensions lazy
1037	# property, the sliceList is ordered according to
1038	# PhysicalDimensions, not Dimensions. But self._SliceList is
1039	# ordered according to Dimensions, so we have to be careful
1040	# about pulling indices from it when computing the slice to
1041	# write to the contained grid.
1042
1043	slicesToWrite = []
1044	for i, d in enumerate(self.GetLazyPropertyValue('PhysicalDimensions')):
1045	offset = self._SliceList[self.Dimensions.index(d)].start
1046	slicesToWrite.append(slice(sliceList[i].start + offset, sliceList[i].stop + offset, sliceList[i].step))
1047
1048	self._Grid._WriteNumpyArray(slicesToWrite, data)
1049
1050
1051	class MaskedGrid(Grid):
1052	__doc__ = DynamicDocString()
1053
1054	def __init__(self, grid, masks, operators, values, unscaledNoDataValue=None, scaledNoDataValue=None):
1055	# TODO: self.__class__.__doc__.Obj.ValidateMethodInvocation()
1056
1057	# Initialize our properties.
1058
1059	if len(values) < len(operators):
1060	values = values + [None] * (len(operators) - len(values))
1061
1062	self._Grid = grid
1063	self._Masks = masks
1064	self._Operators = operators
1065	self._Values = values
1066	self._UnscaledNoDataValue = unscaledNoDataValue
1067	self._ScaledNoDataValue = scaledNoDataValue
1068	self._MaskOffsets = None
1069
1070	self._DisplayName = _(u'%(dn)s, masked where %(maskExpressions)s') % {u'dn': grid.DisplayName, u'maskExpressions': _(u' or ').join(map(lambda s: self._GetMaskDisplayExpression(*s), zip(masks, operators, values)))}
1071
1072	# Initialize the base class.
1073
1074	super(MaskedGrid, self).__init__(self._Grid.ParentCollection, queryableAttributes=self._Grid._QueryableAttributes, queryableAttributeValues=self._Grid._QueryableAttributeValues)
1075
1076	# Our goal is to imitate the contained grid except with ideal
1077	# values for PhysicalDimensions or PhysicalDimensionsFlipped
1078	# (the contained grid takes care of reordering, if needed) and
1079	# with an UnscaledNoDataValue and ScaledNoDataValue (if the
1080	# contained grid does not have any). This task is complicated
1081	# because we have to expose the same queryable attributes and
1082	# have the same parent collection, and those could be used to
1083	# set those lazy properties. Thus we can't just wait to be
1084	# called at _GetLazyPropertyPhysicalValue to return the values
1085	# because if a queryable attribute sets them, we won't ever be
1086	# called.
1087	#
1088	# To work around this, see if values are available without
1089	# accessing physical storage (i.e. they are set by a queryable
1090	# attribute). If they are, then set our modified values now.
1091	# Otherwise, do nothing; we know that we'll be called at
1092	# _GetLazyPropertyPhysicalValue when they are needed.
1093
1094	if self.HasLazyPropertyValue('PhysicalDimensions', allowPhysicalValue=False):
1095	self.SetLazyPropertyValue('PhysicalDimensions', self._GetLazyPropertyPhysicalValue('PhysicalDimensions'))
1096
1097	if self.HasLazyPropertyValue('PhysicalDimensionsFlipped', allowPhysicalValue=False):
1098	self.SetLazyPropertyValue('PhysicalDimensionsFlipped', self._GetLazyPropertyPhysicalValue('PhysicalDimensionsFlipped'))
1099
1100	if self.HasLazyPropertyValue('UnscaledNoDataValue', allowPhysicalValue=False):
1101	self.SetLazyPropertyValue('UnscaledNoDataValue', self._GetLazyPropertyPhysicalValue('UnscaledNoDataValue'))
1102
1103	if self.HasLazyPropertyValue('ScaledNoDataValue', allowPhysicalValue=False):
1104	self.SetLazyPropertyValue('ScaledNoDataValue', self._GetLazyPropertyPhysicalValue('ScaledNoDataValue'))
1105
1106	def _Close(self):
1107	if hasattr(self, '_Grid') and self._Grid is not None:
1108	self._Grid.Close()
1109	for mask in self._Masks:
1110	mask.Close()
1111	super(MaskedGrid, self)._Close()
1112
1113	def _GetDisplayName(self):
1114	return self._DisplayName
1115
1116	def _GetMaskDisplayExpression(self, mask, op, value):
1117	if op in [u'=', u'==', u'!=', u'<>', u'<', u'<=', u'>', u'>']:
1118	return mask.DisplayName + u' ' + op + u' ' + repr(value)
1119
1120	if op in [u'any', u'all']:
1121	if mask.DataType.endswith('8'):
1122	bits = 8
1123	elif mask.DataType.endswith('16'):
1124	bits = 16
1125	elif mask.DataType.endswith('32'):
1126	bits = 32
1127	else:
1128	bits = 64
1129
1130	binaryReprOfValue = ''.join([str(int(value & 1 << i != 0)) for i in range(bits-1, -1, -1)])
1131
1132	if op == u'any':
1133	return mask.DisplayName + u' & 0b' + binaryReprOfValue + ' != 0'
1134	return mask.DisplayName + u' & 0b' + binaryReprOfValue + ' == 0b' + binaryReprOfValue
1135
1136	return mask.DisplayName + u' ' + op
1137
1138	def _GetLazyPropertyPhysicalValue(self, name):
1139
1140	# If the requested property is PhysicalDimensions or
1141	# PhysicalDimensionsFlipped, return values indicating the
1142	# dimensions in the ideal order. The contained grid takes care
1143	# of reordering, if needed.
1144
1145	if name == 'PhysicalDimensions':
1146	return self.Dimensions
1147
1148	if name == 'PhysicalDimensionsFlipped':
1149	return tuple([False] * len(self.Dimensions))
1150
1151	# If the requested property is the UnscaledNoDataValue or
1152	# ScaledNoDataValue, determine the value we should return.
1153
1154	if name == 'UnscaledNoDataValue':
1155	value = self._Grid.GetLazyPropertyValue(name)
1156	if value is not None:
1157	self._LogDebug(_(u'%(class)s 0x%(id)08X: Using the UnscaledNoDataValue of the contained grid (%(value)s).'), {u'class': self.__class__.__name__, u'id': id(self), u'value': repr(value)})
1158	return value
1159
1160	if self._UnscaledNoDataValue is not None:
1161	self._LogDebug(_(u'%(class)s 0x%(id)08X: Using the UnscaledNoDataValue supplied to the MaskedGrid constructor (%(value)s).'), {u'class': self.__class__.__name__, u'id': id(self), u'value': repr(self._UnscaledNoDataValue)})
1162	return self._UnscaledNoDataValue
1163
1164	value = self._GetNoDataValueForDataType(self.UnscaledDataType)
1165	self._LogDebug(_(u'%(class)s 0x%(id)08X: Using the default UnscaledNoDataValue (%(value)s) for UnscaledDataType %(dt)s. The contained grid does not have an UnscaledNoDataValue, nor was one supplied to the MaskedGrid constructor.'), {u'class': self.__class__.__name__, u'id': id(self), u'value': repr(value), u'dt': self.UnscaledDataType})
1166	return value
1167
1168	if name == 'ScaledNoDataValue':
1169	value = self._Grid.GetLazyPropertyValue(name)
1170	if value is not None:
1171	self._LogDebug(_(u'%(class)s 0x%(id)08X: Using the ScaledNoDataValue of the contained grid (%(value)s).'), {u'class': self.__class__.__name__, u'id': id(self), u'value': repr(value)})
1172	return value
1173
1174	if self._ScaledNoDataValue is not None:
1175	self._LogDebug(_(u'%(class)s 0x%(id)08X: Using the ScaledNoDataValue supplied to the MaskedGrid constructor (%(value)s).'), {u'class': self.__class__.__name__, u'id': id(self), u'value': repr(self._ScaledNoDataValue)})
1176	return self._ScaledNoDataValue
1177
1178	value = self._GetNoDataValueForDataType(self.DataType)
1179	if value is not None:
1180	self._LogDebug(_(u'%(class)s 0x%(id)08X: Using the default ScaledNoDataValue (%(value)s) for DataType %(dt)s. The contained grid does not have an ScaledNoDataValue, nor was one supplied to the MaskedGrid constructor.'), {u'class': self.__class__.__name__, u'id': id(self), u'value': repr(value), u'dt': self.DataType})
1181	return value
1182
1183	# Otherwise just get the unaltered value from the contained
1184	# grid.
1185
1186	return self._Grid.GetLazyPropertyValue(name)
1187
1188	@classmethod
1189	def _GetNoDataValueForDataType(cls, dataType):
1190	if dataType == u'int8':
1191	return -128
1192	if dataType == u'uint8':
1193	return 255
1194	if dataType == u'int16':
1195	return -32768
1196	if dataType == u'uint16':
1197	return 65535
1198	if dataType == u'int32':
1199	return -2147483647 # To improve compatibility with ArcGIS, we use -2147483647 rather than -2147483648
1200	if dataType == u'uint32':
1201	return 4294967295L
1202	if dataType == u'float32':
1203	return -3.4028234663852886e+038 # This is the float64 representation of the float32 -3.40282347e+38
1204	if dataType == u'float64':
1205	return -1.7976931348623157e+308
1206	return None
1207
1208	@classmethod
1209	def _TestCapability(cls, capability):
1210	return cls._Grid._TestCapability(capability)
1211
1212	@classmethod
1213	def _GetSRTypeForSetting(cls):
1214	return cls._Grid._GetSRTypeForSetting()
1215
1216	def _SetSpatialReference(self, sr):
1217	return self._Grid._SetSpatialReference(sr)
1218
1219	def _GetCoords(self, coord, coordNum, slices, sliceDims, fixedIncrementOffset):
1220	return self._Grid._GetCoords(coord, coordNum, slices, sliceDims, fixedIncrementOffset)
1221
1222	def _ReadNumpyArray(self, sliceList):
1223
1224	# If we have not yet validated the masks and computed the
1225	# offsets into them, do it now.
1226
1227	if self._MaskOffsets is None:
1228	self._ValidateMasksAndSetOffsets()
1229
1230	# Get the unscaled data from the contained grid.
1231
1232	data = self._Grid.UnscaledData.__getitem__(tuple(sliceList))
1233
1234	# Apply each of the masks.
1235
1236	for i in range(len(self._Masks)):
1237	if self._Operators[i] in [u'=', u'==', u'!=', u'<>', u'<', u'<=', u'>', u'>', u'any', u'all']:
1238	self._LogDebug(_(u'%(class)s 0x%(id)08X: Masking %(dn2)s where %(dn1)s %(op)s %(value)s.'), {u'class': self.__class__.__name__, u'id': id(self), u'dn2': self._Grid.DisplayName, u'dn1': self._Masks[i].DisplayName, u'op': self._Operators[i], u'value': self._Values[i]})
1239	else:
1240	self._LogDebug(_(u'%(class)s 0x%(id)08X: Masking %(dn2)s where %(dn1)s %(op)s.'), {u'class': self.__class__.__name__, u'id': id(self), u'dn2': self._Grid.DisplayName, u'dn1': self._Masks[i].DisplayName, u'op': self._Operators[i]})
1241
1242	# Create a list of slices for retrieving the mask.
1243
1244	maskSliceList = []
1245	for d in range(len(self.Dimensions)):
1246	if self.Dimensions[d] in self._Masks[i].Dimensions:
1247	offset = self._MaskOffsets[i][self._Masks[i].Dimensions.index(self.Dimensions[d])]
1248	maskSliceList.append(slice(sliceList[d].start + offset, sliceList[d].stop + offset))
1249
1250	# Get the mask.
1251
1252	maskData = self._Masks[i].Data.__getitem__(tuple(maskSliceList))
1253
1254	# Perform the requested test on the mask.
1255
1256	if self._Operators[i] in [u'=', u'==']:
1257	maskResult = maskData == self._Values[i]
1258	elif self._Operators[i] in [u'!=', u'<>']:
1259	maskResult = maskData != self._Values[i]
1260	elif self._Operators[i] == u'<':
1261	maskResult = maskData < self._Values[i]
1262	elif self._Operators[i] == u'<=':
1263	maskResult = maskData <= self._Values[i]
1264	elif self._Operators[i] == u'>':
1265	maskResult = maskData > self._Values[i]
1266	elif self._Operators[i] == u'>=':
1267	maskResult = maskData >= self._Values[i]
1268	elif self._Operators[i] == u'any':
1269	maskResult = maskData & self._Values[i] != 0
1270	elif self._Operators[i] == u'all':
1271	maskResult = maskData & self._Values[i] == self._Values[i]
1272	else:
1273	raise ValueError(_(u'Unknown mask operator "%(op)s".') % {u'op': self._Operators[i]})
1274
1275	# Set all cells where the test is True to
1276	# self.UnscaledNoDataValue. Handle the cases where the
1277	# mask has fewer dimensions than the data (for example,
1278	# when a land mask with dimensions yx is used to mask a
1279	# time series of satellite images with dimensions tyx).
1280
1281	if self.Dimensions == self._Masks[i].Dimensions:
1282	data[maskResult] = self.UnscaledNoDataValue
1283	elif self.Dimensions in [u'zyx', u'tyx'] and self._Masks[i].Dimensions == u'yx' or self.Dimensions == u'tzyx' and self._Masks[i].Dimensions == u'zyx':
1284	data[:, maskResult] = self.UnscaledNoDataValue
1285	elif self.Dimensions == u'tzyx' and self._Masks[i].Dimensions == u'tyx':
1286	data.transpose((1,0,2,3))[:, maskResult] = self.UnscaledNoDataValue
1287	else: # self.Dimensions must be u'tzyx' and self._Masks[i].Dimensions must be u'yx'
1288	data[:, :, maskResult] = self.UnscaledNoDataValue
1289
1290	# Return the data and self.UnscaledNoDataValue.
1291
1292	return data, self.UnscaledNoDataValue
1293
1294	def _ValidateMasksAndSetOffsets(self):
1295	maskOffsets = []
1296
1297	for mask in self._Masks:
1298
1299	# Validate that the mask dimensions are a subset of the
1300	# contained grid's dimensions.
1301
1302	if len(self.Dimensions) < len(mask.Dimensions) or not (mask.Dimensions in [u'yx', u'tzyx'] or mask.Dimensions == u'zyx' and self.Dimensions in [u'zyx', 'tzyx'] or mask.Dimensions == u'tyx' and self.Dimensions in [u'tyx', 'tzyx']):
1303	raise ValueError(_(u'%(dn1)s has dimensions (%(dim1)s) that are incompatible with the dimensions of %(dn2)s (%(dim2)s), so it cannot be used as a mask.') % {u'dn1': mask.DisplayName, u'dn2': self._Grid.DisplayName, u'dim1': mask.Dimensions, u'dim2': self._Grid.Dimensions})
1304
1305	# Validate that the mask uses the same coordinate system
1306	# as the contained grid.
1307
1308	if not self.GetSpatialReference('obj').IsSame(mask.GetSpatialReference('obj')):
1309	raise ValueError(_(u'%(dn1)s has a different coordinate system than %(dn2)s, so it cannot be used as a mask.') % {u'dn1': mask.DisplayName, u'dn2': self._Grid.DisplayName})
1310
1311	# If the contained grid has any dimensions for which the
1312	# coordinates depend on the coordinates of other
1313	# dimensions (e.g. the value of z depends on t, y, and x,
1314	# as is the case with certain ROMS datasets), then we
1315	# require that the shape of the mask exactly match that of
1316	# the contained grid (for the mask's dimensions). We do
1317	# not verify the coordinates of the mask as this could be
1318	# a very time consuming operation.
1319
1320	if self.CoordDependencies != tuple([None] * len(self.Dimensions)):
1321	for i in range(len(mask.Dimensions)):
1322	if mask.Shape[i] != self.Shape[self.Dimensions.index(mask.Dimensions[i])]:
1323	raise ValueError(_(u'Because %(dn2)s has dimensions for which the coordinates depend on other dimensions, the MaskedGrid class can only apply masks that have dimensions with the same length as it. %(dn1)s cannot be used as a mask because the %(dim)s dimension has a different length (the mask has length %(len1)i but the grid has length %(len2)i).') % {u'dn1': mask.DisplayName, u'dn2': self._Grid.DisplayName, u'dim': mask.Dimensions[i], u'len1': mask.Shape[i], u'len2': self.Shape[self.Dimensions.index[mask.Dimensions[i]]]})
1324
1325	maskOffsets.append([0] * len(mask.Dimensions))
1326
1327	# Otherwise (there are no coordinate dependencies), verify
1328	# that the mask encloses the grid and has the same
1329	# coordinates.
1330
1331	else:
1332	import numpy
1333
1334	offsets = []
1335
1336	for i in range(len(mask.Dimensions)):
1337	if mask.Dimensions[i] != 't':
1338	gridCoords = numpy.cast['float32'](self.CenterCoords[mask.Dimensions[i]])
1339	maskCoords = numpy.cast['float32'](mask.CenterCoords[mask.Dimensions[i]])
1340	else:
1341	gridCoords = self.CenterCoords[mask.Dimensions[i]]
1342	maskCoords = mask.CenterCoords[mask.Dimensions[i]]
1343
1344	if gridCoords[0] < maskCoords[0] or gridCoords[-1] > maskCoords[-1]:
1345	raise ValueError(_(u'%(dn1)s does not completely enclose %(dn2)s so it cannot be used as a mask.') % {u'dn1': mask.DisplayName, u'dn2': self._Grid.DisplayName})
1346
1347	try:
1348	offset = maskCoords.searchsorted(gridCoords[0])
1349	if len(maskCoords) - offset < len(gridCoords):
1350	raise ValueError
1351	if (maskCoords[offset:offset+len(gridCoords)] != gridCoords).any():
1352	raise ValueError
1353	except:
1354	raise ValueError(_(u'The %(dim)s coordinates of %(dn1)s do not line up with the %(dim)s coordinates of %(dn2)s, so it cannot be used as a mask.') % {u'dn1': mask.DisplayName, u'dn2': self._Grid.DisplayName, u'dim': mask.Dimensions[i]})
1355
1356	offsets.append(offset)
1357
1358	maskOffsets.append(offsets)
1359
1360	self._MaskOffsets = maskOffsets
1361
1362
1363	class DerivedGrid(Grid):
1364	__doc__ = DynamicDocString()
1365
1366	def __init__(self, grids, expression, displayName, dataType, noDataValue=None, queryableAttributes=None, queryableAttributeValues=None):
1367	# TODO: self.__class__.__doc__.Obj.ValidateMethodInvocation()
1368
1369	# Initialize our properties.
1370
1371	self._Grids = grids
1372	self._Expression = expression
1373	self._DisplayName = displayName
1374	self._ValidatedGrids = False
1375
1376	# Initialize the base class.
1377
1378	super(DerivedGrid, self).__init__(queryableAttributes=queryableAttributes, queryableAttributeValues=queryableAttributeValues)
1379
1380	# Set lazy properties that override those of the contained
1381	# grids.
1382
1383	self.SetLazyPropertyValue('UnscaledDataType', dataType)
1384	self.SetLazyPropertyValue('ScaledDataType', None)
1385	self.SetLazyPropertyValue('UnscaledNoDataValue', noDataValue)
1386	self.SetLazyPropertyValue('ScaledNoDataValue', None)
1387	self.SetLazyPropertyValue('ScalingFunction', None)
1388	self.SetLazyPropertyValue('UnscalingFunction', None)
1389
1390	def _Close(self):
1391	if hasattr(self, '_Grids') and self._Grids is not None:
1392	for grid in self._Grids:
1393	grid.Close()
1394	super(DerivedGrid, self)._Close()
1395
1396	def _GetDisplayName(self):
1397	return self._DisplayName
1398
1399	def _GetLazyPropertyPhysicalValue(self, name):
1400
1401	# If the requested property is PhysicalDimensions or
1402	# PhysicalDimensionsFlipped, return values indicating the
1403	# dimensions in the ideal order. The contained grids take care
1404	# of reordering, if needed.
1405
1406	if name == 'PhysicalDimensions':
1407	return self.Dimensions
1408
1409	if name == 'PhysicalDimensionsFlipped':
1410	return tuple([False] * len(self.Dimensions))
1411
1412	# Otherwise just get the unaltered value from the first
1413	# contained grid.
1414
1415	return self._Grids[0].GetLazyPropertyValue(name)
1416
1417	@classmethod
1418	def _TestCapability(cls, capability):
1419	if capability == 'SetSpatialReference':
1420	return NotImplementedError(_(u'DerivedGrid does not support setting the spatial reference.'))
1421	return cls._Grids[0]._TestCapability(capability)
1422
1423	def _GetCoords(self, coord, coordNum, slices, sliceDims, fixedIncrementOffset):
1424	return self._Grids[0]._GetCoords(coord, coordNum, slices, sliceDims, fixedIncrementOffset)
1425
1426	def _ReadNumpyArray(self, sliceList):
1427
1428	# If we have not yet validated the grids, do so now.
1429
1430	if not self._ValidatedGrids:
1431	if len(self._Grids) > 1:
1432	for grid in self._Grids[1:]:
1433	if grid.Shape != self._Grids[0].Shape:
1434	raise ValueError(_(u'%(dn1)s does not have the same shape as %(dn2)s.') % {u'dn1': grid.DisplayName, u'dn1': self._Grids[0].DisplayName})
1435
1436	# TODO: More validation?
1437
1438	self._ValidatedGrids = True
1439
1440	# Evaluate the expression.
1441
1442	import numpy
1443
1444	data = eval(self._Expression)
1445
1446	# If a no data value was provided, mask any cells that have no
1447	# data in any of the grids.
1448
1449	if self.NoDataValue is not None:
1450	for grid in self._Grids:
1451	if grid.NoDataValue is not None:
1452	data[grid.Data.__getitem__(tuple(sliceList)) == grid.NoDataValue] = self.NoDataValue
1453
1454	# Return the data and no data value.
1455
1456	return data, self.NoDataValue
1457
1458
1459	class RotatedGlobalGrid(Grid):
1460	__doc__ = DynamicDocString()
1461
1462	def __init__(self, grid, rotationOffset, rotationUnits=u'Map units'):
1463	self.__class__.__doc__.Obj.ValidateMethodInvocation()
1464
1465	# TODO: Validate that the grid is global?
1466
1467	# Validate that the grid has a constant x increment.
1468
1469	if grid.CoordDependencies[-1] is not None:
1470	raise ValueError(_(u'The provided grid, %(dn)s, does not have a constant x increment. The current implementation of RotatedGlobalGrid only supports grids with constant x increments.') % {u'dn': grid.DisplayName})
1471
1472	# Initialize our properties.
1473
1474	self._Grid = grid
1475	self._XRotationType = rotationUnits
1476
1477	if self._XRotationType == u'cells':
1478	self._XRotationInCells = int(round(rotationOffset))
1479	self._XRotationInMapUnits = self._XRotationInCells * self._Grid.CoordIncrements[-1]
1480	else:
1481	self._XRotationInCells = int(round(rotationOffset / self._Grid.CoordIncrements[-1]))
1482	self._XRotationInMapUnits = self._XRotationInCells * self._Grid.CoordIncrements[-1]
1483
1484	self._DisplayName = _(u'%(dn)s, rotated about the polar axis by %(cells)i cells (%(mapUnits)g map units)') % {u'dn': self._Grid.DisplayName, u'cells': self._XRotationInCells, u'mapUnits': self._XRotationInMapUnits}
1485
1486	# Determine lazy properties that we can cheaply calculate now.
1487	# Transposing and flipping the raw data is handled by the
1488	# contained grid, so we use idealized PhysicalDimensions and
1489	# set PhysicalDimensionsFlipped to False for all dimensions.
1490
1491	lazyPropertyValues = {'PhysicalDimensions': self._Grid.Dimensions,
1492	'PhysicalDimensionsFlipped': tuple([False] * len(self._Grid.Dimensions))}
1493
1494	# Determine whether the contained grid uses a geographic
1495	# coordinate system or a projected coordinate system. If it is
1496	# geographic, we will use the same coordinate system as that
1497	# grid but a different x corner coordinate. To do this, we
1498	# have to override the CornerCoords lazy property. This task
1499	# is complicated because we have to expose the same queryable
1500	# attributes and have the same parent collection, and those
1501	# could be used to set the corner coordinates. Thus we can't
1502	# just wait to be called at _GetLazyPropertyPhysicalValue to
1503	# return CornerCoords because if a queryable attribute sets
1504	# it, we won't ever be called.
1505	#
1506	# To work around this, see if CornerCoords available from the
1507	# contained grid without accessing physical storage (i.e. it
1508	# is already cached by the contained grid or is set by a
1509	# queryable attribute of it or its parents). If it is, then
1510	# set our modified value now. Otherwise, do nothing; we know
1511	# that we'll be called at _GetLazyPropertyPhysicalValue when
1512	# it is needed, and we can retrieve and modify the values from
1513	# the contained grid at that point.
1514
1515	try:
1516	oldCornerCoords = self._Grid.GetLazyPropertyValue('CornerCoords', allowPhysicalValue=False)
1517	sr = self._Grid.GetSpatialReference('obj').Clone()
1518
1519	if sr.IsGeographic() or sr.GetNormProjParm('central_meridian', 123456.) == 123456.:
1520	lazyPropertyValues['SpatialReference'] = sr
1521	if oldCornerCoords is not None:
1522	lazyPropertyValues['CornerCoords'] = tuple(list(oldCornerCoords[:-1]) + [oldCornerCoords[-1] + self._XRotationInMapUnits])
1523
1524	# If the contained grid uses a projected coordinate
1525	# system, we will use the same CornerCoords but change the
1526	# central_meridian to accomplish the rotation.
1527
1528	else:
1529	srAtPrimeMeridian = sr.Clone()
1530	srAtPrimeMeridian.SetNormProjParm('central_meridian', 0.)
1531	srGeographic = self._osr().SpatialReference()
1532	srGeographic.CopyGeogCSFrom(sr)
1533	transformer = self._osr().CoordinateTransformation(srAtPrimeMeridian, srGeographic)
1534	offsetToCentralMeridian = transformer.TransformPoint(self._XRotationInMapUnits, 0.)[0]
1535
1536	sr.SetNormProjParm('central_meridian', sr.GetNormProjParm('central_meridian') + offsetToCentralMeridian)
1537	lazyPropertyValues['SpatialReference'] = sr
1538	if oldCornerCoords is not None:
1539	lazyPropertyValues['CornerCoords'] = oldCornerCoords
1540
1541	except:
1542	self._gdal().ErrorReset()
1543	raise
1544
1545	# Initialize the base class.
1546
1547	super(RotatedGlobalGrid, self).__init__(self._Grid.ParentCollection, queryableAttributes=self._Grid._QueryableAttributes, queryableAttributeValues=self._Grid._QueryableAttributeValues, lazyPropertyValues=lazyPropertyValues)
1548
1549	def _Close(self):
1550	if hasattr(self, '_Grid') and self._Grid is not None:
1551	self._Grid.Close()
1552	super(RotatedGlobalGrid, self)._Close()
1553
1554	def _GetDisplayName(self):
1555	return self._DisplayName
1556
1557	def _GetLazyPropertyPhysicalValue(self, name):
1558	if name == 'CornerCoords':
1559	oldCornerCoords = self._Grid.GetLazyPropertyValue('CornerCoords')
1560	sr = self._Grid.GetSpatialReference('obj')
1561	if sr.IsGeographic() or sr.GetNormProjParm('central_meridian', 123456.) == 123456.:
1562	return tuple(list(oldCornerCoords[:-1]) + [oldCornerCoords[-1] + self._XRotationInMapUnits])
1563	return oldCornerCoords
1564
1565	return self._Grid.GetLazyPropertyValue(name)
1566
1567	@classmethod
1568	def _TestCapability(cls, capability):
1569	return cls._Grid._TestCapability(capability) # TODO: Reject setting the spatial reference
1570
1571	def _GetCoords(self, coord, coordNum, slices, sliceDims, fixedIncrementOffset):
1572	return self._Grid._GetCoords(coord, coordNum, slices, sliceDims, fixedIncrementOffset)
1573
1574	def _ReadNumpyArray(self, sliceList):
1575
1576	# Convert the indices of the requested slice to account for
1577	# the rotation.
1578
1579	xShape = self.Shape[-1]
1580	rotationInCells = self._XRotationInCells % xShape
1581	xStart = (sliceList[-1].start + rotationInCells) % xShape
1582	xStop = (sliceList[-1].stop - 1 + rotationInCells) % xShape + 1
1583
1584	# If the stop index is greater than the start index, it means
1585	# the requested slice does not straddle the left and right
1586	# edges of the contained grid, and we can retrieve the data
1587	# with a single request.
1588
1589	if xStop > xStart:
1590	data = self._Grid.UnscaledData.__getitem__(tuple(list(sliceList[:-1]) + [slice(xStart, xStop)]))
1591
1592	# Otherwise (the requested slice straddles), retrieve the two
1593	# slabs of data and concatenate them together.
1594
1595	else:
1596	import numpy
1597	data = numpy.concatenate((self._Grid.UnscaledData.__getitem__(tuple(list(sliceList[:-1]) + [slice(xStart, xShape)])),
1598	self._Grid.UnscaledData.__getitem__(tuple(list(sliceList[:-1]) + [slice(0, xStop)]))),
1599	axis=len(self.Dimensions) - 1)
1600
1601	# Return the data and self.UnscaledNoDataValue.
1602
1603	return data, self.UnscaledNoDataValue
1604
1605
1606	class MemoryCachedGrid(Grid):
1607	__doc__ = DynamicDocString()
1608
1609	def __init__(self, grid, maxCacheSize=None, xMinBlockSize=None, yMinBlockSize=None, zMinBlockSize=None, tMinBlockSize=None):
1610	# TODO: self.__class__.__doc__.Obj.ValidateMethodInvocation()
1611
1612	# Validate that the caller provided positive values for the
1613	# max cache size and all required dimensions. If not, report a
1614	# debug message (no need for a warning) and set everything to
1615	# zero.
1616
1617	if maxCacheSize is None or maxCacheSize <= 0 or \
1618	xMinBlockSize is None or xMinBlockSize <= 0 or \
1619	yMinBlockSize is None or yMinBlockSize <= 0 or \
1620	u'z' in grid.Dimensions and (zMinBlockSize is None or zMinBlockSize <= 0) or \
1621	u't' in grid.Dimensions and (tMinBlockSize is None or tMinBlockSize <= 0) :
1622	maxCacheSize = 0
1623	self._LogDebug(_(u'MemoryCachedGrid 0x%(id)08X: Wrapping %(dn)s. No caching will be done because one of the required parameters was missing or zero.'), {u'id': id(self), u'dn': grid.DisplayName})
1624	else:
1625	self._LogDebug(_(u'MemoryCachedGrid 0x%(id)08X: Wrapping %(dn)s. maxCacheSize=%(maxCacheSize)s, xMinBlockSize=%(xMinBlockSize)s, yMinBlockSize=%(yMinBlockSize)s, zMinBlockSize=%(zMinBlockSize)s, tMinBlockSize=%(tMinBlockSize)s.'), {u'id': id(self), u'dn': grid.DisplayName, u'maxCacheSize': repr(maxCacheSize), u'xMinBlockSize': repr(xMinBlockSize), u'yMinBlockSize': repr(yMinBlockSize), u'zMinBlockSize': repr(zMinBlockSize), u'tMinBlockSize': repr(tMinBlockSize)})
1626
1627	# Initialize our remaining properties.
1628
1629	self._Grid = grid
1630	self._MaxCacheSize = maxCacheSize
1631
1632	if maxCacheSize > 0:
1633	if grid.Dimensions == u'yx':
1634	self._MinBlockSize = (yMinBlockSize, xMinBlockSize)
1635	elif grid.Dimensions == u'zyx':
1636	self._MinBlockSize = (zMinBlockSize, yMinBlockSize, xMinBlockSize)
1637	elif grid.Dimensions == u'tyx':
1638	self._MinBlockSize = (tMinBlockSize, yMinBlockSize, xMinBlockSize)
1639	else:
1640	self._MinBlockSize = (tMinBlockSize, zMinBlockSize, yMinBlockSize, xMinBlockSize)
1641	else:
1642	self._MinBlockSize = None
1643
1644	self._Cache = []
1645	self._CacheSize = 0
1646
1647	# Initialize the base class.
1648
1649	super(MemoryCachedGrid, self).__init__(self._Grid.ParentCollection, queryableAttributes=self._Grid._QueryableAttributes, queryableAttributeValues=self._Grid._QueryableAttributeValues)
1650
1651	def _Close(self):
1652	if hasattr(self, '_Grid') and self._Grid is not None:
1653	self._Grid.Close()
1654	super(MemoryCachedGrid, self)._Close()
1655
1656	def _GetDisplayName(self):
1657	return self._Grid.DisplayName
1658
1659	def _GetLazyPropertyPhysicalValue(self, name):
1660	return self._Grid.GetLazyPropertyValue(name)
1661
1662	@classmethod
1663	def _TestCapability(cls, capability):
1664	return cls._Grid._TestCapability(capability)
1665
1666	@classmethod
1667	def _GetSRTypeForSetting(cls):
1668	return cls._Grid._GetSRTypeForSetting()
1669
1670	def _SetSpatialReference(self, sr):
1671	return self._Grid._SetSpatialReference(sr)
1672
1673	def _GetCoords(self, coord, coordNum, slices, sliceDims, fixedIncrementOffset):
1674	return self._Grid._GetCoords(coord, coordNum, slices, sliceDims, fixedIncrementOffset)
1675
1676	def _ReadNumpyArray(self, sliceList):
1677
1678	# If the max cache size is zero, just forward the call to the
1679	# contained grid.
1680
1681	if self._MaxCacheSize <= 0:
1682	return self._Grid._ReadNumpyArray(sliceList)
1683
1684	self._LogDebug(_(u'MemoryCachedGrid 0x%(id)08X: Servicing request for slice [%(slice)s].'), {u'id': id(self), u'slice': ','.join(map(lambda s: str(s.start) + ':' + str(s.stop), sliceList))})
1685
1686	# Allocate the numpy array that we will return.
1687
1688	import numpy
1689	data = numpy.zeros(map(lambda s: s.stop-s.start, sliceList), dtype=str(self.UnscaledDataType))
1690
1691	# Try to fill the array completely from the cache.
1692
1693	remainingRegion = []
1694	remainingRegion.extend(sliceList)
1695
1696	if self._CacheSize > 0:
1697	remainingRegion = self._FillArrayFromCache(sliceList, remainingRegion, data)
1698
1699	# If we could not fill the entire array, retrieve and cache
1700	# the region that was not filled in, and fill it in.
1701
1702	if remainingRegion is not None:
1703	regionToRead = []
1704
1705	for i in range(len(sliceList)):
1706	cellsRemaining = remainingRegion[i].stop - remainingRegion[i].start
1707	if cellsRemaining >= self._MinBlockSize[i]:
1708	regionToRead.append(remainingRegion[i])
1709	else:
1710	maxIndex = self.Shape[self.Dimensions.index(self.GetLazyPropertyValue('PhysicalDimensions')[i])]
1711	regionToRead.append(slice(max(0, remainingRegion[i].start - int(math.floor((self._MinBlockSize[i] - cellsRemaining) / 2.))), min(maxIndex, int(math.ceil(remainingRegion[i].stop + (self._MinBlockSize[i] - cellsRemaining) / 2.)))))
1712
1713	self._Cache.insert(0, [regionToRead, self._Grid._ReadNumpyArray(regionToRead)[0]])
1714
1715	bytesAdded = regionToRead[0].stop - regionToRead[0].start
1716	for i in range(1, len(regionToRead)):
1717	bytesAdded *= regionToRead[i].stop - regionToRead[i].start
1718	self._CacheSize += bytesAdded * self._Cache[0][1].dtype.itemsize
1719
1720	self._LogDebug(_(u'MemoryCachedGrid 0x%(id)08X: Added slice [%(slice)s] to the in-memory cache. New cache size = %(size)i bytes.'), {u'id': id(self), u'slice': ','.join(map(lambda s: str(s.start) + ':' + str(s.stop), regionToRead)), u'size': self._CacheSize})
1721
1722	self._FillArrayFromCache(sliceList, remainingRegion, data)
1723
1724	# If the cache is now larger than the maximum allowed
1725	# size, trim it the oldest blocks off the end.
1726
1727	while self._CacheSize > self._MaxCacheSize:
1728	slices, block = self._Cache.pop()
1729
1730	bytesDeleted = slices[0].stop - slices[0].start
1731	for i in range(1, len(slices)):
1732	bytesDeleted *= slices[i].stop - slices[i].start
1733	self._CacheSize -= bytesDeleted * block.dtype.itemsize
1734
1735	self._LogDebug(_(u'MemoryCachedGrid 0x%(id)08X: Deleted slice [%(slice)s] from the in-memory cache. New cache size = %(size)i bytes.'), {u'id': id(self), u'slice': ','.join(map(lambda s: str(s.start) + ':' + str(s.stop), slices)), u'size': self._CacheSize})
1736
1737	# Return successfully.
1738
1739	return data, self.UnscaledNoDataValue
1740
1741	def _FillArrayFromCache(self, sliceList, remainingRegion, data):
1742
1743	# The cache contains zero or more slices, either rectangles
1744	# (if 2D), rectangular parallelepipeds (if 3D), or
1745	# hyperrectangles (if 4D). Ideally we would be able to
1746	# determine, for any collection of slices, whether they
1747	# collectively completely overlap the slice requested by the
1748	# caller and if so, service the request from the slices. If
1749	# there was only partial overlap, we'd "fill in the holes"
1750	# with multiple read requests to the underlying grid.
1751	#
1752	# A completely general algorithm that optimally handles all
1753	# cases appears too complicated to be worth implementing.
1754	# Instead we implement a much simpler algorithm that should
1755	# handle many common cases. This algorithm repeatedly iterates
1756	# through the cached slices, looking for slices that overlap
1757	# the caller's slice such that region that remains
1758	# unoverlapped after each iteration is still rectangular. The
1759	# defining characteristic of candidate slices is that their
1760	# extent exceeds that of the remaining region on all but one
1761	# side, in which case the algorithm continues, or on all
1762	# sides, in which case the caller's slice has been completely
1763	# covered by the cached slices.
1764
1765	i = 0
1766	firstUnusedSlice = 0
1767	while i < len(self._Cache):
1768
1769	# Count how many dimensions of this cached slice fully
1770	# overlap (i.e. enclose) or partially overlap the
1771	# remaining region of this dimension.
1772
1773	dimsWithFullOverlap = []
1774	dimsWithPartialOverlap = []
1775	for j in range(len(self.Dimensions)):
1776	if self._Cache[i][0][j].start <= remainingRegion[j].start and self._Cache[i][0][j].stop >= remainingRegion[j].stop:
1777	dimsWithFullOverlap.append(j)
1778	elif self._Cache[i][0][j].start <= remainingRegion[j].start and self._Cache[i][0][j].stop > remainingRegion[j].start or self._Cache[i][0][j].stop >= remainingRegion[j].stop and self._Cache[i][0][j].start < remainingRegion[j].stop:
1779	dimsWithPartialOverlap.append(j)
1780
1781	#self._LogDebug(_(u'MemoryCachedGrid 0x%(id)08X: i=%(i)i, self._Cache[i][0]=%(cachedSlices)s, dimsWithFullOverlap=%(dimsWithFullOverlap)s, dimsWithPartialOverlap=%(dimsWithPartialOverlap)s'), {u'id': id(self), u'i': i, u'cachedSlices': repr(self._Cache[i][0]), u'dimsWithFullOverlap': repr(dimsWithFullOverlap), u'dimsWithPartialOverlap': repr(dimsWithPartialOverlap)})
1782
1783	# If all of the dimensions of this cached slice fully
1784	# overlap the remaining region, it means the remaining
1785	# region is a subset of (i.e. fully contained by) this
1786	# cached slice.
1787
1788	if len(dimsWithFullOverlap) == len(self.Dimensions):
1789	self._LogDebug(_(u'MemoryCachedGrid 0x%(id)08X: Reading [%(slice1)s] from cached slice [%(slice2)s].'), {u'id': id(self), u'slice1': ','.join(map(lambda s: str(s.start) + ':' + str(s.stop), remainingRegion)), u'slice2': ','.join(map(lambda s: str(s.start) + ':' + str(s.stop), self._Cache[i][0]))})
1790
1791	data.__setitem__(map(lambda s: slice(s[0].start-s[1].start, s[0].stop-s[1].start), zip(remainingRegion, sliceList)), self._Cache[i][1].__getitem__(map(lambda s: slice(s[0].start-s[1].start, s[0].stop-s[1].start), zip(remainingRegion, self._Cache[i][0]))))
1792
1793	self._Cache.insert(0, self._Cache.pop(i))
1794	remainingRegion = None
1795	break
1796
1797	# If all but one of the dimensions of this cached slice
1798	# fully overlap the remaining region and the remaining one
1799	# partially overlaps, it means we can use the slice to
1800	# fill in a block of the rectangular region such that the
1801	# part that remains is a smaller rectangle (2D) /
1802	# rectangular parallelepiped (3D) / hyperrectangle (4D).
1803
1804	if len(dimsWithFullOverlap) == len(self.Dimensions) - 1 and len(dimsWithPartialOverlap) == 1:
1805	cachedRegion = []
1806	for j in range(len(self.Dimensions)):
1807	if j in dimsWithFullOverlap:
1808	cachedRegion.append(remainingRegion[j])
1809	elif self._Cache[i][0][j].start > remainingRegion[j].start:
1810	cachedRegion.append(slice(self._Cache[i][0][j].start, remainingRegion[j].stop))
1811	remainingRegion[j] = slice(remainingRegion[j].start, self._Cache[i][0][j].start)
1812	else:
1813	cachedRegion.append(slice(remainingRegion[j].start, self._Cache[i][0][j].stop))
1814	remainingRegion[j] = slice(self._Cache[i][0][j].stop, remainingRegion[j].stop)
1815
1816	self._LogDebug(_(u'MemoryCachedGrid 0x%(id)08X: Reading [%(slice1)s] from cached slice [%(slice2)s].'), {u'id': id(self), u'slice1': ','.join(map(lambda s: str(s.start) + ':' + str(s.stop), cachedRegion)), u'slice2': ','.join(map(lambda s: str(s.start) + ':' + str(s.stop), self._Cache[i][0]))})
1817
1818	data.__setitem__(map(lambda s: slice(s[0].start-s[1].start, s[0].stop-s[1].start), zip(cachedRegion, sliceList)), self._Cache[i][1].__getitem__(map(lambda s: slice(s[0].start-s[1].start, s[0].stop-s[1].start), zip(cachedRegion, self._Cache[i][0]))))
1819
1820	self._Cache.insert(0, self._Cache.pop(i))
1821	firstUnusedSlice += 1
1822	i = firstUnusedSlice
1823	continue
1824
1825	# Go on to the next slice.
1826
1827	i += 1
1828
1829	# Return the list of slices that define the region not
1830	# fulfilled from the cache, if any.
1831
1832	return remainingRegion
1833
1834	def _WriteNumpyArray(self, sliceList, data):
1835
1836	# If the caller writes anything, invalidate the whole cache.
1837	# This is less than optimal, but we're not really intended to
1838	# be used by callers that need to both write and read.
1839
1840	self._Cache = []
1841	return self._Grid._WriteNumpyArray(sliceList, data)
1842
1843
1844	class AggregateGrid(Grid):
1845	__doc__ = DynamicDocString()
1846
1847	def _GetReportProgress(self):
1848	return self._ReportProgress
1849
1850	def _SetReportProgress(self, value):
1851	# TOOO: Validation
1852	self._ReportProgress = value
1853
1854	ReportProgress = property(_GetReportProgress, _SetReportProgress, doc=DynamicDocString())
1855
1856	def __init__(self, grids, statistic, displayName, reportProgress=False, parentCollection=None, queryableAttributes=None, queryableAttributeValues=None):
1857	# TODO: self.__class__.__doc__.Obj.ValidateMethodInvocation()
1858
1859	# Initialize our properties.
1860
1861	self._Grids = []
1862	self._Grids.extend(grids)
1863	self._Statistic = statistic
1864	self._DisplayName = displayName
1865	self._ReportProgress = reportProgress
1866	self._ReportedStartMessage = False
1867
1868	# Initialize the base class.
1869
1870	super(AggregateGrid, self).__init__(parentCollection=parentCollection, queryableAttributes=queryableAttributes, queryableAttributeValues=queryableAttributeValues)
1871
1872	def _Close(self):
1873	if hasattr(self, '_Grids') and self._Grids is not None:
1874	for grid in self._Grids:
1875	grid.Close()
1876	super(AggregateGrid, self)._Close()
1877
1878	def _GetDisplayName(self):
1879	return self._DisplayName
1880
1881	def _GetLazyPropertyPhysicalValue(self, name):
1882
1883	# If the caller requested PhysicalDimensions or
1884	# PhysicalDimensionsFlipped, return idealized values, as the
1885	# transposing and flipping is handled by the contained grid.
1886
1887	if name == 'PhysicalDimensions':
1888	return self._Grids[0].Dimensions
1889
1890	if name == 'PhysicalDimensionsFlipped':
1891	return tuple([False] * len(self._Grids[0].Dimensions))
1892
1893	# If the caller requested the UnscaledDataType or
1894	# UnscaledNoDataValue, calculate and return them.
1895
1896	if name == 'UnscaledDataType':
1897
1898	# For COUNT, return the smallest integer type that is
1899	# guaranteed hold the largest possible value.
1900
1901	if self._Statistic == u'count':
1902	if len(self._Grids) < 256:
1903	return u'uint8'
1904	if len(self._Grids) < 65536:
1905	return u'uint16'
1906	return u'int32'
1907
1908	# For MIN, MAX, and RANGE, return the data type of the
1909	# first grid. This is guaranteed to work if all of the
1910	# grids have this type. If they do not, it may still work
1911	# but there is the possibility of overflow or underflow.
1912	# If it happens, we will detect it and fail.
1913
1914	if self._Statistic in [u'minimum', u'maximum', u'range']:
1915	return self._Grids[0].DataType
1916
1917	# For all others, return float32, unless the first grid is
1918	# float64, in which case return float64.
1919
1920	if self._Grids[0].DataType == u'float64':
1921	return u'float64'
1922	return u'float32'
1923
1924	if name == 'UnscaledNoDataValue':
1925
1926	# For COUNT, return 0.
1927
1928	if self._Statistic == u'count':
1929	if self.DataType in [u'float32', u'float64']:
1930	return 0.
1931	return 0
1932
1933	# For MIN, MAX, and RANGE, return the NoData value of the
1934	# first grid.
1935
1936	if self._Statistic in [u'minimum', u'maximum', u'range']:
1937	return self._Grids[0].NoDataValue
1938
1939	# For all others, return the smallest representable
1940	# floating point number:
1941	#
1942	# >>> float(numpy.finfo('float32').min)
1943	# -3.4028234663852886e+038
1944	# >>> float(numpy.finfo('float64').min)
1945	# -1.7976931348623157e+308
1946
1947	import numpy
1948
1949	if self.DataType == u'float32':
1950	return float(numpy.finfo('float32').min)
1951	return float(numpy.finfo('float64').min)
1952
1953	# If the caller requested one of the properties related to
1954	# scaling, return None.
1955
1956	if name in ['ScaledDataType', 'ScaledNoDataValue', 'ScalingFunction', 'UnscalingFunction']:
1957	return None
1958
1959	# Otherwise use the value of the property from the first grid.
1960
1961	return self._Grids[0].GetLazyPropertyValue(name)
1962
1963	@classmethod
1964	def _TestCapability(cls, capability):
1965	return cls._Grids[0]._TestCapability(capability) # TODO: Reject writing
1966
1967	def _GetCoords(self, coord, coordNum, slices, sliceDims, fixedIncrementOffset):
1968	return self._Grids[0]._GetCoords(coord, coordNum, slices, sliceDims, fixedIncrementOffset)
1969
1970	def _ReadNumpyArray(self, sliceList):
1971	import numpy
1972
1973	# Initialize arrays needed to perform the aggregation.
1974
1975	shape = [s.stop - s.start for s in sliceList]
1976	count = numpy.zeros(shape, dtype='int32')
1977
1978	if self._Statistic in [u'minimum', u'range']:
1979	smallest = numpy.zeros(shape, dtype=str(self.DataType))
1980	if self.DataType[0] == 'f':
1981	smallest += numpy.finfo(str(self.DataType)).max
1982	else:
1983	smallest += numpy.iinfo(str(self.DataType)).max
1984
1985	if self._Statistic in [u'maximum', u'range']:
1986	largest = numpy.zeros(shape, dtype=str(self.DataType))
1987	if self.DataType[0] == 'f':
1988	largest += numpy.finfo(str(self.DataType)).min
1989	else:
1990	largest += numpy.iinfo(str(self.DataType)).min
1991
1992	if self._Statistic in [u'sum', u'mean']:
1993	sum_ = numpy.zeros(shape, dtype=str(self.DataType))
1994
1995	if self._Statistic == u'standard deviation':
1996	mean = numpy.zeros(shape, dtype=str(self.DataType))
1997	M2 = numpy.zeros(shape, dtype=str(self.DataType))
1998
1999	# If we were asked to report progress and the caller is
2000	# requesting the entire grid, create a progress reporter. If
2001	# they're only requesting part of the grid, do not create a
2002	# progress reporter, even if they asked for one, because it is
2003	# highly likely that they'll just ask for another part of it
2004	# in a moment, which would cause us to report progress again,
2005	# and again...
2006
2007	if self._ReportProgress:
2008	if not self._ReportedStartMessage:
2009	self._LogInfo(_(u'Aggregating %(grids)i grids to create the %(dn)s.') % {u'grids': len(self._Grids), u'dn': self._DisplayName})
2010
2011	if tuple(shape) == self.Shape:
2012	from GeoEco.Logging import ProgressReporter
2013	progressReporter = ProgressReporter(progressMessage1=_(u'Still aggregating: %(elapsed)s elapsed, %(opsCompleted)i grids aggregated, %(perOp)s per grid, %(opsRemaining)i remaining, estimated completion time: %(etc)s.'),
2014	completionMessage=_(u'Aggregation complete: %(elapsed)s elapsed, %(opsCompleted)i grids aggregated, %(perOp)s per grid.'),
2015	abortedMessage=_(u'Aggregation stopped before all grids were aggregated: %(elapsed)s elapsed, %(opsCompleted)i grids aggregated, %(perOp)s per grid, %(opsIncomplete)i not aggregated.'),
2016	loggingChannel=u'GeoEco.Datasets',
2017	arcGISProgressorLabel=_(u'Aggregating %(count)i grids') % {u'count': len(self._Grids)} )
2018	progressReporter.Start(len(self._Grids))
2019	else:
2020	self._ReportedStartMessage = True
2021	progressReporter = None
2022	else:
2023	progressReporter = None
2024
2025	# Iterate through the grids, aggregating their values into the
2026	# appropriate arrays we allocated.
2027
2028	try:
2029	for i in range(len(self._Grids)):
2030
2031	# First validate that this grid has the same dimensions
2032	# and shape as the first grid and a compatible data type.
2033
2034	if self._Grids[i].Dimensions != self._Grids[0].Dimensions:
2035	raise ValueError(_(u'Cannot aggregate %(dn1)s and %(dn2)s because they have different dimensions. The first grid has dimensions %(dim1)s and the second has dimensions %(dim2)s.' ) % {u'dn1': self._Grids[0].DisplayName, u'dn2': self._Grids[1].DisplayName, u'dim1': self._Grids[0].Dimensions, u'dim2': self._Grids[1].Dimensions})
2036
2037	if self._Grids[i].Shape != self._Grids[0].Shape:
2038	raise ValueError(_(u'Cannot aggregate %(dn1)s and %(dn2)s because they have different shapes. The first grid has the shape %(shape1)s and the second has the shape %(shape2)s.' ) % {u'dn1': self._Grids[0].DisplayName, u'dn2': self._Grids[1].DisplayName, u'shape1': self._Grids[0].Shape, u'shape2': self._Grids[1].Shape})
2039
2040	if self._Statistic in [u'minimum', u'maximum', u'range'] and not numpy.can_cast(str(self._Grids[i].DataType), str(self._Grids[0].DataType)):
2041	raise ValueError(_(u'Cannot aggregate %(dn1)s and %(dn2)s because the second grid\'s data type, %(dt2)s, cannot be fully represented by the first grid\'s data type, %(dt1)s. Recreate the grids with the same data type and try again.' ) % {u'dn1': self._Grids[0].DisplayName, u'dn2': self._Grids[1].DisplayName, u'dt1': self._Grids[i].DataType, u'dt2': self._Grids[i].DataType})
2042
2043	# Now aggregate this grid.
2044
2045	data = self._Grids[i].Data.__getitem__(tuple(sliceList))
2046
2047	if self._Grids[i].NoDataValue is not None:
2048	hasData = data != self._Grids[i].NoDataValue
2049	else:
2050	hasData = numpy.ones(data.shape, dtype=bool)
2051
2052	count += numpy.cast['int32'](hasData)
2053
2054	if self._Statistic in [u'minimum', u'range']:
2055	isSmaller = numpy.logical_and(hasData, data < smallest)
2056	smallest[isSmaller] = data[isSmaller]
2057
2058	if self._Statistic in [u'maximum', u'range']:
2059	isLarger = numpy.logical_and(hasData, data > largest)
2060	largest[isLarger] = data[isLarger]
2061
2062	if self._Statistic in [u'sum', u'mean']:
2063	sum_[hasData] += data[hasData]
2064
2065	# To avoid having to make two passes through the grids to
2066	# calculate the standard deviation, we use an online
2067	# algorithm. See "On-line algorithm" in
2068	# http://en.wikipedia.org/wiki/Algorithms_for_calculating_variance.
2069
2070	if self._Statistic == u'standard deviation':
2071	delta = data[hasData] - mean[hasData]
2072	mean[hasData] += delta/count[hasData]
2073	M2[hasData] += delta*(data[hasData] - mean[hasData])
2074
2075	# Report progress.
2076
2077	if progressReporter is not None:
2078	progressReporter.ReportProgress()
2079
2080	except:
2081	if progressReporter is not None:
2082	progressReporter.Stop()
2083	raise
2084
2085	# Calculate and return the result.
2086
2087	if self._Statistic == u'count':
2088	count[count == 0] = self.NoDataValue
2089	return count, self.NoDataValue
2090
2091	if self._Statistic == u'minimum':
2092	if self.NoDataValue is not None:
2093	smallest[count == 0] = self.NoDataValue
2094	return smallest, self.NoDataValue
2095
2096	if self._Statistic == u'maximum':
2097	if self.NoDataValue is not None:
2098	largest[count == 0] = self.NoDataValue
2099	return largest, self.NoDataValue
2100
2101	if self._Statistic == u'range':
2102	result = largest - smallest
2103	if self.NoDataValue is not None:
2104	result[count == 0] = self.NoDataValue
2105	return result, self.NoDataValue
2106
2107	if self._Statistic == u'sum':
2108	sum_[count == 0] = self.NoDataValue
2109	return sum_, self.NoDataValue
2110
2111	if self._Statistic == u'mean':
2112	mean = numpy.zeros(shape, dtype=str(self.DataType))
2113	mean += self.NoDataValue
2114	hasData = count >= 1
2115	mean[hasData] = sum_[hasData] / count[hasData]
2116	return mean, self.NoDataValue
2117
2118	stdev = numpy.zeros(shape, dtype=str(self.DataType))
2119	stdev += self.NoDataValue
2120	hasData = count >= 2
2121	stdev[hasData] = (M2[hasData]/(count[hasData] - 1))**0.5
2122	return stdev, self.NoDataValue
2123
2124
2125	class ClimatologicalGridCollection(DatasetCollection):
2126	__doc__ = DynamicDocString()
2127
2128	def _GetStatistic(self):
2129	return self._Statistic
2130
2131	Statistic = property(_GetStatistic, doc=DynamicDocString())
2132
2133	def _GetBinType(self):
2134	return self._BinType
2135
2136	BinType = property(_GetBinType, doc=DynamicDocString())
2137
2138	def _GetBinDuration(self):
2139	return self._BinDuration
2140
2141	BinDuration = property(_GetBinDuration, doc=DynamicDocString())
2142
2143	def _GetStartDayOfYear(self):
2144	return self._StartDayOfYear
2145
2146	StartDayOfYear = property(_GetStartDayOfYear, doc=DynamicDocString())
2147
2148	def _GetReportProgress(self):
2149	return self._ReportProgress
2150
2151	def _SetReportProgress(self, value):
2152	# TOOO: Validation
2153	self._ReportProgress = value
2154
2155	ReportProgress = property(_GetReportProgress, _SetReportProgress, doc=DynamicDocString())
2156
2157	def __init__(self, grid, statistic, binType, binDuration, startDayOfYear=1, reportProgress=False):
2158	# TODO: self.__class__.__doc__.Obj.ValidateMethodInvocation()
2159
2160	# Validate that the grid has a t dimension and that it can be
2161	# used to produce the desired bin type.
2162
2163	if u't' not in grid.Dimensions:
2164	raise ValueError(_(u'%(dn)s does not have a time dimension. In order to create a climatology, it must have a time dimension.') % {u'dn': grid.DisplayName})
2165
2166	# If the binDuration is monthly, validate that the
2167	# startDayOfYear is not the 28th, 29th, 30th, or 31st day of
2168	# the start month. To keep our calculations simple, we do not
2169	# support those days as starting days.
2170
2171	startDateTime = datetime.datetime(2001, 1, 1) + datetime.timedelta(startDayOfYear - 1) # Use a non-leap year for this calculation
2172	if startDateTime.day > 28:
2173	raise NotImplementedError(_(u'The Start Day Of Year parameter specifies that the first climatological grid should start on %(day)i %(month)s. This is not supported. The grid may not start later than day 28 of any month. We apologize for the inconvenience.') % {u'day': startDateTime.day, u'month': self._Unicode(startDateTime.strftime('%B'))})
2174
2175	# If the grid has a time increment that is larger than the bin
2176	# length, warn the user that this may produce undesirable
2177	# output.
2178
2179	if binType != u'cumulative' and grid.CoordIncrements[0] is not None and grid.TIncrementUnit in ['second', 'minute', 'hour', 'day', 'month', 'year']:
2180	if grid.TIncrementUnit == u'second':
2181	gridTIncrement = grid.CoordIncrements[0]
2182	elif grid.TIncrementUnit == u'minute':
2183	gridTIncrement = grid.CoordIncrements[0] * 60
2184	elif grid.TIncrementUnit == u'hour':
2185	gridTIncrement = grid.CoordIncrements[0] * 60 * 60
2186	elif grid.TIncrementUnit == u'day':
2187	gridTIncrement = grid.CoordIncrements[0] * 60 * 60 * 24
2188	elif grid.TIncrementUnit == u'month':
2189	gridTIncrement = grid.CoordIncrements[0] * 60 * 60 * 24 * 30
2190	else:
2191	gridTIncrement = grid.CoordIncrements[0] * 60 * 60 * 24 * 365
2192
2193	if binType == u'daily':
2194	binLength = binDuration * 60 * 60 * 24
2195	else:
2196	binLength = binDuration * 60 * 60 * 24 * 30 # Monthly
2197
2198	if gridTIncrement > binLength:
2199	if binType == u'monthly':
2200	if binDuration == 1:
2201	binTypeName = u'month'
2202	else:
2203	binTypeName = u'months'
2204	else:
2205	if binDuration == 1:
2206	binTypeName = u'day'
2207	else:
2208	binTypeName = u'days'
2209	self._LogWarning(_(u'The time increment of %(dn)s is longer than the requested climatology bin length (%(gridTIncrement)g %(gridTIncrementUnit)s(s) vs. %(binDuration)i %(binTypeName)s). This may produce unsatisfactory results. The time slices are placed into bins according to their center time coordinate, and are not split between bins. This can lead to problems if the bin is short compared to the time slice. For example, if the bin is 1 day long while the time slice is 1 month long, the entire time slice will be placed in the bin for the day that occurs in the middle of that month. To avoid this kind of problem, increase the bin duration so that the bins are longer than the time slices (preferably a lot longer).') % {u'dn': grid.DisplayName, u'gridTIncrement': grid.CoordIncrements[0], u'gridTIncrementUnit': grid.TIncrementUnit, u'binDuration': binDuration, u'binTypeName': binTypeName})
2210
2211	# Initialize our properties.
2212
2213	self._Grid = grid
2214	self._Statistic = statistic
2215	self._BinType = binType
2216	self._BinDuration = binDuration
2217	self._StartDayOfYear = startDayOfYear
2218	self._ReportProgress = reportProgress
2219	self._AggregateGrids = None
2220
2221	if self._BinType == u'cumulative':
2222	self._DisplayName = _(u'cumulative climatological %(statistic)s of the %(dn)s') % {u'dn': grid.DisplayName, u'statistic': self._Statistic}
2223	elif self._BinType == u'monthly':
2224	if self._BinDuration == 1:
2225	self._DisplayName = _(u'monthly climatological %(statistic)s of the %(dn)s') % {u'dn': grid.DisplayName, u'statistic': self._Statistic}
2226	else:
2227	self._DisplayName = _(u'%(binDuration)i-month climatological %(statistic)s of the %(dn)s') % {u'dn': grid.DisplayName, u'binDuration': self._BinDuration, u'statistic': self._Statistic}
2228	else:
2229	if self._BinDuration == 1:
2230	self._DisplayName = _(u'daily climatological %(statistic)s of the %(dn)s') % {u'dn': grid.DisplayName, u'statistic': self._Statistic}
2231	else:
2232	self._DisplayName = _(u'%(binDuration)i-day climatological %(statistic)s of the %(dn)s') % {u'dn': grid.DisplayName, u'binDuration': self._BinDuration, u'statistic': self._Statistic}
2233
2234	# Create the queryable attributes appropriate for the bin
2235	# type.
2236
2237	qa = [QueryableAttribute(u'ClimatologyBinType', u'Climatology bin type', UnicodeStringTypeMetadata()),
2238	QueryableAttribute(u'ClimatologyBinName', u'Climatology bin name', UnicodeStringTypeMetadata()),
2239	QueryableAttribute(u'Statistic', u'Statistic', UnicodeStringTypeMetadata())]
2240
2241	qav = {u'Statistic': self._Statistic.lower()}
2242
2243	if binType == u'cumulative':
2244	qav[u'ClimatologyBinType'] = u'Cumulative'
2245
2246	if binType == u'monthly':
2247	qa.append(QueryableAttribute(u'FirstMonth', u'First month of this bin', IntegerTypeMetadata(minValue=1, maxValue=12)))
2248	qa.append(QueryableAttribute(u'DayOfFirstMonth', u'Day of the first month', IntegerTypeMetadata(minValue=1, maxValue=31)))
2249	qa.append(QueryableAttribute(u'LastMonth', u'Last month of this bin', IntegerTypeMetadata(minValue=1, maxValue=12)))
2250	qa.append(QueryableAttribute(u'DayOfLastMonth', u'Day of the last month', IntegerTypeMetadata(minValue=1, maxValue=31)))
2251	if self._BinDuration == 1:
2252	qav[u'ClimatologyBinType'] = u'Monthly'
2253	else:
2254	qav[u'ClimatologyBinType'] = u'%imonth' % self._BinDuration
2255
2256	elif binType == u'daily':
2257	qa.append(QueryableAttribute(u'FirstDay', u'First day of the year of this bin', IntegerTypeMetadata(minValue=1, maxValue=366)))
2258	qa.append(QueryableAttribute(u'LastDay', u'Last day of the year of this bin', IntegerTypeMetadata(minValue=1, maxValue=366)))
2259	if self._BinDuration == 1:
2260	qav[u'ClimatologyBinType'] = u'Daily'
2261	else:
2262	qav[u'ClimatologyBinType'] = u'%iday' % self._BinDuration
2263
2264	# Copy the queryable attributes for the grid.
2265
2266	gridQA = grid.GetAllQueryableAttributes()
2267	if gridQA is not None and len(gridQA) > 0:
2268	for a in gridQA:
2269	if grid.ParentCollection is None or grid.ParentCollection.GetQueryableAttribute(a.Name) is None:
2270	qa.append(a)
2271	qav[a.Name] = grid.GetQueryableAttributeValue(a.Name)
2272
2273	# Initialize the base class.
2274
2275	super(ClimatologicalGridCollection, self).__init__(parentCollection=grid.ParentCollection, queryableAttributes=tuple(qa), queryableAttributeValues=qav)
2276
2277	def _Close(self):
2278	if hasattr(self, '_Grid') and self._Grid is not None:
2279	self._Grid.Close()
2280	super(ClimatologicalGridCollection, self)._Close()
2281
2282	def _GetDisplayName(self):
2283	return self._DisplayName
2284
2285	def _QueryDatasets(self, parsedExpression, progressReporter, options, parentAttrValues):
2286
2287	# If we have not created the aggregate grids, do it now.
2288
2289	if self._AggregateGrids is None:
2290	self._AggregateGrids = []
2291
2292	# If the bin type is 'cumulative', create a single
2293	# AggregateGrid for all of the time slices.
2294
2295	if self._BinType == u'cumulative':
2296	self._AggregateGrids.append(AggregateGrid(GridSliceCollection(self._Grid, zQAName=None).QueryDatasets(reportProgress=False), self._Statistic, self._DisplayName, self._ReportProgress, parentCollection=self, queryableAttributeValues={u'ClimatologyBinName': 'cumulative'}))
2297
2298	# Otherwise, if the bin type is 'monthly':
2299
2300	elif self._BinType == u'monthly':
2301
2302	# Produce a list of [startMonthAndDay, endMonthAndDay]
2303	# records, one for each bin, sorted by endDay in
2304	# ascending order.
2305
2306	numBins = int(round(12. / self._BinDuration))
2307	binMonthsAndDays = []
2308	startDateTime = datetime.datetime(2001, 1, 1) + datetime.timedelta(self._StartDayOfYear - 1) # Use a non-leap year for this calculation.
2309	startMonth = startDateTime.month
2310
2311	for i in range(numBins):
2312	if i < numBins - 1:
2313	if startMonth + self._BinDuration <= 12:
2314	endMonthAndDay = int((datetime.datetime(2000, startMonth + self._BinDuration, startDateTime.day) - datetime.timedelta(1)).strftime('%m%d')) # Use a leap year, to ensure that the last day of February is 29
2315	else:
2316	endMonthAndDay = int((datetime.datetime(2000, startMonth + self._BinDuration - 12, startDateTime.day) - datetime.timedelta(1)).strftime('%m%d')) # Use a leap year, to ensure that the last day of February is 29
2317	else:
2318	endMonthAndDay = int((startDateTime - datetime.timedelta(1)).strftime('%m%d'))
2319
2320	binMonthsAndDays.append([startMonth*100 + startDateTime.day, endMonthAndDay])
2321
2322	startMonth += self._BinDuration
2323	if startMonth > 12:
2324	startMonth -= 12
2325
2326	binMonthsAndDays.sort(cmp=lambda a, b: cmp(a[1], b[1]))
2327
2328	self._LogDebug(_(u'ClimatologicalGridCollection 0x%(id)08X: Using bin [startMonthAndDay, endMonthAndDay] definitions: %(binMonthsAndDays)s.'), {u'id': id(self), u'binMonthsAndDays': str(binMonthsAndDays)})
2329
2330	# Create GridSlice instances for each time slice of
2331	# our grid and sort them into the appropriate bins.
2332
2333	binEndMonthsAndDays = [x[1] for x in binMonthsAndDays]
2334	bins = [[] for i in range(numBins)]
2335
2336	for i, t in enumerate(self._Grid.CenterCoords['t']):
2337	bin = bisect.bisect_left(binEndMonthsAndDays, t.month*100 + t.day)
2338	if bin >= numBins:
2339	bin = 0
2340	bins[bin].append(GridSlice(self._Grid, tIndex=i, tQACoordType=u'center'))
2341
2342	# Create an AggregateGrid for each non-empty bin.
2343
2344	for bin, grids in enumerate(bins):
2345	if self._BinDuration == 1:
2346	queryableAttributeValues = {u'ClimatologyBinName': u'month%02i' % (binMonthsAndDays[bin][0] / 100)}
2347	displayName = _(u'%(month)s climatological %(statistic)s of the %(dn)s') % {u'dn': self._Grid.DisplayName, u'month': self._Unicode(datetime.datetime(2001, binMonthsAndDays[bin][0]/100, 1).strftime('%B')), u'statistic': self._Statistic}
2348	else:
2349	queryableAttributeValues = {u'ClimatologyBinName': u'months%02ito%02i' % (binMonthsAndDays[bin][0] / 100, binMonthsAndDays[bin][1] / 100)}
2350	displayName = _(u'months %(firstMonth)i-%(lastMonth)i climatological %(statistic)s of the %(dn)s') % {u'dn': self._Grid.DisplayName, u'firstMonth': binMonthsAndDays[bin][0] / 100, u'lastMonth': binMonthsAndDays[bin][1] / 100, u'statistic': self._Statistic}
2351
2352	if len(grids) > 0:
2353	queryableAttributeValues[u'FirstMonth'] = binMonthsAndDays[bin][0] / 100
2354	queryableAttributeValues[u'DayOfFirstMonth'] = binMonthsAndDays[bin][0] % 100
2355	queryableAttributeValues[u'LastMonth'] = binMonthsAndDays[bin][1] / 100
2356	queryableAttributeValues[u'DayOfLastMonth'] = binMonthsAndDays[bin][1] % 100
2357
2358	self._AggregateGrids.append(AggregateGrid(grids, self._Statistic, displayName, self._ReportProgress, parentCollection=self, queryableAttributeValues=queryableAttributeValues))
2359	else:
2360	self._LogWarning(_(u'No data is available to create the %(dn)s.') % {u'dn': displayName})
2361
2362	# Otherwise the bin type is 'daily':
2363
2364	else:
2365
2366	# Produce a list of [startDay, endDay] records, one
2367	# for each bin, sorted by endDay in ascending order.
2368
2369	numBins = int(round(365. / self._BinDuration))
2370	binDays = []
2371	startDay = self._StartDayOfYear
2372
2373	for i in range(numBins):
2374	if i < numBins - 1:
2375	endDay = startDay + self._BinDuration - 1
2376	else:
2377	endDay = self._StartDayOfYear - 1
2378
2379	if endDay > 365:
2380	endDay -= 365
2381	elif endDay == 0 or endDay == 365:
2382	endDay = 366
2383
2384	binDays.append([startDay, endDay])
2385
2386	startDay += self._BinDuration
2387	if startDay >= 366:
2388	startDay -= 365
2389
2390	binDays.sort(cmp=lambda a, b: cmp(a[1], b[1]))
2391
2392	self._LogDebug(_(u'ClimatologicalGridCollection 0x%(id)08X: Using bin [startDay, endDay] definitions: %(binDays)s.'), {u'id': id(self), u'binDays': str(binDays)})
2393
2394	# Create GridSlice instances for each time slice of
2395	# our grid and sort them into the appropriate bins.
2396
2397	binEndDays = [x[1] for x in binDays]
2398	bins = [[] for i in range(numBins)]
2399
2400	for i, t in enumerate(self._Grid.CenterCoords['t']):
2401	bin = bisect.bisect_left(binEndDays, t.timetuple()[7])
2402	if bin >= numBins:
2403	bin = 0
2404	bins[bin].append(GridSlice(self._Grid, tIndex=i, tQACoordType=u'center'))
2405
2406	# Create an AggregateGrid for each bin.
2407
2408	for bin, grids in enumerate(bins):
2409	if self._BinDuration == 1:
2410	queryableAttributeValues = {u'ClimatologyBinName': u'day%03i' % binDays[bin][0]}
2411	displayName = _(u'day %(day)03i climatological %(statistic)s of the %(dn)s') % {u'dn': self._Grid.DisplayName, u'day': binDays[bin][0], u'statistic': self._Statistic}
2412	else:
2413	queryableAttributeValues = {u'ClimatologyBinName': u'days%03ito%03i' % (binDays[bin][0], binDays[bin][1])}
2414	displayName = _(u'days %(firstDay)03i-%(lastDay)03i climatological %(statistic)s of the %(dn)s') % {u'dn': self._Grid.DisplayName, u'firstDay': binDays[bin][0], u'lastDay': binDays[bin][1], u'statistic': self._Statistic}
2415
2416	if len(grids) > 0:
2417	queryableAttributeValues[u'FirstDay'] = binDays[bin][0]
2418	queryableAttributeValues[u'LastDay'] = binDays[bin][1]
2419
2420	self._AggregateGrids.append(AggregateGrid(grids, self._Statistic, displayName, self._ReportProgress, parentCollection=self, queryableAttributeValues=queryableAttributeValues))
2421	else:
2422	self._LogWarning(_(u'No data is available to create the %(dn)s.') % {u'dn': displayName})
2423
2424	# Execute the parsed expression against each of the aggregate
2425	# grids, returning all of those that match.
2426
2427	results = []
2428
2429	if parsedExpression is None:
2430	results.extend(self._AggregateGrids)
2431	if progressReporter is not None:
2432	progressReporter.ReportProgress(len(self._AggregateGrids))
2433
2434	else:
2435	for aggregatedGrid in self._AggregateGrids:
2436	allAttrValues = {}
2437	for attr in aggregatedGrid.GetAllQueryableAttributes():
2438	allAttrValues[attr] = aggregatedGrid.GetQueryableAttributeValue(attr)
2439
2440	if parsedExpression.eval(allAttrValues):
2441	results.append(aggregatedGrid)
2442	if progressReporter is not None:
2443	progressReporter.ReportProgress()
2444
2445	return results
2446
2447
2448	class SeafloorGrid(Grid):
2449	__doc__ = DynamicDocString()
2450
2451	def __init__(self, grid, queryableAttributes=None, queryableAttributeValues=None):
2452	# TODO: self.__class__.__doc__.Obj.ValidateMethodInvocation()
2453
2454	# Validate that the grid has z coordinates and that those
2455	# coordinates do not depend on any others.
2456
2457	if u'z' not in grid.Dimensions:
2458	raise ValueError(_(u'Cannot create a grid representing the values of the %(dn)s at the seafloor because that grid does not have a depth coordinate.') % {u'dn': grid.DisplayName})
2459
2460	if grid.CoordDependencies[grid.Dimensions.index(u'z')] is not None:
2461	raise ValueError(_(u'Cannot create a grid representing the values of the %(dn)s at the seafloor because the depth cooordinates of that grid depend on the values of other coordinates.') % {u'dn': grid.DisplayName})
2462
2463	# Validate that the grid has an UnscaledNoDataValue. If there
2464	# is no UnscaledNoDataValue, it means we can't tell where
2465	# there is data and where there isn't.
2466
2467	if grid.UnscaledNoDataValue is None:
2468	raise ValueError(_(u'Cannot create a grid representing the values of the %(dn)s at the seafloor because that grid does not have an UnscaledNoDataValue.') % {u'dn': grid.DisplayName})
2469
2470	# Initialize our properties.
2471
2472	self._Grid = grid
2473	self._DisplayName = _(u'seafloor values of the %(dn)s') % {u'dn': grid.DisplayName}
2474
2475	# For our queryable attributes, copy all of those of the grid
2476	# plus those provided by the caller.
2477
2478	qa = []
2479	if grid._QueryableAttributes is not None:
2480	qa.extend(list(grid._QueryableAttributes))
2481	if queryableAttributes is not None:
2482	qa.extend(queryableAttributes)
2483
2484	qav = {}
2485	if grid._QueryableAttributeValues is not None:
2486	qav.update(grid._QueryableAttributeValues)
2487	if queryableAttributeValues is not None:
2488	qav.update(queryableAttributeValues)
2489
2490	# Initialize the base class.
2491
2492	super(SeafloorGrid, self).__init__(parentCollection=grid.ParentCollection, queryableAttributes=tuple(qa), queryableAttributeValues=qav)
2493
2494	def _Close(self):
2495	if hasattr(self, '_Grid') and self._Grid is not None:
2496	self._Grid.Close()
2497	super(SeafloorGrid, self)._Close()
2498
2499	def _GetDisplayName(self):
2500	return self._DisplayName
2501
2502	def _GetLazyPropertyPhysicalValue(self, name):
2503
2504	# If the caller is requesting one of the sequences related to
2505	# the dimensions, get the value from the contained grid but
2506	# remove the z dimension.
2507
2508	if name == 'Dimensions':
2509	if self._Grid.Dimensions == u'tzyx': # It is either tzyx or zyx.
2510	return u'tyx'
2511	return u'yx'
2512
2513	if name == 'Shape':
2514	if self._Grid.Dimensions == u'tzyx':
2515	return (self._Grid.Shape[0], self._Grid.Shape[2], self._Grid.Shape[3])
2516	return (self._Grid.Shape[1], self._Grid.Shape[2])
2517
2518	if name == 'CoordDependencies':
2519	if self._Grid.Dimensions == u'tzyx':
2520	return (self._Grid.CoordDependencies[0], self._Grid.CoordDependencies[2], self._Grid.CoordDependencies[3])
2521	return (self._Grid.CoordDependencies[1], self._Grid.CoordDependencies[2])
2522
2523	if name == 'CoordIncrements':
2524	if self._Grid.Dimensions == u'tzyx':
2525	return (self._Grid.CoordIncrements[0], self._Grid.CoordIncrements[2], self._Grid.CoordIncrements[3])
2526	return (self._Grid.CoordIncrements[1], self._Grid.CoordIncrements[2])
2527
2528	if name == 'PhysicalDimensions':
2529	physicalDimensions = u''
2530	for d in self._Grid.GetLazyPropertyValue('PhysicalDimensions'):
2531	if d != u'z':
2532	physicalDimensions += d
2533	return physicalDimensions
2534
2535	if name == 'PhysicalDimensionsFlipped':
2536	physicalDimensionsFlipped = []
2537	for i, d in enumerate(self._Grid.GetLazyPropertyValue('PhysicalDimensions')):
2538	if d != u'z':
2539	physicalDimensionsFlipped.append(self._Grid.GetLazyPropertyValue('PhysicalDimensionsFlipped')[i])
2540	return tuple(physicalDimensionsFlipped)
2541
2542	# Otherwise get the value from the contained grid.
2543
2544	return self._Grid.GetLazyPropertyValue(name)
2545
2546	@classmethod
2547	def _TestCapability(cls, capability):
2548	return cls._Grid._TestCapability(capability) # TODO: Disallow writing
2549
2550	@classmethod
2551	def _GetSRTypeForSetting(cls):
2552	return cls._Grid._GetSRTypeForSetting()
2553
2554	def _SetSpatialReference(self, sr):
2555	return self._Grid._SetSpatialReference(sr)
2556
2557	def _GetCoords(self, coord, coordNum, slices, sliceDims, fixedIncrementOffset):
2558	return self._Grid._GetCoords(coord, coordNum, slices, sliceDims, fixedIncrementOffset)
2559
2560	def _ReadNumpyArray(self, sliceList):
2561
2562	# Initialize a numpy array to return with the
2563	# UnscaledNoDataValue.
2564
2565	import numpy
2566
2567	shape = [s.stop - s.start for s in sliceList]
2568	result = numpy.zeros(shape, dtype=str(self.UnscaledDataType)) + self.UnscaledNoDataValue
2569	hasNoData = result == self.UnscaledNoDataValue
2570
2571	# Iterate through the depth layers, from deepest to
2572	# shallowest, filling cells of the result array that still
2573	# have no data with values from the depth layer.
2574
2575	zIndex = self._Grid.GetLazyPropertyValue('PhysicalDimensions').index('z')
2576	if self._Grid.GetLazyPropertyValue('PhysicalDimensionsFlipped')[zIndex]:
2577	zLayers = range(self._Grid.Shape[-3])
2578	else:
2579	zLayers = range(self._Grid.Shape[-3] - 1, -1, -1)
2580
2581	for z in zLayers:
2582	sliceToRequest = list() + sliceList
2583	sliceToRequest.insert(zIndex, slice(z, z+1))
2584	data = self._Grid._ReadNumpyArray(sliceToRequest)[0]
2585
2586	sliceToExtract = [slice(None) for i in range(len(sliceList))]
2587	sliceToExtract.insert(zIndex, 0)
2588	data = data.__getitem__(tuple(sliceToExtract))
2589
2590	result[hasNoData] = data[hasNoData]
2591
2592	hasNoData = result == self.UnscaledNoDataValue
2593	if numpy.logical_not(hasNoData).all():
2594	break
2595
2596	# Return the result.
2597
2598	return result, self.UnscaledNoDataValue
2599
2600
2601	class SpatiotemporalAggregateGrid(Grid):
2602	__doc__ = DynamicDocString()
2603
2604	def __init__(self, grid, statistic, xyWindowSize=None, xyIncrement=None, zWindowSize=None, zIncrement=None, tWindowSize=None, tIncrement=None, tUnit=None, tWindowPosition=None):
2605	# TODO: self.__class__.__doc__.Obj.ValidateMethodInvocation()
2606
2607	# TODO: Perform additional validation.
2608
2609	# Fail if the caller's grid has coordinate dependencies for
2610	# any dimension. We only support grids for which there are no
2611	# coordinate dependenices, i.e. the value of coordinate c at
2612	# index i is constant for all values of the other coordinates.
2613	#
2614	# Most datasets adhere to this scheme. One that does not is
2615	# ROMS. It uses a terrain-following system where the z
2616	# coordinate varies with x, y, and t. It is not possible to
2617	# use the SpatiotemporalAggregateGrid with ROMS grids unless
2618	# they have been interpolated to constant z levels.
2619
2620	for i in range(len(grid.Dimensions)):
2621	if grid.CoordDependencies[i] is not None:
2622	raise ValueError(_(u'A spatiotemoral aggregate grid cannot be constructed from %(dn)s, because the values of the %(d)s coordinates of that dataset depend on the values of other coordinates.') % {u'dn': grid.DisplayName, u'd': grid.Dimensions[i]})
2623
2624	# Initialize our properties.
2625
2626	self._Grid = grid
2627	self._Statistic = statistic
2628	self._XYWindowSize = xyWindowSize
2629	self._XYIncrement = xyIncrement
2630
2631	if 'z' in grid.Dimensions:
2632	self._ZWindowSize = zWindowSize
2633	self._ZIncrement = zIncrement
2634	else:
2635	self._ZWindowSize = None
2636	self._ZIncrement = None
2637
2638	if 't' in grid.Dimensions:
2639	self._TWindowSize = tWindowSize
2640	self._TIncrement = tIncrement
2641	self._TUnit = tUnit
2642	self._TWindowPosition = tWindowPosition
2643	else:
2644	self._TWindowSize = None
2645	self._TIncrement = None
2646	self._TUnit = None
2647	self._TWindowPosition = None
2648
2649	self._DisplayName = _(u'spatiotemporal aggregation of %(dn)s') % {u'dn': grid.DisplayName}
2650
2651	# Initialize the base class.
2652
2653	queryableAttributes = tuple(grid.GetAllQueryableAttributes())
2654
2655	queryableAttributeValues = {}
2656	for qa in queryableAttributes:
2657	queryableAttributeValues[qa.Name] = grid.GetQueryableAttributeValue(qa.Name)
2658
2659	super(SpatiotemporalAggregateGrid, self).__init__(parentCollection=grid.ParentCollection, queryableAttributes=queryableAttributes, queryableAttributeValues=queryableAttributeValues)
2660
2661	def _Close(self):
2662	if hasattr(self, '_Grid') and self._Grid is not None:
2663	self._Grid.Close()
2664	super(SpatiotemporalAggregateGrid, self)._Close()
2665
2666	def _GetDisplayName(self):
2667	return self._DisplayName
2668
2669	def _GetLazyPropertyPhysicalValue(self, name):
2670
2671	# If the caller requested PhysicalDimensions or
2672	# PhysicalDimensionsFlipped, return idealized values, as the
2673	# transposing and flipping is handled by the contained grid.
2674
2675	if name == 'PhysicalDimensions':
2676	return self._Grid.Dimensions
2677
2678	if name == 'PhysicalDimensionsFlipped':
2679	return tuple([False] * len(self._Grid.Dimensions))
2680
2681	# If the caller requested CoordDependencies, return None for
2682	# all dimensions because all of our dimensions have constant
2683	# increments, by design.
2684
2685	if name == 'CoordDependencies':
2686	return tuple([None] * len(self._Grid.Dimensions))
2687
2688	# If the caller requested TIncrementUnit and we are not
2689	# aggregating in the t direction, return the contained grid's
2690	# TIncrementUnit. Otherwise return the TIncrementUnit we were
2691	# constructed with.
2692
2693	if name == 'TIncrementUnit':
2694	if self._TUnit is None:
2695	return self._Grid.GetLazyPropertyValue('TIncrementUnit')
2696
2697	return self._TUnit
2698
2699	# If the caller requested TSemiRegularity and we are not
2700	# aggregating in the t direction, return the contained grid's
2701	# TSemiRegularity. Otherwise, if the t window position is not
2702	# 'year-aligned' or it is but the TUnit and TIncrement would
2703	# not result in a truncated or extended window at the end of
2704	# the year, return None. Otherwise return 'annual', indicating
2705	# that the windows are aligned to the first day of the year
2706	# and that the last window may be truncated or extended due to
2707	# remainder days or leap year.
2708
2709	if name == 'TSemiRegularity':
2710	if self._TUnit is None:
2711	return self._Grid.GetLazyPropertyValue('TSemiRegularity')
2712
2713	if self._TWindowPosition != 'year-aligned' or \
2714	self._TUnit == 'year' or \
2715	self._TUnit == 'month' and 12 % self._TIncrement == 0 or \
2716	self._TUnit == 'day' and self._TIncrement == 1 or \
2717	self._TUnit == 'hour' and 24 % self._TIncrement == 0 or \
2718	self._TUnit == 'minute' and (24 * 60) % self._TIncrement == 0 or \
2719	self._TUnit == 'second' and (24 * 60 * 60) % self._TIncrement == 0:
2720	return None
2721
2722	return 'annual'
2723
2724	# If the caller requested TCountPerSemiRegularPeriod, and we
2725	# are not aggregating in the t direction, return the contained
2726	# grid's TCountPerSemiRegularPeriod. Otherwise, if there is no
2727	# TSemiRegularity, return None. Othwerise, compute and return
2728	# the count of time slices per year.
2729
2730	if name == 'TCountPerSemiRegularPeriod':
2731	if self._TUnit is None:
2732	return self._Grid.GetLazyPropertyValue('TCountPerSemiRegularPeriod')
2733
2734	if self.TSemiRegularity is None:
2735	return None
2736
2737	if self._TUnit == u'month':
2738	f, i = math.modf(12. / self._TIncrement)
2739	elif self._TUnit == u'day':
2740	f, i = math.modf(365. / self._TIncrement)
2741	elif self._TUnit == u'hour':
2742	f, i = math.modf(365. * 24. / self._TIncrement)
2743	elif self._TUnit == u'minute':
2744	f, i = math.modf(365. * 24. * 60. / self._TIncrement)
2745	elif self._TUnit == u'second':
2746	f, i = math.modf(365. * 24. * 60. * 60. / self._TIncrement)
2747	else:
2748	raise RuntimeError(_(u'Programming error in this tool: unknown TUnit "%(tunit)s". Please contact the author of this tool for assistance.') % {u'tunit': self._TUnit})
2749
2750	if f != 0:
2751	i += 1
2752
2753	return int(i)
2754
2755	# If the caller requested CoordIncrements, calculate the y and
2756	# x increments by simple multiplication and return the t and z
2757	# increments that were provided to our constructor.
2758
2759	if name == 'CoordIncrements':
2760	coordIncrements = []
2761
2762	if 't' in self._Grid.Dimensions:
2763	if self._TIncrement is not None:
2764	coordIncrements.append(self._TIncrement)
2765	else:
2766	coordIncrements.append(self._Grid.CoordIncrements[0])
2767
2768	if 'z' in self._Grid.Dimensions:
2769	if self._ZIncrement is not None:
2770	coordIncrements.append(self._ZIncrement)
2771	else:
2772	coordIncrements.append(self._Grid.CoordIncrements[-3])
2773
2774	if self._XYIncrement is not None and self._XYIncrement > 1:
2775	tuple(coordIncrements + [self._Grid.CoordIncrements[-2] * self._XYIncrement, self._Grid.CoordIncrements[-1] * self._XYIncrement])
2776
2777	return tuple(coordIncrements + list(self._Grid.CoordIncrements[-2:]))
2778
2779	# For TCornerCoordType, if we are aggregating in the t
2780	# direction, always return 'min'.
2781
2782	if name == 'TCornerCoordType':
2783	if 't' not in self._Grid.Dimensions or self._TWindowSize is None:
2784	return self._Grid.GetLazyPropertyValue('TCornerCoordType')
2785
2786	return 'min'
2787
2788	# If the caller requested CornerCoords, calculate and return
2789	# them.
2790
2791	if name == 'CornerCoords':
2792	cornerCoords = []
2793
2794	# For the t corner coordinate, if we are not aggregating
2795	# in the t direction, just use the corner coordinate of
2796	# the contained grid.
2797
2798	if 't' in self._Grid.Dimensions:
2799	if self._TWindowSize is None:
2800	cornerCoords.append(self._Grid.GetLazyPropertyValue('CornerCoords')[0])
2801
2802	# Otherwise (we are aggregating in the t direction),
2803	# the t corner coordinate is the minimum t coordinate
2804	# of the first aggregation window. That depends on the
2805	# starting date of the contained grid, the t window
2806	# size, and how the caller wants to align the t window
2807	# relative to the center t coordinates of the
2808	# contained grid.
2809
2810	else:
2811	t0 = self._Grid.CenterCoords['t', 0]
2812
2813	if self._TUnit == u'year':
2814	if self._TWindowPosition == u'start':
2815	cornerCoords.append(datetime.datetime(t0.year, 1, 1))
2816	elif self._TWindowPosition in [u'center', u'year-aligned']:
2817	cornerCoords.append(datetime.datetime(t0.year - (self._TWindowSize - 1) / 2, 1, 1))
2818	else:
2819	cornerCoords.append(datetime.datetime(t0.year - (self._TWindowSize - 1), 1, 1))
2820
2821	elif self._TUnit == u'month':
2822	if self._TWindowPosition == u'start':
2823	months = t0.year * 12 + t0.month - 1
2824	elif self._TWindowPosition == u'center':
2825	months = t0.year * 12 + t0.month - (self._TWindowSize - 1) / 2 - 1
2826	elif self._TWindowPosition == u'end':
2827	months = t0.year * 12 + t0.month - self._TWindowSize
2828	else:
2829	months = 1
2830	while months + self._TWindowSize - 1 < t0.month:
2831	months += self._TIncrement
2832	months += t0.year * 12 - 1
2833	cornerCoords.append(datetime.datetime(months / 12, months % 12 + 1, 1))
2834
2835	elif self._TUnit == u'day':
2836	if self._TWindowPosition == u'start':
2837	cornerCoords.append(datetime.datetime(t0.year, t0.month, t0.day))
2838	elif self._TWindowPosition == u'center':
2839	cornerCoords.append(datetime.datetime(t0.year, t0.month, t0.day) - datetime.timedelta(days=(self._TWindowSize - 1) / 2))
2840	elif self._TWindowPosition == u'end':
2841	cornerCoords.append(datetime.datetime(t0.year, t0.month, t0.day) - datetime.timedelta(days=(self._TWindowSize - 1)))
2842	else:
2843	t0DayOfYear = int(t0.strftime('%j'))
2844	dayOfYear = 1
2845	while dayOfYear + self._TWindowSize - 1 < t0DayOfYear and dayOfYear + self._TWindowSize - 1 < 365:
2846	dayOfYear += self._TIncrement
2847	cornerCoords.append(datetime.datetime(t0.year, 1, 1) + datetime.timedelta(days=dayOfYear - 1))
2848
2849	elif self._TUnit == u'hour':
2850	if self._TWindowPosition == u'start':
2851	cornerCoords.append(datetime.datetime(t0.year, t0.month, t0.day, t0.hour))
2852	elif self._TWindowPosition == u'center':
2853	cornerCoords.append(datetime.datetime(t0.year, t0.month, t0.day, t0.hour) - datetime.timedelta(hours=(self._TWindowSize - 1) / 2))
2854	elif self._TWindowPosition == u'end':
2855	cornerCoords.append(datetime.datetime(t0.year, t0.month, t0.day, t0.hour) - datetime.timedelta(hours=(self._TWindowSize - 1)))
2856	else:
2857	t0HourOfYear = int(t0.strftime('%j')) * 24 + t0.hour
2858	hourOfYear = 0
2859	while hourOfYear + self._TWindowSize - 1 < t0HourOfYear and hourOfYear + self._TWindowSize - 1 < 365 * 24 - 1:
2860	hourOfYear += self._TIncrement
2861	cornerCoords.append(datetime.datetime(t0.year, 1, 1) + datetime.timedelta(hours=hourOfYear))
2862
2863	elif self._TUnit == u'minute':
2864	if self._TWindowPosition == u'start':
2865	cornerCoords.append(datetime.datetime(t0.year, t0.month, t0.day, t0.hour, t0.minute))
2866	elif self._TWindowPosition == u'center':
2867	cornerCoords.append(datetime.datetime(t0.year, t0.month, t0.day, t0.hour, t0.minute) - datetime.timedelta(minutes=(self._TWindowSize - 1) / 2))
2868	elif self._TWindowPosition == u'end':
2869	cornerCoords.append(datetime.datetime(t0.year, t0.month, t0.day, t0.hour, t0.minute) - datetime.timedelta(minutes=(self._TWindowSize - 1)))
2870	else:
2871	t0MinuteOfYear = int(t0.strftime('%j')) * 24 * 60 + t0.hour * 60 + t0.minute
2872	minuteOfYear = 0
2873	while minuteOfYear + self._TWindowSize - 1 < t0MinuteOfYear and minuteOfYear + self._TWindowSize - 1 < 365 * 24 * 60 - 1:
2874	minuteOfYear += self._TIncrement
2875	cornerCoords.append(datetime.datetime(t0.year, 1, 1) + datetime.timedelta(minutes=minuteOfYear))
2876
2877	else:
2878	if self._TWindowPosition == u'start':
2879	cornerCoords.append(datetime.datetime(t0.year, t0.month, t0.day, t0.hour, t0.minute, t0.second))
2880	elif self._TWindowPosition == u'center':
2881	cornerCoords.append(datetime.datetime(t0.year, t0.month, t0.day, t0.hour, t0.minute, t0.second) - datetime.timedelta(seconds=(self._TWindowSize - 1) / 2))
2882	elif self._TWindowPosition == u'end':
2883	cornerCoords.append(datetime.datetime(t0.year, t0.month, t0.day, t0.hour, t0.minute, t0.second) - datetime.timedelta(seconds=(self._TWindowSize - 1)))
2884	else:
2885	t0SecondOfYear = int(t0.strftime('%j')) * 24 * 60 * 60 + t0.hour * 60 * 60 + t0.minute * 60 + t0.second
2886	secondOfYear = 0
2887	while secondOfYear + self._TWindowSize - 1 < t0SecondOfYear and secondOfYear + self._TWindowSize - 1 < 365 * 24 * 60 * 60 - 1:
2888	secondOfYear += self._TIncrement
2889	cornerCoords.append(datetime.datetime(t0.year, 1, 1) + datetime.timedelta(seconds=secondOfYear))
2890
2891	# TODO: z
2892
2893	# For the y and x corner coordinates, just use those from
2894	# the contained grid.
2895
2896	return tuple(cornerCoords + list(self._Grid.GetLazyPropertyValue('CornerCoords')[-2:]))
2897
2898	# If the caller requested Shape, calculate and return it. For
2899	# each coordinate, count the number of aggregate cells
2900	# required to reach the last cell of the contained grid.
2901
2902	if name == 'Shape':
2903	shape = []
2904
2905	if 't' in self._Grid.Dimensions:
2906	if self._TWindowSize is None:
2907	shape.append(self._Grid.Shape[0])
2908	else:
2909	tMax = self._GetTCoordsList(0.5, self._Grid.Shape[0])
2910	shape.append(bisect.bisect_right(tMax, self._Grid.CenterCoords['t', -1]) + 1)
2911
2912	if 'z' in self._Grid.Dimensions:
2913	if self._ZWindowSize is None:
2914	shape.append(self._Grid.Shape[-3])
2915	else:
2916	i = 1
2917	increment = self.CoordIncrements[-3]
2918	c = self._Grid.CenterCoords['z', -1]
2919	while i * increment <= c:
2920	i += 1
2921	shape.append(i)
2922
2923	if self._XYWindowSize is None or self._XYWindowSize <= 1:
2924	shape += self._Grid.Shape[-2:]
2925	else:
2926	for coord in ['y', 'x']:
2927	i = 1
2928	increment = self.CoordIncrements[self.Dimensions.index(coord)]
2929	c = self._Grid.CenterCoords[coord, -1]
2930	while i * increment <= c:
2931	i += 1
2932	shape.append(i)
2933
2934	return shape
2935
2936	# If the caller requested the UnscaledDataType or
2937	# UnscaledNoDataValue, calculate and return them.
2938
2939	if name == 'UnscaledDataType':
2940
2941	# For COUNT, return the smallest integer type that is
2942	# guaranteed hold the largest possible value.
2943
2944	if self._Statistic == u'count':
2945
2946	# First, for each dimension, find the maximum number
2947	# of cells that will be binned in that dimension.
2948
2949	maxCells = [0] * len(self.Dimensions)
2950	for i, d in enumerate(self.Dimensions):
2951	maxCoords = self.MaxCoords[d, :]
2952	j = 0
2953	count = 0
2954	for c in self._Grid.CenterCoords[d, :]:
2955	if c < maxCoords[j]:
2956	count += 1
2957	else:
2958	if count > maxCells[i]:
2959	maxCells[i] = count
2960	j += 1
2961	count = 1
2962	if count > maxCells[i]:
2963	maxCells[i] = count
2964
2965	# Now multiply those counts together to get a maximum
2966	# possible count.
2967
2968	maxCount = reduce(lambda x, y: x*y, maxCells)
2969
2970	# Now return the smallest data type that ranges to
2971	# that number.
2972
2973	if maxCount < 256:
2974	return u'uint8'
2975	if maxCount < 65536:
2976	return u'uint16'
2977	return u'int32'
2978
2979	# For MIN, MAX, and RANGE, return the data type of the
2980	# grid.
2981
2982	if self._Statistic in [u'minimum', u'maximum', u'range']:
2983	return self._Grid.DataType
2984
2985	# For all others, return float32, unless the grid is
2986	# float64, in which case return float64.
2987
2988	if self._Grid.DataType == u'float64':
2989	return u'float64'
2990	return u'float32'
2991
2992	if name == 'UnscaledNoDataValue':
2993
2994	# For COUNT, return 0.
2995
2996	if self._Statistic == u'count':
2997	if self.DataType in [u'float32', u'float64']:
2998	return 0.
2999	return 0
3000
3001	# For MIN, MAX, and RANGE, return the NoData value of the
3002	# grid.
3003
3004	if self._Statistic in [u'minimum', u'maximum', u'range']:
3005	return self_Grid.NoDataValue
3006
3007	# For all others, return the smallest representable
3008	# floating point number:
3009	#
3010	# >>> float(numpy.finfo('float32').min)
3011	# -3.4028234663852886e+038
3012	# >>> float(numpy.finfo('float64').min)
3013	# -1.7976931348623157e+308
3014
3015	import numpy
3016
3017	if self.DataType == u'float32':
3018	return float(numpy.finfo('float32').min)
3019	return float(numpy.finfo('float64').min)
3020
3021	# If the caller requested one of the properties related to
3022	# scaling, return None.
3023
3024	if name in ['ScaledDataType', 'ScaledNoDataValue', 'ScalingFunction', 'UnscalingFunction']:
3025	return None
3026
3027	# Otherwise use the value of the property from the grid.
3028
3029	return self._Grid.GetLazyPropertyValue(name)
3030
3031	@classmethod
3032	def _TestCapability(cls, capability):
3033	return cls._Grid._TestCapability(capability) # TODO: Reject writing
3034
3035	def _ReadNumpyArray(self, sliceList):
3036	import numpy
3037
3038
3039	class InpaintedGrid(Grid):
3040	__doc__ = DynamicDocString()
3041
3042	def __init__(self, grid, method=u'Del2b', maxHoleSize=None, xEdgesWrap=False):
3043	# TODO: self.__class__.__doc__.Obj.ValidateMethodInvocation()
3044
3045	# Perform additional validation.
3046
3047	if grid.DataType not in [u'float32', u'float64']:
3048	raise TypeError(_(u'This tool can only fill missing values in floating-point grids. It cannot fill missing values in integer grids.'))
3049
3050	# Initialize our properties.
3051
3052	self._Grid = grid
3053	self._TempDir = None
3054	self._Method = method
3055	self._MaxHoleSize = maxHoleSize
3056	self._XEdgesWrap = xEdgesWrap
3057
3058	if maxHoleSize is None:
3059	self._DisplayName = _(u'%(dn)s with missing cells filled using the %(method)s') % {u'dn': grid.DisplayName, u'method': method}
3060	else:
3061	self._DisplayName = _(u'%(dn)s with clusters of missing cells no larger than %(maxHoleSize)i filled using the %(method)s') % {u'dn': grid.DisplayName, u'maxHoleSize': maxHoleSize, u'method': method}
3062
3063	# Initialize the base class.
3064
3065	queryableAttributes = tuple(grid.GetAllQueryableAttributes())
3066
3067	queryableAttributeValues = {}
3068	for qa in queryableAttributes:
3069	queryableAttributeValues[qa.Name] = grid.GetQueryableAttributeValue(qa.Name)
3070
3071	super(InpaintedGrid, self).__init__(queryableAttributes=queryableAttributes, queryableAttributeValues=queryableAttributeValues)
3072
3073	def _Close(self):
3074	if hasattr(self, '_Grid') and self._Grid is not None:
3075	self._Grid.Close()
3076	self._TempDir = None
3077	super(InpaintedGrid, self)._Close()
3078
3079	def _GetDisplayName(self):
3080	return self._DisplayName
3081
3082	def _GetLazyPropertyPhysicalValue(self, name):
3083
3084	# If the caller requested PhysicalDimensions or
3085	# PhysicalDimensionsFlipped, return idealized values, as any
3086	# transposing and flipping is handled by the contained grid.
3087
3088	if name == 'PhysicalDimensions':
3089	return self._Grid.Dimensions
3090
3091	if name == 'PhysicalDimensionsFlipped':
3092	return tuple([False] * len(self._Grid.Dimensions))
3093
3094	# The contained grid also handles everything related to
3095	# scaling.
3096
3097	if name == 'UnscaledDataType':
3098	return self._Grid.DataType
3099
3100	if name == 'UnscaledNoDataValue':
3101	return self._Grid.NoDataValue
3102
3103	if name in ['ScaledDataType', 'ScaledNoDataValue', 'ScalingFunction', 'UnscalingFunction']:
3104	return None
3105
3106	# Otherwise use the value of the property from the grid.
3107
3108	return self._Grid.GetLazyPropertyValue(name)
3109
3110	@classmethod
3111	def _TestCapability(cls, capability):
3112	return cls._Grid._TestCapability(capability) # TODO: Reject writing
3113
3114	@classmethod
3115	def _GetSRTypeForSetting(cls):
3116	return cls._Grid._GetSRTypeForSetting()
3117
3118	def _SetSpatialReference(self, sr):
3119	return self._Grid._SetSpatialReference(sr)
3120
3121	def _GetCoords(self, coord, coordNum, slices, sliceDims, fixedIncrementOffset):
3122	return self._Grid._GetCoords(coord, coordNum, slices, sliceDims, fixedIncrementOffset)
3123
3124	def _ReadNumpyArray(self, sliceList):
3125
3126	# If we have not created the temporary directory, do it now.
3127
3128	if self._TempDir is None:
3129	self._TempDir = self._CreateTempDirectory()
3130
3131	# Iterate through each requested 2D slice and create an
3132	# inpainted version of it in the temporary directory, if we
3133	# have not done so already. If we catch an exception, call
3134	# _Close() to delete the temporary directory.
3135
3136	try:
3137	if len(self.Dimensions) == 2:
3138	slices = [(sliceList[0], sliceList[1])]
3139	elif len(self.Dimensions) == 2:
3140	slices = [(d, sliceList[1], sliceList[2]) for d in range(sliceList[0])]
3141	else:
3142	slices = [(d1, d2, sliceList[2], sliceList[3]) for d1 in range(sliceList[0]) for d2 in range(sliceList[1])]
3143
3144	for s in slices:
3145	inpaintedFile = os.path.join(self._TempDir, 'slice_%s_inpainted.dat' % '_'.join(map(str, slices[:len(self.Dimensions) - 2])))
3146	if not os.path.exists(inpaintedFile):
3147	self._LogDebug(_(u'%(class)s 0x%(id)08X: Creating %(inpaintedFile)s from slice %(slice)r of %(dn)s.'), {u'class': self.__class__.__name__, u'id': id(self), u'inpaintedFile': inpaintedFile, u'slice': slices, u'dn': self._Grid.DisplayName})
3148
3149	# Write this 2D slice to the temporary directory.
3150
3151	holesFile = os.path.join(self._TempDir, 'slice_%s.dat' % '_'.join(map(str, slices[:len(self.Dimensions) - 2])))
3152
3153	if self.Dimensions == u'yx':
3154	grid = self._Grid
3155	elif self.Dimensions == u'zyx':
3156	grid = GridSlice(self._Grid, zIndex=d)
3157	elif self.Dimensions == u'tyx':
3158	grid = GridSlice(self._Grid, tIndex=d)
3159	else:
3160	grid = GridSlice(self._Grid, tIndex=d1, zIndex=d2)
3161
3162	grid.Data[:].tofile(holesFile)
3163
3164	# Execute InpaintGrid.py to do the inpaint
3165	# operation. This script calls a MATLAB function.
3166	# We prefer to call that function directly right
3167	# here but there is a continuing incompatibility
3168	# between MATLAB DLLs and ArcGIS DLLs (they both
3169	# try to load their own incompatible versions of
3170	# xerces-c_2_7.dll) so we have to do it in a
3171	# separate process.
3172
3173	from GeoEco.DataManagement.Processes import ChildProcess
3174	from GeoEco.Logging import Logger
3175	import win32api
3176
3177	args = [os.path.join(os.path.dirname(__file__), 'InpaintGrid.py'),
3178	holesFile,
3179	inpaintedFile,
3180	str(grid.Shape[0]),
3181	str(grid.Shape[1]),
3182	str(grid.DataType),
3183	str(grid.NoDataValue),
3184	{u'del2a': '1', u'del2b': '0', u'del2c': '2', u'del4': '3', 'spring': '4'}[self._Method.lower()],
3185	str(self._MaxHoleSize),
3186	str(self._XEdgesWrap)]
3187
3188	if self._MaxHoleSize is None:
3189	args[8] = '0'
3190
3191	if self._XEdgesWrap is None:
3192	args[9] = '0'
3193
3194	oldLogInfoAsDebug = Logger.GetLogInfoAsDebug()
3195	Logger.SetLogInfoAsDebug(True)
3196	try:
3197	ChildProcess.ExecuteProgram(win32api.FindExecutable(os.path.join(os.path.dirname(__file__), 'InpaintGrid.py'), os.path.dirname(__file__))[1],
3198	arguments=args,
3199	windowState=u'invisible',
3200	maxRunTime=None)
3201	finally:
3202	Logger.SetLogInfoAsDebug(oldLogInfoAsDebug)
3203	except:
3204	self._Close()
3205	raise
3206
3207	# Allocate an array to return.
3208
3209	import numpy
3210	data = numpy.zeros(map(lambda s: s.stop-s.start, sliceList), dtype=str(self.UnscaledDataType))
3211
3212	# Read the inpainted slices from the temporary directory and
3213	# write them to the array.
3214
3215	for s in slices:
3216	inpaintedFile = os.path.join(self._TempDir, 'slice_%s_inpainted.dat' % '_'.join(map(str, slices[:len(self.Dimensions) - 2])))
3217	data.__setitem__(s, numpy.fromfile(inpaintedFile, str(self.DataType)).reshape(self.Shape[-2], self.Shape[-1]).__getitem__((s[-2], s[-1])))
3218
3219	# Return successfully.
3220
3221	return data, self.NoDataValue
3222
3223
3224	class MosaickedTiles(Grid):
3225	__doc__ = DynamicDocString()
3226
3227	def _GetTiles(self):
3228	return self._Grids
3229
3230	Tiles = property(_GetTiles, doc=DynamicDocString())
3231
3232	def __init__(self, tiles, displayName):
3233	# TODO: self.__class__.__doc__.Obj.ValidateMethodInvocation()
3234
3235	# TODO: Validate that the caller's list of lists of grids has
3236	# the proper structure.
3237
3238	# TODO: Validate that the caller's grids have the same
3239	# characteristics.
3240
3241	# Initialize our properties.
3242
3243	self._Grids = tiles
3244	self._DisplayName = displayName
3245
3246	# Initialize the base class.
3247
3248	queryableAttributes = tuple(tiles[0][0].GetAllQueryableAttributes())
3249
3250	queryableAttributeValues = {}
3251	for qa in queryableAttributes:
3252	queryableAttributeValues[qa.Name] = tiles[0][0].GetQueryableAttributeValue(qa.Name)
3253
3254	super(MosaickedTiles, self).__init__(queryableAttributes=queryableAttributes, queryableAttributeValues=queryableAttributeValues)
3255
3256	def _Close(self):
3257	if hasattr(self, '_Grids') and self._Grids is not None:
3258	for row in self._Grids:
3259	for grid in row:
3260	grid.Close()
3261	super(MosaickedTiles, self)._Close()
3262
3263	def _GetDisplayName(self):
3264	return self._DisplayName
3265
3266	def _GetLazyPropertyPhysicalValue(self, name):
3267
3268	# If the caller requested PhysicalDimensions or
3269	# PhysicalDimensionsFlipped, return idealized values, as any
3270	# transposing and flipping is handled by the contained grids.
3271
3272	if name == 'PhysicalDimensions':
3273	return self._Grids[0][0].Dimensions
3274
3275	if name == 'PhysicalDimensionsFlipped':
3276	return tuple([False] * len(self._Grids[0][0].Dimensions))
3277
3278	# The contained grids also handle everything related to
3279	# scaling.
3280
3281	if name == 'UnscaledDataType':
3282	return self._Grids[0][0].DataType
3283
3284	if name == 'UnscaledNoDataValue':
3285	return self._Grids[0][0].NoDataValue
3286
3287	if name in ['ScaledDataType', 'ScaledNoDataValue', 'ScalingFunction', 'UnscalingFunction']:
3288	return None
3289
3290	# If the caller requested the Shape, calculate it by adding up
3291	# the shapes of the contained grids.
3292
3293	if name == 'Shape':
3294	x = reduce(lambda a, b: a + b, [col.Shape[-1] for col in self._Grids[0]])
3295	y = reduce(lambda a, b: a + b, [row[0].Shape[-2] for row in self._Grids])
3296	return tuple(list(self._Grids[0][0].Shape[:-2]) + [y, x])
3297
3298	# Otherwise use the value of the property from the first
3299	# contained grid. Note: this includes the CornerCoords.
3300	# Because of that, do not change the [0][0] to some other
3301	# indices.
3302
3303	return self._Grids[0][0].GetLazyPropertyValue(name)
3304
3305	@classmethod
3306	def _TestCapability(cls, capability):
3307	return cls._Grid._TestCapability(capability) # TODO: Reject writing
3308
3309	@classmethod
3310	def _GetSRTypeForSetting(cls):
3311	return cls._Grid._GetSRTypeForSetting()
3312
3313	def _SetSpatialReference(self, sr):
3314	return self._Grid._SetSpatialReference(sr)
3315
3316	def _GetCoords(self, coord, coordNum, slices, sliceDims, fixedIncrementOffset):
3317	if coord in 'xy':
3318	raise NotImplementedError(_(u'The MosaickedTiles class does not currently support grids that contain x or y coordinates that have non-constant increments. Please contact the MGET development team for assistance.'))
3319
3320	return self._Grids[0][0]._GetCoords(coord, coordNum, slices, sliceDims, fixedIncrementOffset)
3321
3322	def _ReadNumpyArray(self, sliceList):
3323
3324	# Allocate an array to return.
3325
3326	import numpy
3327	data = numpy.zeros(map(lambda s: s.stop-s.start, sliceList), dtype=str(self.UnscaledDataType))
3328
3329	# Find the indices of the first tile for the requested slice.
3330
3331	xIndexIntoFirstTile = sliceList[-1].start
3332	xFirstTile = 0
3333	while xIndexIntoFirstTile >= self._Grids[0][xFirstTile].Shape[-1]:
3334	xIndexIntoFirstTile -= self._Grids[0][xFirstTile].Shape[-1]
3335	xFirstTile += 1
3336
3337	yIndexIntoFirstTile = sliceList[-2].start
3338	yFirstTile = 0
3339	while yIndexIntoFirstTile >= self._Grids[yFirstTile][0].Shape[-2]:
3340	yIndexIntoFirstTile -= self._Grids[yFirstTile][0].Shape[-2]
3341	yFirstTile += 1
3342
3343	# Fill the array with slices from the tiles.
3344
3345	yTile = yFirstTile
3346	yCellsRead = 0
3347	yCellsRemaining = data.shape[-2]
3348
3349	while yCellsRemaining > 0:
3350	if yTile == yFirstTile:
3351	yStart = yIndexIntoFirstTile
3352	else:
3353	yStart = 0
3354	yStop = min(yStart + yCellsRemaining, self._Grids[yTile][0].Shape[-2])
3355
3356	xTile = xFirstTile
3357	xCellsRead = 0
3358	xCellsRemaining = data.shape[-1]
3359
3360	while xCellsRemaining > 0:
3361	if xTile == xFirstTile:
3362	xStart = xIndexIntoFirstTile
3363	else:
3364	xStart = 0
3365	xStop = min(xStart + xCellsRemaining, self._Grids[0][xTile].Shape[-1])
3366
3367	data.__setitem__(map(lambda s: slice(0, s.stop-s.start), sliceList[:-2]) + [slice(yCellsRead, yCellsRead + yStop - yStart), slice(xCellsRead, xCellsRead + xStop - xStart)],
3368	self._Grids[yTile][xTile].Data.__getitem__(tuple(sliceList[:-2] + [slice(yStart, yStop), slice(xStart, xStop)])))
3369
3370	xTile += 1
3371	xCellsRead += xStop - xStart
3372	xCellsRemaining -= xStop - xStart
3373
3374	yTile += 1
3375	yCellsRead += yStop - yStart
3376	yCellsRemaining -= yStop - yStart
3377
3378	# Return successfully.
3379
3380	return data, self.NoDataValue
3381
3382
3383	###############################################################################
3384	# Metadata: module
3385	###############################################################################
3386
3387	from GeoEco.Metadata import *
3388	from GeoEco.Types import *
3389
3390	AddModuleMetadata(shortDescription=_(u'Classes representing virtual datasets.')) # TODO: Better description
3391
3392	###############################################################################
3393	# Metadata: TimeSeriesGridStack class
3394	###############################################################################
3395
3396	AddClassMetadata(TimeSeriesGridStack,
3397	shortDescription=_(u'TODO: Add description.'))
3398
3399	# TODO: Add metadata
3400
3401	###############################################################################
3402	# Metadata: GridSlice class
3403	###############################################################################
3404
3405	AddClassMetadata(GridSlice,
3406	shortDescription=_(u'TODO: Add description.'))
3407
3408	# TODO: Add metadata
3409
3410	###############################################################################
3411	# Metadata: GridSliceCollection class
3412	###############################################################################
3413
3414	AddClassMetadata(GridSliceCollection,
3415	shortDescription=_(u'TODO: Add description.'))
3416
3417	# TODO: Add metadata
3418
3419	###############################################################################
3420	# Metadata: ClippedGrid class
3421	###############################################################################
3422
3423	AddClassMetadata(ClippedGrid,
3424	shortDescription=_(u'TODO: Add description.'))
3425
3426	# Constructor
3427
3428	_ClipMaxParameterDescription = _(
3429	u"""%(extent)s %(dim)s coordinate or index.
3430
3431	The Clip By parameter determines whether the value is a coordinate or
3432	an index. If a coordinate is provided, it may be an integer or a
3433	floating point number. If an index is provided, it must be an
3434	non-negative integer.
3435
3436	If a value is not provided, the clipped grid will extend to the full
3437	extent in the %(dir)s %(dim)s direction.""")
3438
3439	AddMethodMetadata(ClippedGrid.__init__,
3440	shortDescription=_(u'Constructs a new ClippedGrid instance.'),
3441	isExposedToPythonCallers=True)
3442
3443	AddArgumentMetadata(ClippedGrid.__init__, u'self',
3444	typeMetadata=ClassInstanceTypeMetadata(cls=ClippedGrid),
3445	description=_(u'ClippedGrid instance.'))
3446
3447	AddArgumentMetadata(ClippedGrid.__init__, u'grid',
3448	typeMetadata=ClassInstanceTypeMetadata(cls=Grid),
3449	description=_(u"""Grid to clip."""))
3450
3451	AddArgumentMetadata(ClippedGrid.__init__, u'clipBy',
3452	typeMetadata=UnicodeStringTypeMetadata(allowedValues=[u'Map coordinates', u'Cell indices'], makeLowercase=True),
3453	description=_(
3454	u"""Specifies the type of extent values that will be used to clip the
3455	grid, one of:
3456
3457	* Coordinates - the values are floating-point numbers representing
3458	coordinates in space and time. Typically this means that x and y
3459	coordinates will be in degrees (if a geographic coordinate system is
3460	used) or a linear unit such as meters (if a projected coordinate
3461	system is used), the z coordinate will be in a linear unit such as
3462	meters, and the time coordinate will be a date and time string in a
3463	standard format such as 'YYYY-MM-DD HH:MM:SS' or a Python datetime
3464	instance.
3465
3466	* Indices - the values are integer indices of cells in the grid,
3467	starting at 0 and ending at the maximum index value. For example, if
3468	a 2D grid has 360 columns and 180 rows, the x and y indices will
3469	range from 0 to 359 and 0 to 179, respectively.
3470
3471	If 'Indices' is specified but floating point numbers are provided for
3472	the extent values, the decimal portions of the numbers will be
3473	truncated (e.g. the value 1.7 will be truncated to 1).
3474
3475	If 'Coordinates' is specified and a coordinate value falls within a
3476	cell, rather than exactly on a boundary between two cells, the cell
3477	will be included in the clipped grid (it will not be clipped out).
3478	Because computers cannot represent all floating-point numbers at full
3479	precision, the resulting rounding errors can sometimes produce
3480	unexpected results. For example, consider a grid with a cell size of
3481	0.1. We expect that the cells centered at 0.5 and 1.5 would meet at
3482	coordinate value 0.1, but 0.1 cannot be fully represented by the
3483	computer using standard 64-bit floating-point numbers. The computer
3484	rounds it to 0.10000000000000001. Therefore, if you were to clip the
3485	grid a maximum coordinate of 0.10000000000000001, you would expect
3486	that the cell centered at 1.5 would be included in the resulting grid,
3487	because 0.10000000000000001 falls within that cell. But the computer
3488	actually considers the boundary between the two cells to be at
3489	0.10000000000000001, not 0.1, so that cell would be clipped out."""), # TODO: Add discussion of how small indices correspond to small coordinate values
3490	arcGISDisplayName=_(u'Clip by'))
3491
3492	AddArgumentMetadata(ClippedGrid.__init__, u'xMin',
3493	typeMetadata=FloatTypeMetadata(canBeNone=True),
3494	description=_ClipMaxParameterDescription % {u'extent': _(u'Minimum'), u'dim': u'x', u'dir': _(u'negative')},
3495	arcGISDisplayName=_(u'Minimum x extent'))
3496
3497	AddArgumentMetadata(ClippedGrid.__init__, u'xMax',
3498	typeMetadata=FloatTypeMetadata(canBeNone=True),
3499	description=_ClipMaxParameterDescription % {u'extent': _(u'Maximum'), u'dim': u'x', u'dir': _(u'positive')},
3500	arcGISDisplayName=_(u'Maximum x extent'))
3501
3502	AddArgumentMetadata(ClippedGrid.__init__, u'yMin',
3503	typeMetadata=FloatTypeMetadata(canBeNone=True),
3504	description=_ClipMaxParameterDescription % {u'extent': _(u'Minimum'), u'dim': u'y', u'dir': _(u'negative')},
3505	arcGISDisplayName=_(u'Minimum y extent'))
3506
3507	AddArgumentMetadata(ClippedGrid.__init__, u'yMax',
3508	typeMetadata=FloatTypeMetadata(canBeNone=True),
3509	description=_ClipMaxParameterDescription % {u'extent': _(u'Maximum'), u'dim': u'y', u'dir': _(u'positive')},
3510	arcGISDisplayName=_(u'Maximum y extent'))
3511
3512	AddArgumentMetadata(ClippedGrid.__init__, u'zMin',
3513	typeMetadata=FloatTypeMetadata(canBeNone=True),
3514	description=_ClipMaxParameterDescription % {u'extent': _(u'Minimum'), u'dim': u'z', u'dir': _(u'negative')}, # TODO: What about the sign of z? Does z get smaller or larger with increasing depth?
3515	arcGISDisplayName=_(u'Minimum z extent'))
3516
3517	AddArgumentMetadata(ClippedGrid.__init__, u'zMax',
3518	typeMetadata=FloatTypeMetadata(canBeNone=True),
3519	description=_ClipMaxParameterDescription % {u'extent': _(u'Maximum'), u'dim': u'z', u'dir': _(u'positive')},
3520	arcGISDisplayName=_(u'Maximum z extent'))
3521
3522	AddArgumentMetadata(ClippedGrid.__init__, u'tMin',
3523	typeMetadata=AnyObjectTypeMetadata(canBeNone=True),
3524	description=_(
3525	u"""Minimum t coordinate or index.
3526
3527	The Clip By parameter determines whether the value is a coordinate or
3528	an index. If a coordinate is provided, it may be a date and time
3529	string in a standard format such as 'YYYY-MM-DD HH:MM:SS' or a Python
3530	datetime instance. If an index is provided, it must be an non-negative
3531	integer.
3532
3533	If a value is not provided, the clipped grid will extend to the full
3534	extent in the negative t direction."""),
3535	arcGISDisplayName=_(u'Minimum t extent'))
3536
3537	AddArgumentMetadata(ClippedGrid.__init__, u'tMax',
3538	typeMetadata=AnyObjectTypeMetadata(canBeNone=True),
3539	description=_(
3540	u"""Maximum t coordinate or index.
3541
3542	The Clip By parameter determines whether the value is a coordinate or
3543	an index. If a coordinate is provided, it may be a date and time
3544	string in a standard format such as 'YYYY-MM-DD HH:MM:SS' or a Python
3545	datetime instance. If an index is provided, it must be an non-negative
3546	integer.
3547
3548	If a value is not provided, the clipped grid will extend to the full
3549	extent in the positive t direction."""),
3550	arcGISDisplayName=_(u'Maximum t extent'))
3551
3552	AddResultMetadata(ClippedGrid.__init__, u'clippedGrid',
3553	typeMetadata=ClassInstanceTypeMetadata(cls=Grid),
3554	description=_(u"""Clipped grid."""))
3555
3556	###############################################################################
3557	# Metadata: MaskedGrid class
3558	###############################################################################
3559
3560	AddClassMetadata(MaskedGrid,
3561	shortDescription=_(u'TODO: Add description.'))
3562
3563	# TODO: Add metadata
3564
3565	###############################################################################
3566	# Metadata: DerivedGrid class
3567	###############################################################################
3568
3569	AddClassMetadata(DerivedGrid,
3570	shortDescription=_(u'TODO: Add description.'))
3571
3572	# TODO: Add metadata
3573
3574	###############################################################################
3575	# Metadata: RotatedGlobalGrid class
3576	###############################################################################
3577
3578	AddClassMetadata(RotatedGlobalGrid,
3579	shortDescription=_(u'TODO: Add description.'))
3580
3581	# Constructor
3582
3583	AddMethodMetadata(RotatedGlobalGrid.__init__,
3584	shortDescription=_(u'Constructs a new RotatedGlobalGrid instance.'),
3585	isExposedToPythonCallers=True)
3586
3587	AddArgumentMetadata(RotatedGlobalGrid.__init__, u'self',
3588	typeMetadata=ClassInstanceTypeMetadata(cls=RotatedGlobalGrid),
3589	description=_(u'RotatedGlobalGrid instance.'))
3590
3591	AddArgumentMetadata(RotatedGlobalGrid.__init__, u'grid',
3592	typeMetadata=ClassInstanceTypeMetadata(cls=Grid),
3593	description=_(
3594	u"""Grid to rotate. The grid must have global longitudinal
3595	extent."""))
3596
3597	AddArgumentMetadata(RotatedGlobalGrid.__init__, u'rotationOffset',
3598	typeMetadata=FloatTypeMetadata(),
3599	description=_(
3600	u"""Quantity to rotate the grid by, in the units specified by the
3601	Rotation Units parameter.
3602
3603	Use this parameter to center the grid on a different longitude.
3604	Positive values move the center of the grid to the right (east);
3605	negative values move the center to the west. For example, if the
3606	Rotation Units parameter is 'Cells', the value 10 will cause the
3607	center of the grid to shift 10 cells to the right. Effectively, the 10
3608	left-most cells to be stripped from the left edge and moved to the
3609	right edge."""))
3610
3611	AddArgumentMetadata(RotatedGlobalGrid.__init__, u'rotationUnits',
3612	typeMetadata=UnicodeStringTypeMetadata(allowedValues=[u'Map units', u'Cells'], makeLowercase=True),
3613	description=_(
3614	u"""Specifies the type of values that will be used to rotate the grid,
3615	one of:
3616
3617	* Map units - the unit of rotation will be the same as the linear unit
3618	of the map, typically degrees for data in geographic coordinate
3619	systems and meters for data in projected coordinate systems. Because
3620	the rotation must be in whole cells, the rotation quantity will be
3621	converted to grid cells and rounded to the nearest cell.
3622
3623	* Cells - the unit of rotation will be in grid cells. Because the
3624	rotation must be in whole cells, the rotation quantity will be
3625	rounded to the neareast cell.
3626	"""))
3627
3628	AddResultMetadata(RotatedGlobalGrid.__init__, u'rotatedGrid',
3629	typeMetadata=ClassInstanceTypeMetadata(cls=Grid),
3630	description=_(u"""Rotated grid."""))
3631
3632	###############################################################################
3633	# Metadata: MemoryCachedGrid class
3634	###############################################################################
3635
3636	AddClassMetadata(MemoryCachedGrid,
3637	shortDescription=_(u'TODO: Add description.'))
3638
3639	# TODO: Add metadata
3640
3641	###############################################################################
3642	# Metadata: AggregateGrid class
3643	###############################################################################
3644
3645	AddClassMetadata(AggregateGrid,
3646	shortDescription=_(u'TODO: Add description.'))
3647
3648	# TODO: Add metadata
3649
3650	###############################################################################
3651	# Metadata: ClimatologicalGridCollection class
3652	###############################################################################
3653
3654	AddClassMetadata(ClimatologicalGridCollection,
3655	shortDescription=_(u'TODO: Add description.'))
3656
3657	# TODO: Add metadata
3658
3659	###############################################################################
3660	# Metadata: SeafloorGrid class
3661	###############################################################################
3662
3663	AddClassMetadata(SeafloorGrid,
3664	shortDescription=_(u'TODO: Add description.'))
3665
3666	# TODO: Add metadata
3667
3668	###############################################################################
3669	# Metadata: SpatiotemporalAggregateGrid class
3670	###############################################################################
3671
3672	AddClassMetadata(SpatiotemporalAggregateGrid,
3673	shortDescription=_(u'TODO: Add description.'))
3674
3675	# TODO: Add metadata
3676
3677	###############################################################################
3678	# Metadata: InpaintedGrid class
3679	###############################################################################
3680
3681	AddClassMetadata(InpaintedGrid,
3682	shortDescription=_(u'TODO: Add description.'))
3683
3684	# Constructor
3685
3686	AddMethodMetadata(InpaintedGrid.__init__,
3687	shortDescription=_(u'Constructs a new InpaintedGrid instance.'),
3688	isExposedToPythonCallers=True)
3689
3690	AddArgumentMetadata(InpaintedGrid.__init__, u'self',
3691	typeMetadata=ClassInstanceTypeMetadata(cls=InpaintedGrid),
3692	description=_(u'InpaintedGrid instance.'))
3693
3694	AddArgumentMetadata(InpaintedGrid.__init__, u'grid',
3695	typeMetadata=ClassInstanceTypeMetadata(cls=Grid),
3696	description=_(
3697	u"""Grid to fill. Its ScaledDataType must be float32 or float64."""))
3698
3699	AddArgumentMetadata(InpaintedGrid.__init__, u'method',
3700	typeMetadata=UnicodeStringTypeMetadata(makeLowercase=True, allowedValues=[u'Del2a', u'Del2b', u'Del2c', u'Del4', u'Spring']),
3701	description=_(
3702	u"""Method to use for interpolation and extrapolation of No Data
3703	values. One of:
3704
3705	* Del2a - Performs Laplacian interpolation and linear extrapolation.
3706
3707	* Del2b - Same as Del2a but does not build as large a linear system of
3708	equations. May be faster than Del2a at the cost of some accuracy.
3709
3710	* Del2c - Same as Del2a but solves a direct linear system of equations
3711	for the No Data values. Faster than both Del2a and Del2b but is the
3712	least robust to noise on the boundaries of No Data cells and least
3713	able to interpolate accurately for smooth surfaces.
3714
3715	* Del4 - Same as Del2a but instead of the Laplace operator (also
3716	called the del^2 operator) it uses the biharmoic operator (also
3717	called the del^4 operator). May result in more accurate
3718	interpolations, at some cost in speed.
3719
3720	* Spring - Uses a spring metaphor. Assumes springs (with a nominal
3721	length of zero) connect each cell with every neighbor (horizontally,
3722	vertically and diagonally). Since each cell tries to be like its
3723	neighbors, extrapolation is as a constant function where this is
3724	consistent with the neighboring nodes.
3725	"""))
3726
3727	AddArgumentMetadata(InpaintedGrid.__init__, u'maxHoleSize',
3728	typeMetadata=IntegerTypeMetadata(mustBeGreaterThan=0, canBeNone=True),
3729	description=_(
3730	u"""Maximum size, in cells, that a region of 4-connected No Data cells
3731	may be for it to be filled in.
3732
3733	Use this option to prevent the filling of large No Data regions (e.g.
3734	large clouds in remote sensing images) when you are concerned that
3735	values cannot be accurately guessed for those regions. If this option
3736	is omitted, all regions will be filled, regardless of size."""))
3737
3738	AddArgumentMetadata(InpaintedGrid.__init__, u'xEdgesWrap',
3739	typeMetadata=BooleanTypeMetadata(),
3740	description=_(
3741	u"""If True, the left and right edges of the grid are assumed to be
3742	connected and computations along those edges will consider the values
3743	on the opposite side of the grid."""))
3744
3745	AddResultMetadata(InpaintedGrid.__init__, u'filledGrid',
3746	typeMetadata=ClassInstanceTypeMetadata(cls=Grid),
3747	description=_(u"""Filled grid."""))
3748
3749	# TODO: Add metadata
3750
3751	###############################################################################
3752	# Metadata: MosaickedTiles class
3753	###############################################################################
3754
3755	AddClassMetadata(MosaickedTiles,
3756	shortDescription=_(u'TODO: Add description.'))
3757
3758	# TODO: Add metadata
3759
3760	###############################################################################
3761	# Names exported by this module
3762	###############################################################################
3763
3764	__all__ = ['TimeSeriesGridStack',
3765	'GridSlice',
3766	'GridSliceCollection',
3767	'ClippedGrid',
3768	'MaskedGrid',
3769	'DerivedGrid',
3770	'RotatedGlobalGrid',
3771	'MemoryCachedGrid',
3772	'AggregateGrid',
3773	'ClimatologicalGridCollection',
3774	'SeafloorGrid',
3775	'SpatiotemporalAggregateGrid',
3776	'InpaintedGrid',
3777	'MosaickedTiles']

Note: See TracBrowser for help on using the browser.

Context Navigation

root/MGET/Branches/Jason/PythonPackage/src/GeoEco/Datasets/Virtual.py @ 939

Download in other formats: