A Decompression Pipeline for Accelerating
Out-of-Core Volume Rendering ofTime-Varying Data
Daisuke NAGAYASU a,b Fumihiko INO a,∗
Kenichi HAGIHARA aa Graduate School of Information Science and Technology, Osaka University,

1-3 Machikaneyama, Toyonaka, Osaka 560-8531, Japanb Oracle Corporation Japan, 4-1 Kioi, Chiyoda, Tokyo 102-0094, Japan

Abstract
This paper presents a decompression pipeline capable of accelerating out-of-corevolume rendering of time-varying scalar data. Our pipeline is based on a twostage compression method that cooperatively uses the CPU and GPU (graphics
processing unit) to transfer compressed data entirely from the storage device tothe video memory. This method combines two diﬀerent compression algorithms,
namely packed volume texture compression (PVTC) and Lempel-Ziv-Oberhumer(LZO) compression, allowing us to exploit both temporal and spatial coherence intime-varying data. Furthermore, it achieves fast decompression by taking architectural advantages of each processing unit: a hardware component on the GPU and alarge cache on the CPU, each suited to decompress PVTC and LZO encoded data,
respectively. We also integrate the method with a thread-based pipeline mechanism
to increase the data throughput by overlapping data loading, data decompression,and rendering stages. Our pipelined renderer runs on a quad-core PC and achieves
a video rate of 41 frames per second (fps) in average for 258 × 258 × 208 voxel datawith 150 time steps. It also demonstrates an almost interactive rate of 8 fps for
512 × 512 × 295 voxel data with 411 time steps.
Key words: Volume rendering, time-varying data, pipelined rendering, datacompression, GPU∗ Corresponding author. Tel.: +81 6 6850 6597; fax: +81 6 6850 6599.
Email address: ino@ist.osaka-u.ac.jp (Fumihiko INO).

Preprint submitted to Computers & Graphics

30 March 2008

1

Introduction

Volume rendering [1] of time-varying data is capable of producing animation
sequences that show how the three-dimensional (3-D) structure evolves overtime. This visualization technique plays an increasingly important role for theintuitive understanding of complex time-varying phenomena [2]. For example,
it is useful to interpret simulation results and scanned X-ray images in physicaland life sciences. Thus, real-time rendering of time-varying volume data isrequired to assist scientists eﬀectively in time-series analysis.
One challenging issue in time-varying volume visualization is to develop an
eﬃcient mechanism for handling 4-D data. For example, a time-varying volumeof 512×512×512 voxels with a hundred time steps requires 12.5 GB of memory
if each voxel has 1-byte data. Moreover, recent advances in integrated circuits
allow us to produce large-scale datasets, which could not be entirely stored
in the main memory as well as in the video memory. Due to this increasingdata size, we need fast out-of-core rendering systems, which relax memory
requirements by on-demand loading of data.
Data compression techniques also provide an eﬀective solution to the issue
mentioned above. Actually, such techniques are integrated into many rendering systems [3–10] to reduce data size. Most of them use the graphics processing unit (GPU) [11] for data decompression. These GPU-based methods [3–7]are more eﬃcient compared with CPU-based methods [8,9], which might suﬀer
from low performance due to the raw data transferred from the main memory
to the video memory. In addition, GPU-based methods will be further accelerated by future GPUs, because they are rapidly increasing the performancein terms of ﬂoating point operations and memory bandwidth [12].
Although GPUs are becoming ﬂexible with enhancing programmability, their
special-purpose architecture is still a barrier to achieve both fast decompression and high compression ratio. For example, random accesses are not fastcompared with sequential accesses [13,14], slowing down the processing speed
against complex, irregular computation including branching instructions. Furthermore, GPUs currently do not have large video memory of more than 1.5
GB. Therefore, we need a hybrid method [10] that performs data decompression both on the CPU and GPU in a cooperative manner. This hybrid method,which is our main idea, can oﬄoad GPU cycles to the CPU without any performance loss if the storage device or GPU is kept as a performance bottleneck
after oﬄoading. That is, although this oﬄoading increases the communication amount between the main memory and the video memory, it will notdegrade the rendering performance if the increased amount does not limit thedata throughput. In addition, the oﬄoading should be applied to a sequenceof operations that are suited to the CPU architecture rather than the GPU2architecture.
The objective of our work is to develop an eﬃcient compression mechanism
capable of accelerating out-of-core rendering of time-varying scalar data. To
achieve this, we develop a thread-based decompression pipeline consisting of
two stages, each processed on the CPU and GPU, respectively. We expand
on our preliminary work [10] in this paper, describing how our pipeline mechanism accelerates the hybrid method for higher throughput. The diﬀerence
to the earlier work [10] is the pipeline mechanism that hides decompression
and rendering costs using a multi-core system. Our decompression pipeline
requires a sequence of double-compressed data as the input, each consisting of
three volumes at successive time steps. The double-compressed data is thendecompressed ﬁrst by the CPU, and then by the GPU in a pipeline fashion.
Experimental results are also presented to demonstrate the eﬀectiveness of thepipeline by using a desktop PC equipped with a multi-core CPU.
The rest of the paper is organized as follows. Section 2 introduces related work.
Section 3 gives a brief summary on volume texture compression (VTC) [15]and LZO compression [16], which are the basis of our decompression pipeline.
Section 4 describes the details of the pipeline with the hybrid method. Section
5 shows experimental results evaluating the pipeline in terms of performance,
scalability, and image quality. Section 6 concludes the paper with future work.

2

Related Work

There are many pipelined rendering systems [17–20] that overlap data loading and rendering stages to increase the rendering rate. Since these systems
deal with raw data, data compression techniques are not integrated into them.
On the other hand, many compression-based systems [3–9] are also proposed
to increase the rendering performance. In these systems, only Bernardon etal. [9] uses a pipelined system that overlaps decompression with rendering. To
the best of our knowledge, however, there is no work reporting on a pipelined
rendering system that exploits both of the CPU and GPU for decompression. Thus, it is still not clear how eﬀectively data compression techniques
and pipeline mechanisms contribute to higher performance. Table 1 shows aclassiﬁcation of existing compression-based systems. As we mentioned earlier,they can be classiﬁed into three groups where data is decompressed: a hybrid
method [10]; GPU-based methods [3–7]; and CPU-based methods [8, 9].
Most of prior systems use only the GPU for data decompression. Akiba et
al. [3] propose a combination of two lossless compression techniques that
exploit spatial coherence. Their system ﬁrst applies the wavelet transformation [21] to each volume in order to convert them into a hierarchical represen3tation. This wavelet-based compression allows users to select an appropriate
resolution level according to interactivity and image quality. After this, thesystem applies a packing procedure to a selected subset of the low-resolution
data. This procedure partitions the volume into 3-D blocks, which are then reduced by replacing the same blocks with pointers to the representative block.
On the other hand, data decompression (unpacking) is done at run-time by afragment program running on the GPU. However, this unpacking procedure
could cause random accesses in the worst case, because it exploits the coherence across 3-D blocks. We think that such coherence should be exploited by
the CPU, which has a larger cache with a higher tolerance to lower degrees
of locality. Our system realizes this by performing LZO decompression on theCPU. Furthermore, the system exploits both temporal and spatial coherence,
although it is a lossy compression.
Binotto et al. [4] also realize a lossless compression method. Their method
partitions the volume into 3-D blocks, which are then reduced by replacing
the same blocks with the representative block. In addition to this spatial coherence, temporal coherence is also exploited by replacing blocks between
diﬀerent time steps. Since this replacement requires an exact match between
diﬀerent blocks, their method is eﬀective for highly sparse and temporally
coherent data.
Fout et al. [5] realize lossy compression based on vector quantization. Theirsystem has two modes depending on coherence to be exploited: spatial andtemporal modes. In the spatial mode, their system packs 2 × 2 × 2 neighboringvoxels into a vector, and then applies quantization to the vector in order toapproximate it with a representative vector. Schneider and Westermann [6]also use vector quantization but they exploit both spatial and temporal coherence by sharing representatives between successive time steps. Similar to
Akiba’s packing procedure, vector quantization produces representative vectors and pointers to them. Due to the nature of pointers, this sequence can result in irregular data, which generally decreases compression ratio. Therefore,
when we combine multiple compression methods in a pipeline, such irregular
data should not be generated at the ﬁrst compression stage in order to avoid
disturbing the second compression method. This also implies that complex
decompression algorithms should be executed on the CPU rather than on theGPU.
Lum et al. [7] use the discrete cosine transform (DCT) to achieve lossy compression. Their hardware-assisted solution exploits temporal coherence bytransforming data into a set of coeﬃcients. These coeﬃcients are then quantized to create a more compact representation, allowing us to discard coeﬃcients with higher energy, which are not important in terms of image quality.
Since the discard level can be selected by users, their system is capable ofchoosing the appropriate image quality.
4

