Introdução à Visualização VolumétricaIntrodução à Visualização Volumétrica

porRoberto de Beauclair Seixas, IMPAAnselmo Cardoso de Paiva, UFMA

Marcelo Gattass, PUC-Rio

Visualização VolumétricaVisualizaçã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)Radiografia Padrão (Raios X)

2K x 2K x 2 Bytes4 chapasAngiografia: 40 chapas

Tomografia Computadorizada (TC)Tomografia Computadorizada (TC)

Receptores

Raio X

Exemplos de Imagens de TCExemplos de Imagens de TC

Normalmente 512 x 512 de 1/2 a 2 mm

Ar Gordura Água Rim Pâncreas Músculo Liver OssoEsponjoso

OssoCompacto

-1000 -110 15 0 27 15 35 10 40 10 55 10 200-400 >1000

Números em unidades de Hounsfield (HU)

Exame de 5 a 30 min.

Problemas da TCProblemas 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)Ressonância Magnética Nuclear (MR)

Grande campomagnético

Mesa comtrilhos

Principais Vantagens de MRPrincipais 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 MRDesvantagens 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 MRImagem MR

512 x 512 x 2 Bytes

30 fatias/exame

Ultra-somUltra-som

Imagem de Ultra-somImagem de Ultra-som

512 x 512 x 1 Byte

Medicina NuclearMedicina Nuclear

SPECT (Single Photon Emission Computed Tomography)

PET (Positron Emission Tomography)

256 x 256 x 2 Bytes

Câmera Gama bone scan galium scan

Exames MédicosExames Médicos

SísmicaSísmica

Formação de um traço sísmicoFormação de um traço sísmico

*

Seção SísmicaSeção Sísmica

Reconstrução das FatiasReconstrução das Fatias

Zs

Ys

Xs

volume dasfatias (slices)

Classificação do VoxelClassificação do Voxel

Voxel

densidade

Ar

Gordura

Músculo

Osso

op

acid

ade

Amarelo Vermelho Branco

1.0

Tipos de DadosTipos de Dados

Enumeráveis (material, litologia, ...)

Escalares (temperatura, pressão, ...)

Vetoriais (velocidade, aceleração, ...)

Tensoriais (tensão, deformação, ...)

Estrutura dos DadosEstrutura dos Dados

Grade Cartesiana (i, j, k)

Grade Regular (i*dx, j*dy,k*dz)

Grade Retilínea (x[i],y[j],z[k])

Grade Estruturada (x[i,j,k],y[i,j,k],z[i,j,k])

Grade Não Estruturada {(x[i],y[i],z[i]) ,

e=(v1, v2 , v3)}Grade Estruturada

por blocos

InterpolaçãoInterpolação

pontos dagrade

voxels

valores interpolados

pontos da grade

valor valor

x x

Matriz de voxels – análogo 3D do pixel

– (i, j, k)

– Matriz de células– interpolação trilinear

– imagens mais suaves

Métodos de VisualizaçãoMé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 DiretoRendering Direto

Mapeamento dos elementos de volume direto no espaço da imagem

Apropriados para a visualização de objetos amorfos

Mais lentos

Algoritmo Algoritmo ray-castingray-casting

Raios

Plano da

imagem

Volumede dados

Ordem da imagem. – para cada pixel

• lance um raio e encontre os voxels que são interceptados

– fim para

Lançamento dos RaiosLançamento dos Raios

volume dosraios

volume dosslices

plano devisualização

raio

x

y

z

Partição dos RaiosPartição dos Raios

x

y

z ta

tb

t=min(x,y,z)

amostra

voxel

Interpolação Interpolação

x,i

y,j

z,k

(xa, ya, za)

Dzzk

Dyyj

Dxxi

a

a

a

Dz

Dzzz

Dy

Dyyy

Dx

Dxxx

a

a

a

%

%

%

Interpolação no voxelInterpolação no voxel

1,1,1

1,1,

1,,1

1,,

,1,1

,1,

,,1

,,

1

1

11

1

11

11

