Repository: pyvista/pyvista

Issue:
ImageData.points doesn't handle DirectionMatrix
### Describe the bug, what's wrong, and what you expected.

It seems that the ImageData.points property is not handling the image DirectionMatrix correctly. In fact, it has been ignored.

### Steps to reproduce the bug.

```python

import pyvista as pv
import numpy as np

img = pv.ImageData(dimensions=(2,3,4), spacing=(0.5,1.0,1.5), origin=(0., 0.5, 1.0))
img.SetDirectionMatrix((-1, 0, 0, 0, 1, 0, 0, 0, 1))

img.points
#array([[0. , 0.5, 1. ],
#       [0.5, 0.5, 1. ],
#       [0. , 1.5, 1. ],
#       [0.5, 1.5, 1. ],
#       [0. , 2.5, 1. ],
#       [0.5, 2.5, 1. ],
#       [0. , 0.5, 2.5],
#       [0.5, 0.5, 2.5],
#       [0. , 1.5, 2.5],
#       [0.5, 1.5, 2.5],
#       [0. , 2.5, 2.5],
#       [0.5, 2.5, 2.5],
#       [0. , 0.5, 4. ],
#       [0.5, 0.5, 4. ],
#       [0. , 1.5, 4. ],
#       [0.5, 1.5, 4. ],
#       [0. , 2.5, 4. ],
#       [0.5, 2.5, 4. ],
#       [0. , 0.5, 5.5],
#       [0.5, 0.5, 5.5],
#       [0. , 1.5, 5.5],
#       [0.5, 1.5, 5.5],
#       [0. , 2.5, 5.5],
#       [0.5, 2.5, 5.5]])

np.array([img.GetPoint(i) for i in range(img.number_of_points)])
#array([[ 0. ,  0.5,  1. ],
#       [-0.5,  0.5,  1. ],
#       [ 0. ,  1.5,  1. ],
#       [-0.5,  1.5,  1. ],
#       [ 0. ,  2.5,  1. ],
#       [-0.5,  2.5,  1. ],
#       [ 0. ,  0.5,  2.5],
#       [-0.5,  0.5,  2.5],
#       [ 0. ,  1.5,  2.5],
#       [-0.5,  1.5,  2.5],
#       [ 0. ,  2.5,  2.5],
#       [-0.5,  2.5,  2.5],
#       [ 0. ,  0.5,  4. ],
#       [-0.5,  0.5,  4. ],
#       [ 0. ,  1.5,  4. ],
#       [-0.5,  1.5,  4. ],
#       [ 0. ,  2.5,  4. ],
#       [-0.5,  2.5,  4. ],
#       [ 0. ,  0.5,  5.5],
#       [-0.5,  0.5,  5.5],
#       [ 0. ,  1.5,  5.5],
#       [-0.5,  1.5,  5.5],
#       [ 0. ,  2.5,  5.5],
#       [-0.5,  2.5,  5.5]])

```

### System Information

```shell
Date: Thu Oct 12 14:37:23 2023 CEST

                OS : Linux
            CPU(s) : 8
           Machine : x86_64
      Architecture : 64bit
               RAM : 15.0 GiB
       Environment : Python
       File system : ext4
        GPU Vendor : AMD
      GPU Renderer : RENOIR (renoir, LLVM 15.0.7, DRM 3.49, 6.2.0-34-generic)
       GPU Version : 4.6 (Core Profile) Mesa 23.0.4-0ubuntu1~22.04.1
  MathText Support : False

  Python 3.9.13 (main, May 30 2022, 12:53:10)  [GCC 11.2.0]

           pyvista : 0.42.3
               vtk : 9.2.6
             numpy : 1.23.4
        matplotlib : 3.5.2
            scooby : 0.5.12
             pooch : v1.6.0
            pillow : Module not found
           imageio : 2.19.3
           IPython : 8.4.0
        ipywidgets : 7.7.0
             scipy : 1.10.1
              tqdm : 4.64.0
            meshio : 5.3.4
        jupyterlab : 3.4.2
             trame : Version unknown
      trame_client : Version unknown
      trame_server : 2.9.1
         trame_vtk : Version unknown
     trame_vuetify : Version unknown
      nest_asyncio : Version unknown
```


### Screenshots

_No response_