A CPU-based method is proposed by Strengert et al. [8]. Their system is basedon the wavelet-based compression [22] that exploits spatial coherence at eachtime step. It employs the wavelet transformation to generate a hierarchicaloctree, which is then compressed by an entropy encoding. During rendering,it decompresses data on the CPU and then transfers the raw data to the video
memory through the graphics bus. Therefore, the performance of this method
could be limited by the graphics bus especially if we increase the compression
ratio (for higher performance).
Bernardon et al. [9] also use a CPU-based method but they overlap the decompression time with the rendering time. Their system uses vector quantization
to exploit temporal coherence in unstructured grids. Although this systemdemonstrates successful overlaps, the decompression stage is processed only
on the CPU. Accordingly, the raw data is transferred through the graphics
bus. To avoid this, we think that the decompression stage itself should be
structured into a pipeline.
Finally, a raw renderer is proposed by Furumura et al. [23] to deal with the
problem of increasing data size. They develop a concurrent technique capableof processing simulation and visualization at the same time. This concurrent
technique directly accesses the simulation results on the main memory, eliminating the need to store huge amounts of data. Therefore, data compression
techniques are not needed for visualization. However, the simulation performance can limit the frame rate, because the simulation must be executed atevery frame. In addition, the interactivity in the visualization phase might
be reduced on batch-style high-performance computing systems, which areneeded for the simulation.

3

Preliminaries

3.1 Volume Texture Compression (VTC)

VTC [15] is a lossy compression method that approximates voxel values bylinear interpolation of representative values. This method is supported by
nVIDIA GPUs [11] and is standardized as a part of OpenGL extensions [24].
One remarkable advantage of this method is fast data decompression accelerated by a hardware component on the GPU. Although VTC has four modes interms of the format of input data, we only explain the COMPRESSED RGB
S3TC DXT1 EXT mode in this paper, which compresses 24-bit RGB data.
We use this mode as an underlying method for our pipeline, because it has
the highest compression ratio of 6:1.
5

Fig. 1 shows an overview of VTC. Given a volume in the RGB format, VTC
partitions the volume into 4 × 4 × 1 voxel blocks B1 , B2 , . . . , BN , which arethen applied to an approximation algorithm to generate compressed data.
This algorithm replaces each Bk (1 ≤ k ≤ N) of the blocks, namely 16 voxel
values v1,1 , v1,2 , · · · , v4,4 , with a compressed block Bk containing (1) two representative values Rk = c1 , c2  in 4 bytes and (2) an interpolation table
Tk = {mi,j | 1 ≤ i ≤ 4, 1 ≤ j ≤ 4} in 4 bytes, where mi,j in 2 bits represents the approximation mode for voxel vi,j . Since mi,j has 2-bits information,
four approximation modes are predeﬁned [15, 25]: c1 , c2 , (2c1 + c2 )/3, and(c1 + 2c2 )/3. For example, if the algorithm determines that the value for voxel
vi,j should be approximated by (2c1 + c2 )/3, then it stores mi,j = 2 in table Tk ,according to the predeﬁned scheme. This linear interpolation is independently
applied to every block, so that each block has a diﬀerent table Tk . Therefore,
VTC exploits spatial coherence in a small 4 × 4 × 1 block.
In contrast to the encoding scheme mentioned above, the decoding scheme
returns the approximated value by linear interpolation of two representativevalues c1 and c2 . For example, it returns (2c1 + c2 )/3 as the value of voxel
vi,j if mi,j = 2. Note that this interpolation is quickly done by the on-the-ﬂy
decoder implemented as a hardware component. In addition, the on-the-ﬂy
decompression does not require additional memory to store the decompressed
data, because the entire volume is not decompressed at a time. The GPU also
exploits two parallelisms to further accelerate this decompression: (1) single
instruction multiple data (SIMD) instructions [26] and (2) vector instructions,allowing us to simultaneously process multiple voxels in a volume and multiple RGB channels in a single voxel, respectively. As compared with SIMD
instructions, vector instructions are useful to increase the data granularity,
which contributes to the full utilization of memory bandwidth.

3.2 Lempel-Ziv-Oberhumer (LZO) Compression

LZO [16] is a real-time data decompression/compression library based on alossless compression algorithm. This library uses a dictionary coder [27], whichscans the input data to replace repeating data with references back to theoriginal data. During this scan, a sliding window keeps recent scanned datato perform a replacement operation when a longest match is found in thewindow. Therefore, LZO obtains high compression ratio if the same data sequence frequently appears in the window. Similar to VTC, LZO has multiple
compression modes. In this paper, LZO denotes LZO1X-1, which is regarded
as the fastest mode in most cases.
Note that the LZO library provides fast decompression, because it does not
require complex operations to decode data. It simply replicates the original
6

data to places where the references exist. In addition, this replication can bequickly done if the sliding window is stored entirely in a CPU cache. Actually,
the LZO library has a 64-KB window, which can be stored in an L2 cache.
In addition to this cache optimization, the LZO library has a hash table forsearch acceleration.

4

Decompression Pipeline

We now describe our decompression pipeline designed on the basis of compression methods presented in Section 3. We ﬁrst show the design aspectsof our pipeline, and then describe the hybrid method [10] and the pipeline
mechanism.
Fig. 2 shows an overview of the hybrid method. As we mentioned in Section 1,
the key idea of this method is to perform data decompression both on the CPU
and the GPU. To realize this, we combine two diﬀerent compression methods,
Mc and Mg , each running on the CPU and GPU, respectively. Therefore,
the raw data must be converted into double-compressed data in advance ofrendering. To do this, the hybrid method compresses the raw data ﬁrst by
Mg , and then by Mc . During the visualization phase, this double-compresseddata is decompressed by applying the decompression methods in the opposite
order: ﬁrst by the CPU, and then by the GPU.
We also develop a pipeline mechanism in order to overlap data loading, data
decompression, and rendering stages. This mechanism increases rendering
throughput, because it allows us to process successive data simultaneously at
diﬀerent pipeline stages. In addition, the pipeline mechanism contributes tohide the decompression overhead incurred on the CPU. In our hybrid method,
this overhead is due to the cycles oﬄoaded from the GPU. Therefore, we thinkthat the pipeline mechanism is essential to maximize the eﬀectiveness of thehybrid method.

