Introdução à Visualização Volumétrica
por
Roberto de Beauclair Seixas, IMPA
Anselmo Cardoso de Paiva, UFMA
Marcelo Gattass, PUC-Rio
Visualização Volumétrica

Conjunto de técnicas para visualização de
dados associados a regiões de um volume
– Objetivo : exibição do interior de objetos
volumétricos, a fim de explorar sua estrutura (ou
falta dela) e facilitar sua compreensão
[McCormick, 1987].
Radiografia Padrão (Raios X)
2K x 2K x 2 Bytes
4 chapas
Angiografia: 40 chapas
Tomografia Computadorizada (TC)
Receptores
Raio X
Exemplos de Imagens de TC
Normalmente 512 x
512 de 1/2 a 2 mm
Exame de 5 a 30 min.
Ar
Gordura
Água
Rim
Pâncreas
Músculo
Liver
-1000
-110 15
0
27  15
35 10
40 10
55 10
Números em unidades de Hounsfield (HU)
Osso
Esponjoso
200-400
Osso
Compacto
>1000
Problemas da TC

Pequena resolução temporal para movimento cardíaco;

Presença de artefatos inerentes ao método de aquisição;

Resolução espacial relativamente pequena;

Inabilidade de detecção de doenças em estágios
incipientes que não tenha resultado ainda em
significantes alterações dos coeficientes de densidade
dos tecidos.
Ressonância Magnética Nuclear (MR)
Grande campo
magnético
Mesa com
trilhos
Principais Vantagens de MR

Produzem contraste de tecidos moles superior as outras
modalidades, sem a necessidade de agentes de contraste
externo;

Permitem a detecção de doenças anteriormente ao
aparecimento de grandes mudanças anatômicas ou
fisiológicas;

Fornecem também informação fisiológica e funcional;

As imagens podem ser adquiridas em planos arbitrários,
através de manipulação eletrônica sem necessidade de
mudanças na postura do paciente;

A ausência de radiação ionizante permite a realização de
estudos freqüentes sobre o paciente.
Desvantagens de MR

Dificuldades no estudo de calcificações;

Suscetibilidade a movimentos do paciente durante a
aquisição, por ser um processo lento;

Impossibilidade de aquisição de dados de pacientes em
sistemas artificiais de suporte a vida (UTI);

Inexistência de uma escala de valores absolutos para
um determinado conjunto de dados;

Alto custo.
Imagem MR
512 x 512 x 2 Bytes
30 fatias/exame
Ultra-som
Imagem de Ultra-som
512 x 512 x 1 Byte
Medicina Nuclear
Câmera Gama
SPECT (Single Photon Emission
Computed Tomography)
PET (Positron Emission
Tomography)
256 x 256 x 2 Bytes
bone scan
galium scan
Exames Médicos
Exame
Tamanho
Bytes
MBytes
Imagens
MBytes
da imagem
por pixel
por imagem
por exame
por exame
Raio-X Digital
2K x 2K
2
0.8
5
40
Angiografia
Digital
Tomografia
Computadorizada
Ultra-som
2K x 2K
2
0.8
40
320
512 x 512
2
0.524
25
13
512 x 512
1
0.262
40
10
Medicina Nuclear
256 x 256
2
0.131
6
0.8
Ressonância
Magnética
512 x 512
2
0.524
30
16
Sísmica
Formação de um traço sísmico
*
Seção Sísmica
Reconstrução das Fatias
Ys
Zs
volume das
fatias (slices)
Xs
Classificação do Voxel
opacidade
Voxel
Osso
1.0
Músculo
Gordura
Ar
Amarelo
Vermelho
Branco
densidade
Tipos de Dados

Enumeráveis (material, litologia, ...)

Escalares (temperatura, pressão, ...)

Vetoriais (velocidade, aceleração, ...)

