Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

2d slices in 3d plot #3919

Closed
twmr opened this issue Dec 13, 2014 · 21 comments
Closed

2d slices in 3d plot #3919

twmr opened this issue Dec 13, 2014 · 21 comments

Comments

@twmr
Copy link
Contributor

twmr commented Dec 13, 2014

Is it possible to create plots like the one below with matplotlib? If not, would it be difficult to add support for such plots in the mplot3d backend?

slice

I know that such plots can be created with VTK but VTK is still limited to python 2 👎

@WeatherGod
Copy link
Member

Slices, yes (sort of, with contourf(), not imshow()), but they can not
arbitrarially intersect. Because mplot3d is not a true 3d plotting library,
it can not handle intersecting artists, or even proper layering of separate
artists that have intersecting bounding boxes. If you want intersecting
artists, you would have to contourf() each quadrant separately, and even
then you may get rendering artifacts at different viewing angles.

They haven't ported mayavi yet?! Another project to take a look at for true
3d rendering is glumpy. Inspired from the matplotlib design, it has an
opengl rendering core which allows it to perform 3d plotting.

I hope that helps!

On Sat, Dec 13, 2014 at 2:18 PM, Thomas Hisch notifications@github.com
wrote:

Is it possible to create plots like the one below with matplotlib? If not,
would it be difficult to add support for such plots in the mplot3d backend?

[image: slice]
https://camo.githubusercontent.com/312d2210de78d9b491ddb0b26fc5185976ea0e9b/687474703a2f2f7777772e646172746d6f7574682e6564752f25374572632f636c61737365732f76697369742f74687265655f736c6963652e676966

I know that such plots can be created with VTK but VTK is still limited to
python 2 [image: 👎]


Reply to this email directly or view it on GitHub
#3919.

@twmr
Copy link
Contributor Author

twmr commented Dec 13, 2014

Do you have an example, which shows how to plot data on a plane (x,y,z + color) in a Axes3d ?

mayavi depends on vtk, therefore it is also py2 only.

@WeatherGod
Copy link
Member

You would use the "offset" and the "zdir" arguments to the contourf() call
on a 3d axes:
http://matplotlib.org/mpl_toolkits/mplot3d/tutorial.html#filled-contour-plots

Keep in mind, those examples are surface data rather than a 3D image data,
which is what you would be handling. So, you would pass the coordinate
information for the first two arguments and then the 2D slice of the 3D
image data for the third argument. The zdir argument would indicate the
normal vector for the image slice, and the offset argument would specify
where along that normal vector the image slice should be placed.

I repeat the importance of needing to split up the image slices. You would
need to break up what would be 3 contourf() calls into 12 contourf() (4
quadrants for each slice). Luckily, numpy makes slicing and dicing your
data relatively easy...

On Sat, Dec 13, 2014 at 4:46 PM, Thomas Hisch notifications@github.com
wrote:

Do you have an example, which shows how to plot data on a plane (x,y,z +
color) in a Axes3d ?

mayavi depends on vtk, therefore it is also py2 only.


Reply to this email directly or view it on GitHub
#3919 (comment)
.

@twmr
Copy link
Contributor Author

twmr commented Dec 13, 2014

I see. Does this also work if I rotate the Axes interactively?

@twmr
Copy link
Contributor Author

twmr commented Dec 13, 2014

How do I set the zorder of the planes in the code below?

yplane=0.8
ind = Y[:,0] < yplane
ax.contourf(X[ind,:], Y[ind,:], Z[ind,:], zdir='z', offset=0.3)
ax.contourf(X, Z, Y, zdir='y', offset=yplane)

zorder=xx seems to have no effect.

edit:
https://stackoverflow.com/questions/14824893/how-to-draw-diagrams-like-this/14825951#14825951

@WeatherGod
Copy link
Member

It should work with rotation just fine.

You don't mess around with zorder. mplot3d computes it dynamically to help
compose the 3d scene.

