.

bigarray/__init__.py

@@ -38,8 +38,8 @@ of physical RAM.
 
 from __future__ import print_function
 from wendelin.lib.calc import mul
+from wendelin.lib.xnumpy import _as_strided
 from numpy import ndarray, dtype, sign, newaxis, asarray, argmax, uint8
-from numpy.lib.stride_tricks import DummyArray
 import logging
 
 
@@ -618,29 +618,7 @@ class ArrayRef(object):
         #
         # it is also safe because we checked .shape and .stridev not to escape
         # from bchild data buffer.
-        #
-        # the code below is very close to
-        #
-        #   a = stride_tricks.as_strided(bchild_z0, shape=self.shape, strides=self.stridev)
-        #
-        # but we don't use as_strided() because we also have to change dtype
-        # with shape and strides in one go - else changing dtype after either
-        # via a.dtype = ..., or via a.view(dtype=...) can raise errors like
-        #
-        #   "When changing to a larger dtype, its size must be a
-        #    divisor of the total size in bytes of the last axis
-        #    of the array."
-        aiface = dict(bchild_z0.__array_interface__)
-        aiface['shape']   = tuple(self.shape)
-        aiface['strides'] = tuple(self.stridev)
-        # type: for now we only care that itemsize is the same
-        aiface['typestr'] = '|V%d' % self.dtype.itemsize
-        aiface['descr']   = [('', aiface['typestr'])]
-
-        a = asarray(DummyArray(aiface, base=bchild_z0))
-
-        # restore full dtype - it should not raise here, since itemsize is the same
-        a.dtype = self.dtype
+        a = _as_strided(bchild_z0, tuple(self.shape), tuple(self.stridev), self.dtype)
 
         # restore full array type
         a = a.view(type=self.atype)

lib/xnumpy.py

@@ -23,45 +23,80 @@ from numpy.lib import stride_tricks as npst
 
 # XXX move ArrayRef here.
 
+# _as_strided is similar to numpy.lib.stride_tricks.as_strided, but allows to
+# set all shape/stridev/dtype in one go.
+#
+# It must be used with extreme care, because if there is math error in the
+# arguments, the resulting array can cover wrong memory. Bugs here thus can
+# lead to mysterious crashes.
+def _as_strided(a, shape, stridev, dtype):
+    # the code below is very close to
+    #
+    #   a = stride_tricks.as_strided(a, shape=shape, strides=stridev)
+    #
+    # but we don't use as_strided() because we also have to change dtype
+    # with shape and strides in one go - else changing dtype after either
+    # via a.dtype = ..., or via a.view(dtype=...) can raise errors like
+    #
+    #   "When changing to a larger dtype, its size must be a
+    #    divisor of the total size in bytes of the last axis
+    #    of the array."
+    aiface = dict(a.__array_interface__)
+    aiface['shape']   = shape
+    aiface['strides'] = stridev
+    # type: for now we only care that itemsize is the same
+    aiface['typestr'] = '|V%d' % dtype.itemsize
+    aiface['descr']   = [('', aiface['typestr'])]
+
+    a = np.asarray(npst.DummyArray(aiface, base=a))
+
+    # restore full dtype - it should not raise here, since itemsize is the same
+    a.dtype = dtype
+
+    # XXX restore full array type?
+
+    # we are done
+    return a
+
 # restructure creates view of the array interpreting its minor axis as fully covered by dtype.
 #
-# The minor axis of the array must be C-contiguous and be fully covered by dtype in size.
+# The minor axis of the array must be C-contiguous and be of dtype.itemsize in size.
 #
 # Restructure is similar to arr.view(dtype) + corresponding reshape, but does
 # not have limitations of ndarray.view(). For example:
 #
-# In [1]: a = np.arange(3*3, dtype=np.int32).reshape((3,3))
+#   In [1]: a = np.arange(3*3, dtype=np.int32).reshape((3,3))
 #
-# In [2]: a
-# Out[2]:
-# array([[0, 1, 2],
-#        [3, 4, 5],
-#        [6, 7, 8]], dtype=int32)
+#   In [2]: a
+#   Out[2]:
+#   array([[0, 1, 2],
+#          [3, 4, 5],
+#          [6, 7, 8]], dtype=int32)
 #
-# In [3]: b = a[:2,:2]
+#   In [3]: b = a[:2,:2]
 #
-# In [4]: b
-# Out[4]:
-# array([[0, 1],
-#        [3, 4]], dtype=int32)
+#   In [4]: b
+#   Out[4]:
+#   array([[0, 1],
+#          [3, 4]], dtype=int32)
 #
-# In [5]: dtxy = np.dtype([('x', np.int32), ('y', np.int32)])
+#   In [5]: dtxy = np.dtype([('x', np.int32), ('y', np.int32)])
 #
-# In [6]: dtxy
-# Out[6]: dtype([('x', '<i4'), ('y', '<i4')])
+#   In [6]: dtxy
+#   Out[6]: dtype([('x', '<i4'), ('y', '<i4')])
 
-# In [7]: b.view(dtxy)
-# ---------------------------------------------------------------------------
-# ValueError                                Traceback (most recent call last)
-# <ipython-input-66-af98529aa150> in <module>()
-# ----> 1 b.view(dtxy)
+#   In [7]: b.view(dtxy)
+#   ---------------------------------------------------------------------------
+#   ValueError                                Traceback (most recent call last)
+#   <ipython-input-66-af98529aa150> in <module>()
+#   ----> 1 b.view(dtxy)
 #
-# ValueError: To change to a dtype of a different size, the array must be C-contiguous
+#   ValueError: To change to a dtype of a different size, the array must be C-contiguous
 #
-# In [8]: restructure(b, dtxy)
-# Out[8]: array([(0, 1), (3, 4)], dtype=[('x', '<i4'), ('y', '<i4')])
+#   In [8]: restructure(b, dtxy)
+#   Out[8]: array([(0, 1), (3, 4)], dtype=[('x', '<i4'), ('y', '<i4')])
 #
-# restructure always creates view and never copies data.
+# Restructure always creates view and never copies data.
 def restructure(arr, dtype):
     dtype = np.dtype(dtype) # convenience
 
@@ -84,22 +119,10 @@ def restructure(arr, dtype):
     shape   = arr.shape[:maxis] + arr.shape[maxis+1:]
     stridev = arr.strides[:maxis] + arr.strides[maxis+1:]
 
-    # NOTE cannot use just np.ndarray because if arr is a slice it can give:
+    # NOTE we cannot use just np.ndarray because if arr is a slice it can give:
     #   TypeError: expected a single-segment buffer object
     #return np.ndarray.__new__(type(arr), shape, dtype, buffer(arr), 0, stridev)
-
-    # XXX dup from ArrayRef
-    aiface = dict(arr.__array_interface__)
-    aiface['shape']   = shape
-    aiface['strides'] = stridev
-    # type: for now we only care that itemsize is the same
-    aiface['typestr'] = '|V%d' % dtype.itemsize
-    aiface['descr']   = [('', aiface['typestr'])]
-
-    a = np.asarray(npst.DummyArray(aiface, base=arr))
-
-    # restore full dtype - it should not raise here, since itemsize is the same
-    a.dtype = dtype
+    a = _as_strided(arr, shape, stridev, dtype)
 
     # restore full array type
     a = a.view(type=type(arr))