Realtime Volume Visualization of Time Evolving Nebulae for Immersive Environments

Department of Science and Technology

Institutionen för teknik och naturvetenskap

Linköpings Universitet

Linköpings Universitet

SE-601 74 Norrköping, Sweden

601 74 Norrköping

Examensarbete

LITH-ITN-MT-EX--06/052--SE

Realtime Volume Visualization

of Time Evolving Nebulae for

Immersive Environments

Magnus Bernhardsson

Christian Sörensen

LITH-ITN-MT-EX--06/052--SE

Realtime Volume Visualization

of Time Evolving Nebulae for

Immersive Environments

Examensarbete utfört i medieteknik

vid Linköpings Tekniska Högskola, Campus

Norrköping

Magnus Bernhardsson

Christian Sörensen

Handledare Staffan Klashed

Examinator Anders Ynnerman

Rapporttyp Report category Examensarbete B-uppsats C-uppsats D-uppsats _ ________________ Språk Language Svenska/Swedish Engelska/English _ ________________ Titel Title Författare Author Sammanfattning Abstract ISBN _____________________________________________________ ISRN _________________________________________________________________

Serietitel och serienummer ISSN

Title of series, numbering ___________________________________

Nyckelord

Keyword

Datum

Date

URL för elektronisk version

Avdelning, Institution

Division, Department

Institutionen för teknik och naturvetenskap Department of Science and Technology

2006-11-09

x

x

LITH-ITN-MT-EX--06/052--SE

Realtime Volume Visualization of Time Evolving Nebulae for Immersive Environments

Magnus Bernhardsson, Christian Sörensen

This paper presents the result of our master thesis work in media technology, which is an

implementation of a realtime volume rendering engine capable of visualizing reflection nebulae. The engine is using scattering and emission properties of time evolving volumetric data, created by astrophysicists at the Hayden Planetarium, New York. The visualization is designed towards full dome planetarium displays. Using the latest GPU based techniques for volume rendering and impostor techniques we are able to interactively fly around as well as animate several nebulae embedded in their true astronomical content.

This paper is aimed towards fellow media technology students and people interested in the area of realtime volume renderings. It is written with the presumption that the reader has the same knowledge as we did when starting our thesis work, which implies good confidence in modern 3D computer graphics, some experience in object-oriented programming and basic knowledge of algebra and calculus. It is our hope that this report can be enjoyed even if skipping some of the more technical parts.

Detta dokument hålls tillgängligt på Internet – eller dess framtida ersättare –

under en längre tid från publiceringsdatum under förutsättning att inga

extra-ordinära omständigheter uppstår.

Tillgång till dokumentet innebär tillstånd för var och en att läsa, ladda ner,

skriva ut enstaka kopior för enskilt bruk och att använda det oförändrat för

ickekommersiell forskning och för undervisning. Överföring av upphovsrätten

vid en senare tidpunkt kan inte upphäva detta tillstånd. All annan användning av

dokumentet kräver upphovsmannens medgivande. För att garantera äktheten,

säkerheten och tillgängligheten finns det lösningar av teknisk och administrativ

art.

Upphovsmannens ideella rätt innefattar rätt att bli nämnd som upphovsman i

den omfattning som god sed kräver vid användning av dokumentet på ovan

beskrivna sätt samt skydd mot att dokumentet ändras eller presenteras i sådan

form eller i sådant sammanhang som är kränkande för upphovsmannens litterära

eller konstnärliga anseende eller egenart.

För ytterligare information om Linköping University Electronic Press se

förlagets hemsida

http://www.ep.liu.se/

Copyright

The publishers will keep this document online on the Internet - or its possible

replacement - for a considerable time from the date of publication barring

exceptional circumstances.

The online availability of the document implies a permanent permission for

anyone to read, to download, to print out single copies for your own use and to

use it unchanged for any non-commercial research and educational purpose.

Subsequent transfers of copyright cannot revoke this permission. All other uses

of the document are conditional on the consent of the copyright owner. The

publisher has taken technical and administrative measures to assure authenticity,

security and accessibility.

According to intellectual property law the author has the right to be

mentioned when his/her work is accessed as described above and to be protected

against infringement.

For additional information about the Linköping University Electronic Press

and its procedures for publication and for assurance of document integrity,

please refer to its WWW home page:

http://www.ep.liu.se/

Acknowledgments

We would like to thank all the people at Sciss; Per Hemmingsson, Staffan Klashed, Martin Rasmusson and Johan Öhlund, for all their help and for putting up with our questions. A big thanks goes out to dr. Carter Emmart at AMNH for his inspirational ideas and for giving us an unforgettable time in New York. We would also like to thank prof. Marcus Magnor and his team at the MPI for the feedback on our implementations of their research, and Mordecai-Mark Mac Low at AMNH for providing us with the dataset used in our final visualizations of nebulae.

Finally we would like to thank Janet Asimov, the fund of Anna Whitlock and Anders Ynnerman for providing grants which made our visit to New York possible.

Contents

1 Introduction 1

1.1 Background . . . 1

1.2 Purpose and objectives . . . 1

1.3 Limitations . . . 2

1.4 Outline of report . . . 2

2 Background 5 2.1 Uniview . . . 5

2.2 Nebulae . . . 6

2.2.1 Colors of the Nebulae . . . 8

2.2.2 Creating Nebulae Data for Volume Visualization . . . 8

2.3 Volume rendering . . . 9

2.3.1 Raycasting . . . 10

2.4 GPU-Raycasting . . . 13

2.4.1 The GPU today . . . 13

2.4.2 The Raycasting Technique . . . 14

2.5 Rendering of Emission and Reflection Nebulae . . . 17

2.5.1 Illumination . . . 17 2.5.2 Scattering . . . 17 2.5.3 Observed Illumination . . . 18 2.5.4 Wavelength dependencies . . . 18 2.5.5 Multiscattering . . . 19 3 Implementation 21 3.1 Environment . . . 21 3.2 Design . . . 21 3.2.1 VolumeObject Class . . . 21

3.2.2 VolumeDataHandler and VolumeDataThread Classes . . . . 22

3.2.3 RayCaster and FBOHandler Classes . . . 22

3.3 Data . . . 23 3.3.1 Pre-processing . . . 23 3.3.2 Data management . . . 24 3.3.3 Nebula generator . . . 24 3.4 Rendering . . . 24 ix

3.4.2 Multiscattering . . . 26

3.4.3 Post-processing . . . 26

3.4.4 Animation of volumetric data . . . 27

3.5 Rendering when inside the volume . . . 28

3.6 Impostors . . . 28

3.7 Graphical User Interface . . . 29

4 Results 31 4.1 Visualization . . . 31

4.2 Animation of data . . . 31

4.3 Data format . . . 31

4.4 Interaction with surroundings . . . 32

4.5 Performance . . . 32 4.6 Limitations . . . 32 5 Discussion 35 5.1 Visualization . . . 35 5.2 Data . . . 35 5.3 Performance . . . 36 5.3.1 Preprocessing . . . 36 5.3.2 Rendering . . . 36 5.3.3 Animation . . . 37 5.4 Future Work . . . 37 Bibliography 39 A Figures 41