4.1 Design Aspects
Compression methods Mc and Mg must satisfy the following three requirements to achieve fast volume rendering of time-varying data.
R1. Both Mc and Mg achieve a compression ratio suﬃciently high to pay
decompression cost T2 , where T2 represents the time for data decompression.
R2. Mc and Mg exploit diﬀerent coherences.
7

R3. Mg outputs a sequence of data such that it keeps the coherence that will
be exploited later by Mc .
The ﬁrst requirement R1 is a general constraint to compression methods. Thisrequirement implies that there is a tradeoﬀ relation between compression ratio
and decompression time. From this viewpoint, although a high-compression
method signiﬁcantly reduces the loading time T1 , it is not satisfactory if it
fails to reduce the total time T1 + T2 .
On the other hand, the remaining requirements R2 and R3 are inherent in hybrid (two-stage) compression methods. Requirement R2 indicates that it is not
eﬀective if similar compression methods are repeatedly applied to the volumedata. For example, a combination of a temporal encoder and a spatial encoder
satisﬁes this requirement, because they focus on independent coherences. Requirement R3 comes from the nature of sequential execution of methods Mc
and Mg . For example, Mc and Mg are not a good pair if Mg destroys thecoherence which Mc exploits later in the second compression stage.

4.2 Hybrid Compression

To satisfy requirements R1, R2, and R3, our hybrid method uses packed volume texture compression (PVTC) [10] for Mg and LZO [16] for Mc . PVTC
is a lossy compression method for time-varying volume data containing scalarvalues. It extends VTC [15] to exploit temporal coherence in time-varying
data. Using PVTC and LZO, the hybrid method iteratively compresses threesuccessive volumes from the beginning of volume sequence. In the following
discussion, let Vt denote the volume at time step t, where 0 ≤ t < tmax . Let
Vt (x, y, z) denote the voxel value of point (x, y, z) in volume Vt . Then, thecompression procedure can be written as follows.
(1) t ← 0.
(2) Data packing. Fig. 3 illustrates this step. Given three successive raw
volumes Vt , Vt+1 , and Vt+2 , the method packs them into RGB data toobtain a packed volume Ct/3 . To do this, three scalar values Vt (x, y, z),
Vt+1 (x, y, z), and Vt+2 (x, y, z) in the same location are copied to R, G, andB channels of Ct/3 (x, y, z), respectively. Here, Ct/3 (x, y, z) represents
the voxel value of point (x, y, z) in Ct/3 . Empty values are padded if
tmax mod 3 = 0.
(3) Data compression by VTC. The packed volume Ct/3 is given to a VTC
encoder to obtain a compressed texture Xt/3 .
(4) Data compression by LZO. The compressed texture Xt/3 is given to a
LZO encoder to obtain a double-compressed texture Zt/3 .
(5) t ← t + 3. Goto step (2) if t < tmax .
8The combination of PVTC and LZO satisﬁes requirement R2, because PVTC
and LZO focus on diﬀerent coherences. That is, PVTC exploits temporal andspatial coherence in a small 4 × 4 × 1 block. On the other hand, LZO exploits
spatial coherence across small blocks, because it receives PVTC compressed
textures as inputs. Here, the sliding window in the LZO library covers 8192
VTC compressed blocks at a time, which is equivalent to 512 × 256 × 1 RGB
voxels. Therefore, LZO compression can be regarded as a slice-based compression if a volume consists of 512 × 512 voxel slices. In this sense, LZO andPVTC do not exploit exactly the same coherence, though they exploit spatial
coherence.
With respect to requirement R3, this requirement can be satisﬁed as follows.
• PVTC can keep the coherence between diﬀerent blocks in the input volume,
because it will encode similar blocks into compressed blocks with similar
representative values and interpolation tables. Thus, the coherence between
blocks in the raw volume is also kept in compressed data for later LZO
compression.
• The coherence mentioned above is expressed in the format of compressed
data, because PVTC produces ﬁxed-size data for each block in ﬁrst-in, ﬁrstout (FIFO) order.
From the viewpoint of architectural design, our combination aims at achieving
fast decompression by taking advantage of GPU and CPU architectures. For
example, as we mentioned in Section 3.1, the GPU has hardware acceleration
capabilities which the CPU does not have. On the other hand, the CPU has
a larger cache as compared with the texture cache on the GPU. This larger
cache is suited to duplicating operations in LZO decompression, because these
operations refer to a large area in the sliding window. Actually, the cacheworking set in an nVIDIA GeForce 8800 GTX card is 8 KB per multiprocessor[28], and thus the entire window cannot be stored in the texture cache. Insummary, the CPU cooperates with the GPU by providing a larger cache toaccelerate decompression for large area.

4.3 Pipelined Data Decompression

Once the raw data is converted into double-compressed textures, our pipelined
renderer repeatedly decompresses them during visualization. Note that PVTC
allows us to load data at every three time steps, because each of PVTC compressed textures contains time-series data. The decompression scheme can be
written as follows (see also Fig. 4).
S0. t ← 0.
9

S1. Data loading from the storage device. A double-compressed texture Zt/3
is transferred from the hard disk to the main memory.
S2. Data decompression by LZO. The CPU decompresses Zt/3 to obtain acompressed texture Xt/3 .
S3. Texture-based volume rendering [29, 30] on the GPU.
(a) Texture transfer. Xt/3 is transferred from the main memory to thevideo memory.
(b) Data decompression by PVTC. For all voxels (x, y, z), the approximated value Ṽt (x, y, z) is automatically obtained by simply accessing
the appropriate channel of packed voxel Ct/3 (x, y, z). That is, theVTC decoder returns R, G, or B channel data if t mod 3 = 0, 1,
or 2, respectively. Note that this branching instruction can be eliminated from fragment programs. Instead of fragment processors, theCPU performs this by binding the appropriate fragment programthat refers one of R, G, or B channel data.
(c) Data classiﬁcation by a transfer function. The ﬁnal color and opacityof Ṽt (x, y, z) are determined by a transfer function.
(d) t ← t + 1. Goto step (b) if t mod 3 = 0.
S4. t ← 0 if t ≥ tmax . Goto stage S1.
The above stages S1, S2, and S3 are overlapped by our pipeline mechanism,as shown in Fig. 4. To realize this mechanism by software, we use the POSIX
thread library [31]. Fig. 5 shows a pseudocode of the mechanism. We create
three threads, each responsible for processing one of the three stages. Thismechanism is thread-safe if it satisﬁes the following conditions C1 and C2
during execution.
C1. For all time steps t, where 0 ≤ t < tmax , the volume at time step t isprocessed sequentially from stage S1 to stage S3.
C2. The rendering thread produces images in an ascending order of time stept.
To satisfy condition C1, we create threads such that each of the threads isblocked with pthread cond wait() until it receives a signal from the previous
thread in the pipeline. This previous thread, on the other hand, sends a signal
to the next thread when it ﬁnishes an incoming task. In addition, the rendering thread sends a wake-up signal to the loading thread after completing arendering task.
Condition C2 is essential to guarantee that the stream data cannot pass eachother. If we lack this condition, we cannot build a correct pipeline on threads.
This condition can be satisﬁed by a FIFO policy. To realize this, each thread
has variable “tstep,” which stores the latest time step processed at the corresponding stage. This information is then used to prevent the stream datafrom overtaking each other. That is, a thread is allowed to process the data if
10it has already been processed by the previous threads. Otherwise, the thread
will be blocked with pthread cond wait() placed in a while loop.
Finally, once the time-varying data is sent to the GPU, the ﬁnal image will bequickly generated by the GPU. We employ a texture-based rendering method[29,30], which is fully accelerated by hardware components in the GPU, such as
texture mapping and alpha blending hardware. We currently use 3-D textures
rather than 2-D textures, because 2-D textures consume three times morememory [30] needed for a copy of the dataset for each orientation. Although
this might be a trivial problem for rendering of non-time-varying data, it iscritical for out-of-core (data-intensive) rendering.