Tensoriais (tensão, deformação, ...)
Estrutura dos Dados
Grade Cartesiana
(i, j, k)
Grade Estruturada
(x[i,j,k],y[i,j,k],z[i,j,k])
Grade Regular
(i*dx, j*dy,k*dz)
Grade Estruturada
por blocos
Grade Retilínea
(x[i],y[j],z[k])
Grade Não Estruturada
{(x[i],y[i],z[i]) ,
e=(v1, v2 , v3)}
Interpolação
– Matriz de células
 Matriz de voxels
–
–
análogo 3D do pixel
(i, j, k)
– interpolação trilinear
– imagens mais suaves
voxels
pontos da
grade
pontos da grade
valores interpolados
valor
valor
x
x
Métodos de Visualização

Indiretos: por extração de superfícies implícitas
+
+
–
–

representação por polígonos [Zbuffer]
dados menores
precisa ser refeito quando muda a classificação
dificuldade de modelar objetos amorfos
Diretos: por modelos de iluminação volumétrica
+
+
–
–
geração de imagens diretamente a partir dos dados volumétricos
visualização de múltiplas características, inclusive de dados amorfos
grande volume de dados
não usa (por enquanto) hardware gráfico
Rendering Direto



Mapeamento dos elementos de volume direto
no espaço da imagem
Apropriados para a visualização de objetos
amorfos
Mais lentos
Algoritmo ray-casting
Volume
de dados

Ordem da imagem.
– para cada pixel
•
lance um raio e encontre
os voxels que são
interceptados
– fim para
Lançamento dos Raios
plano de
visualização
y
volume dos
raios
x
z
volume dos
slices
raio
Partição dos Raios
x
ta
z
voxel
tb
y
amostra
t=min(x,y,z)
Interpolação
x,i
(xa, ya, za)
x
i a 
 Dx
y
j a 
 Dy
za 

k
 Dz
z,k
y,j
xa % Dx
x
Dx
ya % Dy
y
Dy
za % Dz
z
Dz
Interpolação no voxel
v( x, y, z )  1  x 1  y 1  z vi , j ,k 
vi +1,j,k +1
vi,j,k +1
vi+1,j,k
y x
vi,j,k
v z
vi+1,j+1,k+1
vi,j +1,k +1
vi+1,j+1,k
vi,j+1,k
 x 1  y 1  z vi 1, j ,k 
1  x  y 1  z vi , j 1,k 
 x  y 1  z vi 1, j 1,k 
1  x 1  y  z vi , j ,k 1 
 x 1  y  z vi 1, j ,k 1 
1  x  y  z vi , j 1,k 1 
 x  y  z vi 1, j 1,k 1
Classificação do Voxel
Voxel
tons e opacidade
cor
Osso
Branco
(255,255,255)
1.0
Músculo
Vermelho
(255,0,0)
Amarelo
(255,255,0)
Gordura
ar
0
amplitude
255
0
velocidade
ou densidade
255
Diferentes Funções de Transferência
Sensibilidade à função de transferência
1
1
0 .8
0 .8
0 .6
0 .6
0 .4
0 .4
0 .2
0 .2
0
0
50
1 00
1 50
2 00
0
0
50
1 00
1 50
2 00
Iluminação de um voxel
C  Iaka  I k (N  L)
luz d
N
Gel A
Gel B
L
N  V ( X )
Estimativa do vetor normal
f  x, y, z 
N
f  x, y, z 






 f xi 1 , y j , zk   f xi 1 , y j , zk  / 2x ,


f xi , y j , zk    f xi , y j 1 , zk   f xi , y j 1 , zk  / 2y ,
 






f
x
,
y
,
z

f
x
,
y
,
z
/
2

z
,
i
j
k 1
i
j
k 1


1a ordem
2a ordem
Interpolações em Visualização Volumétrica
Phong
Debrin
Gouraud
normal
interpolada
normal
do voxel
aC interpolada
aC do voxel
cor interpolada
cor do voxel
Etapa de Composição

Para cada raio:

Processo de acumulação
gera amostras de cor Cl (pi) e opacidades a(pi)
» reamostragem dos dados dos voxels, em k amostras
igualmente espaçadas
I= t Ib +(1-t)I0
I0 = cor do objeto
Ib = cor do fundo
I = cor resultante
t = coeficiente de transparência
Influência de um Voxel
r=0
g=0
b=0
r0
g0
b0
r1
g1
b1
a0
a1
a0
r = a 0 r0
g = a0 g0
b = a0 b0
a  a0
r = a0 r0+1a0  a1 r1
g = a0 g0+1a0  a1 g1
b = a0 b0+1a0  a1 b1
a  a0 +1a0  a1
Composição no raio
a out  a in  (1  a in )a
C out  C
in
 (1  a
in
)a C
ray
Cin
C  C out a out
Cout
ezy= (y)z +(-z)y +cyz
Tripod
z
y
exy= (y)x + (-x)y +cxy
(x0+x,y0 +y,z0 +z)
(x0,y0,z0)
x
exz= (z)x + (-x)z +cxz
exy exz
ezy
+

+
-
-


+
-
Valor inicial e incremento
(1/2,1/2,1/2)
x++
y++
x
exy= (y)x + (-x)y +cxy
(y-x)/2
y
exz= (z)x + (-x)z +cxz
(z-x)/2
z
ezy= (y)z +(-z)y +cyz
(y-z)/2
z++
x
z
y
Caminhamento discreto
Bresenham
Cohen
Efeito da amostragem
Bresenham
Tripod
Partição celular
Volume
de dados
Voxel
Partição
uma amostra no meio da partição
marcador da partição
Partição na grade
nx  (ib  ia )  1
z
n y  ( jb  ja )  1
(ia ,ja ,0) nz  (kb  ka )  1
tzi= ta
tyi
t x  (t xf  t xi ) nx
t y  (t yf  t yi ) n y
tyf=tb tzf
(ib ,Ny ,kb)
t z  (t zf  t zi ) nz
y
Partição celular: algoritmo
Dados: txi,tyi,tzi, txf,tyf,tzf, nx,ny,nz
dtx = txf/nx; dty = tyf/ny; dtz = tzf/nz;
tx=txi; ty=tyi; tz=tzi;
t1 = min(tx,ty,tz) e w é o eixo do mínimo
n=nx+ny+nz;
while ( n > 0 )
tw += dtw;
n--;
t2 = min(tx,ty,tz) e w é o eixo do mínimo
Sample ((t1+t2)/2);
t1=t2;
Otimizações

Velocidade
– Refinamento progressivo
– Terminação adaptativa do raio
– Estruturas Hieráquicas

Qualidade da imagem
– aumento do número de amostra no raio
– lançamento de mais raios
– melhora esquema de interpolação
Refinamento Progresivo
amostragem inicial
primeira subdivisão
segunda subdivisão
subdivisão final
pixels sendo visitados
pixels já visitados
Refinamento Progressivo: Exemplo
2x2
4x4
8x8
16x16
32x32
64x64
128x128
256x256
Algoritmo Ray Casting
 Alto
custo computacional
 Apresenta todo o conjunto de dados
 Facilmente paralelizável
 Utilizado quando se deseja imagens de
alta qualidade [Elvins,1992]
Raio que calcula a cor de um pixel
(i,j)
Imagem
L
ta
(x, y, z)
t I
tb
Volume
t
I = I(x,y,z,-t,) =
radiância ou
intensidade específica
num ponto (x,y,z) ou t
na direção -t,
na freqüência .
Um pouco de física
I(x,n,) = c h   (x,n,)
c = velocidade da luz (2.997925x108 m/s),
h = constante de Planck (6.626x10-34 J.s),
 = freqüência da luz (4.3x1014 a 7.5x1014 Hz ou s-1)
 (x,n,) = densidade de fótons na freqüência 