Chapter 1

Introduction

1.1

Background

During the last years the technique of visualizing volumetric datasets has been extensively researched. The most common area of application is medical visual-ization. When considering the advances in recent graphics programming research, volume rendering can be applied in several other fields as well.

Sciss AB, a software development company out of Norrköping, are currently developing a visualization platform called Uniview for use in immersive environ-ments such as planetariums and other dome displays. The platform is created for the purpose of visualizing and studying astronomical content, for example planets, stars and galaxies. One interesting phenomenon in space which Uniview has not been able to visualize are nebulae, clouds of gas and dust, spanning over several light years. There are many forms of nebulae and their appearance varies from black patches to colorful cloud formations. Since nebulae consist of gas, an appro-priate method for visualizing them could be volume rendering as it describes how light travels through transparent media.

1.2

Purpose and objectives

In the fall of 2005, Sciss AB decided to extend their visualization platform Uniview with the ability to visualize volumetric datasets in real-time. The purpose of this thesis work has been to implement this volume render. Our objectives were to find, implement and evaluate a method for visualizing volumetric data. The volume rendering engine had to reach the following goals:

• Be able to display one channel volume datasets, using transfer functions, at interactive frame rates for different resolutions and quality.

• Be fully integrated into Uniview, embedding it into an astronomical context and creating the possibility for displaying graphics on planetarium domes and other convex displays.

• Handle volumetric datasets of sizes up to 10243_{, using some streaming and}

level of detail techniques if necessary.

• Create and interact with volumetric transfer functions during runtime. • Visualize nebulae, based on the scattering and emissive properties of the gas

of which they consist.

After achieving the first two goals early into the thesis work, we realized that the last three would require too much time to fulfill. Deciding that the most interesting feature for Uniview was the last goal, visualization of nebulae, we agreed to rewrite the objectives into the following, which would suit our interests better.

• Visualize one or several nebulae in real-time based on the emission and scat-tering properties of the gas of which they consist.

• Be able to interact with the nebulae, flying through them and adjusting their physical parameters.

• Visualizing the evolution of a nebula, if we were to be given access to an animated volumetric dataset.

1.3

Limitations

The limitations of this thesis work were primarily defined by the available and sup-ported hardware. Older shader profiles have limited support for some of the more advanced rendering instructions necessary for our implementations. We have made an effort to optimize and support most common hardware used by the customers of Sciss AB.

Our volume renderer would be integrated into an existing scientific visualization framework, and later become a part of a product. This implied several limitations, mostly concerning the amount of shared CPU and GPU resources we could utilize for our algorithms. The structure and layout of our code was also limited by some of the existing designs in Uniview.

Because of the limitations mentioned above, we first gave the esthetical aspects of the visualization a lower priority, focusing on writing efficient algorithms and code. However, when our work was coming to an end and the core algorithms had been implemented, our focus shifted towards creating a visually better result.

1.4

Outline of report

This report is divided into four parts. The first part after this introduction is ex-tensive but necessary for the understanding of the astronomical and mathematical concepts in our thesis. We begin with describing the application Uniview, which has been the basis for our work, followed by a brief explanation of the physics of nebulae and some techniques for generating volumetric nebula data. We give an overview of the graphical processing unit, and introduce the concept of volume

1.4 Outline of report 3

rendering and how it can be performed on the GPU. However, our focus lies on how the volume rendering method can be used for visualizing nebulae. This part requires some knowledge of calculus, but is important for the understanding of the rest of the thesis.

In the next chapter, describing our implementation, we will start off by defining the tools we have used, followed by an in-depth explanation of the implementations of the volume rendering engine in Uniview. We also describe some of the underlying class-structure and coding techniques that we have used.

The two final chapters describe and discuss our results. We present which goals have been achieved and which have not been fulfilled and why. We evaluate the performance of our volume rendering engine and discuss the problems we have run into, before finishing off with an outline of possible future work.

Chapter 2

Background

2.1

Uniview

The beginning of this thesis can be traced back to the year 2001, when a collabo-ration between the American Museum of Natural History (AMNH) and Linköping University was founded. It began when Anders Ynnerman, professor at Linköping University, met doctor Carter Emmart, a director of visualization at AMNH. They decided to start collaborating on their common interest, real time visualization. The same year it was decided that a master student, Staffan Klashed was to be sent to the Hayden Planetarium, a part of AMNH, to perform his thesis work. He cre-ated a visualization of astronomical data over the universe. AMNH provided him with their comprehensive database, The Digital Universe, containing a majority of the known objects in the universe today. The data is collected by astronomers at AMNH, revieving it from several different mappings of the universe. Klashed´s work turned out so well that the following year two other master students, Per Hemmingsson and Martin Rasmusson, continued on his efforts to create an even better visualization platform. This platform was called Uniview.

Their combined work resulted in the forming of the company Sciss AB, located in Norrköping Sweden. Today the company works with the ongoing development of the visualization software Uniview. The application has become a powerful tool for visualizing astronomic objects from the Digital Universe, as well as other data resources.

The most interesting feature of Uniview is that it is especially developed to be able to render on spherical or domed displays, which can be found in many planetariums, hence the collaboration with the Hayden Planetarium. Today, the company is further developing Uniview in many areas, and one is the visualization of volumetric data.

2.2

Nebulae

Many of the pictures taken by the Hubble telescope shows colorful objects, looking like gigantic clouds of glowing gas. These objects are called nebulae and are by many observers the most beautiful sites in the universe. Nebulae are consisting of gas and interstellar dust, remainders from when stars either were created or did collapse. They are categorized as follows:

• Planetary nebulae. If you study these through a telescope you can under-stand where they got their name, they look circular just like a planet. A bright star at its center makes it shine. A planetary nebula which is of-ten regarded as the archetype of planetary nebulae is the Ring Nebula, see figure 2.1.

• Diffuse nebulae can span over several light-years, looking like ragged and sinister clouds. They are also lit up by stars in their inner cavities, but contain far more stars than the planetary nebulae. See figure 2.1.

• Nebulae formed by gas and dust left from when a dying star explodes, so called supernovas. An example is the Crab Nebula, see figure 2.3, already discovered in the 1000th century.

• Dark nebulae or obscuration nebulae contains gas with high density which blocks out all light from stars and therefore appear as dark patches on the sky. As one can imagine these are very hard to spot and study. See figure 2.4 for an image of the Horsehead nebula.

The diffuse nebulae, which are the important ones for this report, are further divided into two different types depending on how they are illuminated:

• Reflective nebulae are the most straightforward type. The reflective gas is simply lit up by stars and the gas has no light of its own, without the stars you would not be able to see it.

• Emission nebulae, have, as the name implies, their own emission aside from the reflectance from the stars. The hot gas glows because of high amounts of energy, coming from the central stars. The energy works as an ionization of the gas, photons knocking electrons out of their orbit, which then emit light when falling back. The ionization is often concentrated to a band or rim, called the ionization front, along the outskirts of the visible part of the nebula and is slowly expanding outwards. A common example of an emission nebula is the Orion Nebula, located in the constellation of Orion.

