# Spherical Volume Rendering

Components:

* UI - linked histogram and transfer function editor
* Algorithm - fast traversal over voxels
* Implementation - which voxels are next
* Application - using it on real data

## User Interface

Typically, volume rendering in yt is conducted in one of three ways:

* Strict integration: $\int_{t_0}^{t_1} f(t) dt$ where the values along a ray are simply summed along the path
* Some type of averaging or functional integration, such as $\mathrm{max}(f(t))$ or a weighted average, such as $\frac{\int_{t_0}^{t_1}w(t)f(t)dt}{\int_{t_0}^{t_1}w(t)dt}$
* An integration of the radiative transfer equation, in the limit of no scattering, where the value $f(t)$ is actually computed as through lookup tables.

Developing the lookup tables, or "transfer functions," is typically conducted by defining values in two dimensions for the red, green, blue and alpha channels.

The yt "transfer function helper" shows distribution of values in the simulation over top of a sample transfer function:

https://yt-project.org/doc/visualizing/transfer_function_helper.html

This often looks something like:


Which produces an image like:


(It is low resolution -- ours will look better!)

What would be a vast improvement would be a way to draw on the transfer function to define its values, with interactive updates of the image. This could be done in the web interface.

An additional component would be to build in methods for rotating the camera position, using something like [threejs](http://threejs.org/).

An example mockup might look like:


`ColorTransferFunction` - which has `r`, `g`, `b`, `a` values, which are essentially wrappers around numpy arrays of fixed length.

To pass back from Javascript, we could accept something that is isomorphic with that, for instance:

`[1.0, 1.0, ...], [0.1, 0.2, ...], [...], [...]`

or something like:

`[ {'r': .., 'b': .., 'g': .., 'a': ..}, ... ]`

`ColorTransferFunction.red.y, ColorTransferFunction.blue.y, ...`

## Algorithm

[Amanatides & Woo](http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.42.3443&rep=rep1&type=pdf) have developed a method for fast voxel traversal in cartesian coordinates.

What we need for this project is a comparative example for spherical or cylindrical geometry.

Essentially, our grids are organized in r, theta and phi. Given a ray's origin $(x,y,z)$ (note that this is functionally equivalent to an $(r, theta, phi)$ coordinate, but we typically will have it specified in cartesian coordinates) and a given set of cells, we want to know which cell will be entered next, with a minimum number of floating point operations.


The cells will be organized such that inside a given collection of cells, they will all have identical edges such as $r_0, r_1, r_2, r_3, ...$ and $\phi_0, \phi_1, \phi_3, ...$ and $\theta_0, \theta_1, \theta_2, ...$. Typically these will all have equal spacing, $\Delta r, \Delta \phi, \Delta \theta$ but this is not universally the case.

If the spacings are identical, we can establish a mapping between integer indices and cell edges:

$$\phi_L = \phi_0 + i \Delta \phi $$ $$\phi_R = \phi_0 + (i + 1) \Delta \phi $$

and similarly for the others, with $j$ and $k$. This establishes a mapping between integer coordinates and cell boundaries. In Cartesian coordinates, the mapping between $i,j,k$ and the coordinates of the cells is straightforward and linearly defined, so the mapping between increments in $i,j,k$ and increments in ray position are linearly defined. However, in spherical and cylindrical coordinates, these mappings are non-linear.

Given our vector direction and our cell spacings, can we *a priori* establish an ordering of cells and determine passage of rays between cells?

## Implementation

Given this mapping, we need to develop a "ray traveral" implementation that takes the values for the grid spacing and ray location and quickly samples and quickly interpolates within a cell.

The intra-cellular interpolation is key, as it requires interpolation in a spherical shell where the values are defined at vertices. What is the fastest way to compute an interpolation kernel and to define an ordering of cell indices?

## Application

This can then be applied to simulations of supernovae, accretion disks, geodynamics and seismology.

## Implementation on GPU

Given our implementation in C/C++, we can then convert this to GLSL shader language.