111),,(

kji

kji

kji

kji

kji

kji

kji

kji

vzyx

vzyx

vzyx

vzyx

vzyx

vzyx

vzyx

vzyxzyxv

vi+1,j+1,k+1

vi,j,k

vi+1,j,k

vi,j+1,k

vi,j,k +1

vi +1,j,k +1

vi,j +1,k +1

vi+1,j+1,k

xy

zv

Classificação do VoxelClassificação do Voxel

Voxel

velocidadeou densidadeamplitude

tons e opacidade

1.0

cor

ar

Gordura

Músculo

Osso

0 255

Branco(255,255,255)

Vermelho (255,0,0)

Amarelo(255,255,0)

0 255

Diferentes Funções de TransferênciaDiferentes Funções de Transferência

Sensibilidade à função de transferênciaSensibilidade à função de transferência

0

0 .2

0 .4

0 .6

0 .8

1

0 50 100 150 2000

0 .2

0 .4

0 .6

0 .8

1

0 50 100 150 200

Iluminação de um voxelIluminação de um voxel

)( LNd

kluz

IakaIC

)(XVN N L

Gel A

Gel B

1a ordem

Estimativa do vetor normalEstimativa do vetor normal

,2/,,,,

,2/,,,,

,2/,,,,

,,

11

11

11

zzyxfzyxf

yzyxfzyxf

xzyxfzyxf

zyxf

kjikji

kjikji

kjikji

kji

zyxf

zyxfN

,,

,,

2a ordem

Interpolações em Visualização VolumétricaInterpolações em Visualização Volumétrica

normalinterpolada

normaldo voxel

cor interpolada

cor do voxel

Phong

Gouraud

C interpolada

C do voxel

Debrin

Etapa de ComposiçãoEtapa de Composição

Para cada raio:gera amostras de cor C(pi) e opacidades (pi)

» reamostragem dos dados dos voxels, em k amostras igualmente espaçadas

Processo de acumulaçãoI= t Ib +(1-t)I0

I0 = cor do objetoIb = cor do fundoI = cor resultantet = coeficiente de transparência

Influência de um VoxelInfluência de um Voxel

r = g = b =

r0

g0

b0

r = r0

g = g0b = b0

r1

g1

b1

r = r0+r1 g = g0+g1b = b0+b1+

Composição no raioComposição no raio

)1( ininout

CininCoutC )1(

outoutCC

Cin

Cout

ray

TripodTripod

(x0,y0,z0)

y

x

z

(x0+x,y0 +y,z0 +z)

exy= (y)x + (-x)y +cxy

e zy=

(y

)z +

(-z

)y +

c yz

exz= (z)x + (-x)z +cxz

exy exz ezy

+

-

+

-

+

-

Valor inicial e incremento Valor inicial e incremento

exy= (y)x + (-x)y +cxy

ezy= (y)z +(-z)y +cyz

exz= (z)x + (-x)z +cxz

(1/2,1/2,1/2)

(y-x)/2

(z-x)/2

(y-z)/2

x++

y

z

y++

x

z

z++

x

y

Caminhamento discretoCaminhamento discreto

Bresenham

Cohen

Efeito da amostragem Efeito da amostragem

Bresenham Tripod

Partição celularPartição celular

Plan

o da

imag

em

Volumede dados

Partição

Voxel uma amostra no meio da partição

marcador da partição

Partição na gradePartição na grade

tzi= ta

tyf=tb y

z

tyi

tzf

(ia ,ja ,0)

(ib ,Ny ,kb)

1)(

1)(

1)(

abz

aby

abx

kkn

jjn

iin

zzizfz

yyiyfy

xxixfx

nttt

nttt

nttt

)(

)(

)(

Partição celular: algoritmoPartição celular: algoritmo

Dados: txi,tyi,tzi, txf,tyf,tzf, nx,ny,nzdtx = 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ínimon=nx+ny+nz; while ( n > 0 )

tw += dtw; n--;t2 = min(tx,ty,tz) e w é o eixo do mínimoSample ((t1+t2)/2);t1=t2;

Otimizações 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 ProgresivoRefinamento Progresivo

subdivisão finalsegunda subdivisão

primeira subdivisãoamostragem inicial

pixels sendo visitados pixels já visitados

Refinamento Progressivo: ExemploRefinamento Progressivo: Exemplo

256x256128x12864x6432x32

16x168x84x42x2

Algoritmo Algoritmo Ray CastingRay 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 pixelRaio que calcula a cor de um pixel

t

Lta

tb

t

(x, y, z)

(i,j)Imagem

Volume

I = I(x,y,z,-t,) = radiância ou intensidade específicanum ponto (x,y,z) ou tna direção -t, na freqüência .

I

Um pouco de físicaUm pouco de física

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

I(x,n,) = c h (x,n,)

Modelo de absorçãoModelo de absorção

n = número de partículas

a = área das partículas = R2

= densidade = n/vol= n/(As s)

t

atdss

etItctIsds

tdI )(

)()()()()(

I

As s

I-

)()(

),()()()()())((

)()(

sas

tIsstIsastIA

asAstI

A

natI

s

s

s

Aproximação para cor de um pixelAproximação para cor de um pixel

b

a

t

at

b

a

t

t

p

dss

p

t

t

ji dtAtetIdtAttcC )(

, )()(

Lta

tb

t(x, y, z)

(i, j)

Imagem

Volume

Ap

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:

Integrando numericamenteIntegrando numericamente

Partição: t0 = ta t1 t2 … tn-1 tn = tb

1

0

)(

,

1

0)(n

k

tt

kkpkji

k

mmm

ettAtIC

b

a

t

at

t

t

p

dss

ji dtAettIC)(

, )(

Estimando a contribuição da partiçãoEstimando a contribuição da partição

kkkkpk CttAtI )(

1

0

1

0

1

0

)(

,

1

0)(n

k

n

k

k

mkk

tt

kkpkjim

k

mmm

eCettAtIC

Ck= intensidade refletida da fonte na partição k através da área Ap

k = (tk) tk= opacidade (núm. de partículas) da partição k

Shading Shading pelo modelo de Phongpelo modelo de Phong

luz

nsdifluzambk VRkLNkCCC )(

N

R

V

Ltk

t

Shading Shading por partículas Kajiyapor partículas Kajiya

luz

luzluz

luzk LVCLVCC ,

V

Ltk

t

Aproximando a absorçãoAproximando a absorção

1

0

1

0,

n

k

k

mkkji

meCC

)1()1(

)1)(1()1(

1

2011

102201100

1

0

1

0,

nnn

k

mm

n

kkkji

C

CCC

CC

)1(62

132

mmm

mme

Compondo de trás para a frenteCompondo de trás para a frente

))))1()(1()(1(((

)1()1()1)(1()1(

11222111000

2011102201100,

nnn

nnnji

CCCC

CCCCC

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 Operadores do OpenGL

)1( srcdstsrcsrcout CCC

glBlendFunc(GL SRC ALPHA, GL ONE MINUS SRC ALPHA)

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 Compondo de frente para atrás

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]);}

)1()1(

)1)(1()1(

2011

102201100,

nnn

ji

C

CCCC

Compondo fatias duas a duasCompondo fatias duas a duas

Cba

ck ,

back

Cfr

ont ,

fr

ont

Cn

ew ,

new

)1(

)1(

frontbackfrontnew

frontbackbackfrontfrontnew CCC

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-amostragemProblemas de re-amostragem

C

00 11

00 11

C

00 11

00 11

0

0 0.5

0.50

0

1.0

1.0

C=1

C=0.75

Re-amostragem da cor ponderadaRe-amostragem da cor ponderada

C

00 11

00 11

00 11

00 11

0

0 0.5

0.50

0

1.0

1.0

C=1

C=1

C

Splatting (Splatting (Westover,1990) 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

SplattingSplatting

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 volumeReconstrução do volume

),,(),,(,,1 1 1

3 kzjyixhkjizyxsignalx y zn

i

n

j

n

kvD

x

yz

(i,j,k)

(x,y,z)

0

,, ),,(),,(),( dttjyixhzyxyxefeito kji

0

),,(),( dttyxhyxfootprint

t

0

0.2

0.4

0.6

0.8

1

),,(),,(),,(,, kzjyixhkjizyxoncontributi kji

FootprintFootprint

footprint

centro da projeção

Grade de espaçamentos diferentesGrade de espaçamentos diferentes

T-1

T-1(x)=x’z

y

x

x’x

Sy

Sx

Shear-Warp (Lacroute, 1995)Shear-Warp (Lacroute, 1995)

1. shear &resample

voxelscanline

2. project& composite

intermediateimage scanline

3. warp &resample

intermediateimage

finalimage

EsquemaEsquema

viewing rays

volume slices

image plane

Mview = Mwarp2D * Mshear3D

shear

warp

viewing rays

image plane

shear

project

intermediate image plane

OtimizaçãoOtimização

voxel scanline

intermediateimage scanline

Resampling and Composition

Voxel with = 0

Voxel with > 0

Opaque pixel

Non opaque pixel

skip work skip work skip

Volume rendering with texturesVolume rendering with textures

3D texture volume rendering3D 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 rendering2D 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 corre-sponding 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.