no ponto x, na dir n e (m-3 sr -1).
Unidade de I = J/s m-2 sr-1 = W/m2 sr
Modelo de absorção
n = número de partículas
I-I
I
a = área das partículas = pR2
r = densidade = n/vol= n/(As s)
s
na
( r ( s ) As s )a
I (t )   I (t )  
I (t )   r ( s )a s I (t )   ( s )s I (t ),
As
As
As
 ( s )  ar ( s )
t
   ( s ) ds
dI (t )
ta
  ( s) I (t )  c(t )  I (t ) e
ds
Aproximação para cor de um pixel
Ap
(i, j)
Imagem
ta
L
(x, y, z)
t
Um modelo simples para determinação
da cor de um pixel considera todos os
fótons emitidos em sua direção no volume
correspondente a seu raio, ou seja:
tb
Volume
tb
tb
ta
ta
Ci , j   c(t ) t Ap dt   I (t )e

t
ta  ( s ) ds
 t  Ap dt
Integrando numericamente
tb
Ci , j   I (t )  t e

t
ta  ( s ) ds
Ap dt
ta
Partição: t0 = ta  t1  t2  …  tn-1  tn = tb
n 1
Ci , j   I (t k ) Ap t k t k e
k 0

k 1
 (tm ) tm
m 0
Estimando a contribuição da partição
I (tk ) Ap tk tk  Cka k
Ck= intensidade refletida da fonte na partição k através da área Ap
ak =  (tk) tk = opacidade (núm. de partículas) da partição k
n 1
Ci , j   I (t k ) Ap t k t k e
k 0

k 1
 (tm ) tm
m 0
n 1
k 1
k 0
m 0
  Cka k  e a m
Shading pelo modelo de Phong
R
V
N
tk
L
t

Ck  Camb   Cluz kdif N  L   k s ( R  V ) n
luz

Shading por partículas Kajiya
V
tk
L
t
Ck   Cluz V , L    Cluz V  L 
luz
luz
Aproximando a absorção
n 1
k 1
k 0
m 0
Ci , j   C ka k  e a m
eam  1  a m 
n 1
k 1
k 0
m 0
a m2 a m3
2

6
   (1  a m )
Ci , j   Cka k  1  a m 
 C0a 0  C1a1 (1  a 0 )  C2a 2 (1  a 0 )(1  a1 )  
 Cn 1a n 1 (1  a 0 )  (1  a n 2 )
