Raycasting SciVis 2008 - Raycasting Ronald Peikert Direct volume - - PDF document

Raycasting

Ronald Peikert SciVis 2008 - Raycasting 3-1

Direct volume rendering

Volume rendering (sometimes called direct volume rendering) stands for methods that generate images directly from 3D scalar data. "Directly" means: no intermediate geometry (such as an isosurface) is generated. Volume rendering techniques

• depend strongly on the grid type
• exist for structured and unstructured grids

Ronald Peikert SciVis 2008 - Raycasting 3-2

• are predominantly applied to uniform grids (3D images).

Direct volume rendering

2D or 3D image data are uniform grids with cell-centered data. Cell-centered data

• are attributed to cells (pixels, voxels) rather than nodes
• can also occur in (finite volume) CFD datasets
• are converted to node data

– by taking the dual grid (easy for uniform grids, n cells → n-1 cells!) – or by interpolating.

Ronald Peikert SciVis 2008 - Raycasting 3-3

Raycasting

Raycasting is historically the first volume rendering technique. It has common with raytracing:

• image-space method: main loop is over pixels of output image
• a view ray per pixel (or per subpixel) is traced backward
• samples are taken along the ray and composited to a single color

Differences are:

• no secondary (reflected, shadow) rays
• transmitted ray is not refracted
• more elaborate compositing functions
• samples are taken at intervals (not at object intersections)

Ronald Peikert SciVis 2008 - Raycasting 3-4

( j )

Raycasting

Sampling interval can be fixed or adjusted to voxels:

uniform sampling voxel-by-voxel traversal

Connectedness of "voxelized" rays:

uniform sampling voxel by voxel traversal (faster)

Ronald Peikert SciVis 2008 - Raycasting 3-5

6-connected (strongest) 18-connected 26-connected (weakest)

Ray templates

A ray template (Yagel 1991) is a voxelized ray which by translating generates all view rays. R t l t d th li b t b i l Ray templates speed up the sampling process, but are obviously restricted to orthographic views. Algorithm: Algorithm:

• Rename volume axes such that z is the one "most orthogonal" to

the image plane. g p

• Create ray template with 3D version of Bresenham algorithm,

giving 26-connected rays which are functional in z coordinate (have exactly one voxel per z-layer)

• Translate ray template in base plane, not in image plane

Ronald Peikert SciVis 2008 - Raycasting 3-6

Ray templates

Incorrect: translated in image plane Correct: translated in base plane base plane

Ronald Peikert SciVis 2008 - Raycasting 3-7

Compositing

Two simple compositing functions can be used for previewing:

• Maximum intensity projection (MIP):

– maximum of sampled vaules – result resembles X-ray image

• Local maximum intensity projection (LMIP):

– first local maximum which is above a prescribed threshold – approximates occlusion – faster & better(!)

LMIP MIP

intensity threshold

Ronald Peikert SciVis 2008 - Raycasting 3-8

camera

Compositing

Comparison of techniques (Y. Sato, dataset of a left kidney): Isosurface vs. raycasting with MIP, LMIP, α-compositing

fast (1 parameter) fast parameter free fast (1 parameter)

• (1 parameter)
• lighting

parameter free full data range noise insensitive

• (1 parameter)

full data range noise insensitive

• full data range

noise insensitive lighting

Ronald Peikert SciVis 2008 - Raycasting 3-9

• cclusion

(transparency)

• (occlusion)
• cclusion

transparency

α-compositing

Assume that each sample on a view ray has color and opacity:

( ) ( )

, , , ,

N N

C C α α

• [

] [ ]

3

0,1 , 0,1

i i

C α ∈ ∈

where the 0th sample is next to the camera and the Nth one is a (fully opaque) background sample:

background

( , , ) 1

N N

C r g b α = =

α-compositing can be defined recursively: Let denote the composite color of samples f,f+1,…,b

b f

C

Recursion formula for back-to-front compositing:

b b b b

C C α =

Ronald Peikert SciVis 2008 - Raycasting 3-10

( )

1

1

b b f f f f f

C C C α α

+

= + −

α-compositing

The first few generations, written with transparency

b

C C α = 1

i i

T α = −

1 1 1 1 b b b b b b b b b b b

C C C C C T C C C T C T T α α α α α α

− − − −

= = + + +

2 2 2 1 1 2 1 2 3 3 3 2 2 3 1 1 2 3 1 2 3 b b b b b b b b b b b b b b b b b b b b b b b b b b

C C C T C T T C C C T C T T C T T T α α α α α α α

− − − − − − − − − − − − − − − − − − − − −

= + + = + + +

reveal the closed formula for α-compositing:

1 i b 1 i b b f i i j i f j f

C C T α

− = =

= ∑

∏

Ronald Peikert SciVis 2008 - Raycasting 3-11

α-compositing

front-to-back compositing can be derived from the closed formula: Let denote the composite transparency of samples f,f+1,…,b

b f

T

b b f j j f

T T

=

= ∏

Then the simultaneous recursion for front-to-back compositing is:

f f f f

C C α =

1 1 1

1

f f f b b b f f b b f

T C C C T α α