Discussion:
pyvista-bot: <!-- NETLIFY DEPLOY COMMENT GENERATED BY ACTIONS_NETLIFY - APP ID SHA256: d3574f413fe51079937885070831e027ca5b10f7facb9e609494e0a481a31f79 -->
🚀 Deployed on https://66d9273dfd91209d11f45500--pyvista-dev.netlify.app
codecov[bot]: ## [Codecov](https://app.codecov.io/gh/pyvista/pyvista/pull/6612?dropdown=coverage&src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pyvista) Report
All modified and coverable lines are covered by tests :white_check_mark:
> Project coverage is 97.14%. Comparing base [(`6647524`)](https://app.codecov.io/gh/pyvista/pyvista/commit/664752494d7ba82a1d79fae72c2a7e1515c95e7c?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pyvista) to head [(`1e4f3c1`)](https://app.codecov.io/gh/pyvista/pyvista/commit/1e4f3c159512fc9d492b1cd011a2625073f33b9d?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pyvista).

<details><summary>Additional details and impacted files</summary>


```diff
@@            Coverage Diff             @@
##             main    #6612      +/-   ##
==========================================
- Coverage   97.37%   97.14%   -0.24%     
==========================================
  Files         143      143              
  Lines       27786    27798      +12     
==========================================
- Hits        27058    27004      -54     
- Misses        728      794      +66     
```

</details>
tkoyama010: pre-commit.ci autofix
user27182: pre-commit.ci autofix

Patch:
diff --git a/pyvista/core/grid.py b/pyvista/core/grid.py
index ddd3c2fe..a65b6e51 100644
--- a/pyvista/core/grid.py
+++ b/pyvista/core/grid.py
@@ -16,6 +16,8 @@ import pyvista
 if TYPE_CHECKING:
     from pyvista.core._typing_core import NumpyArray
 
+from pyvista.core import _validation
+
 from . import _vtk_core as _vtk
 from .dataset import DataSet
 from .filters import ImageDataFilters
@@ -24,6 +26,7 @@ from .filters import _get_output
 from .utilities.arrays import array_from_vtkmatrix
 from .utilities.arrays import convert_array
 from .utilities.arrays import raise_has_duplicates
+from .utilities.arrays import vtkmatrix_from_array
 from .utilities.misc import abstract_class
 
 
@@ -670,7 +673,14 @@ class ImageData(Grid, ImageDataFilters, _vtk.vtkImageData):
         y = np.insert(np.cumsum(np.full(ny, dy)), 0, 0.0) + oy
         z = np.insert(np.cumsum(np.full(nz, dz)), 0, 0.0) + oz
         xx, yy, zz = np.meshgrid(x, y, z, indexing='ij')
-        return np.c_[xx.ravel(order='F'), yy.ravel(order='F'), zz.ravel(order='F')]
+        points = np.c_[xx.ravel(order='F'), yy.ravel(order='F'), zz.ravel(order='F')]
+
+        direction = self.direction_matrix
+        if not np.array_equal(direction, np.eye(3)):
+            return (
+                pyvista.Transform().rotate(direction, point=self.origin).apply(points, copy=False)
+            )
+        return points
 
     @points.setter
     def points(self, points):  # numpydoc ignore=PR01
@@ -889,6 +899,24 @@ class ImageData(Grid, ImageDataFilters, _vtk.vtkImageData):
         """Cast to a rectangular grid and then convert to tetrahedra."""
         return self.cast_to_rectilinear_grid().to_tetrahedra(*args, **kwargs)
 
+    @property
+    def direction_matrix(self):
+        """Set or get the direction matrix.
+
+        The direction matrix is a 3x3 matrix which controls the orientation of the
+        image data.
+
+        Returns
+        -------
+        np.ndarray
+            Direction matrix as a 3x3 NumPy array.
+        """
+        return array_from_vtkmatrix(self.GetDirectionMatrix())
+
+    @direction_matrix.setter
+    def direction_matrix(self, matrix):  # numpydoc ignore: GL08
+        self.SetDirectionMatrix(vtkmatrix_from_array(_validation.validate_transform3x3(matrix)))
+
     @property
     def index_to_physical_matrix(self) -> NumpyArray[float]:
         """Get 4x4 matrix to convert coordinates from index space (ijk) to physical space (xyz).


Test patch:
diff --git a/tests/core/test_grid.py b/tests/core/test_grid.py
index 556c4eab..44620a51 100644
--- a/tests/core/test_grid.py
+++ b/tests/core/test_grid.py
@@ -1093,13 +1093,45 @@ def test_grid_points():
     )
 
 
+def test_imagedata_direction_matrix():
+    # Create image data with a single voxel cell
+    image = pv.ImageData(dimensions=(2, 2, 2))
+    assert image.n_points == 8
+    assert image.n_cells == 1
+
+    initial_bounds = (0.0, 1.0, 0.0, 1.0, 0.0, 1.0)
+    assert image.bounds == initial_bounds
+
+    # Test set/get
+    expected_matrix = pv.Transform().rotate_vector((1, 2, 3), 30).matrix[:3, :3]
+    image.direction_matrix = expected_matrix
+    assert np.array_equal(image.direction_matrix, expected_matrix)
+
+    # Test bounds using a transformed reference box
+    box = pv.Box(bounds=initial_bounds)
+    box.transform(image.index_to_physical_matrix)
+    expected_bounds = box.bounds
+    assert np.allclose(image.bounds, expected_bounds)
+
+    # Check that filters make use of the direction matrix internally
+    image['data'] = np.ones((image.n_points,))
+    filtered = image.threshold()
+    assert filtered.bounds == expected_bounds
+
+    # Check that points make use of the direction matrix internally
+    poly_points = pv.PolyData(image.points)
+    assert np.allclose(poly_points.bounds, expected_bounds)
+
+
 def test_imagedata_index_to_physical_matrix():
     # Create image with arbitrary translation (origin) and rotation (direction)
     image = pv.ImageData()
     vector = (1, 2, 3)
+    rotation = pv.Transform().rotate_vector(vector, 30).matrix[:3, :3]
     image.origin = vector
+    image.direction_matrix = rotation
 
-    expected_transform = pv.Transform().translate(vector)
+    expected_transform = pv.Transform().rotate(rotation).translate(vector)
     ijk_to_xyz = image.index_to_physical_matrix
     assert np.allclose(ijk_to_xyz, expected_transform.matrix)