selections.py 14 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440
  1. # This file is part of h5py, a Python interface to the HDF5 library.
  2. #
  3. # http://www.h5py.org
  4. #
  5. # Copyright 2008-2013 Andrew Collette and contributors
  6. #
  7. # License: Standard 3-clause BSD; see "license.txt" for full license terms
  8. # and contributor agreement.
  9. """
  10. High-level access to HDF5 dataspace selections
  11. """
  12. import numpy as np
  13. from .base import product
  14. from .. import h5s, h5r, _selector
  15. def select(shape, args, dataset=None):
  16. """ High-level routine to generate a selection from arbitrary arguments
  17. to __getitem__. The arguments should be the following:
  18. shape
  19. Shape of the "source" dataspace.
  20. args
  21. Either a single argument or a tuple of arguments. See below for
  22. supported classes of argument.
  23. dataset
  24. A h5py.Dataset instance representing the source dataset.
  25. Argument classes:
  26. Single Selection instance
  27. Returns the argument.
  28. numpy.ndarray
  29. Must be a boolean mask. Returns a PointSelection instance.
  30. RegionReference
  31. Returns a Selection instance.
  32. Indices, slices, ellipses, MultiBlockSlices only
  33. Returns a SimpleSelection instance
  34. Indices, slices, ellipses, lists or boolean index arrays
  35. Returns a FancySelection instance.
  36. """
  37. if not isinstance(args, tuple):
  38. args = (args,)
  39. # "Special" indexing objects
  40. if len(args) == 1:
  41. arg = args[0]
  42. if isinstance(arg, Selection):
  43. if arg.shape != shape:
  44. raise TypeError("Mismatched selection shape")
  45. return arg
  46. elif isinstance(arg, np.ndarray) and arg.dtype.kind == 'b':
  47. if arg.shape == shape:
  48. return PointSelection.from_mask(arg)
  49. # Allow 1D boolean array on the 1st dim
  50. elif arg.shape != shape[:1]:
  51. raise TypeError("Boolean indexing array has incompatible shape")
  52. elif isinstance(arg, h5r.RegionReference):
  53. if dataset is None:
  54. raise TypeError("Cannot apply a region reference without a dataset")
  55. sid = h5r.get_region(arg, dataset.id)
  56. if shape != sid.shape:
  57. raise TypeError("Reference shape does not match dataset shape")
  58. return Selection(shape, spaceid=sid)
  59. if dataset is not None:
  60. selector = dataset._selector
  61. else:
  62. space = h5s.create_simple(shape)
  63. selector = _selector.Selector(space)
  64. return selector.make_selection(args)
  65. class Selection:
  66. """
  67. Base class for HDF5 dataspace selections. Subclasses support the
  68. "selection protocol", which means they have at least the following
  69. members:
  70. __init__(shape) => Create a new selection on "shape"-tuple
  71. __getitem__(args) => Perform a selection with the range specified.
  72. What args are allowed depends on the
  73. particular subclass in use.
  74. id (read-only) => h5py.h5s.SpaceID instance
  75. shape (read-only) => The shape of the dataspace.
  76. mshape (read-only) => The shape of the selection region.
  77. Not guaranteed to fit within "shape", although
  78. the total number of points is less than
  79. product(shape).
  80. nselect (read-only) => Number of selected points. Always equal to
  81. product(mshape).
  82. broadcast(target_shape) => Return an iterable which yields dataspaces
  83. for read, based on target_shape.
  84. The base class represents "unshaped" selections (1-D).
  85. """
  86. def __init__(self, shape, spaceid=None):
  87. """ Create a selection. Shape may be None if spaceid is given. """
  88. if spaceid is not None:
  89. self._id = spaceid
  90. self._shape = spaceid.shape
  91. else:
  92. shape = tuple(shape)
  93. self._shape = shape
  94. self._id = h5s.create_simple(shape, (h5s.UNLIMITED,)*len(shape))
  95. self._id.select_all()
  96. @property
  97. def id(self):
  98. """ SpaceID instance """
  99. return self._id
  100. @property
  101. def shape(self):
  102. """ Shape of whole dataspace """
  103. return self._shape
  104. @property
  105. def nselect(self):
  106. """ Number of elements currently selected """
  107. return self._id.get_select_npoints()
  108. @property
  109. def mshape(self):
  110. """ Shape of selection (always 1-D for this class) """
  111. return (self.nselect,)
  112. @property
  113. def array_shape(self):
  114. """Shape of array to read/write (always 1-D for this class)"""
  115. return self.mshape
  116. # expand_shape and broadcast only really make sense for SimpleSelection
  117. def expand_shape(self, source_shape):
  118. if product(source_shape) != self.nselect:
  119. raise TypeError("Broadcasting is not supported for point-wise selections")
  120. return source_shape
  121. def broadcast(self, source_shape):
  122. """ Get an iterable for broadcasting """
  123. if product(source_shape) != self.nselect:
  124. raise TypeError("Broadcasting is not supported for point-wise selections")
  125. yield self._id
  126. def __getitem__(self, args):
  127. raise NotImplementedError("This class does not support indexing")
  128. class PointSelection(Selection):
  129. """
  130. Represents a point-wise selection. You can supply sequences of
  131. points to the three methods append(), prepend() and set(), or
  132. instantiate it with a single boolean array using from_mask().
  133. """
  134. def __init__(self, shape, spaceid=None, points=None):
  135. super().__init__(shape, spaceid)
  136. if points is not None:
  137. self._perform_selection(points, h5s.SELECT_SET)
  138. def _perform_selection(self, points, op):
  139. """ Internal method which actually performs the selection """
  140. points = np.atleast_2d(np.asarray(points, order='C', dtype='u8'))
  141. if self._id.get_select_type() != h5s.SEL_POINTS:
  142. op = h5s.SELECT_SET
  143. if len(points) == 0:
  144. self._id.select_none()
  145. else:
  146. self._id.select_elements(points, op)
  147. @classmethod
  148. def from_mask(cls, mask, spaceid=None):
  149. """Create a point-wise selection from a NumPy boolean array """
  150. if not (isinstance(mask, np.ndarray) and mask.dtype.kind == 'b'):
  151. raise TypeError("PointSelection.from_mask only works with bool arrays")
  152. points = np.transpose(mask.nonzero())
  153. return cls(mask.shape, spaceid, points=points)
  154. def append(self, points):
  155. """ Add the sequence of points to the end of the current selection """
  156. self._perform_selection(points, h5s.SELECT_APPEND)
  157. def prepend(self, points):
  158. """ Add the sequence of points to the beginning of the current selection """
  159. self._perform_selection(points, h5s.SELECT_PREPEND)
  160. def set(self, points):
  161. """ Replace the current selection with the given sequence of points"""
  162. self._perform_selection(points, h5s.SELECT_SET)
  163. class SimpleSelection(Selection):
  164. """ A single "rectangular" (regular) selection composed of only slices
  165. and integer arguments. Can participate in broadcasting.
  166. """
  167. @property
  168. def mshape(self):
  169. """ Shape of current selection """
  170. return self._sel[1]
  171. @property
  172. def array_shape(self):
  173. scalar = self._sel[3]
  174. return tuple(x for x, s in zip(self.mshape, scalar, strict=True) if not s)
  175. def __init__(self, shape, spaceid=None, hyperslab=None):
  176. super().__init__(shape, spaceid)
  177. if hyperslab is not None:
  178. self._sel = hyperslab
  179. else:
  180. # No hyperslab specified - select all
  181. rank = len(self.shape)
  182. self._sel = ((0,)*rank, self.shape, (1,)*rank, (False,)*rank)
  183. def expand_shape(self, source_shape):
  184. """Match the dimensions of an array to be broadcast to the selection
  185. The returned shape describes an array of the same size as the input
  186. shape, but its dimensions
  187. E.g. with a dataset shape (10, 5, 4, 2), writing like this::
  188. ds[..., 0] = np.ones((5, 4))
  189. The source shape (5, 4) will expand to (1, 5, 4, 1).
  190. Then the broadcast method below repeats that chunk 10
  191. times to write to an effective shape of (10, 5, 4, 1).
  192. """
  193. start, count, step, scalar = self._sel
  194. rank = len(count)
  195. remaining_src_dims = list(source_shape)
  196. eshape = []
  197. for idx in range(1, rank + 1):
  198. if len(remaining_src_dims) == 0 or scalar[-idx]: # Skip scalar axes
  199. eshape.append(1)
  200. else:
  201. t = remaining_src_dims.pop()
  202. if t == 1 or count[-idx] == t:
  203. eshape.append(t)
  204. else:
  205. raise TypeError("Can't broadcast %s -> %s" % (source_shape, self.array_shape)) # array shape
  206. if any([n > 1 for n in remaining_src_dims]):
  207. # All dimensions from target_shape should either have been popped
  208. # to match the selection shape, or be 1.
  209. raise TypeError("Can't broadcast %s -> %s" % (source_shape, self.array_shape)) # array shape
  210. # We have built eshape backwards, so now reverse it
  211. return tuple(eshape[::-1])
  212. def broadcast(self, source_shape):
  213. """ Return an iterator over target dataspaces for broadcasting.
  214. Follows the standard NumPy broadcasting rules against the current
  215. selection shape (self.mshape).
  216. """
  217. if self.shape == ():
  218. if product(source_shape) != 1:
  219. raise TypeError("Can't broadcast %s to scalar" % source_shape)
  220. self._id.select_all()
  221. yield self._id
  222. return
  223. start, count, step, scalar = self._sel
  224. rank = len(count)
  225. tshape = self.expand_shape(source_shape)
  226. # Avoid ZeroDivisionError below (after the shape checks in expand_source)
  227. if any(d == 0 for d in count):
  228. return
  229. chunks = tuple(x//y for x, y in zip(count, tshape, strict=True))
  230. nchunks = product(chunks)
  231. if nchunks == 1:
  232. yield self._id
  233. else:
  234. sid = self._id.copy()
  235. sid.select_hyperslab((0,)*rank, tshape, step)
  236. for idx in range(nchunks):
  237. offset = tuple(x*y*z + s for x, y, z, s in zip(np.unravel_index(idx, chunks), tshape, step, start, strict=True))
  238. sid.offset_simple(offset)
  239. yield sid
  240. class FancySelection(Selection):
  241. """
  242. Implements advanced NumPy-style selection operations in addition to
  243. the standard slice-and-int behavior.
  244. Indexing arguments may be ints, slices, lists of indices, or
  245. per-axis (1D) boolean arrays.
  246. Broadcasting is not supported for these selections.
  247. """
  248. @property
  249. def mshape(self):
  250. return self._mshape
  251. @property
  252. def array_shape(self):
  253. return self._array_shape
  254. def __init__(self, shape, spaceid=None, mshape=None, array_shape=None):
  255. super().__init__(shape, spaceid)
  256. if mshape is None:
  257. mshape = self.shape
  258. if array_shape is None:
  259. array_shape = mshape
  260. self._mshape = mshape
  261. self._array_shape = array_shape
  262. def expand_shape(self, source_shape):
  263. if not source_shape == self.array_shape:
  264. raise TypeError("Broadcasting is not supported for complex selections")
  265. return source_shape
  266. def broadcast(self, source_shape):
  267. if not source_shape == self.array_shape:
  268. raise TypeError("Broadcasting is not supported for complex selections")
  269. yield self._id
  270. def guess_shape(sid):
  271. """ Given a dataspace, try to deduce the shape of the selection.
  272. Returns one of:
  273. * A tuple with the selection shape, same length as the dataspace
  274. * A 1D selection shape for point-based and multiple-hyperslab selections
  275. * None, for unselected scalars and for NULL dataspaces
  276. """
  277. sel_class = sid.get_simple_extent_type() # Dataspace class
  278. sel_type = sid.get_select_type() # Flavor of selection in use
  279. if sel_class == h5s.NULL:
  280. # NULL dataspaces don't support selections
  281. return None
  282. elif sel_class == h5s.SCALAR:
  283. # NumPy has no way of expressing empty 0-rank selections, so we use None
  284. if sel_type == h5s.SEL_NONE: return None
  285. if sel_type == h5s.SEL_ALL: return tuple()
  286. elif sel_class != h5s.SIMPLE:
  287. raise TypeError("Unrecognized dataspace class %s" % sel_class)
  288. # We have a "simple" (rank >= 1) dataspace
  289. N = sid.get_select_npoints()
  290. rank = len(sid.shape)
  291. if sel_type == h5s.SEL_NONE:
  292. return (0,)*rank
  293. elif sel_type == h5s.SEL_ALL:
  294. return sid.shape
  295. elif sel_type == h5s.SEL_POINTS:
  296. # Like NumPy, point-based selections yield 1D arrays regardless of
  297. # the dataspace rank
  298. return (N,)
  299. elif sel_type != h5s.SEL_HYPERSLABS:
  300. raise TypeError("Unrecognized selection method %s" % sel_type)
  301. # We have a hyperslab-based selection
  302. if N == 0:
  303. return (0,)*rank
  304. bottomcorner, topcorner = (np.array(x) for x in sid.get_select_bounds())
  305. # Shape of full selection box
  306. boxshape = topcorner - bottomcorner + np.ones((rank,))
  307. def get_n_axis(sid, axis):
  308. """ Determine the number of elements selected along a particular axis.
  309. To do this, we "mask off" the axis by making a hyperslab selection
  310. which leaves only the first point along the axis. For a 2D dataset
  311. with selection box shape (X, Y), for axis 1, this would leave a
  312. selection of shape (X, 1). We count the number of points N_leftover
  313. remaining in the selection and compute the axis selection length by
  314. N_axis = N/N_leftover.
  315. """
  316. if(boxshape[axis]) == 1:
  317. return 1
  318. start = bottomcorner.copy()
  319. start[axis] += 1
  320. count = boxshape.copy()
  321. count[axis] -= 1
  322. # Throw away all points along this axis
  323. masked_sid = sid.copy()
  324. masked_sid.select_hyperslab(tuple(start), tuple(count), op=h5s.SELECT_NOTB)
  325. N_leftover = masked_sid.get_select_npoints()
  326. return N//N_leftover
  327. shape = tuple(get_n_axis(sid, x) for x in range(rank))
  328. if product(shape) != N:
  329. # This means multiple hyperslab selections are in effect,
  330. # so we fall back to a 1D shape
  331. return (N,)
  332. return shape