5

Experimental Results

In order to evaluate the proposed pipeline in terms of performance, scalability,
and image quality, we implemented it using the C++ language, the OpenGL
library [24], the Cg toolkit [32], and the POSIX thread library [31]. For experiments, we use a desktop PC equipped with 2 GB of main memory and a
320GB SATA disk device. This PC has a Core 2 Quad CPU running at 2.66
GHz clock speed and an nVIDIA Quadro FX 4600 card having 768 MB ofvideo memory. The graphics card is connected to a PCI Express 16X bus. The
experiments are performed using Windows XP with driver version 162.55.
Table 2 summarizes six datasets D1–D6 used for experiments. Datasets D4,
D5, and D6 are high-resolution versions of datasets D1, D2, and D3, respectively. See also Fig. 6 for their visualization results. Datasets D1 and D2 [33]are ﬂuid dynamics data showing a turbulent jet and a turbulent vortex ﬂow,
respectively. Both datasets have low temporal coherence due to moving objects. However, dataset D1 has high spatial coherence at any time step, because it consists of many transparent voxels. Dataset D3 shows a sequenceof lung deformations. It represents a deformation process of nonrigid registration [34, 35] that aligns a lung dataset to another dataset obtained at adiﬀerent time. These datasets are slightly diﬀerent from each other due tobreathing. Since the deformations are small, this dataset has high temporal
coherence. All datasets have voxels of 1-byte scalar data.
Table 2 also shows the compression ratio for each dataset. As we mentioned in
Section 4.2, PVTC is a lossy compression based on VTC, which has the compression ratio of 6:1. On the other hand, the combination of PVTC and LZO
compression contributes to further reduction of the data size: the highest ratioof 12:1 and the lowest ratio of 1.1:1, each obtained by LZO compression. Thus,
the total compression ratio ranges from 6.5:1 to 71.7:1 after LZO compression.
However, LZO compression is not eﬀective for dataset D2. The coherence in11this dataset is not high, so that the match found in the sliding window is
not long enough to obtain a high compression ratio. Actually, the compression ratio is at most 1.1:1 even if we apply LZO compression directly to D2,
without PVTC. This indicates that the data originally lacks spatial coherence,
thus, we think that the compression ratio of 6.5:1 is not due to PVTC, whichkeeps coherence that will be exploited by LZO compression. We also thinkthat hardware supported compression methods such as PVTC, even though
they are lossy, are essential to achieve fast rendering for low-coherence data.

5.1 Rendering Performance

To evaluate the performance gain of the proposed pipeline, we now compare
it with three methods: (1) a straightforward method that transfers raw datafrom the storage device to the video memory; (2) a GPU-based method thatuses PVTC on the GPU but no compression method on the CPU; and (3) ahybrid method [10] that is not structured in a pipeline.
Fig. 7 shows average frame rates for all methods. Datasets are rendered ona screen with an appropriate size: a 256 × 256 pixel screen for D1 and D2; a
512 × 512 pixel screen for D3, D4, and D5; and a 1024 × 1024 pixel screen
for D6. The viewing direction is initially set to z-axis direction, and then it
is rotated 2 degrees around x- and y-axes when the time step is updated. For
every measurement, we rebooted the machine so that disk cashes are cleaned
before measurement. Thus, the frame rates in our experiments are mainly
determined by the compression ratio and the volume size rather than thenumber tmax of time steps (see also Section 5.2).
Our pipelined hybrid method achieves a video rate of 41 frames per second(fps) for dataset D4: 258 × 258 × 208 voxel data with 150 time steps. It also
demonstrates an almost interactive rate of 8 fps for the largest dataset D6:
512 × 512 × 295 voxel data with 411 time steps. As compared with prior
results [3, 10, 36], these performance results are obtained using more recent
hardware with more CPUs, but we think that our results are reasonable inperformance. For example, Akiba et al. [3] achieves 1.1 fps for 256×256×1024
voxel data with 400 time steps. Fout et al. [36] reports approximately 6 fps
for non-time-varying data of 512 × 512 × 512 voxels.
For all datasets, our pipeline mechanism achieves 33–122% improvement overthe hybrid method, which achieves 3.3–7.4 times higher rendering performance
than the raw renderer. The pipeline also achieves 62–201% improvement overthe GPU-based method. One interesting point here is that the hybrid method
fails to improve the performance for dataset D2, as compared with the GPUbased method. This means that LZO compression fails to obtain higher com12

pression ratio due to the lack of higher coherence. Actually, Table 2 shows thatthe total compression ratio of 6.5:1 is mostly obtained by PVTC. Despite thisfact, the pipeline mechanism increases the rendering rate from 66 fps to 111
fps. Therefore, we think that the pipeline mechanism is essential to maximize
the timing advantage of the hybrid method.
Fig. 8 shows the average execution time T and its breakdown normalized to asingle time-step data: (1) average time T1 spent for data loading at stage S1;(2) average time T2 spent for LZO decompression at stage S2; (3) average timeT3 spent for texture transfers at stage S3(a); and (4) average time T4 spent
for texture-based rendering at stage S3 except texture transfers. Note that theexecution time T for the pipelined hybrid method does not equal to the sumof breakdown T1 . . . T4 , because each breakdown is measured concurrently on
the corresponding thread.
We ﬁrst compare the pipelined hybrid method with the hybrid method toanalyze the eﬀects of the pipeline mechanism. In the hybrid method, overheads
T1 . . . T4 account for most of the entire time T . In contrast, the pipelined hybrid
method overlaps these four overheads, reducing execution time T . Actually,
we can see that time T is smaller than the total of time T1 . . . T4 in Fig. 8.
Thus, this overlap explains the 33–122% improvement mentioned before. In
many cases, the pipeline mechanism reduces the entire time T roughly totime T1 , because data loading operations at stage S1 is usually a performancebottleneck in the pipeline. Therefore, the decompression cost T2 is usuallyhidden by time T1 in our pipeline. Thus, the pipeline mechanism is capable ofhiding the overheads of the hybrid method by overlapping overhead T2 withother overheads. It also should be noted that the thread-related overheads
are much smaller than the entire time T . The overheads can be estimated bycomparing time T1 + T2 + T3 + T4 between the pipelined hybrid method andthe hybrid method.
We next analyze the eﬀects of data compression methods. In Fig. 8, we cansee that the GPU-based method achieves 3.0–3.8 times higher performance ascompared with the straightforward method. This improvement is achieved by
the reduction of data size, which decreases loading time T1 . Due to the same
reason, the hybrid method also achieves 1.2–2.0 times higher performance
than the GPU-based method, for all datasets except D2. For dataset D2,
the hybrid method fails to show the performance gain over the GPU-based
method, due to the lower compression ratio. Therefore, the hybrid method
slightly reduces loading time T1 , but it has longer execution time T due toadditional decompression overhead T2 .
Fig. 8 also indicates that most of execution time T is spent for loading operations if we do not use compression methods. Therefore, overlapping theloading stage with other stages is not so eﬀective in this case. Thus, we must13

