yt development: star particle rendering, simple merger trees and documentation

This is the first of a new series of “what’s up with yt” blog posts I’m going to be writing. By keeping this log, I hope that maybe some things that would otherwise get lost in the version control changesets will get brought to greater light. This covers the time period of the first couple weeks in January.

Star Particle Rendering

On the mailing list, the question of adding star particles to a volume rendering was raised. The question was mostly related to the idea of adding an annotation — if you have a star particle somewhere in the simulation, for instance in galactic star formation, you likely want to annotate that somehow. After a bit of discussion, the question was raised about using star particles in the volume rendering as sources of emission; could star particles contribute to the overall image? I consulted with a couple people and then added a first-pass implementation to add them in to the brick casting. Each star particle is defined as a Gaussian, and all star particles within the 99th percentile of the center of a cell are sampled and added on to the emission of a cell during ray-casting.

This can be an expensive operation, and it comes with a few important caveats. The first is that right now, you are not necessarily going to sample the centroid of the Gaussian, and so this simply may not work correctly unless you have many, many star particles clustered together. I hope to eliminate this problem by moving to an integration method that is aware of the underlying function (erf / gaussian) rather than a pure sampling method and rectangular integration. The second caveat is that this cannot, in its current form, be used for radiation transport, as the method is not conservative of flux from star particles: even relative luminosities from star particles cannot be trusted. This is currently a pure-visualization algorithm. The solution to the second issue is not clear to me yet.

To use this, you will need to re-Cythonize your yt/utilities/amr_utils.pyx file. Then rebuild your extensions. As of right now, there are no real convenience functions for rendering star particles, so the following relatively verbose code has to be used to add them appropriately to the volume renderer. This assumed you have created a camera object.

import yt.utilities.amr_utils as au

si = (dd["creation_time"] > 0)
x = dd["particle_position_x"][si][::10]
y = dd["particle_position_y"][si][::10]
z = dd["particle_position_z"][si][::10]
age = dd["ParticleAge"][si]
age = (age - age.min())/(age.max() - age.min())
cc = iw.map_to_colors(age, "hot").squeeze()[:,:3].astype(
    "float64") / 255.0

skc = au.star_kdtree_container()
skc.add_points(x, y, z, cc)
skc.coeff = 1e-5
skc.sigma = 0.5 * pf.h.get_smallest_dx()
cam.volume.initialize_source()
for b in cam.volume.bricks: b.set_star_tree(skc)

Now when you render, it should add the star particles. In the future, once this feature leaves beta, this will be much easier to do. For a complete working example, see http://paste.enzotools.org/show/1482/ .

Here’s an example image of just star particles being rendered. It’s fuzzy and doesn’t look like much yet, but it shows the basic concepts.

Star Rendering

Simple Merger Tree

Enzo, as of 2.0, can now output Friends Of Friends catalogs on the fly during a simulation, along with particle IDs of particles belonging to those halos. Identifying mergers between the halos found with this method was not implemented, so last summer I wrote a small script to do so. I’ve added it to yt/analysis_modules/halo_merger_tree/enzofof_merger_tree.py. You can access it by amods.halo_merger_tree.EnzoFofMergerTree. This function accepts two ID numbers and loads the corresponding catalogs. It then calculates the parentage and childhood percentages of each pair of halos that it can identify.

This means that you get back a set of tuples, describing the percentage of particles from one halo that went to make up a child halo, as well as the percentage of that child halo’s particles that came from a particular parent halo. JohnW then extended this to output a graphviz plot of how the halos merged and formed over time. He also added the ability to avoid halos that are too small and a couple other feature additions.

It’s a very simple system, but it works well for the EnzoFOF calculations. In the future, I see it acting as a supplement to StephenS’s more complicated and capable merger tree identification system.

Documentation

This week my yt development was mostly focused on getting the yt documentation brought up to the current state. To that end, I re-organized it to fit into a better conceptual flow, added in the orientation section I wrote last December, and asked for feedback from a couple people. JeffO in particular gave extremely helpful comments about the structure and flow. I think that overall it’s much better. StephenS also updated his docs to reflect the new import scheme.

The documentation is now de-coupled from the version control repository. This should help keep them up to date, as we no longer have to execute a build and check that build in to update the online documentation. With this change, I think we’re now ready for release, so hopefully sometime today that will go out.

Maestro

ChrisM and JeffO committed and pushed changes to add MAESTRO support to yt. This is very exciting, as it’s another BoxLib code that seems to be supported without a lot of engine rewriting. Hopefully this will be useful to a number of people.

kD-Tree Rendering Improvements

Hi all,

Just sharing a video here that showcases some improvements I've made to the kD-tree rendering that will be making its way to yt for the 2.0 release.

(download)

Just to be clear this is showing the rendering of a cosmology simulation with a 64^3 root grid + 6 AMR levels in real time on 8 processors. The script is run in parallel, with the root processor displaying the results once each frame is finished. The viewpoint is being randomized, showing the power of a kD-tree homogenization that allows a fast back-to-front sorting algorithm for each brick. The big key here is that each processor keeps the data associated with its volume in memory so that a new viewpoint doesn't require additional file I/O.

To help get an idea of what the load balancing is doing, I figured out a way to plot the outline of the bricks with a color corresponding to each processor. This is using the breadth-first load balancing where the top N subtrees are distributed across the N processors. Some of the colors overlap in an odd way because of the order in which they are shown but you can get the general idea.

