examples/01-filter/compute-volume.ipynb

{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "%matplotlib inline"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Volumetric Analysis {#volumetric_example}\n===================\n\nCalculate mass properties such as the volume or area of datasets\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import numpy as np\n\nfrom pyvista import examples"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Computing mass properties such as the volume or area of datasets in\nPyVista is quite easy using the\n`pyvista.DataSetFilters.compute_cell_sizes`{.interpreted-text\nrole=\"func\"} filter and the `pyvista.DataSet.volume`{.interpreted-text\nrole=\"attr\"} property on all PyVista meshes.\n\nLet\\'s get started with a simple gridded mesh:\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# Load a simple example mesh\ndataset = examples.load_uniform()\ndataset.set_active_scalars(\"Spatial Cell Data\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "We can then calculate the volume of every cell in the array using the\n`.compute_cell_sizes` filter which will add arrays to the cell data of\nthe mesh core the volume and area by default.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# Compute volumes and areas\nsized = dataset.compute_cell_sizes()\n\n# Grab volumes for all cells in the mesh\ncell_volumes = sized.cell_data[\"Volume\"]"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "We can also compute the total volume of the mesh using the `.volume`\nproperty:\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# Compute the total volume of the mesh\nvolume = dataset.volume"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "But what if we have a dataset that we threshold with two volumetric\nbodies left over in one dataset? Take this for example:\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "threshed = dataset.threshold_percent([0.15, 0.50], invert=True)\nthreshed.plot(show_grid=True, cpos=[-2, 5, 3])"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "We could then assign a classification array for the two bodies, compute\nthe cell sizes, then extract the volumes of each body. Note that there\nis a simpler implementation of this below in\n`split_vol`{.interpreted-text role=\"ref\"}.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# Create a classifying array to ID each body\nrng = dataset.get_data_range()\ncval = ((rng[1] - rng[0]) * 0.20) + rng[0]\nclassifier = threshed.cell_data[\"Spatial Cell Data\"] > cval\n\n# Compute cell volumes\nsizes = threshed.compute_cell_sizes()\nvolumes = sizes.cell_data[\"Volume\"]\n\n# Split volumes based on classifier and get the volumes\nidx = np.argwhere(classifier)\nhvol = np.sum(volumes[idx])\nidx = np.argwhere(~classifier)\nlvol = np.sum(volumes[idx])\n\nprint(f\"Low grade volume: {lvol}\")\nprint(f\"High grade volume: {hvol}\")\nprint(f\"Original volume: {dataset.volume}\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Or better yet, you could simply extract the largest volume from your\ndataset directly by passing `'largest'` to the `connectivity` and\nspecifying the scalar range of interest.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# Grab the largest connected volume within a scalar range\nscalar_range = [0, 77]  # Range corresponding to bottom 15% of values\nlargest = threshed.connectivity('largest', scalar_range=scalar_range)\n\n# Get volume as numeric value\nlarge_volume = largest.volume\n\n# Display it\nlargest.plot(show_grid=True, cpos=[-2, 5, 3])"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "------------------------------------------------------------------------\n\nSplitting Volumes {#split_vol}\n=================\n\nWhat if instead, we wanted to split all the different connected bodies /\nvolumes in a dataset like the one above? We could use the\n`pyvista.DataSetFilters.split_bodies`{.interpreted-text role=\"func\"}\nfilter to extract all the different connected volumes in a dataset into\nblocks in a `pyvista.MultiBlock`{.interpreted-text role=\"class\"}\ndataset. For example, lets split the thresholded volume in the example\nabove:\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# Load a simple example mesh\ndataset = examples.load_uniform()\ndataset.set_active_scalars(\"Spatial Cell Data\")\nthreshed = dataset.threshold_percent([0.15, 0.50], invert=True)\n\nbodies = threshed.split_bodies()\n\nfor i, body in enumerate(bodies):\n    print(f\"Body {i} volume: {body.volume:.3f}\")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "bodies.plot(show_grid=True, multi_colors=True, cpos=[-2, 5, 3])"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "------------------------------------------------------------------------\n\nA Real Dataset\n==============\n\nHere is a realistic training dataset of fluvial channels in the\nsubsurface. This will threshold the channels from the dataset then\nseparate each significantly large body and compute the volumes for each.\n\nLoad up the data and threshold the channels:\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "data = examples.load_channels()\nchannels = data.threshold([0.9, 1.1])"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Now extract all the different bodies and compute their volumes:\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "bodies = channels.split_bodies()\n# Now remove all bodies with a small volume\nfor key in bodies.keys():\n    b = bodies[key]\n    vol = b.volume\n    if vol < 1000.0:\n        del bodies[key]\n        continue\n    # Now lets add a volume array to all blocks\n    b.cell_data[\"TOTAL VOLUME\"] = np.full(b.n_cells, vol)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Print out the volumes for each body:\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "for i, body in enumerate(bodies):\n    print(f\"Body {i:02d} volume: {body.volume:.3f}\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "And visualize all the different volumes:\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "bodies.plot(scalars=\"TOTAL VOLUME\", cmap=\"viridis\", show_grid=True)"
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.12.1"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}