Search code examples
python-3.xpyvista

Pyvista surface plot?


I need a way to make a 3-dimensional surface plot using millions of datapoints, so I began checking into pyvista which is supposed to do this well. However, pyvista is a bit difficult for me to grasp.

I have x,y,z data where x is time, y is different measurements, and z is the values for those measurements. All I want is for pyvista to show me a surface plot with this information.

For example, if I use this array in matplotlib or other libraries with surface plots:

X = np.array([1,2,3,4,5,6,7,8,9])
Y = np.array([1,2,3,4,5,6,7,8,9])
X, Y = np.meshgrid(X, Y)
Z = X*Y

I get this output:

enter image description here

But if I use the same data on any of the pyvista plots, I get something like this:

import sys

# Setting the Qt bindings for QtPy
import os
os.environ["QT_API"] = "pyqt5"

from qtpy import QtWidgets
from qtpy.QtWidgets import QMainWindow

import numpy as np

import pyvista as pv
from pyvistaqt import QtInteractor

import pandas as pd

class MainWindow(QMainWindow):

    def __init__(self, parent=None, show=True):
        QtWidgets.QMainWindow.__init__(self, parent)

        # create the frame
        self.frame = QtWidgets.QFrame()
        vlayout = QtWidgets.QVBoxLayout()

        # add the pyvista interactor object
        self.plotter = QtInteractor(self.frame)
        vlayout.addWidget(self.plotter.interactor)

        self.frame.setLayout(vlayout)
        self.setCentralWidget(self.frame)

        # simple menu to demo functions
        mainMenu = self.menuBar()
        fileMenu = mainMenu.addMenu('File')
        exitButton = QtWidgets.QAction('Exit', self)
        exitButton.setShortcut('Ctrl+Q')
        exitButton.triggered.connect(self.close)
        fileMenu.addAction(exitButton)

        # allow adding a sphere
        meshMenu = mainMenu.addMenu('Mesh')
        self.add_sphere_action = QtWidgets.QAction('Add Sphere', self)
        self.add_sphere_action.triggered.connect(self.add_sphere)
        meshMenu.addAction(self.add_sphere_action)

        x = np.array([9,8,7,6,5,4,3,2,1])
        y = np.array([9,8,7,6,5,4,3,2,1])
        x, y = np.meshgrid(x, y)
        z = x*y
    
        # z[z < -10] = np.nan  # get rid of missing data. pyvista needs you to do this

        i_res = 2  # display every nth point
        j_res = 2  # display every nth point
        self.grid = pv.StructuredGrid(x[::i_res, ::j_res], y[::i_res, ::j_res], z[::i_res, ::j_res])

        self.z = z
        self.x = x
        self.y = y

        self.plotter.add_mesh(self.grid, scalars=self.grid.points[:, 2], lighting=True, specular=0.5, smooth_shading=True,
                              show_scalar_bar=True)

        if show:
            self.show()

    def add_sphere(self): #changing resolution, not adding a sphere
        i_res = 5  # display every nth point
        j_res = 5  # display every nth point
        self.grid = pv.StructuredGrid(self.x[::i_res, ::j_res], self.y[::i_res, ::j_res], self.z[::i_res, ::j_res])
        self.plotter.update()

    
if __name__ == '__main__':
    app = QtWidgets.QApplication(sys.argv)
    window = MainWindow()
    sys.exit(app.exec())

enter image description here

import pyvista as pv
import numpy as np

# Define a simple Gaussian surface
x = np.array([1,2,3,4,5,6,7,8,9])
y = np.array([1,2,3,4,5,6,7,8,9])
x, y = np.meshgrid(x, y)
z = x*y

# Get the points as a 2D NumPy array (N by 3)
points = np.c_[x.reshape(-1), y.reshape(-1), z.reshape(-1)]
points[0:5, :]

# simply pass the numpy points to the PolyData constructor
cloud = pv.PolyData(points)
cloud.plot(point_size=15)

enter image description here

I managed to get "something" different using this bit of code:

import pandas as pd
import pyvista as pv
import numpy as np

# Load Excel sheet using Pandas
# Note - you may need to `pip install xlrd`
# x = np.array([1,2,3,4,5,6,7,8,9])
# y = np.array([1,2,3,4,5,6,7,8,9])

x = np.array([[1],[2],[3],[4],[5],[6],[7],[8],[9]])
y = np.array([[1],[2],[3],[4],[5],[6],[7],[8],[9]])
# # x, y = np.meshgrid(x, y)
z = x*y

coords = np.hstack((x,y,z))

# Make the structured surface manually
structured = pv.StructuredGrid()
# Set coordinates
structured.points = coords
# Set the dimensions of the structured grid
structured.dimensions = [1, 1, 9]

# Apply an Elevation filter
elevation = structured.elevation()
elevation.plot(show_edges=True, show_grid=True, notebook=False)

enter image description here

But it only provides a single string of data. I haven't been able to get anything else work properly.

Does anyone know why the x,y,z data is doing weird things in pyvista and how I can provide just a normal surface plot? It would be much appreciated, as I am pretty stumped.


Solution

  • Your first version is correct.

    PyVista has excellent documentation, part of which is an extensive collection of examples. You need the one that's called Creating a Structured Surface. This ends up being pretty much the same code as what you originally showed:

    import pyvista as pv
    import numpy as np
    
    # Define a simple linear surface
    x = np.array([1,2,3,4,5,6,7,8,9])
    y = np.array([1,2,3,4,5,6,7,8,9])
    x, y = np.meshgrid(x, y)
    z = x*y
    
    # Create and plot structured grid
    grid = pv.StructuredGrid(x, y, z)
    plotter = pv.Plotter()
    plotter.add_mesh(grid, scalars=grid.points[:, -1], show_edges=True,
                     scalar_bar_args={'vertical': True})
    plotter.show_grid()
    plotter.show()
    

    Here is the (correct!) output: prettier plot of what OP had

    The reason why this looks different is that matplotlib isn't a 3d visualization tool (in fact its 3d tooling infamously uses a 2d renderer that leads to weird quirks). PyVista on the other hand is designed to visualize spatially referenced data. If your x goes from 1 to 9 and your z goes from 1 to 81 then why would it squash the z axis? What PyVista shows is the truth if you set a 1:1:1 aspect ratio along each coordinate axis.

    If you don't want this, you can mess with scaling yourself:

    import pyvista as pv
    import numpy as np
    
    # Define a simple linear surface
    x = np.array([1,2,3,4,5,6,7,8,9])
    y = np.array([1,2,3,4,5,6,7,8,9])
    x, y = np.meshgrid(x, y)
    z = x*y
    
    # Create and plot structured grid
    grid = pv.StructuredGrid(x, y, z)
    plotter = pv.Plotter()
    plotter.add_mesh(grid, scalars=grid.points[:, -1], show_edges=True,
                     scalar_bar_args={'vertical': True})
    plotter.show_grid()
    # scale plot to enforce 1:1:1 aspect ratio
    plotter.set_scale(xscale=1, yscale=x.ptp()/y.ptp(), zscale=x.ptp()/z.ptp())
    plotter.show()
    

    scaled image where the plot fits in a cube

    If you want PyVista to lie about your data, you have to tell it to do so.