+ + +

= − = +

Advantage of front-to-back compositing: early ray termination when

( )

1 1

1

b b f b f

T T α

+ +

= −

Ronald Peikert SciVis 2008 - Raycasting 3-12

g g y y composite transparency falls below a threshold.

The emission-absorption model

How realistic is α-compositing? The emission-absorption model (Sabella 1988) yields a basic l d i ti volume rendering equation

( ) ( )

( )

x b

x x dx

L d

τ

′

′′ ′′ −∫

′ ′

∫

The equation describes the radiance (power per unit area per solid

( ) ( )

x

x

L x x e dx ε ′ ′ = ∫

The equation describes the radiance (power per unit area per solid angle [W/m2/sr]) arriving along a ray at the position x on this ray. The emission function ε(x) describes the photons "emitted" by the ( ) p y volume along the ray. The absorption function τ(x) is the probability that a photon traveling

Ronald Peikert SciVis 2008 - Raycasting 3-13

• ver a unit distance is lost by absorption.

The emission-absorption model

The emission-absorption model is based on Boltzmann's transport equation in statistical physics, but completely ignores scattering. I l d l ( ) i ti ti f ti h i b th In more general models τ(x) is an extinction function having both an absorption term and a scattering term. Instead, in the emission-absorption model:

• incident scattering is modeled by the emission function

incident scattering is modeled by the emission function

• loss by scattering can be thought to be part of the absorption.

Ronald Peikert SciVis 2008 - Raycasting 3-14

The emission-absorption model

Discrete version of emission-absorption model

1

1

i j x

i n n τ

−

− Δ −

∑

f

( )

1

j j j

x i n n x i i i i j

L x xe x e

τ τ

ε ε

=

Δ − Δ = = =

∑ = Δ = Δ

∑ ∑ ∏

matches the α-compositing formula

( )

1

1

i b b f i i j

C C α α

−

= −

∑ ∏

and gives interpretations of "opacity" and "color":

( )

f i i j i f j f = =

∑ ∏

1

i x

e τ

− Δ

The product is called a premultiplied or associated color

1

i x

i

e τ α

Δ

= −

i i i

C x α ε = Δ C C α =

• Ronald Peikert

SciVis 2008 - Raycasting 3-15

The product is called a premultiplied or associated color.

i i i

C C α =

Transfer functions

Transfer functions map raw voxel data to opacities and colors as needed for the α-compositing. Inputs of TF (one or more):

• voxel value

di t it d

( )

∇

( )

s x

• gradient magnitude
• higher derivatives of

Opacity transfer function

( )

s ∇ x

( )

s x

( ) ( )

( ) s s α ∇ x x

Opacity transfer function Color transfer function

• r premultiplied:

( ) ( )

( , , ) s s α ∇ x x

• ( )

( )

( , , ) C s s ∇ x x

• ( )

( )

( ) C s s ∇ x x

• – or premultiplied:

In general TF don't depend on spatial location, exception: for focus+context techniques

( ) ( )

( , , ) C s s ∇ x x

• Ronald Peikert

SciVis 2008 - Raycasting 3-16

exception: for focus+context techniques

Transfer functions

By choosing different opacity transfer functions different types of applications can be achieved. Examples: Examples:

• pacity
• pacity
• pacity

voxel value voxel value gradient magnitude g g

standard application isosurface 3D edge detector

Ronald Peikert SciVis 2008 - Raycasting 3-17

Transfer functions

Example of a bivariate (=2D) transfer function:

• pacity

gradient magnitude

( )

s ∇ x

3

α

1

α

2

α

3

f

2

f

1

f

Ronald Peikert SciVis 2008 - Raycasting 3-18

voxel value ( )

s x

Transfer functions

Example: bivariate transfer function for isosurface of constant "thickness".

• pacity

gradient magnitude

( )

s ∇ x f

Ronald Peikert SciVis 2008 - Raycasting 3-19

voxel value

( )

s x

Transfer functions

The color transfer function allows to make a simple classification. Example:

color (RGB) fat air bone tissue voxel value ( )

s x

Ronald Peikert SciVis 2008 - Raycasting 3-20

Transfer functions

Better (but more expensive) classification than with a transfer function is obtained by segmentation (typically slice by slice, semi-automatic). (typically slice by slice, semi automatic). Pre-classified volume data. Example: "virtual frog" dataset (Lawrence Berkeley Labs).

Ronald Peikert SciVis 2008 - Raycasting 3-21

Transfer functions

Volume rendering of segmented volume data (from VTK book).

Ronald Peikert SciVis 2008 - Raycasting 3-22

Transfer functions

Example: classifications with different transfer functions (McGill Univ.) Example: temperature data for it d l TF

Ronald Peikert SciVis 2008 - Raycasting 3-23

• pacity and color TF (Nissan Res.Ctr.)

Transfer functions

In pre-classification, the voxels can also be lit:

• The gradient is perpendicular to the local isosurface. It can be

used as a normal vector for a Phong lighting (without rendering the isosurface itself).

• Reflection coefficients can be assigned by a separate transfer

f ti (" t i l " i t d f l l ) function ("materials" instead of colors only).

• The diffuse lighting can be applied to the entire volume dataset

as a pre processing since it is independent of the viewing as a pre-processing since it is independent of the viewing direction.

Ronald Peikert SciVis 2008 - Raycasting 3-24

Pre- vs. post-classification

For quality reasons, current volume rendering implementations

• ften use post-classification.

Pre-Classification: 1. Transfer functions are applied to voxels 2. Results are interpolated to sample locations. Post-Classification: Post-Classification: 1. Raw data are interpolated to sample locations. 2. Transfer functions are applied to sampled data.

Ronald Peikert SciVis 2008 - Raycasting 3-25

Pre- vs. post-classification

Pre-classification:

• can be done as

pre processing

l s(x) x α(s(x))

pre-processing, e.g. segmentation, diffuse lighting

voxel sample

Post-classification:

• Interpolation is in

sample inter- polated

• Interpolation is in

the correct space

• Additional samples can

be interpolated on the fly to improve output quality

view ray voxel value transfer function

Ronald Peikert SciVis 2008 - Raycasting 3-26

Pre- vs. post-classification

128 slices pre- classification 128 slices post- classification 128 slices preintegrated 284 slices post- p g p classification

Ronald Peikert SciVis 2008 - Raycasting 3-27

Image credit: K. Engel, U. Stuttgart

Preintegration

Idea (Engel 2001):

• simulate infinitely many interpolated samples between two

successive samples s = s(x ) and s = s(x ) successive samples si = s(xi ) and si+1 = s(xi+1 )

• assuming:

– field s(x) varies linearly between samples ( ) y p – transfer functions don't depend on derivatives

Ronald Peikert SciVis 2008 - Raycasting 3-28

Preintegration

The discrete formula for opacity at a sample was

1

i

x i

e

τ

α

−

= −

• The continuous version, for a sample interval [xi, xi+1], is

i

( )

( )

1 xi +

∫

Assuming now s(x) to be linear between samples we get

( )

( )

1

xi

s x dx i

e

τ

α

− ∫

= −

Assuming now s(x) to be linear between samples, we get with

( )

1 1

1

si i i s

d s ds s s τ

+ +

− −

∫

1 i i

d

+

= − x x

which is called a preintegrated opacity transfer function.

1

1

i i si

i

e α

+

= −

1 i i +

Ronald Peikert SciVis 2008 - Raycasting 3-29

Preintegration

The integral

1 1

( ) ( ) ( )

i i i

s s s s

s ds s ds s ds τ τ τ

+ +

= −

∫ ∫ ∫

can be evaluated by two lookups in a precomputed table of

i

s s

( )

s

s ds τ ′ ′

∫

Alternatively, the preintegrated opacity TF could be tabulated for all possible combinations of (si, si+1, d), especially if the p ( i,

i+1,

), p y sampling distance d is chosen constant.

Ronald Peikert SciVis 2008 - Raycasting 3-30

Preintegration

The composite color of the same interval

( )

( )

( )

1

' '

( )

x i xi

x s x dx

C s x e dx

τ

ε

+

−∫

∫

simplifies for linear s(x) to

( )

( )

i i

i x

C s x e dx ε = ∫

( )

( )

1 1

' ' 1

s i i i si i

d s s ds s s i i i s

d C s e ds s s

τ

ε

+ +

− − +

∫ = −

∫

which is a preintegrated color transfer function. Again, it can be tabulated for all combinations of (si, si+1, d),

• r it can be approximately calculated from tabluated values of

1

i

i i s +

• r it can be approximately calculated from tabluated values of

and (moving the exponential term out

( )

s

s ds τ ′ ′

∫

( )

s

s ds ε ′ ′

∫

Ronald Peikert SciVis 2008 - Raycasting 3-31

• f the integral).

∫ ∫

Raycasting hardware

Today's GPUs are well suited for volume rendering. Object-space methods are more straight-forward, but also some raycasting implementations on GPUs have been done. raycasting implementations on GPUs have been done. In contrast, two examples of special raycasting hardware:

• VIZARD (1997, U. Tübingen):

– FPGA-based system – long pre-processing – 10 Hz for 2563 volume – perspective view raycasting – perspective view raycasting

Ronald Peikert SciVis 2008 - Raycasting 3-32

Raycasting hardware

volumePRO (1999, MERL): – PCI board – no pre-processing no pre processing – 30 Hz for 2563 volume – orthographic raycasting ray templates – ray templates – shear-warp method – z-supersampling

Ronald Peikert SciVis 2008 - Raycasting 3-33

Example: with / without supersampling

Raycasting hardware

ASIC chip of VolumePRO vg500 board has 4 parallel pipelines with stages: g

• interpolation
• gradient estimation

l ifi ti / h di

• classification/shading
• compositing

Image is rendered in a single pass through the volume. Off-chip memory for slices (shown in blue)

Ronald Peikert SciVis 2008 - Raycasting 3-34

Raycasting hardware

Some special features:

• cropping
• clipping planes

some possible combinations of cropping and clipping planes

Ronald Peikert SciVis 2008 - Raycasting 3-35