There are a few more improvements on the way such as parallel kD-tree construction which should lower the overhead for this method by quite a lot, so keep an eye out!

Figuring Out Stereo Volume Rendering

Last week I was approached by a friend and collaborator to prepare some large volume renderings using the software volume renderer in yt. In the past we've successfully made very, very large image renderings using yt -- Sam's even made one at 8192^2, although at extremely high resolution like that sometimes the lack of fidelity in the underlying volume renderer shows up; sometimes even artifacts in the AMR grid boundaries, but that's less common. Making the very large volume renderings isn't too bad -- it scales roughly with the number of pixels, but we can dispatch many frames to be rendered at once on a cluster.

There are a couple other, more important things to consider when making the big volume renderings. For starters, the entire structure of volume rendering in yt was not really created to generate a series of images -- only a single image. The idea was that you would prepare a specific image, make it, and move on. However, for this project, I want to do a zoomin, or possibly a more complicated camera path.

Additionally, one of the first things that we did with the volume rendering was silly: we applied no normalization to the output images. That was a mistake, I see now. Part of the reason for this was uncertainty in the correct normalization -- the bias that the user wanted to apply may not be the natural bias from the image. But more than that, because the rendering algorithm itself was some what holistically settled upon (the original implementation, which we used for shell-style renderings, was not a "correct" implementation of alpha blending) a natural mechanism for scaling did not immediately present itself. One likely exists, possibly dependent on the field of view, I simply do not yet know it. This will have to be rectified, because the mechanism used for scaling a set of images will have to be different than the mechanism for scaling an image in isolation, or else frames will jump in brightness during the movie's course.

The final thing that I wanted to change was to add support for stereo rendering. Rather than repeat any of the amazing discussion from Paul Bourke's website, I'll simply direct you there. Everything you ever wanted to know about stereo rendering. (When I was a first year grad student, we actually bought a copy of his site to use locally -- it was our way of showing support for him putting it online, and it also came with a bunch of source code for example applications.) I first attempted to apply the correct method for stereo, where the view direction is parallel and the total view frustum is shifted.

This did not work. In fact, it made me realize that all this time, the yt method for volume rendering is in fact ... not really a volume rendering method inasmuch as it is a planar ray-casting method. Typically when doing volume rendering, there's a perspective applied to the image: the rays all emanate from a single place, creating a frustum. But for yt, we actually set up a single plane of vectors at the back of the volume and advance that forward across the image. This is good and bad; it's good in the sense that it's more clear precisely what is going on. But it's bad in the sense that correct stereo is more difficult. (Of course, on Bourke's page he has a workaround that may work for this, but I have not yet attempted it.) Here's a rough depiction of the different between the two methods.

The upshot is that stereo doesn't seem to work unless you go with the "toe-in" method that can cause eyestrain after a long time and shows visible parallax at the edges. I'm not sure if this is going to be a problem, but because I am not right now eager to rewrite the rendering backend, this is the way it is for the moment.

To set up the stereo rendering, I separated out the rendering mechanism from the objects to be rendered. Previously, there was a single VolumeRendering object that you could create, raycast through, then discard. I created a new camera object that accepted a homogenized volume and would call "traverse" on that volume, feeding a back and a front point. The Volume is then responsible for passing off fixed-resolution grids to the camera, which accumulates an image buffer by calling the ray traversal functions. The front and back points are essentially the only thing needed to know this order, but the camera also stores its three orientation vectors and its position that describe it in 3D space. By separating out these two conceptual objects, we undo some of the "single, carefully constructed image" bias that was in the original volume renderer. (And, we open ourselves up to being able to use the Camera with a hardware volume renderer, should that day ever come.)

So now we have a camera, and it makes images like this:

It's a little dim, but that's a task for another day. The next step is taking that perspective and turning it into a set of stereo images. To do that, I added a new class called StereoPairCamera. It accepts a Camera object and turns it into two camera objects, where the final interocular distance is calculated relative to the image plane width. As I mentioned above, this only operates via toe-in stereo, so it does this in the simplest manner possible: it moves each of the Left and Right cameras by half of the interocular distance away from the original location, and then recalculates a normal vector to point back at the original center. Now we can generate left and right images:

Unfortunately, on my laptop (which is my primary/exclusive work computer) I don't have the ability to view these pairs. To get around that, I wrote a simple stereo pair image viewer in OpenGL and imposed upon my friends that do have a stereo viz wall to test it out -- and after some fiddling with the interocular distance, we got what appeared to be workable stereo pairs.

The full code for generating camera paths as well as stereo pairs is already in yt, but the documentation is still being written; I might also clean up the interface a bit. Additionally, at some point in the future, the issue of toe-in stereo versus correct parallel-frustum stereo will need to be dealt with; the last thing I really want to do is force people to only use a bad method for generating stereo pairs. Hopefully that is something that can be dealt with at a later time.

Thanks to the wealth of resources out there for making this a relatively easy task: the aforementioned Paul Bourke website on stereo pairs, the PyOpenGL and PIL teams for making the image pair viewer easy, and everyone else whose work I've built on to make things like this.

Filed under: viz

Star Particle Rendering

Simple Merger Tree

Documentation

Maestro