reduce loading time T1 by using compression methods in order to maximize theeﬀectiveness of the pipeline mechanism. Note that loading time T1 is mainly
limited by the seek time of the storage device if we render small datasets, suchas D1 and D2. Actually, the seek time of our storage device (Western Digital
WD3200AAKS) is 8.9 ms, which is close to loading time T1 (approximately
10 ms) in Fig. 8.

We next analyze the behavior of the pipelined hybrid method. Fig. 9 shows atimeline chart explaining how the pipeline renders dataset D5 on the multicore PC. This chart is generated from a trace ﬁle obtained from an execution of
an instrumented program. Although the pipeline runs sequentially for the ﬁrst
three time steps, it is fully parallelized after time step 6. Thus, three threads
are running simultaneously on a multi-core CPU, showing an eﬀective utilization of computing resources. We think that such a parallel implementation isimportant to realize a full utilization in future many-core machines.

Finally, Fig. 10 presents the time-varying properties of the pipelined hybrid
method in terms of the compression ratio and the rendering performance.
It shows that the ﬁrst time step is rendered at approximately 10 fps, due
to the lack of pipeline eﬀects (latency hiding). Since our pipeline transfers acollection of three time steps at a time, the rendering performance decreases
when the rendering thread loads a compressed texture from the main memory.
Fig. 10 shows that this periodical behavior actually occurs at every three time
steps. To resolve this problem, we must split the rendering thread into two
diﬀerent threads: one for texture transfer at stage S3(a); and the other for theremaining stages S3(b)–(d). This modiﬁcation will improve the performanceon a quad-core machine, because it allows us to prefetch a compressed texture
while rendering the current volume. Thus, we can overlap time T3 with T4 .

With respect to the compression ratio, it shows an irregular behavior in Fig.
10(a), but the rendering performance shows a periodical behavior. This meansthat the pipeline mechanism successfully hides the irregularity at the loading
and decoding stages. However, our pipeline cannot hide the overhead incurredat the rendering stage, namely the last stage in the pipeline. For example, we
observe that the cgGLBindProgram() function called for program switchingby the rendering thread spends longer time every after two successive calls of
the function. This explains why the rendering performance decreases at time
step t, as compared to time step t − 1, where t mod 3 = 2. Due to the same
reason, the entire performance can vary according to the viewing point anddirection.
145.2 Scalability Analysis

We ﬁrst analyze the scalability of our system in terms of the capacity. Letn be the volume size. The system assumes that the video memory is large
enough to store a compressed texture Xt/3 , where 0 ≤ t < tmax . Since Xt/3
is generated from three raw volumes Vt , Vt+1 , and Vt+2 with a compression
ratio of 6:1, a compressed texture requires n3 /2 memory space. Therefore, weneed at least 512 MB of video memory to render 1024 × 1024 × 1024 voxel
data. In other words, 512 MB of video memory is not suﬃcient to render suchdata, because it lacks the memory space for frame buﬀer. Thus, we think thatthe maximum cubic size is limited by 1024 × 1024 × 1024 voxels on current
graphics cards.
We next analyze how the rendering performance is aﬀected by the multilevel storage hierarchy consisting of the disk device, main memory, and video
memory. Let σ be the compression ratio obtained by the hybrid method. Lets = n3 /2τ be the ﬁle size of a double-compressed data, where τ = σ/6 represents the compression ratio obtained by LZO. Let F be the average frame
rate for a single time step. Since the bottleneck stage in the pipeline limits
the frame rate, the frame rate is given by:
F = 1/max(T1 , T2 , T3 + T4 ).

(1)

Each of times T1 . . . T4 can be simply modeled as follows: the loading thread
receives s/3 data in average from the disk device; the decoding thread receives
the s/3 data and then sends sτ /3 data; the rendering thread receives the sτ /3data and then loads O(n3 ) single-precision data to produce the image. We
assume here that the rendering cost is proportional to the number of voxels
in a raw volume. Therefore, we obtain
T1 = s/3B1 ,T2 = s(1 + τ )/3B2 ,T3 = sτ /3B3 ,T4 = 2Ksτ /B4 ,

(2)
(3)
(4)
(5)

where K represents the coeﬃcient of the rendering cost and B1 . . . B4 represent
the throughput of data loading, that of LZO decoding, that of texture transfer,and that of rendering, respectively. Equations (2)–(5) indicate that the volume
size n does not determine the bottleneck stage in the pipeline. The bottleneckstage is determined by the compression ratio τ of LZO and by the throughputs
B1 . . . B4 .
The throughputs can be estimated by using the measured values of compres15

sion ratio and execution time presented in Table 2 and Fig. 8, respectively.
According to this estimation, the lowest throughputs we obtained are as follows: B1 = 13; B2 = 452; B3 = 205; and B4 = 1352K. Using these measured
throughputs and Equations (2)–(5), we obtain the following relations:T1 > T2 ≥ T3 + T4 ,T1 ≥ T3 + T4 > T2 ,T3 + T4 > T1 > T2 ,T3 + T4 > T2 ≥ T1 ,if τ ≤ 0.3,
if 0.3 < τ ≤ 8.3,
if 8.3 < τ < 33.8,
otherwise.

(6)
(7)
(8)
(9)

In summary, the loading thread or the rendering thread can be a performancebottleneck in our desktop PC. The rendering thread determines the frame rate
F when the compression ratio σ is higher than 49.8:1. Otherwise, the loading
thread determines the frame rate F .
Actually, Figs. 8(c) and (f) show that the rendering thread becomes a performance bottleneck when viewing datasets D3 (σ = 67.0) or D6 (σ = 71.7). In
particular, this bottleneck is critical if we render higher resolution data such as
512×512×512 voxel data. Therefore, we think that 512×512×512 voxel datais the maximum volume size that can be interactively rendered by our current
system. We also think that this size will be increased by future GPUs, which
achieve shorter rendering time T4 and transfer time T3 by a higher memory
bandwidth and by a new graphics bus, respectively.

5.3 Image Quality

