epochpic · ShaunD137 · Jul 28, 2026 · Jul 10, 2026 · Jul 11, 2026 · Jul 11, 2026
diff --git a/docs/animation.md b/docs/animation.md
@@ -35,9 +35,10 @@ resolved data has 2 dimensions.
 
 | Dimensions |      Plotting function        |        Notes          |
 | ---------- | ----------------------------- | --------------------- |
-|    `2`     | <inv:#xarray.plot.line>       |                       |
-|    `3`     | <inv:#xarray.plot.pcolormesh> |                       |
-|    `>3`    | <inv:#xarray.plot.hist>       | Not fully implemented |
+|     2      | <inv:#xarray.plot.line>       |                       |
+|     3      | <inv:#xarray.plot.pcolormesh> |                       |
+|     4      | <inv:#mpl_toolkits.mplot3d.axes3d.Axes3D.voxels>| WIP |
+|     >4     | <inv:#xarray.plot.hist>       | Not fully implemented |
 
 
 ### 1D simulation
@@ -103,6 +104,10 @@ anim.show()
 
 ### 3D simulation
 
+```{warning}
+[Voxel](./plotting.md#voxel-plots) animations are not currently implemented. 
+```
+
 Opening a 3D simulation as a multi-file dataset and plotting it will
 return a <inv:#xarray.plot.hist>. However, this may not be
 desirable. We can plot a 3D simulation along a certain plane in the

diff --git a/docs/figures/voxel_plot_high_res.png b/docs/figures/voxel_plot_high_res.png
diff --git a/docs/index.md b/docs/index.md
@@ -21,6 +21,7 @@ why_sdf_xarray.md
 loading_data.md
 understanding_datasets.md
 unit_conversion.md
+plotting.md
 animation.md
 ```
 

diff --git a/docs/plotting.md b/docs/plotting.md
@@ -0,0 +1,95 @@
+---
+file_format: mystnb
+kernelspec:
+  name: python3
+---
+
+# Plotting
+
+You can plot datasets using
+[`xarray.DataArray.epoch.plot`](project:#sdf_xarray.dataarray_accessor.EpochAccessor.plot).
+This is a custom <project:#sdf_xarray> plotting routine that builds on top of
+<inv:#xarray.DataArray.plot>, so you keep the familiar <inv:#xarray> plotting
+behaviour while using <project:#sdf_xarray> conveniences (see 
+[here](project:#sdf_xarray.dataarray_accessor.EpochAccessor.plot) for details).
+Under the hood, plotting is still handled by <inv:#matplotlib>, which means you
+can use the full <inv:#matplotlib> API to customise your figure.
+
+| Dimensions | Base plotting function        |        Notes          |
+| ---------- | ----------------------------- | --------------------- |
+|     1      | <inv:#xarray.plot.line>       |                       |
+|     2      | <inv:#xarray.plot.pcolormesh> |                       |
+|     3      | <inv:#mpl_toolkits.mplot3d.axes3d.Axes3D.voxels>| Custom functionality built for <project:#sdf_xarray> |
+|     >3     | <inv:#xarray.plot.hist>       | Not fully implemented |
+
+## Line plots
+
+One dimensional data will be plotted as a line. Multiple lines can be easily
+plotted at the same time.
+
+```{code-cell} ipython3
+import sdf_xarray as sdfxr
+import matplotlib.pyplot as plt
+ds = sdfxr.open_dataset("tutorial_dataset_1d/0010.sdf")
+da_1 = ds["Derived_Number_Density_Electron"]
+da_2 = ds["Derived_Number_Density_Ion"]
+da_1.epoch.plot(label = "Electron")
+da_2.epoch.plot(label = "Ion")
+plt.legend()
+plt.show()
+```
+
+## Mesh plots
+
+Two dimensional data will be plotted as a mesh.
+
+```{code-cell} ipython3
+ds = sdfxr.open_dataset("tutorial_dataset_2d/0010.sdf")
+da = ds["Derived_Number_Density_Electron"]
+da.epoch.plot()
+plt.show()
+```
+
+## Voxel plots
+
+Three dimensional data will be plotted using voxels. This behaviour is not native
+to <inv:#xarray>, it has been custom built specifically for this package.
+
+```python
+ds = sdfxr.open_dataset("tutorial_dataset_3d/0000.sdf")
+da = ds["Derived_Number_Density"]
+da.epoch.plot(vmin = 1e27)
+plt.show()
+```
+
+![voxel_plot_high_res](./figures/voxel_plot_high_res.png)
+
+```{warning}
+Voxel plots can be extremely computationally expensive and may take longer than
+expected to plot.
+```
+
+Because voxel plots are very expensive, even with relatively small data arrays,
+plotting can be sped up by using [`xarray.DataArray.epoch.resize`](project:#sdf_xarray.dataarray_accessor.EpochAccessor.resize),
+which uses interpolation to reduce the resolution.
+
+```{code-cell} ipython3
+ds = sdfxr.open_dataset("tutorial_dataset_3d/0005.sdf")
+da = ds["Derived_Number_Density"]
+da_resized = da.epoch.resize((20, 20, 20))
+da_resized.epoch.plot(vmin = 1e27)
+plt.show()
+```
+
+## Histograms
+
+When the data array has four or more dimensions (or is specified by the user), it will
+be plotted as a histogram.
+
+```{code-cell} ipython3
+ds = sdfxr.open_dataset("tutorial_dataset_3d/0005.sdf")
+da = ds["Derived_Number_Density"]
+
+da.epoch.plot(hist = True)
+plt.show()
+```
diff --git a/docs/understanding_datasets.md b/docs/understanding_datasets.md
@@ -29,23 +29,10 @@ ds = sdfxr.open_mfdataset("tutorial_dataset_1d/*.sdf")
 ds["Electric_Field_Ex"]
 ```
 
