1{
2 "cells": [
3  {
4   "cell_type": "markdown",
5   "metadata": {},
6   "source": [
7    "Here, we explain how to use TransferFunctionHelper to visualize and interpret yt volume rendering transfer functions.  Creating a custom transfer function is a process that usually involves some trial-and-error. TransferFunctionHelper is a utility class designed to help you visualize the probability density functions of yt fields that you might want to volume render.  This makes it easier to choose a nice transfer function that highlights interesting physical regimes.\n",
8    "\n",
9    "First, we set up our namespace and define a convenience function to display volume renderings inline in the notebook.  Using `%matplotlib inline` makes it so matplotlib plots display inline in the notebook."
10   ]
11  },
12  {
13   "cell_type": "code",
14   "execution_count": null,
15   "metadata": {
16    "collapsed": false
17   },
18   "outputs": [],
19   "source": [
20    "import yt\n",
21    "import numpy as np\n",
22    "from IPython.core.display import Image\n",
23    "from yt.visualization.volume_rendering.transfer_function_helper import TransferFunctionHelper\n",
24    "from yt.visualization.volume_rendering.render_source import VolumeSource\n",
25    "\n",
26    "def showme(im):\n",
27    "    # screen out NaNs\n",
28    "    im[im != im] = 0.0\n",
29    "    \n",
30    "    # Create an RGBA bitmap to display\n",
31    "    imb = yt.write_bitmap(im, None)\n",
32    "    return Image(imb)"
33   ]
34  },
35  {
36   "cell_type": "markdown",
37   "metadata": {},
38   "source": [
39    "Next, we load up a low resolution Enzo cosmological simulation."
40   ]
41  },
42  {
43   "cell_type": "code",
44   "execution_count": null,
45   "metadata": {
46    "collapsed": false
47   },
48   "outputs": [],
49   "source": [
50    "ds = yt.load('Enzo_64/DD0043/data0043')"
51   ]
52  },
53  {
54   "cell_type": "markdown",
55   "metadata": {},
56   "source": [
57    "Now that we have the dataset loaded, let's create a `TransferFunctionHelper` to visualize the dataset and transfer function we'd like to use."
58   ]
59  },
60  {
61   "cell_type": "code",
62   "execution_count": null,
63   "metadata": {
64    "collapsed": false
65   },
66   "outputs": [],
67   "source": [
68    "tfh = TransferFunctionHelper(ds)"
69   ]
70  },
71  {
72   "cell_type": "markdown",
73   "metadata": {},
74   "source": [
75    "`TransferFunctionHelpler` will intelligently choose transfer function bounds based on the data values.  Use the `plot()` method to take a look at the transfer function."
76   ]
77  },
78  {
79   "cell_type": "code",
80   "execution_count": null,
81   "metadata": {
82    "collapsed": false
83   },
84   "outputs": [],
85   "source": [
86    "# Build a transfer function that is a multivariate gaussian in temperature\n",
87    "tfh = TransferFunctionHelper(ds)\n",
88    "tfh.set_field(('gas', 'temperature'))\n",
89    "tfh.set_log(True)\n",
90    "tfh.set_bounds()\n",
91    "tfh.build_transfer_function()\n",
92    "tfh.tf.add_layers(5)\n",
93    "tfh.plot()"
94   ]
95  },
96  {
97   "cell_type": "markdown",
98   "metadata": {},
99   "source": [
100    "Let's also look at the probability density function of the `mass` field as a function of `temperature`.  This might give us an idea where there is a lot of structure. "
101   ]
102  },
103  {
104   "cell_type": "code",
105   "execution_count": null,
106   "metadata": {
107    "collapsed": false
108   },
109   "outputs": [],
110   "source": [
111    "tfh.plot(profile_field=('gas', 'mass'))"
112   ]
113  },
114  {
115   "cell_type": "markdown",
116   "metadata": {},
117   "source": [
118    "It looks like most of the gas is hot but there is still a lot of low-density cool gas.  Let's construct a transfer function that highlights both the rarefied hot gas and the dense cool gas simultaneously."
119   ]
120  },
121  {
122   "cell_type": "code",
123   "execution_count": null,
124   "metadata": {
125    "collapsed": false
126   },
127   "outputs": [],
128   "source": [
129    "tfh = TransferFunctionHelper(ds)\n",
130    "tfh.set_field(('gas', 'temperature'))\n",
131    "tfh.set_bounds()\n",
132    "tfh.set_log(True)\n",
133    "tfh.build_transfer_function()\n",
134    "tfh.tf.add_layers(8, w=0.01, mi=4.0, ma=8.0, col_bounds=[4.,8.], alpha=np.logspace(-1,2,7), colormap='RdBu_r')\n",
135    "tfh.tf.map_to_colormap(6.0, 8.0, colormap='Reds')\n",
136    "tfh.tf.map_to_colormap(-1.0, 6.0, colormap='Blues_r')\n",
137    "\n",
138    "tfh.plot(profile_field=('gas', 'mass'))"
139   ]
140  },
141  {
142   "cell_type": "markdown",
143   "metadata": {},
144   "source": [
145    "Let's take a look at the volume rendering. First use the helper function to create a default rendering, then we override this with the transfer function we just created."
146   ]
147  },
148  {
149   "cell_type": "code",
150   "execution_count": null,
151   "metadata": {
152    "collapsed": false
153   },
154   "outputs": [],
155   "source": [
156    "im, sc = yt.volume_render(ds, [('gas', 'temperature')])\n",
157    "\n",
158    "source = sc.get_source()\n",
159    "source.set_transfer_function(tfh.tf)\n",
160    "im2 = sc.render()\n",
161    "\n",
162    "showme(im2[:,:,:3])"
163   ]
164  },
165  {
166   "cell_type": "markdown",
167   "metadata": {},
168   "source": [
169    "That looks okay, but the red gas (associated with temperatures between 1e6 and 1e8 K) is a bit hard to see in the image. To fix this, we can make that gas contribute a larger alpha value to the image by using the ``scale`` keyword argument in ``map_to_colormap``."
170   ]
171  },
172  {
173   "cell_type": "code",
174   "execution_count": null,
175   "metadata": {
176    "collapsed": false
177   },
178   "outputs": [],
179   "source": [
180    "tfh2 = TransferFunctionHelper(ds)\n",
181    "tfh2.set_field(('gas', 'temperature'))\n",
182    "tfh2.set_bounds()\n",
183    "tfh2.set_log(True)\n",
184    "tfh2.build_transfer_function()\n",
185    "tfh2.tf.add_layers(8, w=0.01, mi=4.0, ma=8.0, col_bounds=[4.,8.], alpha=np.logspace(-1,2,7), colormap='RdBu_r')\n",
186    "tfh2.tf.map_to_colormap(6.0, 8.0, colormap='Reds', scale=5.0)\n",
187    "tfh2.tf.map_to_colormap(-1.0, 6.0, colormap='Blues_r', scale=1.0)\n",
188    "\n",
189    "tfh2.plot(profile_field=('gas', 'mass'))"
190   ]
191  },
192  {
193   "cell_type": "markdown",
194   "metadata": {},
195   "source": [
196    "Note that the height of the red portion of the transfer function has increased by a factor of 5.0. If we use this transfer function to make the final image:"
197   ]
198  },
199  {
200   "cell_type": "code",
201   "execution_count": null,
202   "metadata": {
203    "collapsed": false
204   },
205   "outputs": [],
206   "source": [
207    "source.set_transfer_function(tfh2.tf)\n",
208    "im3 = sc.render()\n",
209    "\n",
210    "showme(im3[:,:,:3])"
211   ]
212  },
213  {
214   "cell_type": "markdown",
215   "metadata": {},
216   "source": [
217    "The red gas is now much more prominent in the image. We can clearly see that the hot gas is mostly associated with bound structures while the cool gas is associated with low-density voids."
218   ]
219  }
220 ],
221 "metadata": {
222  "kernelspec": {
223   "display_name": "Python 3",
224   "language": "python",
225   "name": "python3"
226  },
227  "language_info": {
228   "codemirror_mode": {
229    "name": "ipython",
230    "version": 3
231   },
232   "file_extension": ".py",
233   "mimetype": "text/x-python",
234   "name": "python",
235   "nbconvert_exporter": "python",
236   "pygments_lexer": "ipython3",
237   "version": "3.5.1"
238  }
239 },
240 "nbformat": 4,
241 "nbformat_minor": 0
242}
243