In order to evaluate the quality of rendered images, we measure PSNR (peaksignal-to-noise ratio) values and analyze the error distribution. A higher PSNR
value here represents higher quality with less noise. In general, PSNR valuesof at least 30 dB are desired for rendered images. As a rule of thumb, a
PSNR value over 40 dB indicates that the image is very close to the original
and a value between 30 to 40 dB means that the distortion is visible but
acceptable [37].
Fig. 11 shows PSNR values measured for each dataset. For all datasets, weobserve PSNR values of more than 33 dB. These results are competitive toprior work [7, 36] that achieves similar quality ranging from 28 to 49 dB.
Fig. 12 shows the breakdown of pixel errors in produced images. It represents
how exactly the hybrid method produces the original pixels. In most cases,
more than 95% of pixels have an error within pixel value 3, as compared withthe exact image produced from raw data. We think that this distribution is
acceptable because pixel values are in the range [0,255].
16

Fig. 13 shows two pairs of images, each produced from the raw data and fromthe PVTC compressed data. Fig. 14 shows subtraction images produced fromeach pair. By comparing them in zoomed views, we can see some artifacts
in Fig. 13(b) and Fig. 13(d), especially where neighboring pixels do not have
similar values. More precisely, we can see a glow eﬀect or rim in the lower rightof the zoomed subimage in Fig. 13(b), where this is not a case in Fig. 13(a).
This artifact could mislead the user into a diﬀerent interpretation of simulation
results when the original data is not available for comparison. Nevertheless,these artifacts do not signiﬁcantly change the overall appearance of produced
images. Thus, we think that the image quality of PVTC is tolerable for timeseries analysis.

6Conclusion

We have presented a decompression pipeline for accelerating volume renderingof time-varying scalar data. The proposed pipeline is based on a hybrid compression method and a thread-based pipeline mechanism. The hybrid method
uses a combination of two diﬀerent compression methods, namely PVTC and
LZO compression, each running on the GPU and on the CPU. Our combination is designed to exploit both temporal and spatial coherence in time-varying
data. It also makes a good use of hardware components on the GPU and alarge cache on the CPU. The pipeline mechanism overlaps data loading, data
decompression, and rendering stages, in order to increase the rendering rate.
Experimental results show that our pipelined method achieves a video rateof 41 fps in average for 258 × 258 × 208 voxel data with 150 time steps. It
also demonstrates an almost interactive rate of 8 fps for 512 × 512 × 295
voxel data with 411 time steps. Our pipeline mechanism achieves 33–122%
improvement over the nonpipelined hybrid method, which achieves 3.3–7.4
times higher rendering performance than the raw renderer. We also ﬁnd thatmost of the execution time is spent for loading operations if we do not use
compression methods. Therefore, we think that loading time must be reducedby compression methods in order to maximize the performance beneﬁt of thepipeline mechanism. Furthermore, the pipeline mechanism is also useful to hide
the decompression overhead inherent in multi-stage compression methods.
One future work is to develop an adaptive framework for visualization of higher
resolution data. Although our compression pipeline provides an eﬀective solution to visualization of large time-step data, it assumes that the capacityof video memory is large enough to store the entire of a compressed texture.
A scalable, parallelization scheme [38] will be useful to reduce the memory
requirement for a graphics card. We also think that a multiresolution data
structure [39] contributes to progressive visualization that increases the inter17

activity for large volume.

Acknowledgments

This work was partly supported by MEXT Grant-in-Aid for Scientiﬁc Research
for Young Researchers (19700061), JSPS Grant-in-Aid for Scientiﬁc Research(B)(2)(18300009), and the Global COE Program “in silico medicine” at Osaka
University. We are also grateful to the anonymous reviewers for their valuablecomments.

References[1] M. Levoy, Display of surfaces from volume data, IEEE Computer Graphics andApplications 8 (3) (1988) 29–37.
[2] K.-L. Ma, Visualizing time-varying volume data, Computing in Science andEngineering 5 (2) (2003) 34–42.
[3] H. Akiba, K.-L. Ma, J. Clyne, End-to-end data reduction and hardware
accelerated rendering techniques for visualizing time-varying non-uniform gridvolume data, in: Proc. 4th Int’l Workshop Volume Graphics (VG’05), 2005, pp.
31–39.
[4] A. P. Binotto, J. L. Comba, C. M. Freitas, Real-time volume rendering of timevarying data using a fragment-shader compression approach, in: Proc. 6th IEEE
Symp. Parallel and Large-Data Visualization and Graphics (PVG’03), 2003, pp.
69–76.
[5] N. Fout, K.-L. Ma, J. Ahrens, Time-varying, multivariate volume data
reduction, in: Proc. 20th ACM Symp. Applied Computing (SAC’05), 2005, pp.
1224–1230.
[6] J. Schneider, R. Westermann, Compression domain volume rendering, in: Proc.
14th IEEE Visualization Conf. (VIS’03), 2003, pp. 293–300.
[7] E. B. Lum, K.-L. Ma, J. Clyne, A hardware-assisted scalable solution forinteractive volume rendering of time-varying data, IEEE Trans. Visualization
and Computer Graphics 8 (3) (2002) 286–301.
[8] M. Strengert, M. Magallón, D. Weiskopf, S. Guthe, T. Ertl, Large volume
visualization of compressed time-dependent datasets on GPU clusters, Parallel
Computing 31 (2) (2005) 205–219.
[9] F. F. Bernardon, S. P. Callahan, J. L. Comba, C. T. Silva, An adaptive
framework for visualizing unstructured grids with time-varying scalar ﬁelds,
Parallel Computing 33 (6) (2007) 391–405.

18

[10] D. Nagayasu, F. Ino, K. Hagihara, Two-stage compression for fast volume
rendering of time-varying scalar data, in: Proc. 4th Int’l Conf. Computer
Graphics and Interactive Techniques in Australasia and Southeast Asia
(GRAPHITE’06), 2006, pp. 275–284.
[11] J. Montrym, H. Moreton, The GeForce 6800, IEEE Micro 25 (2) (2005) 41–51.
[12] J. D. Owens, D. Luebke, N. Govindaraju, M. Harris, J. Krüger, A. E. Lefohn,
T. J. Purcell, A survey of general-purpose computation on graphics hardware,
Computer Graphics Forum 26 (1) (2007) 80–113.
[13] I. Buck, K. Fatahalian, P. Hanrahan, GPUBench: Evaluating GPU performancefor numerical and scientiﬁc applications, in: Proc. 1st ACM Workshop GeneralPurpose Computing on Graphics Processors (GP2 ’04), 2004, p. C-20.
[14] B. He, N. K. Govindaraju, Q. Luo, B. Smith, Eﬃcient gather and scatter
operations on graphics processors, in: Proc. High Performance Networking andComputing Conf. (SC’07), 2007, 12 pages (CD-ROM).
[15] OpenGL
Extension
Registry,
Gl nv texture
compression vtc, http://oss.sgi.com/projects/ogl-sample/registry/NV/
texture compression vtc.txt (2004).
[16] M.F.X.J. Oberhumer, LZO real-time data compression library, http://www.
oberhumer.com/opensource/lzo/ (Oct. 2005).
[17] H. Yu, K.-L. Ma, A study of I/O methods for parallel visualization of large-scaledata, Parallel Computing 31 (2) (2005) 167–183.
[18] T. Chiueh, K.-L. Ma, A parallel pipelined renderer for time-varying volume
data, in: Proc. 2nd Int’l Symp. Parallel Architectures, Algorithms and Networks(I-SPAN’97), 1997, pp. 9–15.
[19] W. Bethel, B. Tierney, J. Lee, D. Gunter, S. Lau, Using high-speed WANs andnetwork data caches to enable remote and distributed visualization, in: Proc.
High Performance Networking and Computing Conf. (SC’00), 2000.
[20] X. Cavin, C. Mion, A. Filbois, COTS cluster-based sort-last rendering:
Performance evaluation and pipelined implementation, in: Proc. 16th IEEE
Visualization Conf. (VIS’05), 2005, 8 pages (CD-ROM).
[21] J. Clyne, The multiresolution toolkit: Progressive access for regular gridded
data, in: Proc. 3rd IASTED Int’l Conf. Visualization, Imaging, and Image
Processing (VIIP’03), 2003, pp. 152–157.
[22] S. Guthe, M. Wand, J. Gonser, W. Straßer, Interactive rendering of large volumedata sets, in: Proc. 13th IEEE Visualization Conf. (VIS’02), 2002, pp. 53–60.
[23] T. Furumura, L. Chen, Parallel simulation of strong ground motions during
recent and historical damaging earthquakes in tokyo, japan, Parallel Computing
31 (2) (2005) 149–165.