-## Plotting
+## Visualisation
 
-You can plot datasets using
-[`xarray.DataArray.epoch.plot`](project:#sdf_xarray.dataarray_accessor.EpochAccessor.plot).
-This is a custom <project:#sdf_xarray> plotting routine that builds on top of
-<inv:#xarray.DataArray.plot>, so you keep the familiar <inv:#xarray> plotting
-behaviour while using <project:#sdf_xarray> conveniences (see 
-[here](project:#sdf_xarray.dataarray_accessor.EpochAccessor.plot) for details).
-Under the hood, plotting is still handled by <inv:#matplotlib>, which means you
-can use the full <inv:#matplotlib> API to customise your figure.
-
-```{code-cell} ipython3
-# This is discretized in both space and time
-ds["Electric_Field_Ex"].epoch.plot()
-plt.title("Electric field along the x-axis")
-plt.show()
-```
+`sdf-xarray` has built-in plotting and animation functionality, more details can
+be found on [](./plotting.md) and [](./animation.md).
 
 ## Dimension slicing
 
@@ -56,6 +43,7 @@ for easy indexing. To quickly determine the number of time steps available,
 you can check the size of the time dimension.
 
 ```{code-cell} ipython3
+ds = sdfxr.open_mfdataset("tutorial_dataset_1d/*.sdf")
 # This corresponds to the number of individual SDF files loaded
 print(f"There are a total of {ds['time'].size} time steps")
 
@@ -66,15 +54,15 @@ print(f"The time at the 20th simulation step is {sim_time:.2e} s")
 
 You can select and extract a single simulation snapshot using the integer
 index of the time step with the <inv:#xarray.Dataset.isel> function. This can be
-done by passsing the index to the `time` parameter (e.g., `time=0` for
+done by passing the index to the `time` parameter (e.g., `time=0` for
 the first snapshot).
 
 ```{code-cell} ipython3
 ds["Electric_Field_Ex"].isel(time=20)
 ```
 
 We can also use the <inv:#xarray.Dataset.sel> function if you wish to pass a
-value intead of an index.
+value instead of an index.
 
 ```{tip}
 If you know roughly what time you wish to select but not the exact value
@@ -122,6 +110,51 @@ print(f"Total particle energy absorbed: {ds["Total_Particle_Energy_in_Simulation
 print(f"The laser absorption fraction: {ds["Laser_Absorption_Fraction_in_Simulation"][-1].values:.1f} %")
 ```
 
+## Limit
+
+The `.epoch` accessor has the ability to "limit" data arrays. This is a helper
+function which applies <inv:#xarray.DataArray.where> to the upper and lower bounds
+of each dimension. You can use `None` to preserve the existing limit.
+
+```{code-cell} ipython3
+ds = sdfxr.open_dataset("tutorial_dataset_2d/0020.sdf")
+da = ds["Derived_Number_Density_Electron"]
+da.epoch.plot()
+print(f"Orignal shape: {da.shape}")
+plt.show()
+```
+
+```{code-cell} ipython3
+da_limit = da.epoch.limit(((-2e-6, 2e-6), (None, 0)))
+da_limit.epoch.plot()
+print(f"New shape: {da_limit.shape}")
+plt.show()
+```
+
+## Resize
+
+The data arrays can be interpolated with <project:#sdf_xarray.dataarray_accessor.EpochAccessor.resize>
+to either increase or decrease the resolution while mantaining the general structure.
+Decreasing the size may be useful to download large data from HPCs to local machines
+or ploting with [](./plotting.md#voxel-plots). Be aware that the original
+data will not be preserved.
+
+```{code-cell} ipython3
+ds = sdfxr.open_dataset("tutorial_dataset_2d/0020.sdf")
+da = ds["Derived_Number_Density_Electron"]
+
+da.epoch.plot()
+print(f"Orignal shape: {da.shape}")
+plt.show()
+```
+
+```{code-cell} ipython3
+da_resized = da.epoch.resize((50, 50))
+da_resized.epoch.plot()
+print(f"New shape: {da_resized.shape}")
+plt.show()
+```
+
 ## Visualisation on HPC Machines
 
 In many cases you will be running EPOCH simulations via a HPC cluster and your

diff --git a/docs/why_sdf_xarray.md b/docs/why_sdf_xarray.md
@@ -7,7 +7,7 @@ There are several benefits to using this package over the [`sdf_helper`](https:/
 - [data to numpy and pandas](https://docs.xarray.dev/en/stable/user-guide/pandas.html) - Data can be easily converted to [numpy](https://www.numpy.org/) or [pandas](https://pandas.pydata.org/) if you prefer to work with raw data.
 - [data labelling](project:understanding_datasets.md) - All SDF files come with a bunch of information about them such as the `dimensions`, `units` and `name` of the variable. In <inv:#xarray> these are added to each variable making it much easier to read.
 - [lazy/partial loading](https://docs.xarray.dev/en/stable/internals/internal-design.html#lazy-loading) - This reduces the RAM requirements when loading large SDF files as we only load the actual array values when they are needed.
-- [plotting](project:understanding_datasets.md#plotting) - You can easily plot variables across several dimensions including time using built in <inv:#matplotlib> support.
+- [plotting](project:plotting.md) - You can easily plot variables across several dimensions including time using built in <inv:#matplotlib> support.
 - [animating](project:animation.md) - Visualise your data across time or slice through another dimension and make beautiful GIFs for your presentations.
 - [deck loading](loading-input-deck) - <project:#sdf_xarray> automatically loads in the associated `input.deck` utilised when creating the simulation so that you can access setup variables.
 - [easier unit conversion](project:unit_conversion.md) - Convert dimensions and variable units to other types with the help of [`pint`](https://pint.readthedocs.io/en/stable).
diff --git a/pyproject.toml b/pyproject.toml
@@ -69,6 +69,7 @@ docs = [
   "sphinx-design>=0.7.0",
   "pooch>=1.8.2",
   "tqdm",
+  "scipy",
 ]
 build = ["cibuildwheel[uv]"]
 lint = ["ruff"]
@@ -80,6 +81,7 @@ test = [
   "tqdm",
   "pint",
   "pint-xarray",
+  "scipy",
 ]
 
 [project.entry-points."xarray.backends"]
@@ -133,6 +135,10 @@ ignore = [
   "PLR2004", # magic-comparison
   "B9",      # flake8-bugbear opinionated warnings
   "PLR0913", # remove maximum number of arguments in a function
+  "PLR0917", # remove maximum number of positional arguments in a function
+  "BLE001",  # allow capturing of base exception
 ]
-exclude = ["**/_version.py"]
-# Auto-generated file
+exclude = ["**/_version.py"] # Auto-generated file
+
+[tool.pytest.ini_options]
+filterwarnings = ['ignore::UserWarning']