{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Multi-dimensional datasets\n",
    "\n",
    "This is the continuation of [1-D datasets and tables](introduction.ipynb).\n",
    "\n",
    "## Creation, slicing, and visualization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import scipp as sc"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To create variables with more than one dimension we specify a list of dimension labels and provide data with a corresponding shape.\n",
    "When inserted into a dataset it is important to note that while the dimensions extents have to match, individual variables may have transposed memory layout."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "d = sc.Dataset(\n",
    "    data={\n",
    "        'alice': sc.Variable(dims=['z', 'y', 'x'], values=np.random.rand(10, 10, 10),\n",
    "                             variances=0.1*np.random.rand(10, 10, 10)),\n",
    "        'bob': sc.Variable(dims=['x', 'z'], values=np.arange(0.0, 10.0, 0.1).reshape(10, 10),\n",
    "                           variances=0.1*np.random.rand(10, 10))\n",
    "    },\n",
    "    coords={\n",
    "        'x': sc.arange('x', 11.0, unit=sc.units.m),\n",
    "        'y': sc.arange('y', 11.0, unit=sc.units.m),\n",
    "        'z': sc.arange('z', 11.0, unit=sc.units.m)})"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that in this example the coordinates are exceeding the shape of the data by 1.\n",
    " This means that the coordinates represent bin edges:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.show(d)\n",
    "d"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To slice in multiple dimensions, we can simply chain the slicing notation used previously for 1D data.\n",
    "This gives us a number of different options for visualizing our data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.table(d['x', 5]['z', 2])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can plot and item of a `Dataset` using:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.plot(d['bob'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plotting a 3-dimensional data cube will show a 2D image with a slider to navigate through the third dimension (note that interactive sliders will only appear in a Jupyter notebook and not in the documentation pages):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.plot(d['alice'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, by extracting a 1D variable, we obtain a 1D plot:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.plot(d['x', 8]['y', 2])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that this is now plotted as a histogram since the coordinate in the dataset represents bin edges, in contrast to the 1D data plotted in [1-D datasets and tables](introduction.ipynb).\n",
    "\n",
    "Operations automatically broadcast based on dimension labels.\n",
    "That is, if one of the operands lacks one (or multiple) dimensions that the other operands have, the operand is considered constant along those dimensions.\n",
    "Its values are implicitly \"duplicated\" so the shape matches.\n",
    "In contrast to `numpy` or `MATLAB` there is no need to keep track of dimension order:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "d['alice'] -= d['bob']\n",
    "d['alice'] -= d['alice']['y', 5]\n",
    "sc.plot(d['alice']['x', 4])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercise 1\n",
    "\n",
    "Remove the X and Z surface layer of the volume, i.e., remove the first and last slice in each of the dimensions `'x'` and `'z'`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Solution 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "d = d['x', 1:-1]['z', 1:-1].copy()\n",
    "d"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note the important call to `copy()`.\n",
    "If we omit it, `d` will just be a multi-dimensional slice of the larger volume (which is kept alive), wasting memory and preventing further modification, such as insertion of other variables.\n",
    "\n",
    "Note also that if we had also sliced `'y'` the result would not contain the data for `'bob'` since this item does not depend on Y."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## More advanced operations with multi-dimensional datasets\n",
    "Operations like `concatenate` and `merge` work just like with one-dimensional datasets."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercise 2\n",
    "- Try to concatenate the dataset with itself along the X dimensions. Why does this fail?\n",
    "- Make a copy of the dataset, add an offset to the X coordinate to fix the issue, and try to concatenate again."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Solution 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "try:\n",
    "    d = sc.concatenate(d, d, 'x')\n",
    "except RuntimeError:\n",
    "    print(\"Failed as expected!\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With a data extent of, e.g. `8` in this case, bin edges have extent `9`.\n",
    "Naive concatenation would thus lead a new data extent of `16` and a coordinate extent of `18`, which is meaningless and thus prevented.\n",
    "In this `concatenate` merges the last edge of the first input with the first edge of the second input, if compatible."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "offset = d.copy()\n",
    "offset.coords['x'] += sc.scalar(8.0, unit=sc.units.m)\n",
    "combined = sc.concatenate(d, offset, 'x')\n",
    "sc.plot(combined['alice'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Another available operation is `rebin`.\n",
    " This is only for count-data or count-density-data, so we have to set an appropriate unit first:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "d['alice'].unit = sc.units.counts\n",
    "d['bob'].unit = sc.units.counts"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Before `rebin` we have the following:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.to_html(d)\n",
    "sc.plot(d['z', 0])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We rebin onto a coarser grid, in this case combining two neightboring bins:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "new_x = sc.array(dims=['x'], unit=sc.units.m, values=d.coords['x'].values[::2])\n",
    "d = sc.rebin(d, 'x', new_x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The result looks as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.to_html(d)\n",
    "sc.plot(d['z', 0])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Interaction with `numpy`\n",
    "\n",
    "Variable in a dataset are exposed in a `numpy`-compatible buffer format, so we can directly hand them to `numpy` functions:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "d['alice'].values = np.sin(d['alice'].values)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In contrast to the 1-D case considered earlier, the `values` are now a multi-dimensional array:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "d['alice'].values"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercise 3\n",
    " 1. Use `sc.mean` to compute the mean of the data for Alice along the Z dimension.\n",
    " 2. Do the same with `numpy`, what are the complications you encounter, that are not present when using the dataset?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Solution 3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "help(sc.mean)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mean = sc.mean(d['alice'], 'z')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "When using `numpy` to compute the mean:\n",
    "- We must remember (or lookup) which dimension corresponds to the Z dimensions.\n",
    "- We need a separate call for values and variances.\n",
    "- We need to manually scale the variance with the inverse of the number of data points to get the variance of the mean (standard deviation of the mean scales with `1/sqrt(N)`)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np_values = np.mean(d['alice'].values, axis=0)\n",
    "np_variances = np.mean(d['alice'].variances, axis=0)\n",
    "np_variances /= d['alice'].shape[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(mean.values)\n",
    "print(mean.variances)\n",
    "print(np_values)\n",
    "print(np_variances)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.7.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}