19

[24] D. Shreiner, M. Woo, J. Neider, T. Davis, OpenGL Programming Guide, 5th
Edition, Addison-Wesley, Reading, MA, 2005.
[25] OpenGL
Extension
Registry, Gl ext texture compression dxt1, http://oss.sgi.com/projects/
ogl-sample/registry/EXT/texture compression dxt1.txt (2004).
[26] A. Grama, A. Gupta, G. Karypis, V. Kumar, Introduction to Parallel
Computing, 2nd Edition, Addison-Wesley, Reading, MA, 2003.
[27] J. Ziv, A. Lempel, A universal algorithm for sequential data compression, IEEE
Trans. Information Theory 23 (3) (1977) 337–343.
[28] nVIDIA Corporation, CUDA Programming Guide Version 1.1, http://developer.nvidia.com/cuda/ (Nov. 2007).
[29] B. Cabral, N. Cam, J. Foran, Accelerated volume rendering and tomographic
reconstruction using texture mapping hardware, in: Proc. 4th Symp. Volume
Visualization (VVS’94), 1994, pp. 91–98.
[30] M. Hadwiger, J. M. Kniss, K. Engel, C. Rezk-Salama, High-quality volume
graphics on consumer PC hardware, in: SIGGRAPH 2002, Course Notes 42,
2002.
[31] B. Nichols, B. Buttlar, J. P. Farrell, Pthreads Programming, O’Reilly &Associates, Newton, MA, 1996.
[32] W. R. Mark, R. S. Glanville, K. Akeley, M. J. Kilgard, Cg: A system forprogramming graphics hardware in a C-like language, ACM Trans. Graphics
22 (3) (2003) 896–897.
[33] K.-L.
Ma,
Time-Varying∼
Volume Data Repository, http://www.cs.ucdavis.edu/ ma/ITR/tvdr.html(2003).
[34] J. V. Hajnal, D. L. Hill, D. J. Hawkes (Eds.), Medical Image Registration, CRC
Press, Boca Raton, FL, 2001.
[35] F. Ino, K. Ooyama, K. Hagihara, A data distributed parallel algorithm for
nonrigid image registration, Parallel Computing 31 (1) (2005) 19–43.
[36] N. Fout, K.-L. Ma, Transform coding for hardware-accelerated volume
rendering, IEEE Trans. Visualization and Computer Graphics 13 (6) (2007)
1600–1607.
[37] Y. Wang, J. Ostermann, Y.-Q. Zhang, Video Processing and Communications,
Prentice Hall PTR, Upper Saddle River, NJ, 2001.
[38] H. Childs, M. Duchaineau, K.-L. Ma, A scalable, hybrid scheme for volume
rendering massive data sets, in: Proc. 6th Eurographics Symp. Parallel Graphics
and Visualization (EGPGV’06), 2006, pp. 153–162.
[39] E. LaMar, B. Harmann, K. I. Joy, Multiresolution techniques for interactive
texture-based volume visualization, in: Proc. 10th IEEE Visualization Conf.
(VIS’99), 1999, pp. 355–361.

20

Compressedvoxel block B’1

...

...

4 x 4 x 1 voxel block B1
...

...

Approximation

BN

Volume in RGB format

B’N
Compressed RGB texture

2 bits per voxel representing
interpolation modes

24 bits per voxel
v1,1

m1,1

c1
+

v4,4

...

...

...

...
Partition

Approximation

Bk: voxel block(48 bytes)

c2

m4,4

Tk:
Rk:
Interpolation Representative RGB
table (4 bytes)
colors (4 bytes)

B’k: Compressed voxel block (8 bytes)

Fig. 1. Volume texture compression (VTC).

Vt: Raw volume~

Vt: Approximated volume

GPU

Mg: Compression method for the GPU (PVTC)

Xt: Compressed texture

CPU

Mc: Compression method for the CPU (LZO)

Zt: Double-compressed texture
Preprocessing phase

Storage

Visualization phaseFig. 2. Overview of hybrid compression method. This method uses two diﬀerent
compression methods, Mc and Mg , each running on the CPU and the GPU, respectively. Note that raw data must be compressed in advance of visualization.

21

Compressed voxel block eachcontaining time-series data

Vt+2(x,y,z)

Block 1

Block 2

R

G B

R

G B

t

t+1 t+2t

t+1 t+2

...

...

B

Vt+1(x,y,z)

G
R

...

VTC

C t/3 (x,y,z)

...

...

Vt (x,y,z)Compressed texture X t/3

Data packing

Fig. 3. Packed volume texture compression (PVTC). Time-series scalar voxels inthe same location are packed into an RGB voxel. The packed RGB data is thencompressed using VTC, which generates a sequence of compressed blocks, eachcontaining time-series data.

DoubleStoragecompressed data

Loading thread

Z t/3

Z t/3 +1

Data

Data

loadingloading

Z t/3

CPU Decoding thread

...

LZO

LZO

decoding

decoding

PVTC compressed data X t/3
Rendering thread

X t/3

GPU

VolumeVolumerendering

rendering

PVTC

PVTC

decoding

decodingTimestep tTimestep t+3

...
...

Time

Fig. 4. Decompression pipeline implemented using multiple threads running on amulti-core CPU. The pipeline consists of three stages, each responsible for dataloading, LZO decoding, and volume rendering.

22

// send a signal to the decoding thread
pthread_mutex_lock(&ld->mutex);ld->tstep += 3;// update time step
pthread_cond_signal(&ld->ready);
pthread_mutex_unlock(&ld->mutex);// receive a signal from the rendering thread
pthread_mutex_lock(&rd->mutex);
while (&ld->tstep > &rd->tstep+lookahead) {pthread_cond_wait(&rd->ready,&rd->mutex);}
pthread_mutex_unlock(&rd->mutex);

Rendering thread// receive a signal from the loading thread
pthread_mutex_lock(&ld->mutex);
while (&ld->tstep <= &dc->tstep) {pthread_cond_wait(&ld->ready,&ld->mutex);}
pthread_mutex_unlock(&ld->mutex);// send a signal to the loading thread
pthread_mutex_lock(&rd->mutex);rd->tstep++;// update time step
pthread_cond_signal(&rd->ready);
pthread_mutex_unlock(&rd->mutex);// Stage S2
LZO decoding(&dc->tstep);// send a signal to the rendering thread
pthread_mutex_lock(&dc->mutex);dc->tstep += 3;// update time step
pthread_cond_signal(&dc->ready);
pthread_mutex_unlock(&dc->mutex);