If the gas in a nebula starts to contract and gets denser, the pressure and temperature slowly starts rising. If the temperature gets high enough, the gas ignites, and in a chain reaction the surrounding gas forms a big globe of exploding gas, giving birth to a new star. These newborn stars are extremely hot and therefore the nebulae surrounding these birthplaces are the brightest ones.

2.2 Nebulae 7

Figure 2.1. Planetary nebula (The Ring)

Figure 2.2. Diffuse nebula (Orion

Neb-ula)

Figure 2.3. Supernova (Crab Nebula)

Figure 2.4. Dark nebula (Horsehead

2.2.1

Colors of the Nebulae

The colorful images seen in photos, for example from the Hubble telescope, can sometimes be deceiving. Most of the light sent out from the nebulae is actually infrared or ultraviolet and thus invisible to the human eye. In order to be able to visualize this infrared light sent out from dying stars or the ultraviolet light from newly born stars, these different wavelengths are translated into colors in the visible range of the spectrum. Looking at the three images taken of the Helix Nebula (figure 2.5) different types of color interpretations are clearly visible. Often the colors are chosen to bring out and enhance the level of detail and features of the nebula that the astronomers study. For example, in the image of the Eagle Nebula (figure A.1), the colors represent the appearance of sulphur ions, hydrogen and doubly-ionized oxygen in the red, green and blue channel respectively.

This fact has one major implication on our thesis work, implementing a physi-cally correct optical model of a nebula will not result in renderings looking like the color-enhanced images of nebulae taken by the Hubble telescope. When looking at a reflection nebula, the only thing that the naked eye would see is the effect of the light scattering in the dust. The scattering has the result that a reflection nebula gets a faint blue color. As blue light is of a shorter wavelength, it is more prone to forward scattering compared to the red light, which is cast back into the rear of the nebula.

Figure 2.5. Three different images taken of the Helix Nebula (NGC 7293).

2.2.2

Creating Nebulae Data for Volume Visualization

Since nebulae essentially are big clouds of gas with large cavities in the middle, the process of creating volume data would be to create clouds with some modeling constraints mimicking the solar winds from the stars, which blow the gas away from the center. However, this has shown to be more difficult than expected. There are some recent and impressive techniques for easily creating clouds for real-time applications [13], see figure 2.6. These approaches are very popular and widely used in the computer games of today, but they do not create the three-dimensional datasets that we need.

A better approach is to generate three-dimensional noise, for example Perlin-based [11], and tweak it into looking like nebulae, which has been done by Magnor

2.3 Volume rendering 9

Figure 2.6. Procedural clouds. Figure 2.7. Perlin noise generated

neb-ula.

Figure 2.8. Simulated nebula.

et al. [8], and this produced some really appealing results, see figure 2.7. Pro-cedural noise is however difficult to easily tweak, and a better way would be to use image based modeling, creating 3D data from analyzing a set of 2D images. Some attempts of this has been made on galaxies [5], but the 2D images used are of course limited to only one viewpoint and this is not enough when using image based modeling.

To create realistic volume datasets, the best way is to use real physical simula-tions of how the gas in nebulae behaves. Astronomers at AMNH have made some good efforts at this, creating a highly advanced simulation of how a high concen-tration of gas expands away from a central, high energy star, and also calculating at what rate different species are being ionized. Fortunately, astrophysicist Mac Low at AMNH provided us with samples from one of their latest simulations which we could use in our work, see figure 2.8.

2.3

Volume rendering

Visualizing the exact inner structures of nebulae, as well as imitate their complex lighting conditions, requires a powerful rendering method which is easy to tweak

and can handle transparent objects. One approach of doing just that is volume rendering.

Volume renderings are today used for the visualization of medical data, but could be applied in many other areas, such as scientific simulations in metrology, geology and of course astronomy. These simulations produce larger sets of sampled values in a 3D structure, called voxels or volume data elements. This is the same type of dataset that we received from AMNH, created by a scientific simulation.

Often, an approach that renders the data in a simplified way, using slices or other 2D projection techniques, is adopted instead of using a volume rendering approach. The reasons why volume rendering is not used are issues with poor interactivity and quality, but during the last years this has started to change.

The most common techniques for advanced volume renderings are raycasting and proxy geometry. A lot of recent volume rendering research has contributed in several improvements to both these techniques. These improvements have made it possible to render high-quality volumes at interactive frame rates [6].

2.3.1

Raycasting

Raycasting is a rendering method using the principle of how rays of light, which are being cast from an observer, travels through volumetric data. These rays sample the dataset and transform it into a two dimensional projection suitable for displaying on a screen.

For each pixel in the image that is to be rendered, a ray is cast out from the position of the observer. All these rays sample the volume dataset at equal distances, see figure 2.9. The rays will practically never sample the dataset at the exact positions of the voxels, hence a interpolation method is needed. Often a trilinear interpolation is used, taking into account neighbouring voxels. These sampled voxel values build up to the final pixel values as the rays traverse through the whole dataset. When adding new values to the final image, the raycasting method is using the volume rendering integral 2.1. This equation is based on the physics of light and the law of energy conservation. It is the basis for all volume rendering methods, describing how light travels through transparent media, taking both absorption and scattering into account.

I(D) = I0e− RD 0 τ (t)dt+ D Z 0 L(s)τ (s)e−RsDτ (t)dtds (2.1)

We will not go into the derivation of this algorithm, but dividing it into parts, reveals a logical structure. The total amount of outgoing intensity I(D), going through a volume of thickness D, equals the sum of the initial intensity I0, and

the emission of each voxel L(s). The initial intensity is attenuated through the volume, with an attenuation factor τ at each point through the volume. The emission of each voxel is both attenuated by the voxels own attenuation factor, and by the voxels lying on the remaining path between the voxel and the viewer. The method to solve this equation during rendering is to use something called

2.3 Volume rendering 11

Figure 2.9. Sampling the volume.

alpha blending. When traversing the volume data in a front-to-back order, the blending at each step through the volume follows the two equations below.

Couti= Couti−1+ 1 − αouti−1
αsrciCsrci (2.2)

αouti = αouti−1+ 1 − αouti−1
αsrci (2.3)

The equations describe the resulting color Couti after each ray step i. It is the

sum of the color from the last step Couti−1, and the color at the current position

along the ray Csrci, attenuated by the voxel itself αsrci, and the accumulated

attenuation from the viewer to the voxel 1 − αouti−1. The same relationship holds

for the opacity values 2.3.

The distance between the centres of two voxels is the sampling rate of the dataset. The Nyqvist theorem states that the sampling rate need to be at least twice the highest frequency in the signal we wish to sample, if it is to be reproduced entirely. However, in many cases that is impossible and would require enormous amounts of data, thus making sampling aliasing artefacts unavoidable in volume renderings. Also, when tracing the rays through the dataset during rendering, they need to sample data at, for example 512 points to fully reproduce a 2563

dataset. This is easier to achieve and there even exist optimizations making a double sampling rate superfluous [2].

An alternative to a full raycasting solution is iso-surfaces. Here, a surface is constructed from terminating the rays when they find a voxel with a value lying above a certain threshold. These surfaces can be lit, shaded and are a powerful way of visualizing isolated objects in a volumetric dataset.

Proxy Geometry

The easiest way of implementing interactive volume renderings is to use a tech-nique called proxy geometry. This techtech-nique consists of rendering geometry in the

Figure 2.10. Parallel projection. Figure 2.11. Perspective projection.

position of the volume, in the simplest case parallel planes, and using these planes to render the volumetric data.

By storing the volumetric data in a 3D texture, the texture can be directly applied to the proxy planes. Using the built in trilinear interpolation on the graphic cards of today, each pixel fragment gets a correct interpolated texture value. By finally blending the planes/fragments with the graphics API blending equations 2.4 2.5, the volumetric rendering integral is automatically calculated.

Changing the standard blending equation, several different renderings can be achieved. An example is Maximum Intensity Projection (MIP) that shows only the sample along each ray through the volume with the highest value, thus imitating an x-ray image.

Cdst= (1 − αsrc) Cdst+ αsrcCsrc (2.4)

αdst= (1 − αsrc) αdst+ αsrcαsrc (2.5)

There are two ways of setting up the proxy plane geometry, object-aligned or view-aligned, both having advantages and disadvantages. The easiest method is to render the parallel planes object-aligned. This geometry is static, meaning that the planes are always rendered in one position, not depending on the position of the camera. This will only work until the camera is moved too far to the side of the volume. When this happen the planes have to be repositioned in order to not become parallell towards the camera, and therefore invisible.

A method that produces renderings with higher quality, but involving more computations, are the view-aligned technique. In this method the planes are al-ways drawn parallel to the viewer. Keeping them this way will alal-ways produce a consistent volume rendering, which is not the case when using object-aligned planes since the rearrangement of planes often introduces a noticeable transition. There are however still some problems connected to the use of proxy planes. Looking at images 2.10 and 2.11 the difference between a parallel and a perspective projection of the rays, can be seen. A perspective projection does not keep the sampling distances equal, introducing noticeable distortions at the edges of the volume when using a wide field of view.

To solve this problem, a more complex proxy geometry than planes can be used. One example is spherical proxy-shells, where the proxies are drawn as sections of

2.4 GPU-Raycasting 13

spheres growing bigger the farther away they get from the viewer, see figure 2.12. The disadvantage of using spherical shells for proxy geometry is that the shells includes more computations to create than planes, but in a dome environment, demanding an extreme field of view, it is the only solution if using proxy geometry.

Figure 2.12. Proxy Shells.

2.4

GPU-Raycasting

There are advantages and disadvantages to both methods. A software approach gives the programmer more control, easy data access and optimization algorithms such as early ray termination are easier to achieve than in a hardware approach. The accelerated method takes advantage of the parallelism of modern graphics cards, having pipelines for simultaneous processing of several vertices and frag-ments. This translates into processing all the rays cast through the volume, at the same time, which leads to increased performance.

While there are a few graphic cards that are specialized for accelerated ray-casting in development, up to a couple of years ago there was no alternative to the software-based approach. With the technique presented by Krüger and West-ermann [6], raycasting was however made available for advanced real-time volume applications on most of today’s hardware.

There are advantages and disadvantages to both methods. A software approach gives the programmer more control and easy access to the volume data. This makes it easier to implement optimization algorithms, such as early ray termination, than in a hardware approach. The accelerated method takes advantage of the parallelism of modern graphics cards, having pipelines for simultaneous processing of several vertices and fragments. This means that several of the cast rays can be processed at the same time, which leads to increased performance.

2.4.1

The GPU today

As known, the GPU is a dedicated graphics rendering device that can be pro-grammed to perform a wide range of tasks and complex algorithms faster than a typical CPU. When first introduced the GPU only had a simple structure giving the programmer a few instructions to do simple graphical effects. The GPU:s of

today has evolved into a fully programmable environment, capable of doing much more than just graphical tasks. The GPU, working on both pixel and fragment level, have up till now mostly been used for texturing and lighting of surfaces, although often in very advanced applications. Recently many programmers are looking more into GPGPU, General-Purpose GPU computation, and the more advanced programming possibilities has opened up for a wide variety of tasks to be performed on the GPU. Particle systems, fluid dynamics, compression, raycast-ing, collision detection, cloth simulation, are just some examples of areas using the GPU.

As the capabilities of the GPU:s are advancing fast there is a need to categorize the features of the different processors. The earliest generation, the Pixel Shader 1.0, lack many of the features in the newer Pixel Shader 3.0 generation, such as conditionals, better floating point precision, higher number of allowed instructions and the possibility to bind a greater number of textures. This fact often leads to problems with backward compatibility and requires the developer not only to write one version of a shader but several to support older graphic cards.

2.4.2

The Raycasting Technique

As described earlier, the raycasting technique cast rays through the volume dataset, sampling values at equal distances and compositing them using equations derived from the volume rendering integral. To achieve this in a GPU-accelerated version it is needed to perform several rendering passes. Each pass outputs an intermedi-ate texture that is used in the forthcoming passes. Each pass use several advanced shader-programs, and the following sections explain each pass in depth.

Colorcube pass

A colorcube is a common object in computer graphics and consist of a cube with a length, width and height of 1 unit. The cube is colored by the position value at each point, resulting in a cube as in figure 2.13 with RGB-colors varying from (0, 0, 0) to (1, 1, 1). The cube is rendered at the position of where the volume is to be rendered, so it can be thought of as a bounding box for the raycasting. The result describes the points where the viewing rays enter the volume.

Direction pass

In this pass the Colorcube is rendered again, but now with only the back faces rendered and thus representing the exit points of the volume. By giving the result from the previous pass as a texture input to the shader program, it can easily do a texture lookup to fetch this value. By subtracting the color values in this pass with the result from the previous pass each pixel value generates a three dimensional vector. This vector, which is the viewing direction through the volume is normalized and stored in the RGB values in a texture. Besides the direction, the length through the cube is also stored in the alpha channel of the texture before normalization. This is a simple setup to derive the directions and lengths of the rays through the volume. The resulting texture is seen in figure 2.15.

2.4 GPU-Raycasting 15

Figure 2.13. Colorcube. Figure 2.14. Inverted Colorcube.

Figure 2.15. Direction.

Raycasting pass

Using the direction and length data derived in the previous pass, it is now possible to do the actual sampling of the volumetric data. The starting input to the shader is the same as in the first pass, i.e. the front faces of the bounding cube, since it serves as starting positions for the raycasting. An ordinary volume raycaster consists of a loop in the shader stepping along each ray, using the direction from the second pass, and sampling the volumetric data from a 3D texture and compositing the values according to the volume integral (equation 2.1). Using the built-in trilinear interpolation methods in the GPU, this sampling is easily implemented.

When implementing this on older shader models there are some problems with the raycasting pass. Since older shader models do not support more than 96 instructions per pass, it is impossible to step a ray through the volume. The shader programs are actually often quite short, but the compilation unfolds all

loops so they quickly grow in length. The solution to this is to perform several passes and traversing the rays a short length every time. This is the case in the article Acceleration Techniques for GPU-based Volume Rendering [6], we have however adopted the newer technique presented in the article Advanced GPU Raycasting [12] which just consist of a single raycasting pass. Newer shader models support up to 32768 instructions, but even this imposes problems if each ray were to consist of 1024 steps. After the compositing is done, the result is written to a resulting texture, ready to be put on screen.

Optimization Techniques

Several optimization techniques to improve the speed and quality of the above algorithm have been developed. Krüger and Westermann [6], and Scharsach [12] presents three possible methods in recent research. One method is to perform empty space skipping, which consists of refining the entry and exit points of the volume, clearly unnecessary to be described as a simple cube. The surface can be much more complex, not containing the parts of the volume data which does not contribute to the final rendering. Drawing a more advanced geometry instead of the simple Colorcube does not introduce much of a performance loss.

Another optimization is early ray termination, which occur when a ray trav-eling through the volume has reached an accumulated opacity value of one, the ray can be terminated since no more voxels will contribute to its final color. An advanced and efficient technique called pre-integration is presented by Engel [2], which consists of making pre-calculations of the volume integral (equation 2.3.1), between all possible values in the volume mapped to a specific transfer function. When these are made available to the raycaster as an input 2D texture, the al-gorithm can make much less steps along each ray since parts of it are already calculated.

We further discuss the use of these optimizations, and why they did not quite fit our needs, in the implementation chapter 3.4.1.

Rendering when inside the volume

The commonly used method of doing raycasting inside the bounding box of a volume object is described in the article Advanced GPU Raycasting [12]. The problem lies in the generation of the volume front faces, i.e. the Colorcube and therefore the raycasting entry points. When you get too close, these faces get clipped, or not rendered at all when completely inside the bounding box. The method described by Scharsach, if somewhat simplified, consist of rendering a camera-facing plane at the position of the near clipping-plane. This plane is color-coded by position and rendered before the Colorcube. Without writing to the depth-buffer the plane will replace the holes introduced in the front faces of the Colorcube when clipping.

2.5 Rendering of Emission and Reflection Nebulae 17

2.5

Rendering of Emission and Reflection

Nebu-lae

To do a realistic volumetric rendering of reflective nebulae the standard volume rendering technique has to be extended quite a lot. Researchers Magnor et. al. at the Max Planck Institute has made thorough research of how this can be done, presented in the paper Reflection Nebulae Visualization [8] and the thesis work Rendering and Reconstruction of Astronomical Objects [5]. The following sections summarize the most important aspects from these two papers, for our need to visualize nebulae. We start by considering what factors play a role in the light of nebulae. As mentioned earlier, the main reason why we see the nebulae are the very hot and bright stars in the middle lighting up the cloud of gas and dust. This light is both scattered and absorbed by the dust particles of which the cloud consist. This leads to the denser the cloud, the dimmer the escaping light becomes. The amount of extinction depends on both the size of the dust particles and the wavelengths of the light. To describe the scattering in interstellar dust, the main component of nebulae, only two parameters are needed besides the properties of the gas, and they turn out to be nearly constant throughout our galaxy. They are the albedo a, the percentage of incident light that scatters, and the scattering anisotropy g, a description of the angle in which light scatters.

2.5.1

Illumination

Each voxel in the dataset receives a specific amount of radiance from the stars. Using the equation 2.6 this illumination can be calculated. The radiance received at each voxel, Lill, depends on the amount of light the star emits, Φstar, the

distance r to the star and also the extinction or optical depth τopt of the dust,

accumulated between the voxel and the star blocking the light. Note that it is not depending on the viewing angle and is therefore suitable to be precomputed before the actual rendering.

Lill= Φstar 4πr2e −Rr 0 τopt(r 0_{)/r dr}0 (2.6) The optical depth, can further be derived from the ratio between the scattering depth τsct and the albedo, τopt = τsct/a. The scattering depth which describes

the amount of scattering events a photon makes along a path, can be expressed as the ratio between the scattering coefficient σsct and the size l of the voxels in

the dataset, τsct = σsct/l. Finally, the scattering coefficient σsct is in its turn

proportional to the density in the voxel.

2.5.2

Scattering

When this radiance Lill reaches a voxel it scatters anisotropically, and only a

small amount travels in the direction towards the observer. This fraction of light is denoted P (τsct, θ) and is a function of the scattering depth and the scattering

describes the probability p of all incident light with the anisotropy factor g, a constant, to scatter in a given scattering direction θ.

p (θ) = 1 − g

2

2 (1 + g2_{− 2g cos θ)}3/2 (2.7)

As the scattering probability P (τsct, θ) is the same for all volumes, and only

depends on a direction and the dust density, it is possible to precompute this distribution. This can be done by using a Monte Carlo simulation with the Henyey-Greenstein phase function, which is a quite complicated setup and thoroughly explained in the article by Magnor et al. [8].

The final amount of scattered light Lsctviewed by an observer is derived from

the equation 2.8. The ionization of the gas is an extra light contribution Lem

which is specified for each voxel in the dataset. Since it emits from each voxel it is a simple task o just include it in the scattered light.

Lsct= Lill· P (τsct, θ) + Lem (2.8)

An interesting fact is that the only way to really see the color of a nebulae change with viewing-angle is to do a computer simulation, as the distance to visible nebulae makes it impossible to view them from more than one direction.

2.5.3

Observed Illumination

After the scattering occurs, there is still some extinction to take into account, namely the blocking dust that lies in between the viewer and the voxel from where the light travels. To calculate the observed illumination at all points a similar extinction-integral as in equation 2.6 is needed for the scattered light:

Lobs= Lscte− Rl

0τopt(l 0_{)/l dl}0

(2.9) Each voxel in the volume can be assigned an amount of observed illumination, Lobs, and to produce a final image of the nebula each one of these needs to be

added up along the traced rays from the viewing location.

2.5.4

Wavelength dependencies

In the above calculations no consideration has been made regarding the wavelength dependency, other than what the ionization contribute with. This would result in a nebula with only gray hues, and as seen this is not quite the case. The truth is that the scattering coefficient σsct is wavelength dependant, i.e. it varies with

the color of the light. The parameter used in astrophysics to describe the relative extinction between the different bands of wavelengths is called RV, the ratio of

total to selective visual extinction. Using a standard value of RV = 3.1 the relative

factors are ABlue/AGreen= 1.324 and ARed/AGreen= 0.748 [8]. For denser clouds

of gas a higher value of RV is achieved, and thus less light is scattered. These

2.5 Rendering of Emission and Reflection Nebulae 19

this leads to scattered light consisting of blue hues, while the unscattered light is red.

From the earlier equations you can see that a varying σsctcreates

wavelength-dependencies on the local illumination 2.6, the scattering 2.8 and the extinction 2.9. During rendering this has to be taken into account to get a proper visual result.

2.5.5

Multiscattering

So far we have made the assumption that light only scatters at one point, but in reality it is constantly scattering throughout the gas. Creating a real-time applica-tion with proper multiscattering is impossible due to the fact that it would be far too demanding to calculate. But in the article Reflection Nebulae Visualization [8] a method for simulating some of the multiscattering behavior of light is presented. It is based on doing a series of renderings of the nebula in different resolutions and adding the differences between the smaller resolutions and the higher ones to the final image. The method is explained in more detail in the implementation chapter 3.4.2.

Chapter 3

Implementation

3.1

Environment

Uniview has been developed and extended with new features for some years now, and it has grown to be a quite big and intricate application. It is developed using object oriented C++, and utilizes OpenGL for rendering graphics. At the time of coding Uniview also was running with the help of SGI Performer, an API that aids the development of simulations and real-time applications. Fortunately Uniview has a rather nice modular design which made it easy for us to create our volume renderer as a separate module which could be inserted into the main rendering system without too much effort and learning how to use SGI Performer.

When we started to work on this thesis our main goal was to implement a vol-ume rendering solution capable of rendering arbitrary volvol-ume datasets. We made implementations of both a raycaster on the GPU and a proxy geometry solution, but soon realized that the raycaster was the only way to go. The intricate lighting conditions in nebulae require the ability to implement advanced blending equa-tions and adjusted optical models, which the GPU raycaster technique provides. In the following sections we have decided to exclude the implementation of the traditional volume renderers in favor for more extensive outlining of the nebulae renderings.

3.2

Design

The UML diagram seen in figure 3.1 shows the design of the implementation. Since both Sciss and we wanted a versatile solution which easily could be extended, we implemented the volume renderer using an object oriented approach, dividing up the code into different objects depending on function.

3.2.1

VolumeObject Class

The VolumeObject is an object class which was derived from Univiews

RenderOb-ject and makes it possible for Uniview to communicate with the obRenderOb-ject, trigger

Figure 3.1. Simple class diagram of implementation.

rendering and set parameters.

3.2.2

VolumeDataHandler and VolumeDataThread Classes

The class named VolumeDataHandler contains functions that manages the large amount of data required for volume visualization. During startup it reads low-resolution versions of all the volumes and a precalculated scattering table and stores them in the system memory for fast access.

When reading high-resolution data from disk the data handler utilizes the threaded VolumeDataThread object. Using a separate thread for disk access it was made possible to stream new volume data in small chunks from disk during runtime, without any frame rate issues.

3.2.3

RayCaster and FBOHandler Classes

The RayCaster and FBOHandler manages the GPU-raycasting which will be ex-plained in more detail in the next section. The raycasting code was implemented in Nvidia’s CG (C for Graphics) which at the time was the shader language most commonly used in Uniview. The fact that we also had some previous experience with CG made the decision to use CG easy.

To handle all the textures required by the raycasting algorithm to render onto, and read from, we created a FBOHandler class. When rendering to textures instead of the ordinary OpenGL renderbuffer which outputs graphics on the screen, an extension to the standard API is needed. To our knowledge there are two extensions to choose from, Pixel buffers [9] and Framebuffer Objects [10]. We made a decision early on to use the, at that moment, brand new Framebuffer Object extension. There are several advantages in using Framebuffer Objects over Pixel buffers, they are easier to set up, need less memory, and they are much faster to bind and switch between.

3.3 Data 23

3.3

Data

3.3.1

Pre-processing

The nebula data which we used in this thesis work was received from AMNH, see 2.2.2. The data contained several datasets, both over time and containing a large amount of properties, such as density, temperature, pressure and gas compo-sitions. For our purposes we only needed the gas optical density and the ionization levels for the available species, which were oxygen and hydrogen.

The first step was to make the data fit our objectives, and that included using a low-pass filter to get rid of an artifact from the simulation, down sampling the data, and also fading the edges of the data to get rid of the blocky appearance of the volume. The artifact introduced in the simulations consisted of a large gap appearing at a cross-section through the middle of the dataset. It was the result of poor boundary constraints, and the astrophysicists at AMNH could not easily improve the simulation. We therefore chose to apply gaussian low-pass filters to make them less apparent.

For our renderings we needed to calculate the radiance received at each voxel from the illuminating star. In order to calculate this radiance we implemented a pre-processing application that calculated the amount of light received from the star at each voxel by solving the equation 2.6, sampling at double frequency in order to fulfill the Nyquist theorem. Unfortunately we made a major mistake not taking into account the wavelength dependency of the local illumination, which we discuss further on. Instead of calculating the proper three channels, we only included one, being identical for the red, green and blue channel. Even so, the calculations are quite time-consuming. When down sampling the data, we only used a simple interpolation method taking into account six neighbours to a voxel. This was not sufficient and introduced some artifacts which made the volume looking like it was sampled too low. We removed these artifacts using a low-pass filter. The remaining data to be extracted from the simulation dataset was the ionization levels. To receive the proper values of these, we made use of the fact that the ionization is relative the square of the density [7].

The only step left was to down sample the data for use in the multiscattering algorithm 2.5.5 and in our level of detail system. Using our preprocessing applica-tion all data was saved in a binary raw format with increasing resoluapplica-tions of the data cubes from sizes of 163_{up to 256}3_{. The final data consisted of cubes}

contain-ing four channels; the optical density, the voxel illumination and ionization states for hydrogen and oxygen. A correctly implemented dataset would have included six channels, where the voxel illumination would have introduced two extra ones. Besides of the data preprocessing we also needed to calculate the scattering probability distribution. As explained in section 2.5 this distribution can be pre-generated. Using MATLAB from MathWorks, specialized on easy implementation of mathematical calculations, we created an application that generated this dis-tribution from the Monte Carlo calculations explained by Magnor et al. [8] and exported it as a raw data file. From this data we created a 2D lookup-texture over the scattering probability used in the raycasting pass. The coordinates in the

texture corresponds to the inputs to the Henyey-Greenstein phase function 2.7; the deflecting angle θ and the cosine of the incident scattering angle g.

3.3.2

Data management

Since the data format that we used was a collection of data-cubes of increasing resolution, the files grew large very fast depending on the resolution. A volume of size 1283 _{with four channels using 16-bit (2 byte) floats has for example a}

size of 16MB. At a doubled resolution, 2563_{, the volume data is instead 134MB}

large. Fortunately the quality was not really improved when choosing the higher resolution so we decided to limit ourselves to 1283 _{cubes, but still, this amount of}

data implied problems.

When this data had been read into the system memory by the threaded reader, a 3D texture had to be created and uploaded to the GPU memory. This step could not be done inside our data loader as its thread does not have access to the OpenGL context created on the main Uniview thread. Creating a large texture and trans-ferring it to the GPU memory on the main thread however caused the application to freeze for several seconds. We solved this problem by pre-allocating memory on the graphics memory for the 3D texture needed for the volume rendering, and then during runtime uploading smaller portions of the new texture one chunk at a time using the convenient OpenGL function glTexSubImage.

Using this technique the VolumeDataHandler object could load a high-resolution volume dataset without any noticeable performance drops when the viewer is close up to a nebula, and when positioned further away use lower resolution versions. This implementation made it possible to manage rendering of several volumes at the same time. Although, by using the pre-allocated GPU memory we were limited to only render one of the nebulae at a time at full resolution.

3.3.3

Nebula generator

We did a serious attempt to develop a volume nebula generator. Using a modeling interface similar to a 3D modeling software we defined areas where a nebula should evolve and where the stars inside should be situated. Creating a realistic gas cloud in a procedural way was however more difficult than we expected. Trying both Perlin noise [11] and wavelet noise [1] with poor result we gave up the idea of ourselves creating nebulae. It can however be done, which Magnor shows in his article Reflection Nebula Visualization [8], but developing such a generator was too big of a project.

3.4

Rendering

3.4.1

GPU-Raycasting

The GPU-raycasting technique which we have adopted for our main volume ren-dering is a modified version of the now commonly used method explained in sec-tion 2.4. Although, we experienced errors in precision of the vector describing the

3.4 Rendering 25

viewing-direction through the volume, when we were using the proposed method. These errors were most likely introduced in the texture lookup. So instead of us-ing the texture difference between passes one and two, we used the vector gous-ing from the position of the camera to the front faces, which is easy to calculate in the shader. For calculating the length of the ray through the volume, the texture difference was still needed.

A simple volume raycaster would as explained earlier, basically consist of a loop stepping along each ray and sampling the data from a 3D texture and compositing the values according to the raycasting equation 2.2. The emission and scattering raycaster however is a bit more complicated. At each step along the ray, the density along with the optical density, illumination and emission is sampled for each voxel from the 3D RGBA-texture created in the preprocessing step.

Scattering probabilities for each point along the ray is extracted by calculating the scattering angle between the current position and the star, i.e. in which angle the light has to scatter to travel towards the observer (see figure 3.2). To account for the scattering being different depending on wavelengths the optical density also has to be scaled properly using the extinction coefficients 2.5.4. These resulting scattering probabilities (one for each color channel), illumination color from the star, emission colors and the dusts optical density are then accumulated along the ray using the equation 2.9.

Figure 3.2. Light from illuminating star scattered towards observer.

In earlier versions of the GPU-raycaster we used optimization techniques as pre-integration, early ray-termination and empty space-skipping to optimize the performance. All three of these techniques unfortunately had to be removed from the version with scattering and emission. Empty space skipping was removed because there were no empty spaces in the dataset to skip since all density values have to be included. To use early ray-termination we would have to introduce an if-conditional and a break-statement inside of the shaders main raycasting loop. The overhead that this code introduced resulted in a higher performance loss than the actual gain so it proved to be useless. Even the quite effective pre-integration was removed because there is of course no simple transfer function to

use in the pre-integration, as the voxel color values depends on several factors, such as viewing direction and scattering probabilities, making the ordinary 2D pre-integration texture into a cumbersome high-dimensional texture.

3.4.2

Multiscattering

For multiscattering to occur, we extend the raycasting to implement the method referenced in chapter 2.5.5. The method renders the volume in consecutive passes using down sampled versions of the high resolution data (see figure 3.3). Images from one rendering pass are down sampled and subtracted from the images ren-dered in the next pass, with lower resolution. Finally the difference images are up sampled and added to the full resolution volume rendering to create the mul-tiscattered result. It is important to note that the density values in these lower resolution versions also have to be adjusted for their lower sampling distance as described by Magnor et al. [8]. In our case the multiscattering renderings consisted of four passes, using 1283_{, 64}3_{, 32}3 _{and 16}3 _{sized data, more passes proved to be}

unnecessary.

Figure 3.3. The mulitscattering procedure.

3.4.3

Post-processing

When the volume rendering is done, and the different layers of scattering are added together, the only remaining, but optional, step is to do post-processing of the colors. In order to see all of the hues and shades in the image we perform simple histogram normalization, i.e.

Kout = (Kvol− Kmin)/(Kmax− Kmin) (3.1)

The shading pipeline does not really support this kind of operation without performing several passes, where the volume is rendered in pass one and the mini-mum and maximini-mum values are read from the resulting image in order to use these

3.4 Rendering 27

in a second pass where the end result is achieved. On top of this, a read back from the renderbuffer is really expensive and can not be used in an interactive ap-plication. Our solution was to use a technique called asynchronous read back [9], where you can start several parallel read backs from the same buffer, and thus read the whole buffer faster. Even so, there is a big performance hit and since the color range for the rendering is pretty much consistent from all viewing angles, we decided to only read the buffer for one frame and use these values until the user demands an update.

After the normalization we perform both an opacity correction, based on a user selected value, and an exposure correction according to the HDR technique mentioned in Rendering and Reconstruction of Astronomical Objects [5]. The physics involved in the opacity correction are not really accurate, since it should really consist of tweaking the optical density of the gas and that would require a recalculation of the voxel illumination 2.8. We still decided to use this tool to be able to easily experiment with the visual impact of the nebulae. The exposure correction resembles how the colors of the astronomical pictures are based on what frequencies of light are mapped to which colors.

3.4.4

Animation of volumetric data

The data which is uploaded onto the graphics board for a high-resolution nebula rendering is as said a 16MB large 3D texture. Loading that amount of data to the system memory and thereafter uploading it to the GPU memory takes quite some time. Uploading a new dataset each frame when animating, was at the time of writing impossible, it simply cannot be done on the hardware of today.

Fortunately the data that we received was not suitable to animate in this manner anyway since the simulation only included 20 time steps with quite big time intervals. Our first attempts included replacing the current data on the graphics board by parts. Sending small cubes, 43, at each frame is not very demanding, but the animation did not look very pleasing. The blocky structures were clearly visible and it looked like what it was, replacing data, not animating data.

The approach that we decided to stick to was to simply load all, or whatever amount of data that fits, onto the graphics memory. Having all frames within reach of the GPU it is an easy task to select two frames of 3D-data and do a linear interpolation between these.

F ramei= (1 − αi) × F rame1+ αi× F rame2 (3.2)

It sounds a bit too easy, but the results are really smooth and highly interactive animations. The downsides are of course the initial loading and the data occupying almost all free GPU memory. Our method is easy, but will not be a sufficient solution when the visualizations in Uniview grows and no longer can afford to give up the large amount of GPU memory needed, but this is discussed further in the results section.

3.5

Rendering when inside the volume

When trying to render inside the volume, we used the method explained in chap-ter 2.4.2. Though, this technique did not adopt very well do our needs, as the precision of the resulting directions was far too low, resulting in a jerky volume when the camera was moving. So, instead of using the front faces drawn in the second pass as the origin of the direction rays, we used the position of the camera. The difference between the camera-position and the proposed near-clipping plane is practically none, but the camera position is less prone to precision problems. This gave us high precision ray directions even when moving around inside.

3.6

Impostors

The applications using volume rendering today are often specialized for visualizing objects up-close, for example in medical examinations. However, in Uniview, where you travel throughout the whole known universe, the case is often that the volumes are very far away, almost not visible. For this case, we could implement a powerful performance enhancement. When zooming away from a volume we switch to an impostor-mode, where the volume is rendered at a lower quality. Basically we use three simplifications of the quality

• As mentioned earlier, we use data cubes at a size of 323_{instead of 128}3_{, being}

able to cut down the raycasting loop to four time less steps than before. • When viewing the volume from afar it is usually really small, and thus we

can lower the resolution of the rendering since it is scaled down anyway. Instead of a typical rendering to a 256 by 256 square the impostor rendering is performed to a 64 by 64 square and that leads to four times fewer fragments to process.

• The volume is not rendered every frame. We start by recording the position of the camera and then render the volume, saving it to a texture. This texture is applied to a plane which is always facing the camera, thus creating an impostor-rendering of the true volume. When the camera has moved a certain number of degrees from where the texture was created, the volume is updated. Our tests have shown that a camera-movement of one degree is a good value to use.

Looking at the figure 3.4, the small differences between a full powered volume rendering and the simplified impostor rendering can be seen. When it is time for a switch between the impostor-mode and the true volume rendering we perform a smooth transition by simply rendering both the volume and the impostor and changing their opacity gradually. Since the impostor rendering can be performed at a very low cost, it is possible to render a very large amount of volumes at the same time, this is described further in the result section.

3.7 Graphical User Interface 29

Figure 3.4. Comparison between normal rendering and impostor rendering.

3.7

Graphical User Interface

In order to experiment and observe how different parameters impact the visualiza-tion we implemented a simple graphical user interface from which we can control some basic parameters, see appendix A.2. It is for example possible to control the dust albedo, the emissive colors, some resolution parameters, the exposures in three channels and a non-discrete control of data animation. To aid the user when setting the post-processing parameters, we also implemented an overlay rendering of a simple histogram that makes it easier to choose good minimum and maximum bounds for the gamut.

Chapter 4

Results

4.1

Visualization

Our volume rendering algorithm has successfully been implemented into the as-tronomical visualization software Uniview. It supports both the simulation of scattering and emissive properties of nebulae in real-time. The renderer can fur-thermore animate the data between several time steps using interpolation creating the effect of an evolving nebula. With the use of several optimizations mentioned below, a large number of nebulae can be rendered at the same time, and high-resolution datasets can be interacted with and even when navigated through. We have also implemented a graphical user interface as a part of the volume ren-dering engine, for easy user interaction of renren-dering properties such as physics, appearance, quality properties and animation during runtime.

4.2

Animation of data

At the end of the thesis work we were given access to the latest graphics card from NVIDIA, a Quadro FX5500, featuring a 1024-megabyte texture memory. This made it possible for us to animate up to 20 time steps at a resolution of 128x128x128 at interactive frame rates. Using 20 time steps a nice looking anima-tion ends up somewhere between 20 to 40 seconds, but since we can fade between time steps at any rate we want, this time is fully adjustable. It should be noted that the distance over which the nebula expands are several light years, so a true real-time animation should last several thousands of years.

4.3

Data format

To feed the implementation with data, we developed a new binary data file format. It holds a database containing all the data which our rendering algorithms require and its structure is optimized for fast reading performance. The file structure con-sists of nebula density, star illumination and emission information at each voxel,

stored with 16 bit per channel. For the multiscattering algorithm we require vol-ume data in several resolutions and therefore the file format also supports storing of a volume in more than one resolution. Using our own preprocessing application this database can be created from either raw volume data or data stored in the HDF5-format, widely used in the scientific visualization society.

4.4

Interaction with surroundings

The nebulae can be situated anywhere in space with surrounding stars making the visualization far more realistic than earlier applications. One big drawback with volume rendering is that it is quite performance demanding to render other objects inside of it. Therefore we do not yet have a good solution for rendering the central stars inside the nebulae. There also exist a special case in Uniview, when it comes to the rendering of the surrounding stars. An implementation feature makes them render behind the nebulae even if they should lie in front of them.

4.5

Performance

Using our implemented algorithms we are able to render 1283_{sized volumes with}

several time steps at interactive frame rates. It is possible to easily adjust the performance by setting several quality parameters, for example the number of rays being cast into the volume, or the number of steps they take through the volume. These settings are of course a trade off between quality and performance. To further optimize the performance of the rendering we have implemented several performance boosting techniques such as using an impostor rendering technique when volumes are positioned far away from the camera. Instead of always doing a full volume rendering of the volume, an impostor texture is created and applied to a screen-aligned plane at the volume position. This impostor texture is updated only when the camera angle to the nebula has changed more than one degree. To create a smooth transition between the impostor and normal rendering we fade between the two when closing or pulling away from the volume. Using this impostor technique we are able to render several volumes at the same time with almost no performance loss.

In tables 4.1, 4.2 and 4.3 some performance results are presented. All mea-surements have been recorded during normal usage in Uniview. This implies that the number of milliseconds in the tables does not represent the true performance, but rather the complete visualization pipeline of Uniview. The hardware used was a computer equipped with an AMD Athlon 64 X2 Dual 4200+ processor with 2 Gigabytes of RAM and an NVidia GeForce 7900 GT graphics card.

4.6

Limitations

There are a few limitations when rendering volumetric data using our algorithms. First, our calculations of the local illumination of each voxel in the volume dataset

4.6 Limitations 33

does not take into account the wavelength dependency of the scattering depth τsct. This fact most likely introduces some errors in the lighting calculations of

the nebula and this is really unfortunate, since it does not make our renderings as accurate as the ones we based our work upon [8]. We further discuss this in the next chapter.

Second, there is a limit to the amount of texture memory on the graphics card. Most consumer graphics cards on the market today support a total of 256 megabytes of textures. This restricts the resolution of the volume data and the number of time steps that can be stored on the graphics card at one time during animation. We tried to get around this problem using a streaming method. This solution was however abandoned because of the limited bandwidth and transfer methods between the system memory and texture memory resulting in a stuttering animation. We therefore choose the other alternative, uploading all volume data to the texture memory on the graphics card during initialization. The downside to this is the loading time before the actual animation can begin. On some computers this can take a long time, especially if it is equipped with slow hard drive disks.

Resolution 128Raysteps 64Raysteps remarks 512x512 85 ms 35 ms

484x484 65 ms 33 ms max fps reached 384x384 50 ms 33 ms max fps reached 256x256 33 ms 33 ms max fps reached

Table 4.1. Full volume rendering of 128x128x128 volume at different resolutions and

number of ray steps

Resolution 128Raysteps 64Raysteps remarks 512x512 70 ms 33 ms max fps reached 484x484 60 ms 33 ms max fps reached 384x384 35 ms 33 ms max fps reached 256x256 33 ms 33 ms max fps reached

Table 4.2. Full volume rendering of a animated 64x64x64 volume (23 timesteps) at

different resolutions and number of ray steps

N umberof impostors F rametime remarks 1 33 ms max fps reached 2 33 ms max fps reached 4 33 ms max fps reached 8 33 ms max fps reached