{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### X LINES OF PYTHON\n",
    "\n",
    "# Physical units with `pint`\n",
    "\n",
    "This notebook goes with [a blog post on the same subject](https://agilescientific.com/blog/2019/8/19/x-lines-of-python-physical-units).\n",
    "\n",
    "Have you ever wished you could carry units around with your quantities — and have the computer figure out the best units and multipliers to use?\n",
    "\n",
    "[`pint`](https://pint.readthedocs.io/en/0.9/tutorial.html) is a nince, compact library for doing just this, handling all your [dimensional analysis](https://en.wikipedia.org/wiki/Dimensional_analysis) needs. It can also detect units from strings. We can define our own units, it knows about multipliers (kilo, mega, etc), and it even works with `numpy` and `pandas`.\n",
    "\n",
    "Install `pint` with `pip` or `conda`, e.g.\n",
    "\n",
    "    pip install pint\n",
    "    \n",
    "**NB** If you are running this on Google Colaboratory, you must uncomment these lines (delete the initial `#`) and run this first:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#!pip install pint\n",
    "#!pip install git+https://github.com/hgrecco/pint-pandas#egg=Pint-Pandas-0.1.dev0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To use it in its typical mode, we import the library then instantiate a `UnitRegistry` object. The registry contains lots of physical units. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pint\n",
    "\n",
    "units = pint.UnitRegistry()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'0.9'"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pint.__version__"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Attaching and printing units"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "68 meter"
      ],
      "text/latex": [
       "$68\\ \\mathrm{meter}$"
      ],
      "text/plain": [
       "68 <Unit('meter')>"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "thickness = 68 * units.m\n",
    "thickness"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In a Jupyter Notebook you see a 'pretty' version of the quantity. In the interpreter, you'll see something slightly different (the so-called `repr` of the class):\n",
    "\n",
    "    >>> thickness\n",
    "    <Quantity(68, 'meter')>\n",
    "\n",
    "We can get at the magnitude, the units, and the dimensionality of this quantity:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(68, <Unit('meter')>, <UnitsContainer({'[length]': 1.0})>)"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "thickness.magnitude, thickness.units, thickness.dimensionality"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can also use the following abbreviations for magnitude and units:\n",
    "\n",
    "    thickness.m, thickness.u\n",
    "    \n",
    "For printing, we can use Python's string formatting:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'4624 meter ** 2'"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f'{thickness**2}'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "But `pint` extends the string formatting options to include special options for `Quantity` objects. The most useful option is `P` for 'pretty', but there's also `L` for $\\LaTeX$ and `H` for HTML. Adding a `~` (tilde) before the option tells `pint` to use unit abbreviations instead of the full names:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "4624 meter²\n",
      "4624 m²\n",
      "4624\\ \\mathrm{m}^{2}\n",
      "4624 m<sup>2</sup>\n"
     ]
    }
   ],
   "source": [
    "print(f'{thickness**2:P}')\n",
    "print(f'{thickness**2:~P}')\n",
    "print(f'{thickness**2:~L}')\n",
    "print(f'{thickness**2:~H}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Doing maths\n",
    "\n",
    "If we multiply by a scalar, `pint` produces the result you'd expect:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "136 meter"
      ],
      "text/latex": [
       "$136\\ \\mathrm{meter}$"
      ],
      "text/plain": [
       "136 <Unit('meter')>"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "thickness * 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that you must use units when you need them:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "ename": "DimensionalityError",
     "evalue": "Cannot convert from 'meter' to 'dimensionless'",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mDimensionalityError\u001b[0m                       Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-10-ecd15f6b9f4c>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mthickness\u001b[0m \u001b[0;34m+\u001b[0m \u001b[0;36m10\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      2\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      3\u001b[0m \u001b[0;31m# This is meant to produce an error...\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/anaconda3/envs/xlines/lib/python3.7/site-packages/pint/quantity.py\u001b[0m in \u001b[0;36m__add__\u001b[0;34m(self, other)\u001b[0m\n\u001b[1;32m    752\u001b[0m             \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mto_timedelta\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;34m+\u001b[0m \u001b[0mother\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    753\u001b[0m         \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 754\u001b[0;31m             \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_add_sub\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mother\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0moperator\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0madd\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    755\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    756\u001b[0m     \u001b[0m__radd__\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0m__add__\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/anaconda3/envs/xlines/lib/python3.7/site-packages/pint/quantity.py\u001b[0m in \u001b[0;36mwrapped\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m     73\u001b[0m         \u001b[0;32melif\u001b[0m \u001b[0misinstance\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mother\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mlist\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mand\u001b[0m \u001b[0misinstance\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mother\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtype\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     74\u001b[0m             \u001b[0;32mreturn\u001b[0m \u001b[0mNotImplemented\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 75\u001b[0;31m         \u001b[0mresult\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mf\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m     76\u001b[0m         \u001b[0;32mreturn\u001b[0m \u001b[0mresult\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     77\u001b[0m     \u001b[0;32mreturn\u001b[0m \u001b[0mwrapped\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/anaconda3/envs/xlines/lib/python3.7/site-packages/pint/quantity.py\u001b[0m in \u001b[0;36m_add_sub\u001b[0;34m(self, other, op)\u001b[0m\n\u001b[1;32m    663\u001b[0m                                _to_magnitude(other, self.force_ndarray))\n\u001b[1;32m    664\u001b[0m             \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 665\u001b[0;31m                 \u001b[0;32mraise\u001b[0m \u001b[0mDimensionalityError\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_units\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m'dimensionless'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    666\u001b[0m             \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m__class__\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmagnitude\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0munits\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    667\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mDimensionalityError\u001b[0m: Cannot convert from 'meter' to 'dimensionless'"
     ]
    }
   ],
   "source": [
    "thickness + 10\n",
    "\n",
    "# This is meant to produce an error..."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's try defining an area of $60\\ \\mathrm{km}^2$, then multiplying it by our thickness. To make it more like a hydrocarbon volume, I'll also multiply by net:gross `n2g`, porosity `phi`, and saturation `sat`, all of which are dimensionless:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "285.59999999999997 kilometer<sup>2</sup> meter"
      ],
      "text/latex": [
       "$285.59999999999997\\ \\mathrm{kilometer}^{2} \\cdot \\mathrm{meter}$"
      ],
      "text/plain": [
       "285.59999999999997 <Unit('kilometer ** 2 * meter')>"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "area = 60 * units.km**2\n",
    "n2g = 0.5 * units.dimensionless  # Optional dimensionless 'units'...\n",
    "phi = 0.2                        # ... but you can just do this.\n",
    "sat = 0.7  \n",
    "\n",
    "volume = area * thickness * n2g * phi * sat\n",
    "volume"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can convert to something more compact:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "285599999.99999994 meter<sup>3</sup>"
      ],
      "text/latex": [
       "$285599999.99999994\\ \\mathrm{meter}^{3}$"
      ],
      "text/plain": [
       "285599999.99999994 <Unit('meter ** 3')>"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "volume.to_compact()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Or be completely explicit about the units and multipliers we want:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "285599999.99999994 meter<sup>3</sup>"
      ],
      "text/latex": [
       "$285599999.99999994\\ \\mathrm{meter}^{3}$"
      ],
      "text/plain": [
       "285599999.99999994 <Unit('meter ** 3')>"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "volume.to('m**3')  # Or use m^3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `to_compact()` method can also take units, if you want to be more explicit; it applies multipliers automatically:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "285.6 gigaliter"
      ],
      "text/latex": [
       "$285.6\\ \\mathrm{gigaliter}$"
      ],
      "text/plain": [
       "285.6 <Unit('gigaliter')>"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "volume.to_compact('L')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Oil barrels are already defined (**careful**, they are abbreviated as `oil_bbl` not `bbl` &mdash; that's a 31.5 gallon barrel, about the same as a beer barrel). "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "1.7963699560354092 gigaoil_barrel"
      ],
      "text/latex": [
       "$1.7963699560354092\\ \\mathrm{gigaoil_barrel}$"
      ],
      "text/plain": [
       "1.7963699560354092 <Unit('gigaoil_barrel')>"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "volume.to_compact('oil_barrel')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If we use string formatting (see above), we can get pretty specific:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'The volume is 1.80\\\\ \\\\mathrm{Goil_bbl}'"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f\"The volume is {volume.to_compact('oil_barrel'):~0.2fL}\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Defining new units\n",
    "\n",
    "`pint` defines hundreads of units ([here's the list](https://github.com/hgrecco/pint/blob/master/pint/default_en.txt)), and it knows about tonnes of oil equivalent... but it doesn't know about barrels of oil equivalent ([for more on conversion to BOE](https://en.wikipedia.org/wiki/Barrel_of_oil_equivalent)). So let's define a custom unit, using the USGS's conversion factor:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "units.define('barrel_of_oil_equivalent = 6000 ft**3 = boe')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's suspend reality for a moment and imagine we now want to compute our gross rock volume in BOEs..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "1680978.135942857 barrel_of_oil_equivalent"
      ],
      "text/latex": [
       "$1680978.135942857\\ \\mathrm{barrel_of_oil_equivalent}$"
      ],
      "text/plain": [
       "1680978.135942857 <Unit('barrel_of_oil_equivalent')>"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "volume.to('boe')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "1.680978135942857 megabarrel_of_oil_equivalent"
      ],
      "text/latex": [
       "$1.680978135942857\\ \\mathrm{megabarrel_of_oil_equivalent}$"
      ],
      "text/plain": [
       "1.680978135942857 <Unit('megabarrel_of_oil_equivalent')>"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "volume.to_compact('boe')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Getting units from strings\n",
    "\n",
    "`pint` can also parse strings and attempt to convert them to `Quantity` instances:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "2.34 kilometer"
      ],
      "text/latex": [
       "$2.34\\ \\mathrm{kilometer}$"
      ],
      "text/plain": [
       "2.34 <Unit('kilometer')>"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "units('2.34 km')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This looks useful! Let's try something less nicely formatted."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "2340.0 kilometer"
      ],
      "text/latex": [
       "$2340.0\\ \\mathrm{kilometer}$"
      ],
      "text/plain": [
       "2340.0 <Unit('kilometer')>"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "units('2.34*10^3 km')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "-12000.0 foot"
      ],
      "text/latex": [
       "$-12000.0\\ \\mathrm{foot}$"
      ],
      "text/plain": [
       "-12000.0 <Unit('foot')>"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "units('-12,000.ft')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "3.2 meter"
      ],
      "text/latex": [
       "$3.2\\ \\mathrm{meter}$"
      ],
      "text/plain": [
       "3.2 <Unit('meter')>"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "units('3.2 m')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can also use the `Quantity` constructor, like this:\n",
    "\n",
    "    >>> qty = pint.Quantity\n",
    "    >>> qty('2.34 km')\n",
    "    2.34 kilometer\n",
    "    \n",
    "But the `UnitRegistry` seems to do the same things and might be more convenient."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## `pint` with `uncertainties`\n",
    "\n",
    "Conveniently, `pint` works well with [`uncertainties`](https://pythonhosted.org/uncertainties/). Maybe I'll do an _X lines_ on that package in the future. Install it with `conda` or `pip`, e.g.\n",
    "\n",
    "    pip install uncertainties"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "27.4+/-2.1 gigaoil_barrel"
      ],
      "text/latex": [
       "$27.4+/-2.1\\ \\mathrm{gigaoil_barrel}$"
      ],
      "text/plain": [
       "27.373256472920524+/-2.138535661946916 <Unit('gigaoil_barrel')>"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from uncertainties import ufloat\n",
    "\n",
    "area = ufloat(64, 5) * units.km**2  # 64 +/- 5 km**2\n",
    "(thickness * area).to('Goil_bbl')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## `pint` with `numpy`\n",
    "\n",
    "`pint` works fine with NumPy arrays:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "\\[\\begin{pmatrix}5520000.0 & 6120000.0 & 6375000.0 & 8480000.0\\end{pmatrix} kilogram/(meter<sup>2</sup> second)\\]"
      ],
      "text/latex": [
       "$\\begin{pmatrix}5520000.0 & 6120000.0 & 6375000.0 & 8480000.0\\end{pmatrix}\\ \\frac{\\mathrm{kilogram}}{\\left(\\mathrm{meter}^{2} \\cdot \\mathrm{second}\\right)}$"
      ],
      "text/plain": [
       "array([5520000., 6120000., 6375000., 8480000.]) <Unit('kilogram / meter ** 2 / second')>"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "vp = np.array([2300, 2400, 2550, 3200]) * units.m/units.s\n",
    "rho = np.array([2400, 2550, 2500, 2650]) * units.kg/units.m**3\n",
    "\n",
    "z = vp * rho\n",
    "z"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For some reason, this sometimes doesn't render properly. But we can always do this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[5520000. 6120000. 6375000. 8480000.] kilogram / meter ** 2 / second\n"
     ]
    }
   ],
   "source": [
    "print(z)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As expected, the magnitude of this quantity is just a NumPy array:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([5520000., 6120000., 6375000., 8480000.])"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "z.m"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## `pint` with `pandas`\n",
    "\n",
    "**Note** that this functionality is fairly new and is still settling down. YMMV.\n",
    "\n",
    "To use `pint` (version 0.9 and later) with `pandas` (version 0.24.2 works; 0.25.0 does not work at the time of writing), we must first install `pintpandas`, which must be done from source; [get the code from GitHub](https://github.com/hgrecco/pint-pandas). Here's how I do it:\n",
    "\n",
    "    cd pint-pandas\n",
    "    python setup.py sdist\n",
    "    pip install dist/Pint-Pandas-0.1.dev0.tar.gz\n",
    "    \n",
    "You could also do:\n",
    "\n",
    "    pip install git+https://github.com/hgrecco/pint-pandas#egg=Pint-Pandas-0.1.dev0\n",
    "\n",
    "Once you have done that, the following should evaluate to `True`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pint._HAS_PINTPANDAS"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To use this integration, we pass special `pint` data types to the `pd.Series()` object:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Vp</th>\n",
       "      <th>Vs</th>\n",
       "      <th>rho</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2300.0</td>\n",
       "      <td>1200</td>\n",
       "      <td>2400.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2400.0</td>\n",
       "      <td>1200</td>\n",
       "      <td>2550.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2550.0</td>\n",
       "      <td>1250</td>\n",
       "      <td>2500.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3200.0</td>\n",
       "      <td>1300</td>\n",
       "      <td>2650.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "       Vp    Vs     rho\n",
       "0  2300.0  1200  2400.0\n",
       "1  2400.0  1200  2550.0\n",
       "2  2550.0  1250  2500.0\n",
       "3  3200.0  1300  2650.0"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import pandas as pd\n",
    "\n",
    "df = pd.DataFrame({\n",
    "    \"Vp\": pd.Series(vp.m, dtype=\"pint[m/s]\"),\n",
    "    \"Vs\": pd.Series([1200, 1200, 1250, 1300], dtype=\"pint[m/s]\"),\n",
    "    \"rho\": pd.Series(rho.m, dtype=\"pint[kg/m**3]\"),\n",
    "})\n",
    "\n",
    "df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0     9075366233.766233\n",
       "1          9792000000.0\n",
       "2     10483220521.25506\n",
       "3    12550276023.391813\n",
       "Name: E, dtype: pint[kilogram / meter / second ** 2]"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import bruges as bg\n",
    "\n",
    "df['E'] = bg.rockphysics.moduli.youngs(df.Vp, df.Vs, df.rho)\n",
    "\n",
    "df.E"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can't convert the units of a whole `Series` but we can do one:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "9.075366233766236 gigapascal"
      ],
      "text/latex": [
       "$9.075366233766236\\ \\mathrm{gigapascal}$"
      ],
      "text/plain": [
       "9.075366233766236 <Unit('gigapascal')>"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df.loc[0, 'E'].to('GPa')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So to convert a whole series, we can use `Series.apply()`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0     9.075366233766236 gigapascal\n",
       "1     9.792000000000003 gigapascal\n",
       "2    10.483220521255063 gigapascal\n",
       "3    12.550276023391817 gigapascal\n",
       "Name: E, dtype: object"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df.E.apply(lambda x: x.to('GPa'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Bonus: dataframe display with units\n",
    "\n",
    "We *could* subclass dataframes to tweak their `_repr_html_()` method, which would allow us to make units show up in the Notebook representation of the dataframe..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [],
   "source": [
    "class UnitDataFrame(pd.DataFrame):\n",
    "    def _repr_html_(self):\n",
    "        \"\"\"New repr for Jupyter Notebook.\"\"\"\n",
    "        html = super()._repr_html_()  # Get the old repr string.\n",
    "        units = [''] + [f\"{dtype.units:~H}\" for dtype in self.dtypes]\n",
    "        style = \"text-align: right; color: gray;\"\n",
    "        new = f'<tr style=\"{style}\"><th>' + \"</th><th>\".join(units) + \"</th></tr></thead>\"\n",
    "        return html.replace('</thead>', new)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Vp</th>\n",
       "      <th>Vs</th>\n",
       "      <th>rho</th>\n",
       "    </tr>\n",
       "  <tr style=\"text-align: right; color: gray;\"><th></th><th>m/s</th><th>m/s</th><th>kg/m<sup>3</sup></th></tr></thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2300.0</td>\n",
       "      <td>1200</td>\n",
       "      <td>2400.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2400.0</td>\n",
       "      <td>1200</td>\n",
       "      <td>2550.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2550.0</td>\n",
       "      <td>1250</td>\n",
       "      <td>2500.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3200.0</td>\n",
       "      <td>1300</td>\n",
       "      <td>2650.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "       Vp    Vs     rho\n",
       "0  2300.0  1200  2400.0\n",
       "1  2400.0  1200  2550.0\n",
       "2  2550.0  1250  2500.0\n",
       "3  3200.0  1300  2650.0"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df = UnitDataFrame({\n",
    "    \"Vp\": pd.Series(vp.m, dtype=\"pint[m/s]\"),\n",
    "    \"Vs\": pd.Series([1200, 1200, 1250, 1300], dtype=\"pint[m/s]\"),\n",
    "    \"rho\": pd.Series(rho.m, dtype=\"pint[kg/m**3]\"),\n",
    "})\n",
    "\n",
    "df"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Cute."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "----\n",
    "\n",
    "&copy; Agile Scientific 2019, licensed CC-BY"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "xlines",
   "language": "python",
   "name": "xlines"
  },
  "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.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}