Search code examples
pythonmatplotlibvectormatplotlib-3d

How do I plot a 3D vector field?


The Data

I have a vector field, which is 0 in all components except for the z component. I just have the data for one slice of this field. My goal is to show this slice in a 3D plot.

The slice:

enter image description here

import numpy as np
import matplotlib.pyplot as plt

# DATA FIELD (1,N,M)
dz = np.array([[[0.24884899, 0.24884899, 0.24884899, 0.24884899, 0.24884899,
         0.24884899],
        [0.248849  , 0.248849  , 0.248849  , 0.248849  , 0.248849  ,
         0.248849  ],
        [0.24885767, 0.24885513, 0.24885108, 0.24885113, 0.2488552 ,
         0.24885767],
        [0.2451304 , 0.24563262, 0.24642831, 0.24641793, 0.24561579,
         0.2451304 ],
        [0.0764377 , 0.12581053, 0.09866768, 0.10043774, 0.12461962,
         0.0764377 ],
        [0.03382106, 0.03394624, 0.03414449, 0.03414171, 0.03394174,
         0.03382106]]])

dx = np.zeros(np.shape(dz))
dy = np.zeros(np.shape(dz))

# DATA POINTS (N,)
X = np.array([0. , 0.2, 0.4, 0.6, 0.8, 1. ])
Y = X 
Z = np.array([-500., -360., -220.,  -80.,   60.,  200.])

To create the plot:

# COMPUTE LENGTH OF VECTORS
length = np.sqrt(dx[0]**2+ dz[0]**2)

# PLOT 2D Cross-Section of vector field
fig = plt.figure(dpi=300)
Q = plt.quiver(X, Z, dx[0], dz[0], length, units='xy' ,angles='xy', scale=0.005,  
       pivot = "tail", headaxislength = 5, headlength = 5, cmap='jet')
fig.colorbar(Q)
plt.gca().invert_yaxis()
plt.ylabel("Z")
plt.xlabel("X")
plt.title("2D Cross-Section")

The Problem

So far so good, but when I try to plot it in 3D, it seems that I am doing something wrong, since

  1. the arrows are not visible any more
  2. the plot seems to be rotated
# Create 3D Quiver Plot with color gradient
# Source: https://stackoverflow.com/questions/65254887/how-to-plot-with-matplotlib-a-3d-quiver-plot-with-color-gradient-for-length-giv

from mpl_toolkits.mplot3d import Axes3D # <--- This is important for 3d plotting 

def plot_3d_quiver(x, y, z, u, v, w):
    # COMPUTE LENGTH OF VECTOR -> MAGNITUDE
    c = np.sqrt(np.abs(v)**2 + np.abs(u)**2 + np.abs(w)**2)

    c = (c.ravel() - c.min())/c.ptp()
    # Repeat for each body line and two head lines
    c = np.concatenate((c, np.repeat(c, 2)))
    # Colormap
    c = plt.cm.jet(c)

    fig = plt.figure(dpi = 300)
    ax = fig.gca(projection = '3d')
    ax.quiver(x, y, z, u, v, w, colors = c, length = 1, arrow_length_ratio = 0.5)
    plt.gca().invert_zaxis()
    plt.show()


# Create Mesh !
xi, yi, zi = np.meshgrid(X, Y, Z, indexing = 'xy')

# slice mesh
skip_every = 1
skip_slice = 6
skip3D = (slice(None, None, skip_slice), slice(None, None, skip_every), slice(None, None, skip_every))

plot_3d_quiver(xi[skip3D], yi[skip3D], zi[skip3D], dx, dy, dz)

enter image description here

What can I try next to resolve this?

Edit

As pointed out by @Marc the arrows are there in the 3D plot, but apparently, the arrow size is relative to the size of the axis. See more here: https://github.com/matplotlib/matplotlib/issues/11746

So far there seems to be no other fix than to scale all axis to roughly the same values. So, if -as suggested by @Marc- I divide the z-axis by 1000, I can see the arrows again:

plot_3d_quiver(xi[skip3D], yi[skip3D], zi[skip3D]/1000, dx, dy, dz)

Then we get:

enter image description here

which clearly shows that the orientation is somehow off.

Edit 2

What I would like to achieve:

The data that I have is just for one cross-section (one plane) of the 3D space (-> dimensions (1,N,M)), hence if I would show it in 3D, it should look something like this:

enter image description here

Full Data

In the end, I would like to correctly plot the full vector field. You can download it here: https://filebin.net/csjvs7xde5lfbwek

# %% 
import pickle
import numpy as np
import matplotlib.pyplot as plt

# Import Full Data
with open('full_data.pickle', 'rb') as handle:
    full_data = pickle.load(handle)

# Axis
X = np.linspace(0,1,101)
Y = np.linspace(0,1,10)
Z = np.linspace(-500,200,101)

for cross_section in full_data["cross_sections"].keys():
    
    # extract field components in x, y, and z 
    dx,dy,dz = full_data["cross_sections"][cross_section]
    
    # Make them numpy imediatley
    dx = np.array(dx)
    dy = np.array(dy)
    dz = np.array(dz)
    

    
    # ------------------------
    # plot cross-sections z-x
    # -----------------------
    
    length=np.sqrt(dx**2+dz**2) 

    fig = plt.figure(dpi=300)
    
    skip_every = 10
    skip_slice = 2
    skip1D= slice(None,None,skip_every)
    skip2D=(slice(None,None,skip_every),slice(None,None,skip_every))
    
    # If U and V are 2-D arrays and X and Y are 1-D, 
    # and if len(X) and len(Y) match the column and row dimensions of U, 
    # then X and Y will be expanded with numpy.meshgrid().
    # To make the head a triangle, make headaxislength the same as headlength.
    Q = plt.quiver(X[skip1D], Z[skip1D], dx[skip2D], dz[skip2D],\
                length[skip2D], units='xy' ,angles='xy', scale=0.02,  
                pivot = "tail", headaxislength = 5, headlength = 5, cmap='jet')
    fig.colorbar(Q)
    plt.title("Cross-Section: " + str(cross_section))
    # Invert y axis to have same represntation than MC GRATING
    plt.gca().invert_yaxis()
    
    plt.pause(0.01)
  
# ------------------------
# plot full 3D Vector Field ?? 
# -----------------------

Solution

  • As discussed in the comments:

    Arrowheads

    The arrow heads don't show, as your z-axis data is much larger compared to x/y-axis values. You'd have to rescale or supply parameters to quiver.

    Data Creation

    There is no general good solution on how you should create your data. It very much depends on your use-case. In my opinion, an easily understandable and tuned approach is the following:

    import numpy as np, matplotlib.pyplot as plt, pandas as pd
    from mpl_toolkits.mplot3d import Axes3D
    
    def plot_3d_quiver(x, y, z, u, v, w):
        # COMPUTE LENGTH OF VECTOR -> MAGNITUDE
        c = np.sqrt(np.abs(v)**2 + np.abs(u)**2 + np.abs(w)**2)
    
        c = (c.ravel() - c.min())/c.ptp()
        # Repeat for each body line and two head lines
        c = np.concatenate((c, np.repeat(c, 2)))
        # Colormap
        c = plt.cm.jet(c)
    
        fig = plt.figure(dpi = 300)
        ax = fig.gca(projection = '3d')
        ax.quiver(x, y, z, u, v, w, colors = c, length = .5, arrow_length_ratio = 0.2)
        plt.gca().invert_zaxis()
        plt.show()
    
    axis = np.r_[-1:1:.2]
    X,Y,Z = np.meshgrid(axis,axis,axis)
    
    zeros = np.zeros_like(X)
    U,V,W = (zeros, zeros, zeros+np.linspace(0,-1,len(axis))[:,np.newaxis,np.newaxis])
        
    plot_3d_quiver(X,Y,Z,U,V,W)
    
    

    enter image description here

    edit The problem lays in how you expect the 3D-mesh points to be ordered with respect to your data. You can match your data dz to the meshgrid by switching your third axis with the second:

    plot_3d_quiver(xi[skip3D], yi[skip3D], zi[skip3D], dx, dy, np.moveaxis(dz,2,1))
    

    Or, if you'd use two dimensions only, you can transpose the matrix:

    plot_3d_quiver(xi[0], yi[0], zi[0], dx[0], dy[0], dz[0].T)