Compondo de trás para a frente
Ci , j  C0a 0  C1a1 (1  a 0 )  C2a 2 (1  a 0 )(1  a1 )    Cn1a n1 (1  a 0 )(1  a n2 )
 ((C0a 0  (1  a 0 )(C1a1  (1  a1 )(C2a 2    (1  a n2 )Cn1a n1 )))
cor = c[n-1]*alfa[n-1];
for (i=n-2; i>=0; i--) {
cor = c[i]*alfa[i]+(1-alfa[i])*cor;
}
Operadores do OpenGL
glBlendFunc(GL SRC ALPHA, GL ONE MINUS SRC ALPHA)
Cout  Csrca src  Cdst (1  a src )
cor = c[n-1]*alfa[n-1];
for (i=n-2; i>=0; i--) {
cor = c[i]*alfa[i]+cor*(1-alfa[i]);
}
Compondo de frente para atrás
Ci , j  C0a 0  C1a1 (1  a 0 )  C2a 2 (1  a 0 )(1  a1 )  
 Cn1a n1 (1  a 0 )(1  a n2 )
cor = c[0]*alfa[0];
trp = 1- alfa[0];
for (i=1; i<n && trp>tol; i++) {
cor += trp*(c[i]*alfa[i]);
trp *= (1-alfa[i]);
}
Cnew , anew
Cback , aback
Cfront , afront
Compondo fatias duas a duas
Cnew  C fronta front  Cbacka back (1  a front )
a new  a front  a back (1  a front )
cor = c[0]*alfa[0];
opc = alfa[0];
for (i=1; i<n && opc<1-tol; i++) {
cor += (1-opc)*(c[i]*alfa[i]);
opc += (1-opc)*alfa[i];
}
Problemas de re-amostragem
C
a
C
a
0
0
1
1
0
0
1
1
0
0
0
0.5
1
1.0
1
0
0
0
0.5
1
1.0
1
C=1
a1
C=0.75
a1

Re-amostragem da cor ponderada
a C0
0
1
1
0
0
1
1
0
0
0
0.5
1
1.0
1
0
0
0
0.5
1
1.0
1
a
aC
a
C=1
a1
C=1
a1

Splatting (Westover,1990)
 projeção
de voxel equivale a atirar uma
bola de neve num prato de vidro
–contribuição da neve decresce com distância ao
centro
–região de influência de voxel no espaço da
imagem é constante, a menos de uma translação
Splatting
 Algoritmo:
–transformação: ordem para percorrer o volume
–iluminação: modelo local
–reconstrução: determina a porção da imagem
influenciada por cada voxel
–determinação da visibilidade: compõe as
contribuições acumuladas em um buffer
Reconstrução do volume
nx
ny
t
nz
signal3 D  x, y, z    r (i, j , k ) hv (x  i, y  j , z  k )
i 1 j 1 k 1
(x,y,z)
x
contributioni, j , k ( x, y, z)  r (i, j, k )  h( x  i, y  j, z  k )

efeitoi , j ,k ( x, y )  r ( x, y, z )  h( x  i, y  j , t ) dt
(i,j,k)
z
y
0
1

footprint ( x, y)   h( x, y, t ) dt
0
0.8
0.6
0.4
0.2
0
Footprint
centro da projeção
footprint
Grade de espaçamentos diferentes
T-1
x
z
Sy
y
Sx
x
q
T-1(x)=x’
x’
Shear-Warp (Lacroute, 1995)
intermediate
image scanline
1. shear &
resample
3. warp &
resample
2. project
& composite
voxel
scanline
intermediate
image
final
image
Esquema
viewing rays
shear
shear
project
volume
slices
intermediate
image plane
warp
image plane
Mview = Mwarp2D * Mshear3D
Otimização
voxel scanline
Resampling and Composition
intermediate
image scanline
skip
work skip
work
skip
Voxel with a = 0
Opaque pixel
Voxel with a > 0
Non opaque pixel
Volume rendering with textures
3D texture volume rendering

Load the volume data into a 3D texture. This processing is done
once for a particular data volume.

Choose the number of slices. Usually this matches the texel
dimensions of the volume data cube.

Find the desired viewpoint and view direction

Compute a series of polygons that cut through the data
perpendicular to the direction of view. Use texture coordinate
generation to texture the slice properly with respect to the 3D
texture data.

Use the texture transform matrix to set the desired orientation the
textured images on the slices.

Render each slice as a textured polygon, from back to front.

As the viewpoint and direction of view changes, recompute the
data slice positions and update the texture transformation matrix
as necessary.
2D texture volume rendering







Generate the three sets of 2D textures from the volume data. Each set of 2D
textures is oriented perpendicular to one of volume’s major axes. This
processing is done once for a particular data volume.
Choose the number of slices. Usually this matches the texel dimensions of the
volume data cube.
Find the desired viewpoint and view direction.
Find the set of 2D textures most perpendicular to the direction of view.
Generate data slice polygons parallel to the 2D texture set chosen. Use texture
coordinate generation to texture each slice properly with respect to its corresponding 2D texture in the texture set.
Use the texture transform matrix to set the desired orientation the textured
images on the slices.
Render each slice as a textured polygon, from back to front. A blend operation
is performed at each slice; the type of blend depends on the desired effect.
As the viewpoint and direction of view changes, recompute the data slice
positions and update the texture transformation matrix as necessary. Always
orient the data slices to the 2D texture set that is most closely aligned with it.