On Sat, Dec 13, 2014 at 5:43 PM, Thomas Hisch notifications@github.com
wrote:

How do I set the zorder of the planes in the code below?

yplane=0.8
ind = Y[:,0] < yplane
ax.contourf(X[ind,:], Y[ind,:], Z[ind,:], zdir='z', offset=0.3)
ax.contourf(X, Z, Y, zdir='y', offset=yplane)

zorder=xx seems to have no affect.


Reply to this email directly or view it on GitHub
#3919 (comment)
.

@petehuang
Copy link
Contributor

Seems like we've resolved this @WeatherGod?

@twmr
Copy link
Contributor Author

twmr commented Jan 6, 2017

The last time I tested it, it did not work satisfactorily AFAIR.

@WeatherGod
Copy link
Member

WeatherGod commented Jan 6, 2017 via email

@louity
Copy link

louity commented Jul 12, 2017

As explained in this stackoverflow answer, matplotlib won't be able to make the 3 slices plots on the same figure "until OpenGL support is added to all of the backends".

Indeed, occlusions of the different slices will not be rendered correctly, and you will get this kind of result (I plot the 3d function (x,y,z) -> x^2 + y^2 + z^2 over the domain [-1,1]^3) :
slice_plot_xyz

However, you can plot each slice on a separate figure :
slice_plot_x
slice_plot_y
slice_plot_z

or plot slices that do not occlude each other from the point of view you are looking from.

code :

import matplotlib.pyplot as plt
import numpy as np
from mpl_toolkits.mplot3d import Axes3D