Variables
typedef struct stage_tag {pthread_mutex_t mutex;
pthread_cond_t ready;// ready signal
int
tstep;// time step} stage_t;
stage_t *ld, *dc, *rd;// each for stage S1, S2, S3
const int lookahead = 7;// Stage S3(a)if (rd->tstep % 3 == 0) {glCompressedTexImage3d(&rd->tstep);}
glutPostRedisplay();// receive a signal from the decoding thread
pthread_mutex_lock(&dc->mutex);
while (&dc->tstep <= &rd->tstep) {pthread_cond_wait(&dc->ready,&dc->mutex);}
pthread_mutex_unlock(&dc->mutex);// Stage S3(b) and (c)Texture-based rendering(&rd->tstep);

Idle callback function// Stage S1
Data loading(&ld->tstep);

Decoding thread

Display callback function

Loading thread

Fig. 5. Pseudocode of our thread-based pipeline mechanism. See text for details.
Variable “lookahead” limits the number of stream data in the pipeline. This pipelineis allowed to process nine time steps of the data simultaneously, because each of thethree stages can have the data containing three time steps.

23

(a)

(b)

(c)
Fig. 6. Rendering results obtained from (a) D1: small jet, (b) D2: small vortex, and(c) D3: middle lung.

24

Rendering performance (fps)120
Raw
GPU-based (PVTC)
Hybrid (PVTC+LZO)
Pipelined hybrid

100
80
60
40
20
0
D1: Smalljet

D2: Small D3: Middle D4: Middle D5: Middle D6: Large
vortexlung
jet
vortexlung

Fig. 7. Rendering performance in frames per second (fps). The performance is theaverage value.

25

80

60
50
40

T’4: RenderingT’3: CPU-GPU transfer
T’2: LZO decompression
T’1: Data loading
T: Execution time

50
Execution time (ms)

Execution time (ms)60

T’4: RenderingT’3: CPU-GPU transfer
T’2: LZO decompression
T’1: Data loading
T: Execution time

70

30
20

40
30
20
10

10
0

0
Raw

GPU-based

HybridPipelined hybrid

Raw

GPU-based

(a)
180

300

120
100

T’4: RenderingT’3: CPU-GPU transfer
T’2: LZO decompression
T’1: Data loading
T: Execution time

250
Execution time (ms)

Execution time (ms)

140

80
60
40

200
150
100
50

20
0

0
Raw

GPU-based

HybridPipelined hybrid

Raw

GPU-based

(c)
300

1400

150

Pipelined hybrid

T’4: RenderingT’3: CPU-GPU transfer
T’2: LZO decompression
T’1: Data loading
T: Execution time

1200
Execution time (ms)200

Hybrid(d)T’4: RenderingT’3: CPU-GPU transfer
T’2: LZO decompression
T’1: Data loading
T: Execution time

250
Execution time (ms)Pipelined hybrid(b)T’4: RenderingT’3: CPU-GPU transfer
T’2: LZO decompression
T’1: Data loading
T: Execution time

160

Hybrid

100

1000

50

800
600
400
200

0

0
Raw

GPU-based

HybridPipelined hybrid

Raw

GPU-based(e)HybridPipelined hybrid(f)Fig. 8. Breakdown analysis of execution time required for rendering of a volume.
Each of subﬁgures (a)–(f) corresponds to datasets D1–D6, respectively, showing the
average time spent for rendering of a single time-step data.

Time step t=0

t=3

t=6

t=9

t=12

t=15

Loading thread

Decoding thread

Rendering thread
Time (ms)
0

50

100

150

200

250

300

350

400

450

500

Fig. 9. Timeline chart showing how dataset D5 is rendered by the pipelined hybrid
method. A rectangle represents a busy state of the corresponding thread.

26

30
25

40

20
30
15
20

10
Rendering performance

10

5

Compression ratio

0
20

40

60
80
100
t: Time step

120

12.6

60

12.4

50

12.2

40
12
30
11.8

20
Rendering performance
10

11.4
0

140

11.6

Compression ratio

0

0
0

70

σ: Compression ratio

50

Rendering performance (fps)35

σ: Compression ratio

Rendering performance (fps)

60

20

(a)

40
60
t: Time step

80

(b)

Fig. 10. Compression ratio and rendering performance for (a) datasets D4 and (b)D5 with diﬀerent time steps.

60
47.3

PSNR (dB)

50
40

46.1

45.8

47.8

36.1

33.7

30
20
10
0
D1: Smalljet

D2: SmallvortexD3: Middle D4: Middle D5: Middle
lungjet
vortex

D6: LargelungFig. 11. Image quality of rendering results. PSNR values are averaged over all time
steps.

27

100

Ratio (%)

8060
Error >= 5
Error = 4
Error = 3
Error = 2
Error = 1
Error = 0

40

20

0
D1: Smalljet

D2: SmallvortexD3: Middle D4: Middle D5: Middle
lungvortexjet

D6: LargelungFig. 12. Breakdowns of pixel errors in rendered images. Pixel values are in the range[0, 255].

(a)

(b)

(c)

(d)

Fig. 13. Comparison of rendering results using datasets D4 at time step 20 and D5
at time step 65. The left-hand side is the result from raw data while the right-handside is that from PVTC-compressed data.

28

(a)

(b)

Fig. 14. Subtraction images for (a) datasets D4 and (b) D5, each produced from apair of rendering results shown in Fig. 13. Subtracted pixels are rendered on thegray campus.

29

Table 1
Comparison of rendering systems. Compression ratio σ is given by σ = D1 /D2 ,where D1 and D2 represent raw size and reduced size, respectively. Note that thistable does not show a fair comparison, because compression ratio and image quality
are not measured under the same condition, such as datasets, machines, and transfer
functions.
StrategyHybridGPU

CPU

SystemPipelined

Lossless

Supported

Ratio

Image quality

decompression

compression

coherence

σ

(dB)

No

Temporal and spatial

6–72

PSNR=33–48

1.5–50—This paper

Yes

[10]

No

[3]

No

Yes

Spatial

[4]

No

Yes

Temporal and spatial

7–154

—

[5]

No

No

Temporal or spatial

3.2–3.4

PSNR=32–71

[6]

No

No

Temporal and spatial

8–60

SNR=10–26

[7]

No

No

Temporal

2–8

PSNR=28–49

[8]

No

Yes

Spatial

—

—

[9]

Yes

No

Temporal

7.1–34.7

SNR=5–43

Table 2
Datasets used for experiments.
DatasetVolume size

Time Raw ﬁle size per Compression ratio(voxel)step time step (MB) PVTC PVTC+LZO

D1: small jet

129 × 129 × 104 150

1.7

11.5

D2: small vortex

128 × 128 × 128

99

2.0

6.5

D3: middle lung

256 × 256 × 148 411

9.3

D4: middle jet

258 × 258 × 208 150

13.2

21.5

99

16.0

12.0

512 × 512 × 295 411

73.8

71.7

D5: middle vortex 256 × 256 × 256
D6: large lung

30

6

67.0