def plot_3D_array_slices(array):
    min_val = array.min()
    max_val = array.max()
    n_x, n_y, n_z = array.shape
    colormap = plt.cm.YlOrRd

    x_cut = array[n_x//2,:,:]
    Y, Z = np.mgrid[0:n_y, 0:n_z]
    X = n_x//2 * np.ones((n_y, n_z))
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.plot_surface(X, Y, Z, rstride=1, cstride=1, facecolors=colormap((x_cut-min_val)/(max_val-min_val)), shade=False)
    ax.set_title("x slice")

    y_cut = array[:,n_y//2,:]
    X, Z = np.mgrid[0:n_x, 0:n_z]
    Y = n_y//2 * np.ones((n_x, n_z))
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.plot_surface(X, Y, Z, rstride=1, cstride=1, facecolors=colormap((y_cut-min_val)/(max_val-min_val)), shade=False)
    ax.set_title("y slice")

    z_cut = array[:,:,n_z//2]
    X, Y = np.mgrid[0:n_x, 0:n_y]
    Z = n_z//2 * np.ones((n_x, n_y))
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.plot_surface(X, Y, Z, rstride=1, cstride=1, facecolors=colormap((z_cut-min_val)/(max_val-min_val)), shade=False)
    ax.set_title("z slice")

    plt.show()


n_pts = 100
r_square = (np.mgrid[-1:1:1j*n_pts, -1:1:1j*n_pts, -1:1:1j*n_pts]**2).sum(0)
plot_3D_array_slices(r_square)

@louity
Copy link

louity commented Jul 12, 2017

You might actually overcome this issue using matplotlib Poly3DCollection, see example, but as you can read in the doc :"this class (Poly3DCollection) does a bit of magic with the _facecolors and _edgecolors properties."... Let me know if you manage to do it :-).

@louity
Copy link

louity commented Jul 12, 2017

There is an example of solution using plotly.

@tacaswell
Copy link
Member

Also see https://stackoverflow.com/questions/14824893/how-to-draw-intersecting-planes/14825951#14825951 . Manually chopping your slices up can also be made to work.

@jiadongdan
Copy link

Is it type of plot still possible? As of matplotlib v3, it seems there no proper way to handle this

@WeatherGod
Copy link
Member

WeatherGod commented Jul 30, 2019 via email

@github-actions
Copy link

github-actions bot commented Mar 9, 2023

This issue has been marked "inactive" because it has been 365 days since the last comment. If this issue is still present in recent Matplotlib releases, or the feature request is still wanted, please leave a comment and this label will be removed. If there are no updates in another 30 days, this issue will be automatically closed, but you are free to re-open or create a new issue if needed. We value issue reports, and this procedure is meant to help us resurface and prioritize issues that have not been addressed yet, not make them disappear. Thanks for your help!

@github-actions github-actions bot added the status: inactive Marked by the “Stale” Github Action label Mar 9, 2023
@github-actions github-actions bot closed this as not planned Won't fix, can't repro, duplicate, stale Apr 9, 2023
@KrisSatie
Copy link

Is it type of plot still possible? As of matplotlib v3, it seems there no proper way to handle this

@rcomer rcomer added the status: closed as inactive Issues closed by the "Stale" Github Action. Please comment on any you think should still be open. label May 30, 2023
@dongjiaping
Copy link

yes, we need a better way to solve make 3 slices plots on the same figure.

@scottshambaugh scottshambaugh added status: won't or can't fix and removed status: inactive Marked by the “Stale” Github Action status: closed as inactive Issues closed by the "Stale" Github Action. Please comment on any you think should still be open. labels May 16, 2024
@scottshambaugh
Copy link
Contributor

Relabeling the reason for closure as "can't fix", as this is a renderer issue. However, plotting slices can be done by plotting each quadrant of the slice around the intersection separately. Here's an example (which could probably be made more compact but demonstrates the concept):

import matplotlib.pyplot as plt
import numpy as np

def plot_quadrants(ax, slice_data, fixed_coord, coord_val, colormap, shape, min_val, max_val):
    half_shape = (shape[1] // 2, shape[2] // 2)
    quadrants = [
        slice_data[:half_shape[0], :half_shape[1]],
        slice_data[:half_shape[0], half_shape[1]:],
        slice_data[half_shape[0]:, :half_shape[1]],
        slice_data[half_shape[0]:, half_shape[1]:]
    ]
    
    for i, quad in enumerate(quadrants):
        if fixed_coord == 'x':
            Y, Z = np.mgrid[0:half_shape[0], 0:half_shape[1]]
            X = coord_val * np.ones_like(Y)
            Y_offset = (i // 2) * half_shape[0]
            Z_offset = (i % 2) * half_shape[1]
            ax.plot_surface(X, Y + Y_offset, Z + Z_offset, rstride=1, cstride=1, facecolors=colormap((quad-min_val)/(max_val-min_val)), shade=False)
        elif fixed_coord == 'y':
            X, Z = np.mgrid[0:half_shape[0], 0:half_shape[1]]
            Y = coord_val * np.ones_like(X)
            X_offset = (i // 2) * half_shape[0]
            Z_offset = (i % 2) * half_shape[1]
            ax.plot_surface(X + X_offset, Y, Z + Z_offset, rstride=1, cstride=1, facecolors=colormap((quad-min_val)/(max_val-min_val)), shade=False)
        elif fixed_coord == 'z':
            X, Y = np.mgrid[0:half_shape[0], 0:half_shape[1]]
            Z = coord_val * np.ones_like(X)
            X_offset = (i // 2) * half_shape[0]
            Y_offset = (i % 2) * half_shape[1]
            ax.plot_surface(X + X_offset, Y + Y_offset, Z, rstride=1, cstride=1, facecolors=colormap((quad-min_val)/(max_val-min_val)), shade=False)

def plot_3D_array_slices(array):
    colormap = plt.cm.YlOrRd
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')

    min_val = array.min()
    max_val = array.max()
    
    x_slice = array[array.shape[0]//2, :, :]
    y_slice = array[:, array.shape[1]//2, :]
    z_slice = array[:, :, array.shape[2]//2]
    
    plot_quadrants(ax, x_slice, 'x', array.shape[0]//2, colormap, array.shape, min_val, max_val)
    plot_quadrants(ax, y_slice, 'y', array.shape[1]//2, colormap, array.shape, min_val, max_val)
    plot_quadrants(ax, z_slice, 'z', array.shape[2]//2, colormap, array.shape, min_val, max_val)
    
    plt.show()

n_pts = 100
r_square = (np.mgrid[-1:1:1j*n_pts, -1:1:1j*n_pts, -1:1:1j*n_pts]**2).sum(0)
plot_3D_array_slices(r_square)

image

@timhoffm
Copy link
Member

Thanks @scottshambaugh! I think this would be good as a gallery example. There were some indices mixed up (noticed when using different shapes for the dimensions). Below is a corrected version, which also simplifies the plot_quadrants interface:

import matplotlib.pyplot as plt
import numpy as np


def plot_quadrants(ax, array, fixed_coord, cmap):
    """For a given 3d *array* plot a plane with *fixed_coord*, using four individual quadrants."""
    nx, ny, nz = array.shape
    index = {
        'x': (nx//2, slice(None), slice(None)),
        'y': (slice(None), ny//2, slice(None)),
        'z': (slice(None), slice(None), nz//2),
    }[fixed_coord]
    plane_data = array[index]

    n0, n1 = plane_data.shape
    quadrants = [
        plane_data[:n0//2, :n1//2],
        plane_data[:n0//2, n1//2:],
        plane_data[n0//2:, :n1//2],
        plane_data[n0//2:, n1//2:]
    ]

    min_val = array.min()
    max_val = array.max()
    
    cmap = plt.get_cmap(cmap)
    
    for i, quadrant in enumerate(quadrants):
        facecolors = cmap((quadrant - min_val) / (max_val-min_val))
        if fixed_coord == 'x':
            Y, Z = np.mgrid[0:ny//2, 0:nz//2]
            X = nx//2 * np.ones_like(Y)
            Y_offset = (i // 2) * ny//2
            Z_offset = (i % 2) * nz//2
            ax.plot_surface(X, Y + Y_offset, Z + Z_offset, rstride=1, cstride=1, facecolors=facecolors, shade=False)
        elif fixed_coord == 'y':
            X, Z = np.mgrid[0:nx//2, 0:nz//2]
            Y = ny//2 * np.ones_like(X)
            X_offset = (i // 2) * nx//2
            Z_offset = (i % 2) * nz//2
            ax.plot_surface(X + X_offset, Y, Z + Z_offset, rstride=1, cstride=1, facecolors=facecolors, shade=False)
        elif fixed_coord == 'z':
            X, Y = np.mgrid[0:nx//2, 0:ny//2]
            Z = nz//2 * np.ones_like(X)
            X_offset = (i // 2) * nx//2
            Y_offset = (i % 2) * ny//2
            ax.plot_surface(X + X_offset, Y + Y_offset, Z, rstride=1, cstride=1, facecolors=facecolors, shade=False)


def figure_3D_array_slices(array, cmap=None):
    """Plot a 3d array using three intersecting centered planes."""
    fig = plt.figure()
    ax = fig.add_subplot(projection='3d')
    ax.set_box_aspect(array.shape)
    plot_quadrants(ax, array, 'x', cmap=cmap)
    plot_quadrants(ax, array, 'y', cmap=cmap)
    plot_quadrants(ax, array, 'z', cmap=cmap)
    return fig, ax


nx, ny, nz = 70, 100, 50
r_square = (np.mgrid[-1:1:1j*nx, -1:1:1j*ny, -1:1:1j*nz]**2).sum(0)

figure_3D_array_slices(r_square, cmap='viridis_r')
plt.show()

grafik

@scottshambaugh
Copy link
Contributor

scottshambaugh commented May 16, 2024

That is cleaner, full disclosure I had chatgpt do my example so I'm not surprised it messed up indexing! Agreed it'd make a nice gallery example, I can whip that up.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests