Dissertação

dissertacao.allyson.cabral.pdf
Documento PDF (9.5MB)
                    Visualização da Curvatura de
Objetos Implícitos em um Sistema
Extensível

Allyson Ney Teodosio Cabral

Maceió
Fevereiro de 2010

Allyson Ney Teodosio Cabral

Visualização da Curvatura de
Objetos Implı́citos em um
Sistema Extensı́vel

Dissertação de Mestrado na área de concentração de
Computação Gráfica submetida em 11 de fevereiro
de 2010 à Banca Examinadora, designada pelo
Colegiado do Programa de Pós-Graduação em
Matemática da Universidade Federal de Alagoas,
como parte dos requisitos necessários à obtenção
do grau de mestre em Matemática.

Orientador:

Vinı́cius Mello

Maceió, Fevereiro de 2010

Catalogação na fonte
Universidade Federal de Alagoas
Biblioteca Central
Divisão de Tratamento Técnico
Bibliotecária Responsável: Maria Auxiliadora Gonçalves da Cunha
C117v

Cabral, Allyson Ney Teodosio.
Visualização da curvatura de objetos implícitos em um sistema extensível. /
Allyson Ney Teodosio Cabral, 2010.
51 f. : il., grafs.
Orientador: Vinícius Moreira Mello.
Dissertação (mestrado em Matemática) – Universidade Federal de
Alagoas. Instituto de Matemática, 2010.
Bibliografia: f. 50-51.

1. Curvatura. 2. Curvatura – Visualização volumétrica . 3. B-Spline – Funções .
4. GPU – Linguagem de computação. 5. GLSL – Linguagem de computação.
6. Interpolação tricúbica. I. Título.
CDU: 514.772.2

Resumo
Neste trabalho, estudaremos a visualização da curvatura de superfı́cies definidas implicitamente por funções do tipo f : [0, 1]3 → [0, 1], usando a técnica de lançamento de raios
(ray casting). Como em geral conhecemos apenas valores amostrados de f , estudaremos
um método de interpolação tricúbica, a fim de calcular as derivadas de segunda ordem
precisamente.
A implementação computacional deste trabalho foi desenvolvida na forma de módulos
do framework de visualização e processamento de imagens Voreen, o qual se beneficia do
poder de processamento das placas gráficas atuais para acelerar o processo de visualização.
Palavras-Chaves: Curvatura, Visualização Volumétrica, Interpolação Tricúbica, BSpline, GPU, GLSL, Objetos Implı́citos.

i

Abstract
In this work we study the curvature visualization problem on surfaces implicitly defined
by functions f : [0, 1]3 → [0, 1], using the ray casting technique. As we usually know only
sampled values of f , we study the tricubic interpolation method to compute second order
derivatives accurately.
This work’s implementation was designed as modules to the framework for volume rendering and image processing named Voreen, that uses the processing capability of graphics
cards to improve the rendering tasks.
Keywords: Curvature, Volume Rendering, Tricubic Interpolation, B-Spline, GPU,
GLSL, Implicit Objects.

ii

Agradecimentos
Agradeço a Deus, a razão de tudo em minha vida. Aos meus pais Cleide e Paulo por
ajudarem a tornar meus sonhos possı́veis e por estarem sempre ao meu lado em todos
os momentos. Ao meu orientador, Vı́nicius, pela paciência e dedicação durante o desenvolvimento desse trabalho. Ao professor Adán Corcho pelo carinho na condução da
Coordenação do programa de Pós-Graduação. Aos professores Adelailson Peixoto e Ralph
Teixeira pelas valiosas contribuições para esta dissertação. Aos amigos da PUC-Rio, professor Thomas Lewiner, Thales, Clarissa, Maria, Allan, Renata e Alex pela receptividade
e valiosas sugestões. Aos amigos de Laboratório, Douglas, Marcelo, Michel, Fábio Bóia,
Ailton, Ruff e Alisson pela importante companhia nos últimos dois anos. Aos amigos da
turma de mestrado, Alex, Ana, Adalgisa, Isnaldo, Kenerson, Fábio Henrique, Isadora,
Gregório, Rodrigo, Natália e Viviane. Um agradecimento especial aos amigos Douglas
Cedrim e Fabiane Queiroz pelas correções e ajuda na elaboração de imagens. À FAPEAL
pelo apoio financeiro.

iii

Sumário
1 Introdução
1.1 Objetivo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.2 Estrutura do Trabalho . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1
2
2

2 Visualização de Objetos Implı́citos
2.1 Objetos Implı́citos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2 Visualização de Objetos Implı́citos . . . . . . . . . . . . . . . . . . . . . .
2.3 Visualização Volumétrica . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4
4
6
9

3 Interpolação e Cálculo de Curvatura
3.1 Interpolação Cúbica . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.2 Implementação da Inversão . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.3 Derivadas das Funções Spline . . . . . . . . . . . . . . . . . . . . . . . . .
3.4 Avaliação Rápida . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.5 Caso Tridimensional . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.6 Curvatura . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

13
14
16
18
19
20
21

4 Visualização de Objetos Implı́citos Baseada em GPU
4.1 GPU . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2 OpenGL Shading Language . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2.1 Vertex Processor . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2.2 Geometry Processor . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2.3 Fragment Processor . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.3 GPGPU e CUDA . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.4 Voreen: Modularidade . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.4.1 Classes de Processadores . . . . . . . . . . . . . . . . . . . . . . . .
4.4.2 Criando Processadores . . . . . . . . . . . . . . . . . . . . . . . . .

25
25
27
28
29
30
32
32
34
35

5 Módulos
5.1 Módulo de Funções . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5.2 Módulo de Interpolação Tricúbica . . . . . . . . . . . . . . . . . . . . . . .
5.3 Módulos de Cálculos em GLSL (interpolação tricúbica e curvatura) . . . .
5.4 Módulo de Renderização de Isosuperfı́cies . . . . . . . . . . . . . . . . . . .
5.5 Módulo de Renderização de Volumes . . . . . . . . . . . . . . . . . . . . .

38
38
39
40
43
44

6 Resultados

46

7 Considerações Finais

50

iv

Lista de Figuras
1.1

2.1
2.2

Resultados obtidos através do módulo de visualização de curvatura implementado. Nas imagens, azul indica curvaturas negativas, verde curvaturas
próximas de zero e vermelho curvaturas positivas. . . . . . . . . . . . . . .

3

2.7
2.8

Visualização dos tecidos ósseo e cutâneo (obtida de [7]). . . . . . . . . . . . 5
Descrição do algoritmo de traçado de raios. Os quadrados vermelhos representam os pixels cujos raios intersectam o volume. . . . . . . . . . . . . 7
Componentes do modelo de iluminação . . . . . . . . . . . . . . . . . . . . 8
Resultado do algoritmo de ray casting . . . . . . . . . . . . . . . . . . . . 9
Efeito de transparência aplicado ao volume possibilitando a visualização de
diferentes isosuperfı́cies. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
Descrição do algoritmo de visualização volumétrica. Para cada pixel um
raio atravessa o volume acumulando os valores de cor e opacidade encontrados. 10
Modelo de absorção. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11
Algoritmo de visualização aplicado a um volume. . . . . . . . . . . . . . . 12

3.1
3.2
3.3
3.4
3.5
3.6

Funções B-spline. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15
Gráfico do sinal h[k]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
Interpolação Cúbica. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
Derivadas da função spline cúbica. . . . . . . . . . . . . . . . . . . . . . . 18
Funções B-spline bidimensionais e suas derivadas parciais de primeira ordem. 21
Funções B-spline bidimensionais e suas derivadas parciais de segunda ordem. 22

2.3
2.4
2.5
2.6

4.1

Gráficos de comparação entre peformances da CPU e da GPU. O primeiro
gráfico mostra como a GPU tem obtido um crescimento muito maior da sua
capacidade de processamento em relação à CPU enquanto o segundo traz
uma comparação do crescimento na capacidade de transferência de dados
em GB/s (dados obtidos de [17]). . . . . . . . . . . . . . . . . . . . . . . .
4.2 Diferença de arquitetura entre as GPUs e as CPUs, em verde as unidades
lógicas e aritméticas que são responsáveis pelo processamento dos dados e
em laranja as unidades de memória. . . . . . . . . . . . . . . . . . . . . . .
4.3 Pipeline de renderização. . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.4 Fluxo de dados no pipeline de renderização . . . . . . . . . . . . . . . . . .
4.5 Etapa do pipeline de renderização relacionada ao vertex shader . . . . . . .
4.6 Etapa do pipeline de renderização relacionada ao geometry shader . . . . .
4.7 Etapa do pipeline de renderização relacionada ao fragment shader . . . . .
4.8 Resultado da execução dos shaders . . . . . . . . . . . . . . . . . . . . . .
4.9 Visão Geral do Voreen . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.10 Rede de processadores. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
v

26

27
28
28
29
30
31
31
33
33

LISTA DE FIGURAS

vi

4.11 Texturas com os pontos de entrada e saı́da do volume em forma de cubo. . 34
4.12 Propriedade camera compartilhadas entre componentes. . . . . . . . . . . . 35
5.1
5.2
5.3
5.4
5.5
5.6
5.7
5.8

6.1

6.2
6.3
6.4
6.5
6.6

FunctionSelector . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Imagens Geradas com o processador FunctionSelector . . . . . . . . . .
CubicSplineProcessor . . . . . . . . . . . . . . . . . . . . . . . . . . . .
CurvatureRaycaster . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Resultado obtido com o CurvatureRaycaster e suas propriedades . . . . .
Editor para funções de transferência . . . . . . . . . . . . . . . . . . . . . .
CurvatureRaycaster . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Resultado obtido com o CurvatureVolumeRaycaster mostrando diversas
porções do volume com efeito de transparência e levando em consideração
os valores de curvatura para o cálculo das cores. . . . . . . . . . . . . . . .
Resultados obtidos utilizando o mapa de cores JET para visualização de
valores de curvatura. Vermelho indica curvaturas positivas, azul curvaturas
negativas, verde valores próximos de zero, e em preto está delimitada a
curva de nı́vel zero. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Resultados obtidos utilizando o volume smile. . . . . . . . . . . . . . . . .
Resultados obtidos utilizando o volume cow. . . . . . . . . . . . . . . . . .
Erros que ocorrem utilizando diferentes formas de representação de valores.
Visualização da curvatura Gaussiana em volumes diferentes. . . . . . . . .
Diferentes nı́veis do volume cubo. . . . . . . . . . . . . . . . . . . . . . . .

38
39
39
43
44
44
45

45

46
47
48
48
49
49

Capı́tulo 1
Introdução
Os métodos de visualização de volumes consistem em um conjunto de técnicas cujo objetivo é a geração de imagens a partir de dados trimensionais. As abordagens tradicionais
[12] realizam a extração de malhas poligonais a fim de fazer uso delas no processo de
visualização. Tais métodos são denominados métodos indiretos. Nesses casos, a malha
construı́da pode ser utilizada para realizar processamentos posteriores. Entretanto, esses métodos possibilitam apenas a visualização de uma parte do volume por vez, mais
especificamente a sua superfı́cie de bordo. Essa limitação esconde informações que estão
no interior do volume e que são de grande importância para, por exemplo, a geração de
imagens médicas. Como alternativa para solução deste problema, foram desenvolvidos
métodos que eliminam a etapa de geração de malhas e que tornaram possı́vel visualizar
em uma mesma imagem porções diferentes do volume. Essas técnicas, conhecidas como
métodos diretos, em geral, utilizam algoritmos de lançamento de raios que percorrem o
volume de dados em um processo de integração de valores que resulta nas cores finais dos
pixels da imagem. Exemplos clássicos de algoritmos que utilizam o lançamento de raios
são encontrados em [5], [10] e [3].
O conjunto de dados de entrada para algoritmos de renderização de volumes consiste
em um conjunto de superfı́cies de nı́vel descritas por uma aplicação f : [0, 1]3 → [0, 1].
Esses algoritmos realizam uma série de operações sobre os pontos da superfı́cie, como o
cálculo do gradiente, para as quais algumas condições de regularidades são necessárias. Na
maioria dos casos é suficiente que f seja continuamente diferenciável, entretanto existem
aplicações que necessitam do cálculo de derivadas de segunda ordem, a exemplo das aplicações que utilizam o cálculo de curvatura [9] [8]. Para esses casos podem ser utilizadas
técnicas de interpolção tricúbica que garantem a exatidão das derivadas até a segunda
ordem [22].
Uma desvantagem das técnicas de visualização direta é que para cada imagem do
volume gerada todo o conjunto de dados precisa ser percorrido. Entretanto, como esses
algoritmos são divididos em etapas independentes, vários trabalhos passaram a utilizar a

1

1.1. OBJETIVO

2

programação paralela em placas gráficas para acelerar a obtenção dos resultados. Dois
exemplos desses trabalhos são [25] e [7]. Este último foi um dos muitos trabalhos recentes
resultantes de um projeto que deu origem a uma ferramenta de visualização de volumes
chamada Voreen que utiliza a GPU para realizar a visualização em tempo real em um
sistema modular e interativo que permite a associação de várias técnicas com esse intuito.

1.1

Objetivo

O objetivo desse trabalho é fazer um estudo sobre os problemas envolvidos na visualização
de volumes utilizando informações de curvatura e foi baseado no trabalho desenvolvido
por Kindlmann et al. em [9]. Kindlmann utilizou os passos descritos por Möller et. al [16]
para cálculo de curvaturas em superfı́cies para a aplicação na visualização de volumes.
A idéia inicial para o trabalho foi otimizar esse cálculo implementando o trabalho de
Kindlmann para a execução em placas gráficas. Durante as pesquisas nos deparamos com
o framework Voreen que, além de ser uma ferramenta muito poderosa para a geração de
imagens de volumes, possui código fonte aberto estando disponı́vel para colaborações.
Como contribuição foram desenvolvidos módulos até então inexistentes no Voreen. O
primeiro deles permite a visualização de volumes levando em consideração informações
de curvatura que utiliza os passos descritos por Kindlmann. O segundo módulo permite
a avaliação de funções amostradas através de interpolação tricúbica utilizando o método
descrito em [22]. Em ambos foi utilizada a alta capacidade dos múltiplos processadores
existentes nas placas gráficas para paralelizá-los e acelerar os cálculos envolvidos. Além
desses, foram implementados módulos que realizam a visualização de volumes utilizando
o algoritmo de ray casting1 . Estes componentes foram incorporados ao Voreen e estão
descritos no capı́tulo 5.
A imagem abaixo mostra alguns resultados obtidos através dos módulos implementados. Nelas aparece a visualização das curvaturas Gaussiana (K) e média (H) bem como
das curvaturas principais (k1 e k2 ) sobre um volume utilizando o método de interpolação
tricúbica.

1.2

Estrutura do Trabalho

O próximo capı́tulo deste trabalho (capı́tulo 2) apresentará métodos para a visualização
direta de objetos implı́citos que eliminam a necessidade da extração de malhas. No capı́tulo seguinte descreveremos a teoria necessária das técnicas de interpolação e do cálculo
de curvaturas (capı́tulo 3). No capı́tulo 4 é realizada uma descrição de duas ferramentas
1

Deixemos claro que existem diferenças entre os termos ray casting e ray tracing utilizados durante
muito tempo na literatura como sinônimos. Atualmente, faz-se uma distinção clara entre os dois métodos
já que o primeiro nunca utiliza o lançamento recursivo de raios secundários.

1.2. ESTRUTURA DO TRABALHO

(a) K

(b) H

3

(c) k1

(d) k2

Figura 1.1: Resultados obtidos através do módulo de visualização de curvatura implementado. Nas imagens, azul indica curvaturas negativas, verde curvaturas próximas de zero
e vermelho curvaturas positivas.
computacionais utilizadas para o desenvolvimento deste trabalho, a linguagem para programação em placas gráficas GLSL e o framework Voreen deixando para o capı́tulo 5 a
tarefa de descrever os módulos desenvolvidos para o framework. Por fim apresentaremos
os resultados, as conclusões e algumas propostas para trabalhos posteriores (capı́tulos 6
e 7).

Capı́tulo 2
Visualização de Objetos Implı́citos
O processo de modelagem tradicional fornece como resultado uma malha de polı́gonos,
independente do método utilizado para sua geração (superfı́cies de subdivisão, NURBS
etc.). Esse fato é tão comum que os hardwares gráficos (GPUs) produzidos atualmente
estão otimizados para a visualização de malhas de polı́gonos. Nestes casos, objetos sólidos
são representados por sua superfı́cie de bordo (B-Rep), não sendo possı́vel, em um mesmo
passo de renderização1 , visualizar seu interior.
Existem abordagens que utilizam a representação implı́cita de objetos, que consiste
em caracterizá-los como o conjunto solução de uma função. Para visualizar esses objetos,
alguns métodos utilizam a extração de malhas, tais como o marching cubes [12], que é
um algoritmo utilizado para poligonização de superfı́cies definidas implicitamente. Entretanto, há métodos que não necessitam da extração de malhas, além de oferecerem a
possibilidade de visualizar o interior do objeto, o que vem a ser importante em muitas
aplicações, como na geração de imagens médicas, por exemplo. Aliado a isso, a constante evolução e facilidade de programação das placas gráficas permitem que elas sejam
utilizadas também nesses métodos.
Neste capı́tulo, faremos uma breve descrição do método de visualização de objetos
implı́citos baseada em lançamento de raios (ray casting) [5]. Inicialmente, nos preocuparemos apenas em fazer uma descrição conceitual do problema, deixando os aspectos de
implementação para o capı́tulo 4.

2.1

Objetos Implı́citos

Um conceito fundamental para este trabalho é o de objeto implı́cito. De maneira simples,
um objeto implı́cito pode ser entendido como o conjunto solução de uma determinada
função. Mais precisamente:
1

Em diversos momentos utilizaremos a palavra renderização, que é um estrangerismo da palavra
proveniente da lı́ngua inglesa rendering, para indicar o processo de visualização.

4

2.1. OBJETOS IMPLÍCITOS

5

Definição 2.1 Um subconjunto O ⊂ Rn é chamado objeto implı́cito se existe uma função
f : Rn → R e um número real c tal que O = f −1 (c).
A esfera é um exemplo de objeto implı́cito, pois se considerarmos a função g : R3 → R
definida por g(x, y, z) = x2 + y 2 + z 2 e um r ∈ R, o objeto implı́cito S(r) = f −1 (r)
representa a esfera de raio r.
Uma ilustração prática desses conjuntos são os dados obtidos através do processo de
tomografia computadorizada, onde são obtidos valores que indicam, para cada ponto do
espaço, o ı́ndice de absorção de raios X do tecido do corpo humano presente neste ponto.
Processados esses valores, pode-se obter imagens como a figura 2.1. Se considerarmos que
c1 e c2 indicam, respectivamente, os ı́ndices de absorção dos tecidos ósseo e cutâneo, a
imagem da esquerda é uma representação do conjunto f −1 (c1 ) e a da direita do conjunto
f −1 (c2 ).
Como qualquer conjunto do Rn pode ser descrito como um objeto implı́cito 2 , a definição acima torna-se muito geral. Na maioria dos casos essa generalidade é desnecessária e,
além disso, vários problemas podem ser simplificados se condições de regularidade forem
impostas. Nesta seção, enunciaremos alguns resultados com o objetivo de caracterizar
os objetos implı́citos do ponto de vista da geometria diferencial. Para mais detalhes,
consultar [14].
Definição 2.2 Um objeto implı́cito O = f −1 (c) é regular se f é diferenciável e ∇f (p) 6=
0, para todo p ∈ O, onde ∇f (p) indica o gradiente de f no ponto p.
Essa definição pode ser estendida naturalmente impondo-se que f seja de classe C k ,
com k ≥ 1. Em muitas aplicações, é necessário que f seja ao menos de classe C 2 , como
2

Para definir um conjunto A ⊂ Rn como um objeto implı́citito basta tormar a função caracterı́stica
IA : Rn → {0, 1} que tem valor 1 em todos os pontos de A e 0 nos pontos x ∈ Rn − A e fazer A = f −1 (1).

(a) f −1 (c1 )

(b) f −1 (c2 )

Figura 2.1: Visualização dos tecidos ósseo e cutâneo (obtida de [7]).

2.2. VISUALIZAÇÃO DE OBJETOS IMPLÍCITOS

6

veremos adiante. Um número real c é denominado valor regular de f se O = f −1 (c) é um
objeto implı́cito regular.
Dois resultados teóricos são muito importantes no estudo dos objetos implı́citos regulares: o teorema da função implı́cita e o teorema de Sard. O primeiro pode ser utilizado para
mostrar que todo objeto implı́cito regular O = f −1 (c) é localmente o gráfico de uma função de n variáveis sobre a qual pode-se definir um espaço vetorial tangente Tp O = ∇f (p)⊥ .
O segundo assegura que “quase todo”3 c ∈ R é valor regular de f . Na prática, considerando um intervalo [a, b] representado no computador e a função f , se existe algum valor
c ∈ [a, b] tal que f −1 (c) não é regular, pequenas perturbações no valor de c garantem essa
regularidade. Esse fato torna aceitável a suposição, feita pela maioria dos algoritmos de
visualização, de que a superfı́cie a ser visualizada é de fato regular. Descrições formais
dos dois teoremas podem ser encontrados em [11].
Aproveitando a notação definida acima e os resultados mencionados, já que é possı́vel
definir um plano tangente em p, podemos calcular o vetor normal nesse ponto da superfı́cie
da seguinte forma
∇f
(p).
N (p) =
k∇f k

2.2

Visualização de Objetos Implı́citos

Diversos métodos podem ser utilizados para a visualização de objetos implı́citos, dentre
eles estão algoritmos clássicos, como marching cubes [12] que, como foi dito anteriormente,
é um método de poligonização, e o ray casting [5] que utiliza o lançamento de raios para
gerar a visualização do objeto, eliminando a etapa de extração de malhas. Este último
será descrito brevemente mais adiante, para o caso particular de objetos implı́citos.
Um fato comum nestes métodos é que a visualização de um objeto implı́cito depende
de uma função de atributo a : O → A que associa a cada ponto p ∈ O um elemento em um
conjunto de atributos A. O conteúdo desse elemento depende da aplicação em questão,
mas tipicamente é consituı́do pela cor, pelo vetor normal e por propriedades fı́sicas da
superfı́cie no ponto p.
Por exemplo, se uma certa aplicação necessita dos atributos vetor normal e cor em cada
ponto p ∈ O, então a(p) = (N (p), C(p)). Caso o atributo seja constante em Oc = f −1 (c)
podemos usar a notação a(c).
O processo de visualização consiste em dar um sentido visual a essas informações,
unindo a informação geométrica em O e as propriedades do objeto descritas pela função
a. Para tanto, é utilizado um algoritmo de iluminação que simula a interação desses dados
com a luz presente na cena, gerando como resultado uma imagem.
3
“Quase todo” neste caso significa que o conjunto dos valores que não são regulares (valores crı́ticos)
tem medida nula.

2.2. VISUALIZAÇÃO DE OBJETOS IMPLÍCITOS

7

Consideradas as definições acima, podemos agora descrever o algoritmo de ray casting
para objetos implı́citos ilustrado na figura 2.2.

Figura 2.2: Descrição do algoritmo de traçado de raios. Os quadrados vermelhos representam os pixels cujos raios intersectam o volume.
Dado objeto O = f −1 (c), o primeiro passo consiste em considerarmos um grid regular,
que representa a imagem resultante do algoritmo, onde cada elemento desse grid é chamado
de pixel. Para cada pixel Pi,j é traçado um raio r(t) que parte da origem o e passa pelo
ponto Pi,j , ou seja, r(t) = Pi,j + t(Pi,j − o). Conhecendo o raio e a equação da superfı́cie,
podemos calcular a interseção entre o raio r e o objeto O fazendo f (r(t)) = c. Note que
chegamos assim a uma equação na variável t. A primeira interseção de cada um dos raios
com o objeto, que corresponde ao menor valor de t que satisfaz a equação acima, nos
fornecerá pontos pi,j sobre os quais podemos aplicar a função a para obtermos sua cor
original, C(i, j), e o vetor normal, N (i, j).
A etapa final do algoritmo consiste em aplicarmos um modelo de iluminação sobre
os pontos utilizando as informações adquiridas nos passos anteriores. O resultado desse
algoritmo é a cor final de cada pixel. Como ilustração, para este exemplo, utilizaremos
modelo de iluminação de Phong-Blinn [20] que descreveremos a seguir.
A imagem acima ilustra as componentes do modelo de iluminação em questão. Esse
modelo procura, de forma heurı́stica, descrever o processo de reflexão da luz sobre a
superfı́cie. Sendo assim, o primeiro elemento a ser considerado é a reflexão difusa que
é provocada pela absorção e irradiação uniformemente distribuı́da da luz sobre o objeto
iluminado. Esse efeito é modelado pela Lei de Lambert, descrita pela equação
I = kd (N · L),

(2.1)

onde I é a intensidade refletida, kd é o coeficiente de reflexão difusa determinado pelo
material que constitui o objeto, N é o vetor normal à superfı́cie no ponto em questão e L
é o vetor que representa a direção da fonte de luz.
Aliada à componente do parágrafo anterior, para simplificar o modelo, a contribuição

2.2. VISUALIZAÇÃO DE OBJETOS IMPLÍCITOS

8

Figura 2.3: Componentes do modelo de iluminação
da energia luminosa provocada pela reflexão da luz por outros objetos contidos na cena
é considerada em um novo elemento chamado de luz ambiente. Dessa maneira, o modelo
pode ser novamente descrito como
I = Ia ka + kd (N · L),

(2.2)

onde Ia representa a intensidade da luz ambiente e ka é a constante de reflexão ambiente.
Finalmente, podemos considerar a componente especular que simula propriedades de
reflexão da superfı́cie. Essa componente é responsável pelo efeito de brilho intenso em uma
determinada região da superfı́cie conhecido como highlight. Chegamos então à seguinte
expressão:
I = Ia ka + kd (N · L) + ks (N · H)α ,
(2.3)
onde ks é o coeficiente de reflexão especular, α é o expoente de reflexão especular e H é
o vetor que indica a direção de highlight máximo [1], este vetor está situado entre a fonte
de luz e o observador a uma distância angular definida como a metade da existente entre
os dois vetores.
Encontrada a intensidade de luz refletida Ii,j em um determinado ponto pi,j , é calculada
então a cor Ci,j do pixel Pi,j da seguinte maneira:
Ci,j = Ii,j C(i, j).
A figura 2.4 mostra um resultado do ray casting de um objeto utilizando o modelo de
iluminação descrito.

2.3. VISUALIZAÇÃO VOLUMÉTRICA

9

Figura 2.4: Resultado do algoritmo de ray casting

2.3

Visualização Volumétrica

Como já foi mencionado, em muitas aplicações é necessário visualizar o interior de um
sólido, mas a representação por objetos implı́citos nos permite visualizar apenas a superfı́cie deles. Veremos nesta seção como estender o algoritmo descrito na seção anterior para
visualização do interior de um sólido descrito por uma função escalar f : [0, 1]3 → [0, 1]. A
idéia é semelhante a descrita na introdução do capı́tulo. A função f classifica os pontos do
cubo [0, 1]3 atribuindo um valor de densidade a cada um deles. O processo de visualização
que utiliza esses valores de densidade é chamado de Visualização Volumétrica.
Uma maneira de interpretar a visualização volumétrica é considerar que
[0, 1]3 = ∪c∈[0,1] f −1 (c),
ou seja, o cubo [0, 1]3 é a união de todos as superfı́cies de nı́vel da função f . Assim,
podemos aplicar as técnicas de visualização de objetos implı́citos para cada um dos nı́veis
c ∈ [0, 1], integrando-os de alguma forma. É possı́vel combinar os atributos contidos em
cada um desses objetos de maneira a visualizar partes diferentes do volume simultaneamente através de um efeito de transparência ilustrado na figura 2.5.

Figura 2.5: Efeito de transparência aplicado ao volume possibilitando a visualização de
diferentes isosuperfı́cies.
O segredo para realizar essa integração é considerar que existe um atributo ω, chamado

2.3. VISUALIZAÇÃO VOLUMÉTRICA

10

opacidade, associado a cada ponto p ∈ [0, 1]3 . Este atributo pode ser constante em cada
superfı́cie Oc = f −1 (c).
Durante o ray casting, um mesmo raio pode intersectar diferentes superfı́cies de nı́vel,
com diferentes valores de opacidade e cor. No algoritmo descrito anteriormente, considerávamos apenas a primeira interseção. O problema então é: como compor esses valores
para gerar o resultado final?

Figura 2.6: Descrição do algoritmo de visualização volumétrica. Para cada pixel um raio
atravessa o volume acumulando os valores de cor e opacidade encontrados.
Tomemos como base a ilustração anterior. Como mencionado no inı́cio da seção, uma
série de raios serão traçados através do volume. O primeiro cálculo a ser feito é o da
quantidade de energia luminosa que atravessa o volume na direção de cada um desses
raios e chega ao observador. Assim, dado um pixel Pi,j , sua intensidade de cor é dada pela
seguinte equação:
Z
tA

Ii,j =

e

R
− tt τ (s)ds
A

Iv (t)dt.

(2.4)

tB

Nesta equação, Iv indica a intensidade luminosa em cada um dos pontos considerados e
tA e tB são valores escalares que servem para fixar a posição dos pontos de entrada e saı́da
r(tA ) e r(tB ) do raio r no volume.
Para entender a equação acima, considere que a energia luminosa Iv esteja atravessando
um volume cilı́ndrico de área At e a altura ∆s onde existem n partı́culas de área Ap que
impedem a passagem de luz, conforme ilustra a figura 2.7.
Supondo-se que a espessura ∆s seja suficientemente pequena e que as partı́culas não
obscureçam umas às outras, a redução de energia luminosa pode ser calculada por:
∆Iv = −

nAp
(ρ(s)At ∆s)Ap
Iv = −
Iv = −ρ(s)Ap Iv ∆s = −τ (s)Iv ∆s,
At
At

onde ρ(s) indica a densidade de partı́culas opacas e τ (s) = ρ(s)Ap o coeficiente de oclusão

2.3. VISUALIZAÇÃO VOLUMÉTRICA

11

Figura 2.7: Modelo de absorção.
por unidade de comprimento ao longo do raio.
No limite em que ∆s tende para zero, tem-se
dIv
= −τ (s)Iv ,
ds
que integrada resulta em
It = Iv (t)e

R
− tt τ (s)ds
0

.

(2.5)

Integrando todas as contribuições It ao longo do raio r de tA a tB chega-se a equação
de rendering 2.4.
A integral da equação de rendering de volumes pode ser avaliada usando-se a soma de
Riemann sobre uma partição t0 = tA < t1 < t2 < ... < tn = tB , de onde obtém-se
Ii,j =

n−1
X

e−

Pk−1

l=0 τl ∆tl

Ik ∆tk =

k=0

n−1
X

Ik ∆tk ·

k=0

com


Ik ≡ I

tk + tk+1
2



!
e−τl ∆tl

,

(2.6)

l=0


e τl ≡ τ

k−1
Y

tl + tl+1
2


.

Definindo ωl ≡ 1 − e−τl , de onde τl = − ln(1 − ωl ), como a opacidade da amostra l;
Ck ≡ (Ik /ωk ) ∆tk , a cor da amostra k e ck ≡ Ck ωk , a multiplicação dos dois, obtém-se:
Ii,j ∼
=

n−1
X
k=0

ck ·

k−1
Y

!

(1 − ωk ) .

(2.7)

l=0

Antes de encerrar o capı́tulo, é importante mencionar que existe um tipo de função de
atributo especial que associa para cada ponto p um valor de cor e opacidade e que possibilita a alteração desses valores para efeito de visualização, essas funções são conhecidas
na literatura como funções de transferência.
Para mais informações sobre os aspéctos fı́sicos envolvidos no processo de visualização
de volumes e detalhes sobre as fórmulas que acabamos de descrever, consultar [19].

2.3. VISUALIZAÇÃO VOLUMÉTRICA

Figura 2.8: Algoritmo de visualização aplicado a um volume.

12

Capı́tulo 3
Interpolação e Cálculo de Curvatura
Como foi mencionado no capı́tulo introdutório, o objetivo deste trabalho é o estudo dos
problemas envolvendo o cálculo de curvatura em superfı́ceis representadas implicitamente
e a aplicação dos seus resultados no processo de visualização de volumes. Entretanto,
para que isso seja possı́vel, alguns aspectos devem ser garantidos.
Em muitas aplicações, desejamos obter uma função contı́nua f: R → R a partir de um
sinal discreto g[k], de tal modo que f (k) = g[k] para todo k ∈ Z. Esse problema é denominado problema de interpolação de sinais discretos. Nestes casos, podemos empregar
interpolação linear, de modo que
f (x) = (1 − t)g[i] + tg[i + 1],
onde i = bxc denota a parte inteira e t = {x} a parte fracionária de x. Interpolação linear
é uma solução muito eficiente para o problema de interpolação. Entretanto, em alguns
problemas, como no cálculo de curvaturas, mais uma condição é necessária, a de que f
seja de classe C 2 . Neste caso, deve-se recorrer a outro método de interpolação, pois em
geral as funções obtidas por interpolação linear não são diferenciáveis nos pontos k ∈ Z.
Uma abordagem geral para o problema da interpolação é considerar que a função
f pertence a um certo espaço de funções que possua as propriedades desejadas. Em
aplicações de Computação Gráfica, os espaços mais usados são os espaços das funções
splines de grau n,
S n = {f ∈ C n−1 (R)|f |[k,k+1] é um polinômio de grau n},
ou seja, S n é o espaço das funções de classe C n−1 que são polinomiais por partes nos
intervalos [k, k + 1] da reta. Note que para cada sinal discreto g[k] existe uma única
função em S 1 que satisfaz à condição de interpolação, exatamente a função obtida por
interpolação linear.
Nos modelos tradicionais de visualização espera-se apenas que a superfı́cie seja regular,
13

3.1. INTERPOLAÇÃO CÚBICA

14

em outras palavras, que f definida em R3 seja de classe C 1 , dessa forma podem ser
definidos planos tangentes em cada um dos seus pontos. Esse fato é muito importante
pois garante que vetores normais possam ser calculados sobre esses pontos para serem
utilizados, por exemplo, no processo de iluminação da cena. Para este trabalho mais uma
condição será imposta, a de que f seja de classe C 2 , ou seja, que suas derivadas até a
segunda ordem sejam contı́nuas. Essa garantia é fundamental para o cálculo da curvatura
sobre pontos em superficies descritas por f .
Entretanto, para que f seja de classe C 2 , devemos buscá-la no espaço S 3 , conhecido
como espaço das splines cúbicas. Neste capı́tulo, mostraremos como resolver o problema
de interpolação para splines cúbicas, ou seja, dado um sinal g[k], determinar a única
função f ∈ S 3 tal que f (k) = g[k], para todo k ∈ Z. Além disso, será apresentado um
método para obtenção de valores de curvatura sobre superfı́cies.

3.1

Interpolação Cúbica

Para uma boa compreensão do problema de interpolação cública, é necessária uma prévia
caracterização dos espaços S n . Analisando o método de interpolação linear, é fácil ver
que se f ∈ S 1 então
X
f (x) =
c[k]B1 (x − k),
(3.1)
k∈Z

onde B1 (x) é a função B-spline de grau 1 exibida na figura 3.1(a) e c[k] = g[k]. O prefixo
B antes da palavra spline significa que as translações inteiras da função B1 formam uma
base do espaço vetorial S 1 . Note que como B1 possui suporte compacto, para cada x o
somatório em (3.1) é finito.
É possı́vel mostrar, que para cada espaço S n , com n > 1, existe uma função Bn ,
chamada função B-spline de grau n, tal que o conjunto {Bn (· − k)}k∈Z é uma base do
espaço S n , ou seja, cada função f ∈ S n pode ser escrita da forma
f (x) =

X

c[k]Bn (x − k),

(3.2)

k∈Z

para coeficientes c[k] unicamente determinados por f .
Existem várias fórmulas recursivas para as funções Bn . É possı́vel mostrar, por exemplo, que Bn+1 = Bn ? B0 , onde B0 é a função caracterı́stica do intervalo (− 21 , 12 ], conhecida
como função box, e a operação f ? g é o produto de convolução entre duas funções,
Z ∞
f (t)g(x − t)dt.

f ? g(x) =
−∞

O gráfico das funções Bn para n = 0, · · · , 3 pode ser visto na figura 3.1.

3.1. INTERPOLAÇÃO CÚBICA

15

(a) B0 (x)

(b) B1 (x)

(c) B2 (x)

(d) B3 (x)

Figura 3.1: Funções B-spline.
Aplicando a definição recursiva, concluı́mos que

B3 (x) =



 0



1
· (2 − |x|)3
6
2
− 12 |x|2 · (2 − |x|)
3

, se |x| ≥ 2
, se 1 ≤ |x| < 2
, se |x| < 1

(3.3)

Como o suporte da função B3 tem largura 4, em geral cada ponto x da função f (x) sofre a
influência de quatro funções de base B3 (·−k), mais exatamente para k = i−1, i, i+1, i+2
onde i = bxc. Mas se x é inteiro, apenas três funções de base contribuem para o valor de
f (x), aquelas com k = i − 1, i, i + 1, onde i = x. Portanto, temos que para todo k ∈ Z
f (k) = c[k − 1]B3 (k − (k − 1)) + c[k]B3 (k − k) + c[k + 1](k − (k + 1))
= c[k − 1]B3 (1) + c[k]B3 (0) + c[k + 1]B3 (−1)
4
1
1
=
c[k − 1] + c[k] + c[k + 1].
6
6
6
Voltando ao problema de interpolação, o problema de encontrar a função f ∈ S 3 que
interpola o sinal discreto g[k] resume-se a determinar coeficientes c[k] tais que
4
1
1
g[k] = c[k − 1] + c[k] + c[k + 1]
6
6
6

(3.4)

para todo k ∈ Z. Note que os coeficientes c[k] podem ser interpretados como um sinal
discreto de modo que a equação (3.4) pode ser colocada na forma de uma convolução

3.2. IMPLEMENTAÇÃO DA INVERSÃO

16

Figura 3.2: Gráfico do sinal h[k].
discreta
g[k] = ([

1 4 1
] ? c)[k],
6 6 6

onde
(g ? h)[k] =

X

g[i]h[k − i].

i∈Z

Considerando que o produto de convolução é associativo, temos que, ao menos formalmente,
1 4 1 −1
c[k] = ([
] ? g)[k],
6 6 6
e o problema se reduz a encontrar a inversa com respeito à convolução do sinal discreto
[ 16 46 16 ]. A técnica da Transformada Z pode ser utilizada para esse fim, mas nesse caso
em particular é fácil verificar diretamente que
h[k] =

−6z |k|
z ,
1 − z2

√
com z = −2 + 3, é a inversa de [ 16 46 16 ], basta usar a definição de convolução discreta e
o fato que z 2 + 1 = −4z.

3.2

Implementação da Inversão

Vemos pelo gráfico de h[k] (figura 3.2) que seu valor absoluto cai exponencialmente a
medida que k vai para o infinito. Assim podemos truncar o sinal h[k] para um valor de
k0 suficientemente grande, ou seja, considerar que h[k] = 0 se |k| > k0 , e obter os valores
dos coeficientes c[k] convoluindo g[k] com o sinal truncado. Mas é possı́vel implementar
o cálculo dos coeficientes c[k] mais eficientemente, do seguinte modo: temos que
h[k] =

−6z +
(y [k] + y + [−k] − 1),
2
1−z

3.2. IMPLEMENTAÇÃO DA INVERSÃO

17

(a) g[k]

(b) f (x)

Figura 3.3: Interpolação Cúbica.
onde

(
y + [k] =

0 , se k < 0
.
z k , se k ≥ 0

Como y + [k + 1] = zy + [k], se k > 0, temos o = y + ? g pode ser calculada recursivamente
com regra de atualização
o[k] := g[k] + zo[k − 1].
Por outro lado, como c = h ? g, temos que c e o satisfazem a relação
c[k] = z(c[k + 1] − o[k]).
Assim, o valor dos coeficientes c[k] pode ser calculado em duas etapas: com k variando
de 1 a N , fazemos
o[k] := g[k] + zo[k − 1],
e com k variando de N − 1 a 0, fazemos
c[k] := z(c[k + 1] − o[k]).
Resta apenas determinar os valores iniciais dessas duas recursões, o[0] e c[N ], respectivamente. Mas isso depende de como o sinal é estendido para k < 0 e k > N . Dependendo
da extensão do sinal, teremos condições iniciais distintas. Na extensão por reflexão, na
qual, g[k] = g[−k], se −N ≤ k < 0 e g[k] = g[2N − k], se N < k ≤ 2N , temos das
definições de o e c que
o[0] =

k0
X
k=0

z k g[k] e c[N ] =

−6z
(2o[N ] − g[N ]).
1 − z2

Esse algoritmo é bastante eficiente, pois são efetuadas apenas duas multiplicações e adições
para cada entrada k. Na próxima seção, analisaremos como calcular as derivadas de f até
a segunda ordem.

3.3. DERIVADAS DAS FUNÇÕES SPLINE

18

(b) B30 (x)

(a) B3 (x)

(c) B300 (x)

Figura 3.4: Derivadas da função spline cúbica.

3.3

Derivadas das Funções Spline

Agora que sabemos como avaliar a função f utilizando B-Splines cúbicas, podemos analisar o problema de calcular as derivadas desta função utilizando os coeficientes obtidos.
Encontrar as derivadas da função f se torna uma tarefa simples quando a caracterizamos como
X
c[k]B3 (x − k),
f (x) =
k∈Z

pois, pela linearidade da derivação,
f 0 (x) =

X

c[k]B30 (x − k) e f 00 (x) =

k∈Z

X

c[k]B300 (x − k).

k∈Z

onde B30 e B300 indicam as derivadas de primeira e segunda ordem da função B3 descritas
pelas equações 3.5 e 3.6.


, se |x| ≥ 2
 0


1 2


, se − 2 < x ≤ −1
 2 x + 2x + 2
0
3 2
B3 (x) =
− 2 x − 2x
, se − 1 < x ≤ 0


3 2

x − 2x
, se 0 < x ≤ 1


2

 1 2
− 2 x + 2x − 2 , se 1 < x < 2


0
, se |x| ≥ 2





, se − 2 < x ≤ −1
 x+2
00
B3 (x) =
−3x − 2 , se − 1 < x ≤ 0



3x − 2
, se 0 < x ≤ 1




−x + 2 , se 1 < x < 2
Estas derivadas estão representadas na figura 3.4.

(3.5)

(3.6)

3.4. AVALIAÇÃO RÁPIDA

3.4

19

Avaliação Rápida

Apresentaremos agora uma versão otimizada para o cálculo da convolução vista na seção
3.1. Esse resultado foi obtido por Sigg e Hadwiger em [24] com o intuito de aproveitar
melhor o potencial oferecido pelas placas gráficas. Eles levaram em consideração que o
acesso linear à textura é muito mais rápido em relação a dois acessos de vizinho mais
próximo, mesmo as duas operações acessando o mesmo número de texels1 . O objetivo é
se beneficiar disso durante o cálculo da convolução.
Como o suporte da função B3 tem largura 4 podemos escrever a equação 3.2 da seguinte
forma:
f (x) = Bi−1 · gi−1 + Bi · gi + Bi+1 · gi+1 + Bi+2 · gi+2 .
(3.7)
onde g[i] = gi e B3 (x − i) = Bi .
A idéia básica do algoritmo idealizado por Sigg e Hadwiger é reescrever a equação 3.7
como a soma ponderada de interpolações lineares entre cada par de valores da função.
Como foi visto, a interpolação linear é calculada através de uma combinação convexa da
seguinte forma
f (x) = (1 − t)g[i] + tg[i + 1],
onde i = bxc denota a parte inteira e t = {x} a parte fracionária de x e g[k] = f (k) ∀k ∈ Z.
Podemos reescrever uma combinação linear qualquer a · g[i] + b · g[i + 1] construı́da através
de uma combinação convexa com a e b quaisquer como
(a + b) · f (i + b/(a + b))

(3.8)

sempre que a condição de combinação convexa 0 ≤ (b/(a + b)) ≤ 1 for satisfeita. Assim,
ao invés de realizar duas buscas na memória de textura para recuperar os valores g[i]
e g[i + 1] e uma interpolação linear, faz-se apenas uma busca em i + b/(a + b) e uma
multiplicação por (a + b).
Podemos então reescrever a equação 3.7 como
f (x) = q0 (x) · f (x − h0 (x)) + q1 (x) · f (x + h1 (x)),

(3.9)

onde as novas funções g0 , g1 , h0 e h1 são definidas da seguinte forma
q0 (x) = Bi−1 + Bi

i
+x
h0 (x) = 1 − Bi−1B+B
i

i+2
q1 (x) = Bi+1 + Bi+2 h1 (x) = 1 + Bi+1B+B
.
i+2

Através desse método, a convolução pode ser realizada utilizando apenas dois acessos
lineares à textura e uma interpolação linear, o que é muito mais rápido que quatro acessos
1

Texel é a menor unidade de textura.

3.5. CASO TRIDIMENSIONAL

20

de vizinho mais próximo.

3.5

Caso Tridimensional

Os passos descritos nas seções anteriores para o cálculo de interpolação utilizando funções B-spline cúbicas podem ser facilmente estendidos para o caso trimensional. Neste
contexto, f passa a ser avaliada a partir de um produto tensorial.
No caso geral temos então o seguinte
X

f (x, y, z) =

c[i, j, k]B3 (x − i)B3 (y − j)B3 (z − k).

(3.10)

i,j,k∈Z

Podemos ainda calcular as derivadas parciais de primeira e segunda ordem de f da
seguinte forma
fx (x, y, z) =

X

c[i, j, k]B30 (x − i)B3 (y − j)B3 (z − k),

i,j,k∈Z

fy (x, y, z) =

X

c[i, j, k]B3 (x − i)B30 (y − j)B3 (z − k),

i,j,k∈Z

fz (x, y, z) =

X

c[i, j, k]B3 (x − i)B3 (y − j)B30 (z − k),

i,j,k∈Z

fxx (x, y, z) =

X

c[i, j, k]B300 (x − i)B3 (y − j)B3 (z − k),

i,j,k∈Z

fxy (x, y, z) =

X

c[i, j, k]B30 (x − i)B30 (y − j)B3 (z − k),

i,j,k∈Z

fxz (x, y, z) =

X

c[i, j, k]B30 (x − i)B3 (y − j)B30 (z − k),

i,j,k∈Z

fyy (x, y, z) =

X

c[i, j, k]B3 (x − i)B300 (y − j)B3 (z − k),

i,j,k∈Z

fyz (x, y, z) =

X

c[i, j, k]B3 (x − i)B30 (y − j)B30 (z − k),

i,j,k∈Z

fzz (x, y, z) =

X

c[i, j, k]B3 (x − i)B3 (y − j)B300 (z − k).

i,j,k∈Z

Pela natureza das equações apresentadas, percebe-se facilmente que fxy = fyx , fxz = fzx
e fyz = fzy .
Como pode ser notado, o cálculo desses somatórios possui um alto custo computacional, entretanto, as considerações que acabamos de fazer sobre o método de avaliação
rápida para o caso unidimensional podem ser facilmente estendidas para o cenário tridi-

3.6. CURVATURA

21

mensional. Nesse caso, é possı́vel realizar o somatório de 64 termos descrito na equação
3.10 utilizando apenas 8 acessos à textura.
No caso 3D os pesos e as distâncias entre as amostras podem ser calculados separadamente para cada uma das dimensões. Os seus valores finais são dados por
g→ =
j

Y

→

gjk e h → =
j

X→
ek · hjk

→

onde k denota os eixos e ek os vetores da base.
As figuras 3.5 e 3.6 apresentam respectivamente as derivadas de primeira e segunda
ordem da função f para o caso 2D.

(a) B3 (x, y)

(b) ∂x B3 (x, y)

(c) ∂y B3 (x, y)

Figura 3.5: Funções B-spline bidimensionais e suas derivadas parciais de primeira ordem.

3.6

Curvatura

Nesta seção utilizaremos os conhecimentos adquiridos durante este capı́tulo para entender
como calcular valores de curvatura sobre uma superfı́cie. Para tanto, seguiremos os passos
descritos por Kindlmann et al. em [9] que sustenta seu método no fato de a curvatura
em uma superfı́cie ser definida como a taxa de variação da normal na vizinhança de um
ponto da superfı́cie com respeito a uma variação infinitesimal da posição do ponto.

3.6. CURVATURA

22

(a) ∂xx B3 (x, y)

(b) ∂xy B3 (x, y)

(c) ∂yy B3 (x, y)

Figura 3.6: Funções B-spline bidimensionais e suas derivadas parciais de segunda ordem.
Seja f uma função do cubo [0, 1]3 no intervalo [0, 1]. Convencionemos que os valores
de f aumentam à medida que nos afastamos do centro do cubo (padrão seguido pelos
equipamentos de tomografia computadorizada) e lembremos do capı́tulo 2 que a normal
∂f ∂f
sobre a superfı́cies é definida como N = g/|g| onde g = ∇f = [ ∂f
]. As informações
∂x ∂y ∂z
T
sobre a curvatura estão em ∇N que é em uma matriz 3 × 3 e pode ser expandida da
seguinte forma:
∇N

T

=
=
=
=

 T

 T
∇g
g∇T |g|
g
=−
−
−∇
|g|
|g|
|g|2




1
g∇T (g T g)1/2
1
g∇T (g T g)
−
H−
=−
H−
|g|
|g|
|g|
2|g|(g T g)1/2




1
g(2g T H)
1
gg T
−
H−
=−
1− 2H
|g|
2|g|2
|g|
|g|

1
−
I − N N T H.
|g|

3.6. CURVATURA

23

onde I é a matriz identidade e H a matriz Hessiana definida por

∂ 2 f /∂x2 ∂ 2 f /∂x∂y ∂ 2 f /∂x∂z


H =  ∂ 2 f /∂y∂x ∂ 2 f /∂y 2 ∂ 2 f /∂y∂z  .
∂ 2 f /∂z∂x ∂ 2 f /∂z∂y ∂ 2 f /∂z 2


O produto externo de N com si próprio, N N T , é um operador linear que projeta
vetores sobre o subespaço vetorial unidimensional gerado por N . O operador I − N N T
projeta sobre o complemento ortogonal desse subespaço, que vem a ser o plano tangente
à isosuperfı́cie em p. Fazendo P = I − N N T temos
∇N T = −

1
P H.
|g|

(3.11)

A equação 3.11 fornece uma boa descrição de ∇N T . A matriz Hessiana representa a
maneira como o gradiente g varia de acordo com uma função de mudanças infinitesimais
de posição em R3 . Essa variação em g possui uma componente ao longo de g (o gradiente
pode variar o tamanho) e uma componente no plano tangente (o gradiente pode mudar de
direção). Apenas o último componente é relevante para descrever a curvatura. Esse componente pode ser isolado apenas com uma multiplicação à esquerda por P . Finalmente, o
fator escalar 1/|g| converte mudanças infinitesimais do gradiente g (não normalizado) em
mudanças infinitesimais do vetor normal unitário N .
As matrizes P e H são ambas simétricas, entretanto ∇N T não é simétrica em geral.
Apesar disso, se v está sobre o plano tangente tem-se P v = v e v T P = v T , portanto para
u e v no plano tangente,
v T P Hu = v T Hu = uT Hv = uT P Hv.
Isso siginifica que a restrição de ∇N T = −P H/|g| ao plano tangente é simétrica e,
portanto, existe uma base ortonormal {p1 , p2 } para o plano tangente na qual ∇N T é uma
matriz diagonal 2 × 2. Essa base pode ser facilmente estendida para uma base ortonormal
para todo R3 , {p1 , p2 , N }. Nessa base, a derivada da normal a superfı́cie é dada por



k1 0 σ1


∇N T =  0 k2 σ2  .
0 0 0
A última linha é nula porque nenhuma variação na posição fará a normal mudar de
tamanho. Mudanças na posição do plano tangente ao longo de p1 e p2 provocam mudanças
em N na mesma direção com proporções k1 e k2 respectivamente. Pela definição de
curvatura em superfı́cies [2], p1 e p2 são as direções principais enquanto k1 e k2 são as

3.6. CURVATURA

24

curvaturas principais. A medida que nos movemos na direção da normal, se aproximando
ou afastando da superfı́cie, a normal inclina de acordo com os valores de σ1 e σ2 .
Multiplicar ∇N T por P isola k1 e k2 na base {p1 , p2 , N }, dessa maneira obtemos


 

1 0 0
k1 0 0

 

G = ∇N T P = ∇N T  0 1 0  =  0 k2 0  .
0 0 0
0 0 0

(3.12)

A avaliação da curvatura em questão é baseada em G, que é chamado de tensor
geométrico. O traço da matriz G é k1 + k2 e a norma de Frobenius de G, denotada por
p
p
|G|F , é definida como tr(GGT ) = k12 + k22 .
Sintetizando o que foi visto, enumeremos os passos necessários para calcular curvatura
em pontos arbitrários de um campo escalar.
1. Calcular as derivadas parciais de primeira ordem e obter o gradiente g. Calcular
N = −g/|g| e P = I − N N T .
2. Calcular as derivadas parciais de segunda ordem e obter a Hessiana. Calcular G =
−P HP/|g|.
3. Calcular o traço T e a norma de Frobenius F de G e obter as curvaturas principais
k1 =

T+

√

2F 2 − T 2
T−
e k2 =
2

√

2F 2 − T 2
.
2

Lembramos que a curvatura gaussiana é definida como K = k1 · k2 e a curvatura média
como H = (k1 + k2 )/2. A implementação dos algoritmos vistos aqui será discutida no
capı́tulo 5.

Capı́tulo 4
Visualização de Objetos Implı́citos
Baseada em GPU
Nos últimos anos, acompanhamos um desenvolvimento formidável dos hardwares especializados no processamento gráfico impulsionado pela indústria de jogos e filmes de animação. Essas aplicações exigem cada cada vez mais realismo e, consequentemente, o
processamento de um volume cada vez maior de informações em um curto espaço de
tempo.
Existe um cenário de pesquisa crescente em torno desses equipamentos e para facilitar
a tarefa de programação, antes realizada utilizando linguagem de montagem (assembly),
diversas linguagens de alto nı́vel foram criadas para programação de placas gráficas.
Neste capı́tulo, conheceremos um pouco mais sobre as possibilidades oferecidas pelas placas gráficas e as restrições impostas por sua estrutura, que privilegia o alto grau
de paralelismo no processamento ao custo de uma capacidade menor de armazenamento
de dados. Faremos uma breve descrição de uma das linguagens mais utilizadas na programação desses processadores, a OpenGL Shading Language [21], que foi utilizada no
desenvolvimento deste trabalho. Por fim, será apresentado um framework extensı́vel especializado em aplicações de visualização utilizando hardware gráfico chamado Voreen, para
o qual foram desenvolvidos módulos para análise dos problemas envolvendo o cálculo de
curvatura.

4.1

GPU

A unidade gráfica de processamento (GPU) é um tipo especial de processador criado para
poupar a unidade lógica de processamento (CPU) das tarefas relacionadas à visualização.
Essa unidade está presente em diversos dispostivos, de computadores pessoais a video
games, e por sua especificidade possui uma capacidade de processamento gráfico muito
maior se comparada aos processadores comuns. Além disso, a demanda por esse tipo de
25

4.1. GPU

26

processamento só tem aumentado com a exigência cada vez mais comum de respostas em
tempo real para a geração de gráficos em alta definição.
Com esse objetivo, as placas gráficas têm evoluı́do muito ao longo dos anos para
proporcionar um grau cada vez maior de paralelismo às aplicações com o uso de processamento multithread em múltiplos precessadores com enorme poder de processamento e
memórias de acesso super rápido.

.
Figura 4.1: Gráficos de comparação entre peformances da CPU e da GPU. O primeiro
gráfico mostra como a GPU tem obtido um crescimento muito maior da sua capacidade
de processamento em relação à CPU enquanto o segundo traz uma comparação do crescimento na capacidade de transferência de dados em GB/s (dados obtidos de [17]).
A razão para a discrepância entre as operações com valores em ponto flutuante, pre-

4.2. OPENGL SHADING LANGUAGE

27

sente no gráfico da figura 4.1, entre a CPU e a GPU está na diferença de arquitetura
existente entre elas. Esta última foi projetada para a realização de cálculos intensivos
com alto grau de paralelismo. Por esse motivo, são desenvolvidas com um número muito
maior de transistores dedicados ao processamento de dados em relação à sua capacidade
de armazenamento, como ilustrado na imagem 4.2.

Figura 4.2: Diferença de arquitetura entre as GPUs e as CPUs, em verde as unidades
lógicas e aritméticas que são responsáveis pelo processamento dos dados e em laranja as
unidades de memória.

4.2

OpenGL Shading Language

OpenGL Shading Language, informalmente conhecida como GLSL, é uma linguagem de
alto nı́vel para programação de placas gráficas baseada na sintaxe das linguagens C e
C++. GLSL possui a vantagem de ter sido incorporada ao OpenGL [23], uma API para
processamento gráfico multiplataforma presente na maioria dos equipamentos existentes.
Para entendermos um pouco melhor a estrutura da GLSL analisemos a figura 4.3.
Os elementos que aparecem na imagem descrevem as etapas que constituem o pipeline
de renderização que é responsável por transformar primitivas geométricas em imagens.
Primeiro, as primitivas que representam o modelo a ser visualizado são passadas como
entrada. Em um primeiro passo, os vértices do modelo passam por um processo de
transformação de coordenadas que os transporta do sistema de coordenadas do mundo
para o da câmera. Além disso, se existe uma fonte de luz na cena, uma primeira mudança
nas cores dos vértices é realizada seguindo um modelo de iluminação como o descrito no
capı́tulo 2. Para acelerar o processo de renderização, os vértices situados fora do campo de
visão do observador são eliminados em um processo conhecido como clipping. Esta etapa
do pipeline está representada no diagrama da figura 4.3 como o passo de transformação
nos vértices.
O último estágio do pipeline consiste em unir as primitivas em um processo chamado
de rasterização que indica quais pixels da tela elas cobrem. Um fragmento é gerado para
cada pixel, cada um deles contendo informações dos vértices aos quais está ligado como
cor e coordenadas de textura.

4.2. OPENGL SHADING LANGUAGE

28

Figura 4.3: Pipeline de renderização.
O pipeline de renderização é composto por elementos chamados de processadores.
Atualmente, existem três tipos de processadores, vertex processors, geometry processsors
e fragment processors. Estes processadores estão representados na figura 4.4 que mostra
também que as etapas do processo de visualização são independentes, ou seja, os vértices
podem ser processados ao mesmo tempo em que pixels estão sendo processados. Além
da independência entre as etapas, dentro de um mesmo estágio do pipeline os elementos
podem ser processados simultaneamente por processadores diferentes.

Figura 4.4: Fluxo de dados no pipeline de renderização
Nas próximas linhas descreveremos melhor cada um dos processadores do ponto de
vista da GLSL.

4.2.1

Vertex Processor

O vertex processor é uma unidade programável que processa vértices e os atributos relacionados a eles. Os programas escritos para rodar nesses processadores são conhecidos
como vertex shaders.
Esses elementos podem ser usados com várias finalidades como transformação de vértices, transformação e normalização de vetores normais, geração e transformação de coordenadas de textura etc. Entretanto, apenas um vértice é processado por vez, não sendo

4.2. OPENGL SHADING LANGUAGE

29

Figura 4.5: Etapa do pipeline de renderização relacionada ao vertex shader
possı́vel executar operações que dependam das informações contidas em uma série de vértices simultaneamente. O código fonte a seguir mostra um exemplo simples de vertex
shader.




1
2
3
4
5
6
7

// Vert ex Shader Main
void main ( void ) {
// Pass v e r t e x c o l o r t o n e x t s t a g e
gl FrontColor = gl Color ;
// Transform v e r t e x p o s i t i o n b e f o r e p a s s i n g i t
g l P o s i t i o n = gl ModelViewProjectionMatrix ∗ gl Vertex ;
}





Código Fonte 4.1: Exemplo de Vertex Shader
Assim como na linguagem C, para cada shader uma função principal é obrigatória e a
sua declaração é análoga. Nas linhas seguintes, a cor do vértice é delegada a um próximo
passo e a posição é transformada para o espaço da tela.

4.2.2

Geometry Processor

Um geometry processor é uma unidade responsável pelo processamento de vértices que
fazem parte de uma mesma primitiva, linhas, triângulos etc., e tem como saı́da uma
sequência de vértices formando outras primitivas. Os programas escritos para rodar nesse
tipo de unidade são chamados de geometry shaders.
Continuando nosso exemplo, o código abaixo exemplifica um geometry shader.

1
2
3
4
5
6

#e x t e n s i o n GL EXT geometry shader4 : e n a b l e
void main ( ) {
// I t e r a t e s o v e r a l l v e r t i c e s i n t h e i n p u t p r i m i t i v e
f o r ( i n t i = 0 ; i < g l V e rt i ce s I n ; ++i ) {
// Pass c o l o r and p o s i t i o n t o n e x t s t a g e
gl FrontColor = gl FrontColor [ i ] ;




4.2. OPENGL SHADING LANGUAGE

gl Position = gl PositionIn [ i ] ;
// Done with t h i s v e r t e x
EmitVertex ( ) ;

7
8
9

}
// Done with t h e i n p u t p r i m i t i v e
EndPrimitive ( ) ;

10
11
12
13

30

}





Código Fonte 4.2: Exemplo de Geometry Shader
A primeira linha do código, habilita a extensão geometry shader definida pelo shader
model 4.0. Essa diretiva é sempre necessária ao utilizar um geometry shader. A função
principal possui um loop cujo fator de iteração é o número de vértices da primitiva.

Figura 4.6: Etapa do pipeline de renderização relacionada ao geometry shader

4.2.3

Fragment Processor

Um fragment processor é uma unidade que opera sobre fragmentos e os dados associados
a eles. Os programas escritos para esses processadores são chamados de fragment shaders.
Cada iteração recebe a cor e a posição do vértice correspondente e associa esses valores à
cor e à posição dos vértices de saı́da. A função EmitVertex, como o próprio nome indica,
emite o vértice para saı́da, enquanto a função EndPrimitive encerra a primitiva.
Um fragment shader não pode alterar a posição de um fragmento. O acesso a fragmentos vizinhos também não é permitido. Os valores calculados nesses shaders são utilizados
para atualizar dados na memória de um framebuffer ou de uma textura.
Esse tipo de processador é utilizado para realizar operações gráficas como as de interpolação, acesso a textura, soma de cores etc. As linhas de código a seguir trazem um
exemplo simples de fragment shader.

1
2

// Fragment Shader Main
void main ( void ) {




4.2. OPENGL SHADING LANGUAGE

31

Figura 4.7: Etapa do pipeline de renderização relacionada ao fragment shader
// Pass f r a g m e n t c o l o r
gl FragColor = gl Color ;

3
4
5

}





Código Fonte 4.3: Exemplo de Fragment Shader
Nesse exemplo, o identificador gl_Color refere-se a um tipo diferente do mostrado no
vertex shader. No primeiro caso, gl_Color foi associado à saı́da do pipeline como um
atributo do vértice. No caso do fragment shader ele foi definido dentro do pipeline como
um valor interpolado.
A imagem 4.8 mostra a saı́da do pipeline descrito acima aplicado ao volume smile,
onde o resultado é apenas a silhueta do objeto visto que, para efeito de simplificação,
nenhum cálculo de iluminação foi definido.
Mais detalhes e alguns exemplos básicos da utilização da linguagem GLSL podem ser
encontrados em [13].

Figura 4.8: Resultado da execução dos shaders

4.3. GPGPU E CUDA

4.3

32

GPGPU e CUDA

O crescimento vertiginoso dos hardwares gráficos provocou o surgimento de um outro tipo
de demanda para as GPUs, agora os desenvolvedores de softwares empregam o seu alto
poder de processamento em atividades diversas que lidam com a manipulação de grandes
volumes de dados. Com o objetivo de aproveitar todo o potencial oferecido pelos múltiplos
processadores das placas gráficas algumas ferramentas foram criadas em uma vertente de
pesquisa denominada GPGPU [6] (sigla para General-Purpose Computation on Graphics
Hardware) que vem crescendo em várias áreas como simulação fı́sica e processamento de
sinais.
Foi com esse pensamento que a NVIDIA [18] apresentou em novembro de 2006 uma
nova arquitetura e linguagem para programação paralela em GPU chamada CUDA [17].
O objetivo dessa ferramenta é desenvolver aplicações que se adaptem bem à GPU, ou
seja, que possam ser divididas em núcleos de processamento independentes. Assim como
GLSL, CUDA foi desenvolvida para programadores familiarizados com a linguagem C e
em breve dará suporte a outras linguagens como C++ e FORTRAN.
A linguagem CUDA foi brevemente testada durante o desenvolvimento deste trabalho,
mas, para os fins em que foi empregada, não ofereceu nenhuma vantagem que justificasse
seu uso em um primeiro momento. Apesar disso, existe um grande potencial para a
utilização dessa ferramenta em algoritmos baseados em lançamento de raios, já que na
maioria dos casos o processamento em cada pixel é independente. Portanto, a possibilidade
de utilizarmos a linguagem para o processamento de curvatura em superfı́cies, bem como
para outras finalidades, não foi descartada.

4.4

Voreen: Modularidade

Esta seção é dedicada à descrição do sistema utilizado nesse trabalho para desenvolver
os módulos de visualização utilizados no estudo dos problemas envolvidos no cálculo de
curvatura. O Voreen (contração de Volume Rendering Engine) é um sistema open source
desenvolvido em C++ sob os termos da GPL 2.0 [4]. O objetivo dessa ferramenta é fornecer o acesso a algoritmos de rendering de volumes e processamento de imagens baseados
em GPU que permitem a visualização interativa de conjuntos de dados volumétricos de
maneira flexı́vel. Essa flexibilidade é resultante da caracterı́stica modular do framework,
que permite a associação de várias técnicas de visualização através de estruturas chamadas
de processadores (processors). A figura 4.9 mostra uma visão geral do sistema.
O rendering de volumes é realizado através da criação de redes (networks) que são
constituı́das de processadores, cada um deles responsável por uma tarefa especı́fica. O
usuário pode facilmente estender esses elementos para, por exemplo, integrar novos algoritmos de processamento de volumes.

4.4. VOREEN: MODULARIDADE

33

Figura 4.9: Visão Geral do Voreen

Figura 4.10: Rede de processadores.
Os dados são trocados entre os diferentes processadores através de portas, inport para
receber dados e outport para fornecê-los. Os tipos de porta indicam as conexões possı́veis
de entrada/saı́da entre os processadores, apenas portas de dados compatı́veis podem ser
conectadas. Na imagem 4.10, os quadrados coloridos representam os diferentes tipos de

4.4. VOREEN: MODULARIDADE

34

portas interligados pelas setas em preto.
O Voreen foi construı́do seguindo o paradigma de orientação a objetos (OO) onde cada
objeto possui sua própria maneira de controlar e manipular dados. Esse paradigma não
se aplica à arquitetura do OpenGL que tem seu estado determinado por variáveis globais
e onde os dados podem ser acessados em qualquer parte do programa. Sendo assim, os
desenvolvedores precisam ser cuidadosos ao escolher entre seguir os princı́pios de OO ou
obter alta peformance em processamento gráfico.

4.4.1

Classes de Processadores

Nos próximos parágrafos serão apresentadas algumas classes importantes de processadores do Voreen. No capı́tulo 4 serão descritos os principais processadores implementados
durante o desenvolvimento deste trabalho como partes constituintes dos módulos criados
para estudo de curvatura e métodos de interpolação.
O processador VolumeSource permite ao usuário criar seu próprio conjunto de volumes
de dados incluindo arquivos de tipos suportados pelo Voreen [26]. A partir desse conjunto,
volumes podem ser facilmente selecionados para visualização.
Um EntryExitPoints cria imagens contendo a posição de pontos de entrada e saı́da
armazenados como cores. Esses pontos, que posteriormente serão armazenados como
texturas, serão utilizados no processo de visualização direta como dados de entrada do
algoritmo. Para volumes em forma de cubo, os dados armazenados nessas texturas formam
um cubo de cores RGB (figura 4.11).

Figura 4.11: Texturas com os pontos de entrada e saı́da do volume em forma de cubo.
Finalmente, os processadores do tipo RayCaster recebem como entrada os pontos de
entrada e saı́da, assim como o volume de dados. A partir desse momento, o algoritmo de
ray casting é executado em um fragment program na GPU usando a linguagem GLSL.
Cada um dos processadores mencionados possui uma série de propriedades que podem
ser alteradas para gerar diferentes tipos de visualização em tempo real. Essas propriedades
podem ter seus valores compartilhados entre processadores para que a alteração em seu
estado tenha efeito em todos que a compartilham. A figura 4.12 mostra dois processadores
que compartilham as mesmas propriedades de câmera.

4.4. VOREEN: MODULARIDADE

35

Figura 4.12: Propriedade camera compartilhadas entre componentes.
Existem muitos outros processadores que podem ser associados aos mencionados para
criar diferentes resultados. Optou-se por omitir esses processadores pois eles não são o
foco deste trabalho. Para mais informações consultar [15] que possui uma descrição mais
detalhada do framework ou o site do projeto que possui entre outras informações, o código
fonte do sistema para análise e possı́veis colaborações [26].

4.4.2

Criando Processadores

Nesta seção será apresentado um exemplo simples de criação de componente para o Voreen.
Este exemplo será importante para o entendimento do capı́tulo 5 que apresentará os
módulos desenvolvidos para esse trabalho e agregados ao Voreen.
O processador em questão terá como tarefa converter uma imagem em RGB para uma
imagem em escala de cinza. Este processador está implementado no Voreen com o nome
de Grayscale relacionado ao grupo de compenentes para processamento de imagens.

1
2
3
4
5
6
7
8
9
10
11
12

Grayscale : : Grayscale ( )
: ImageProcessor ( " pp_grayscale " ) , // l o a d s f r a g m e n t s h a d e r p p g r e y s c a l e . f r a g
s a t u r a t i o n ( " saturation " , " Saturation " , 0 . 0 f ) ,
i n p o r t ( Port : : INPORT, " inport " ) ,
o u t p o r t ( Port : : OUTPORT, " outport " )
{
// r e g i s t e r p r o p e r t i e s and p o r t s :
addProperty ( s a t u r a t i o n ) ;
addPort ( i n p o r t ) ;
addPort ( o u t p o r t ) ;
}




4.4. VOREEN: MODULARIDADE

13
14
15
16
17
18

19
20
21
22
23
24
25
26
27
28
29
30
31
32

36

void Grayscale : : process ( ) {
// a c t i v a t e and c l e a r o utput r e n d e r t a r g e t
outport . activateTarget () ;
glClear (GL COLOR BUFFER BIT | GL DEPTH BUFFER BIT) ;
// bind r e s u l t from p r e v i o u s p r o c e s s o r
i n p o r t . b i n d T e x t u r e s ( tm . getGLTexUnit ( shadeTexUnit ) , tm . getGLTexUnit ( depthTexUnit )
);
// a c t i v a t e s h a d e r and s e t u n i f o r m s :
program −>a c t i v a t e ( ) ;
s e t G l o b a l S h a d e r P a r a m e t e r s ( program ) ;
program −>s e t U n i f o r m ( " shadeTex_ " , tm . getTexUnit ( shadeTexUnit ) ) ;
program −>s e t U n i f o r m ( " depthTex_ " , tm . getTexUnit ( depthTexUnit ) ) ;
program −>s e t U n i f o r m ( " saturation_ " , s a t u r a t i o n . g e t ( ) ) ;
// r e n d e r s c r e e n a l i g n e d quad :
glDepthFunc (GL ALWAYS) ;
renderQuad ( ) ;
glDepthFunc ( GL LESS ) ;
program −>d e a c t i v a t e ( ) ;
g l A c t i v e T e x t u r e (GL TEXTURE0) ;
LGL ERROR;
}





Código Fonte 4.4: Este exemplo mostra a criação de um processador para criar imagens
em escala de cinza.
No código fonte apresentado acima iniciaremos explicando o que acontece no método
process que deve ser implmententado por todas as classes que herdam de Processor,
esse é o casso da classe GrayScale. A linha 15 ativa a textura que será passada como
saı́da do processador. A linha 18 registra as texturas que serão passadas como dado de
entrada para o fragment shader pp_grayscale na linha 22 através da chamada ao método
setUniform. As linhas 20 e 29 alteram o estado do objeto program_ ativando e desativando o processamento do código do fragment shader em questão. O primeiro método
apresentado no código é o construtor do processador que inicializa alguns atributos do objeto e as linhas 9 e 10 criam as portas que servirão de comunicação entre os processadores
para troca de informações, uma porta de entrada para a imagem original e uma porta de
saı́da para a imagem em tons de cinza.

1
2
3
4

uniform SAMPLER2D TYPE shadeTex ;
uniform SAMPLER2D TYPE depthTex ;
uniform TEXTURE PARAMETERS t e x t u r e P a r a m e t e r s ;
uniform f l o a t s a t u r a t i o n ;

5
6
7
8
9
10

vec4 toGrayScale ( i n vec4 c o l o r , i n f l o a t s a t u r a t i o n ) {
float brightness = (0.30 ∗ color . r ) + (0.59 ∗ color . g) + (0.11 ∗ color . b) ;
vec4 g r a y s c a l e = vec4 ( b r i g h t n e s s , b r i g h t n e s s , b r i g h t n e s s , c o l o r . a ) ;
return mix ( g r a y s c a l e , c o l o r , s a t u r a t i o n ) ;
}

11
12
13
14

void main ( ) {
v e c 2 p = g l F r a g C o o r d . xy ∗ screenDimRCP ;
gl FragColor = toGrayScale ( t e x t u r e L o o k u p 2 D n o r m a l i z e d ( shadeTex , t e x t u r e P a r a m e t e r s , p
) , saturation ) ;




4.4. VOREEN: MODULARIDADE

gl FragDepth = t e x t u r e L o o k u p 2 D n o r m a l i z e d ( depthTex , t e x t u r e P a r a m e t e r s , p ) . z ;

15
16

37

}





Código Fonte 4.5: Fragment shader associado ao processador Grayscale.
No código do fragment shader acima as primeiras linhas são declarações de variáveis
que armazenarão as informações de textura e uma váriavel float que guardará um valor
de saturação. Todos esses dados são provenientes do processador Grayscale e serão
apenas processados pelo shader sem sofrer nenhuma alteração.
A linha 13 recupera a coordenada do pixel que está sendo processado enquanto a linha
14 faz a chamada ao método que irá converter a cor original do pixel em tons de cinza e
passar essa informação adiante.

Capı́tulo 5
Módulos
Para este trabalho foram desenvolvidos alguns módulos empregados no cálculo de curvatura e na avaliação de funções utilizando interpolação tricúbica. Esses módulos foram
incorporados ao framework Voreen na forma de processadores. O objetivo desse capı́tulo
é descrever os aspectos mais relevantes da construção desses processadores e apresentar
os algoritmos desenvolvidos utilizando a teoria dos capı́tulos iniciais.

5.1

Módulo de Funções

Este processador foi criado para eliminar a necessidade de leitura de arquivos com os dados
volumétricos. Não há nenhuma entrada para esse componente e a sua saı́da é uma textura
com um volume de dados indicado pelo usuário através de uma de suas propriedades.
Os volumes são gerados através de funções que descrevem o volume implicitamente. A
imagem 5.1 mostra o processador em destaque com as suas propriedades que permitem
escolher o tipo utilizado para armazenar os valores (inteiro de 8 e 16 bits e ponto flutuante)
e a resolução do volume (entre 24 e 28 ). Já a imagem 5.1 traz alguns exemplos de volumes
gerados com o processador ao qual demos o nome de FunctionSelector.

Figura 5.1: FunctionSelector

38

5.2. MÓDULO DE INTERPOLAÇÃO TRICÚBICA

(a) Cubo

(b) Toro

(c) Smile

(d) Esfera

39

Figura 5.2: Imagens Geradas com o processador FunctionSelector

5.2

Módulo de Interpolação Tricúbica

O módulo de interpolação tricúbica encontra-se implementado na classe TricubicSplineProcessor, esta classe implementa o processador CubicSplineProcessor (em destaque
na figura 5.3). A entrada desse processador é um volume de dados como o descrito nos
capı́tulos anteriores e saı́da é um outro volume com as mesmas caracterı́sticas. Entretanto,
os valores dos voxels são convertidos em coeficientes B-Splines que serão utilizados para
a avaliação dos valores da função nos pontos do volume utilizando interpolação tricúbica.

Figura 5.3: CubicSplineProcessor

5.3. MÓDULOS DE CÁLCULOS EM GLSL (INTERPOLAÇÃO TRICÚBICA E
CURVATURA)

40

A única propriedade do processador habilita ou desabilita o processamento do volume.

5.3

Módulos de Cálculos em GLSL (interpolação tricúbica e curvatura)

Esse módulos equivalem aos shaders criados para o cálculo de curvatura e para avaliação
de funções utilizando interpolação tricúbica. Utilizaremos o código abaixo para descrever
o que acontece no módulo de curvatura.

1
2
3

f l o a t m a t r i x t r a c e ( i n mat3 m) {
return m[ 0 ] [ 0 ] +m[ 1 ] [ 1 ] +m [ 2 ] [ 2 ] ;
}

4
5
6
7
8
9
10

f l o a t f r o b e n i u s n o r m ( i n mat3 m) {
f l o a t x = m[ 0 ] [ 0 ] ∗m[ 0 ] [ 0 ] +m[ 0 ] [ 1 ] ∗m[ 0 ] [ 1 ] +m[ 0 ] [ 2 ] ∗m[ 0 ] [ 2 ] +
m[ 1 ] [ 0 ] ∗m[ 1 ] [ 0 ] +m[ 1 ] [ 1 ] ∗m[ 1 ] [ 1 ] +m[ 1 ] [ 2 ] ∗m[ 1 ] [ 2 ] +
m[ 2 ] [ 0 ] ∗m[ 2 ] [ 0 ] +m[ 2 ] [ 1 ] ∗m[ 2 ] [ 1 ] +m[ 2 ] [ 2 ] ∗m [ 2 ] [ 2 ] ;
return s q r t ( x ) ;
}

11
12
13

14
15

mat3 c a l c n n T ( i n v e c 3 n ) {
mat3 nnT = mat3 ( n . x∗n . x , n . y∗n . x , n . z ∗n . x , n . x∗n . y , n . y∗n . y , n . z ∗n . y , n . x∗n . z , n .
y∗n . z , n . z ∗n . z ) ;
return nnT ;
}

16
17
18
19
20

mat3 c a l c P ( i n mat3 nnT ) {
mat3 I = mat3 ( 1 . 0 ) ;
return I−nnT ;
}

21
22
23
24
25
26
27

// n g − g r a d i e n t norm
mat3 c a l c G ( i n mat3 P , i n mat3 H, i n f l o a t n g ) {
mat3 G = mat3 ( 0 . 0 ) ;
G = (P∗H∗P) ∗ ( 1 / n g ) ;
return G;
}

28
29
30
31
32
33
34
35
36
37
38
39
40
41
42

v e c 2 computeCurvature ( i n v e c 3 grad , i n mat3 h e s s i a n ) {
vec3 g
= grad ;
f l o a t n g = s q r t ( dot ( g , g ) ) ;
g
= normalize ( g ) ;
vec3 n
= v e c 3 (−g . x , −g . y , −g . z ) ;
mat3 nnT = c a l c n n T ( n ) ;
mat3 P
= c a l c P ( nnT ) ;
mat3 G
= c a l c G (P , h e s s i a n , n g ) ;
f l o a t t r a c e = m a t r i x t r a c e (G) ;
f l o a t f norm = f r o b e n i u s n o r m (G) ;
f l o a t k1 = ( t r a c e + s q r t (max ( 0 . 0 , 2∗ f norm ∗ f norm−t r a c e ∗ t r a c e ) ) ) / 2 . f ;
f l o a t k2 = ( t r a c e − s q r t (max ( 0 . 0 , 2∗ f norm ∗ f norm−t r a c e ∗ t r a c e ) ) ) / 2 . f ;
return v e c 2 ( k1 , k2 ) ;
}




5.3. MÓDULOS DE CÁLCULOS EM GLSL (INTERPOLAÇÃO TRICÚBICA E
CURVATURA)

41




Código Fonte 5.1: Módulo de cálculo de curvaturas.
Entre as linhas 1 e 10 estão duas funções extremamente simples, a primeira calcula o
traço de uma matriz e o segundo a norma de Frobenius descritas no capı́tulo 3.
A linha 12 corresponde a declaração da função que realiza o cálculo resultante do
produto entre as matrizes linha e coluna que representam o vetor normal N e a sua
transposta N T . O resultado é uma matriz 3 × 3 que em GLSL corresponde ao tipo mat3.
Na linha 23 está definido o cálculo da matriz G onde a função possui como parâmetros
de entrada as matrizes P e H e a norma do gradiente.
A função computeCurvature realiza o cálculo das curvaturas principais e recebe como
dados de entrada a matriz Hessiana e o gradiente. Este trecho do código apenas faz
chamadas as funções definidas acima e nas linhas 39 e 40 calcula as curvaturas principais
utilizando as fórmulas definidas no capı́tulo 3 onde
√
√
T − 2F 2 − T 2
2F 2 − T 2
e k2 =
.
k1 =
2
2
O código abaixo refere-se a implementação da avaliação de funções utilizando b-splines.
Foram selecionados dois métodos que merecem destaque. No primeiro deles interpolate_tricubic o loop principal representa o somátorio
T+

f (x, y, z) =

X

c[i, j, k]B3 (x − i)B3 (y − j)B3 (z − k),

i,j,k∈Z

definido no capı́tulo 3. O segundo, interpolate_tricubic_fast, implementa o método
rápido de avaliação que faz apenas 8 consultas a memória de textura representados no
código pela função texture3D.

1

2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18

f l o a t i n t e r p o l a t e t r i c u b i c ( sampler3D tex , VOLUME PARAMETERS volumeParameters , v e c 3 c o o r d )
{
v e c 3 n=volumeParameters . d a t a s e t D i m e n s i o n s ;
c o o r d = n∗ c o o r d ; // − 0 . 5 ;
// s h i f t t h e c o o r d i n a t e from [ 0 , e x t e n t ] t o [ − 0 . 5 , e x t e n t − 0 . 5 ]
vec3 c o o r d g r i d = coord − 0 . 5 ;
vec3 index = f l o o r ( c o o r d g r i d ) ;
vec3 f r a c t i o n = c o o r d g r i d − index ;
f l o a t r e s u l t =0;
f o r ( f l o a t z=−1; z < 2 . 5 ; z++) // r a n g e [ −1 , 2 ] {
f l o a t b z = b s p l i n e ( z−f r a c t i o n . z ) ;
f o r ( f l o a t y=−1; y < 2 . 5 ; y++)
{
f l o a t b y = b s p l i n e ( y−f r a c t i o n . y ) ;
f o r ( f l o a t x=−1; x < 2 . 5 ; x++)
{
f l o a t b x = b s p l i n e ( x−f r a c t i o n . x ) ;
vec3 p
= ( i n d e x+v e c 3 ( x , y , z ) +0.5) /n ;
f l o a t c = t e x t u r e 3 D ( tex , p ) . a ;
result
+= c ∗ b x ∗ b y ∗ b z ;
}




5.3. MÓDULOS DE CÁLCULOS EM GLSL (INTERPOLAÇÃO TRICÚBICA E
CURVATURA)
}
}
return r e s u l t ;

19
20
21
22

42

}

23
24

25
26
27
28
29
30
31
32
33
34
35

36

37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56

f l o a t i n t e r p o l a t e t r i c u b i c f a s t ( sampler3D tex , VOLUME PARAMETERS volumeParameters , v e c 3
coord ) {
v e c 3 n=volumeParameters . d a t a s e t D i m e n s i o n s ;
c o o r d = n∗ c o o r d ; // − 0 . 5 ;
// s h i f t t h e c o o r d i n a t e from [ 0 , e x t e n t ] t o [ − 0 . 5 , e x t e n t − 0 . 5 ]
vec3 c o o r d g r i d = coord − 0 . 5 ;
vec3 index = f l o o r ( c o o r d g r i d ) ;
vec3 f r a c t i o n = c o o r d g r i d − index ;
v e c 3 w0 , w1 , w2 , w3 ;
b s p l i n e w e i g h t s ( f r a c t i o n , w0 , w1 , w2 , w3 ) ;
v e c 3 g0 = w0 + w1 ;
v e c 3 g1 = w2 + w3 ;
v e c 3 h0 = ( w1 / g0 ) − 0 . 5 + i n d e x ;
// h0 = w1/ g0 − 1 , move from [ − 0 . 5 , e x t e n t
−0.5] to [ 0 , extent ]
v e c 3 h1 = ( w3 / g1 ) + 1 . 5 + i n d e x ;
// h1 = w3/ g1 + 1 , move from [ − 0 . 5 , e x t e n t
−0.5] to [ 0 , extent ]
h0 = ( h0 ) /n ;
h1 = ( h1 ) /n ;
// f e t c h t h e e i g h t l i n e a r i n t e r p o l a t i o n s
// w e i g h t i n g and f e t c h i n g i s i n t e r l e a v e d f o r p e r f o r m a n c e and s t a b i l i t y r e a s o n s
f l o a t t e x 0 0 0 = t e x t u r e 3 D ( tex , v e c 3 ( h0 . x , h0 . y , h0 . z ) ) . a ;
f l o a t t e x 1 0 0 = t e x t u r e 3 D ( tex , v e c 3 ( h1 . x , h0 . y , h0 . z ) ) . a ;
t e x 0 0 0 = mix ( tex100 , tex000 , g0 . x ) ;
// weigh a l o n g t h e x−d i r e c t i o n
f l o a t t e x 0 1 0 = t e x t u r e 3 D ( tex , v e c 3 ( h0 . x , h1 . y , h0 . z ) ) . a ;
f l o a t t e x 1 1 0 = t e x t u r e 3 D ( tex , v e c 3 ( h1 . x , h1 . y , h0 . z ) ) . a ;
t e x 0 1 0 = mix ( tex110 , tex010 , g0 . x ) ;
// weigh a l o n g t h e x−d i r e c t i o n
t e x 0 0 0 = mix ( tex010 , tex000 , g0 . y ) ;
// weigh a l o n g t h e y−d i r e c t i o n
f l o a t t e x 0 0 1 = t e x t u r e 3 D ( tex , v e c 3 ( h0 . x , h0 . y , h1 . z ) ) . a ;
f l o a t t e x 1 0 1 = t e x t u r e 3 D ( tex , v e c 3 ( h1 . x , h0 . y , h1 . z ) ) . a ;
t e x 0 0 1 = mix ( tex101 , tex001 , g0 . x ) ;
// weigh a l o n g t h e x−d i r e c t i o n
f l o a t t e x 0 1 1 = t e x t u r e 3 D ( tex , v e c 3 ( h0 . x , h1 . y , h1 . z ) ) . a ;
f l o a t t e x 1 1 1 = t e x t u r e 3 D ( tex , v e c 3 ( h1 . x , h1 . y , h1 . z ) ) . a ;
t e x 0 1 1 = mix ( tex111 , tex011 , g0 . x ) ;
// weigh a l o n g t h e x−d i r e c t i o n
t e x 0 0 1 = mix ( tex011 , tex001 , g0 . y ) ;
// weigh a l o n g t h e y−d i r e c t i o n
return mix ( tex001 , tex000 , g0 . z ) ;
// weigh a l o n g t h e z−d i r e c t i o n
}



Código Fonte 5.2: Módulo de avaliação de funções através de interpolação tricúbica.
Antes de encerrar a seção, existe um fato importante que deve ser mencionado. Durante o processo de visualização, a função de transferência utiliza valores no intervalo
[0, 1] para acessar os valores de cor e opacidade associados aos pontos. Por esse motivo, é
necessária a normalização dos valores de curvatura encontrados já que eles possivelmente
não se encaixarão nesse intervalo. Dois tipos de mormalização foram implementados.
O primeiro método para a normalização é linear e utiliza os valores de curvatura
positiva máxima (k max = max{k, 0}) e o valor de curvatura negativa mı́nima (k min =
min{k, 0}) para o cálculo, onde a curvatura k varia entre todos os valores de curvatura



5.4. MÓDULO DE RENDERIZAÇÃO DE ISOSUPERFÍCIES

43

encontrados. O valor normalizado kn de k é então dado por
kn = k/min k se k < 0 e k = k/max k se k ≥ 0.
O segundo método adota uma estratégia não linear. O valor de kn é dado utilizando
a função exponencial. Assim,
kn = (ek·t )/2 se k < 0 e kn = 1 − (e−k·t )/2 caso contrário.
O valor t é uma constante de atenuação.

5.4

Módulo de Renderização de Isosuperfı́cies

O módulo de renderização de isosuperfı́cies está implementado no processador CurvatureRaycaster e utiliza o algoritmo de ray casting descrito na seção 2.2 do capı́tulo 2
para realizar a visualização de isosuperfı́cies. A navegação entre as diferentes isosuperfı́cies do volume é realizada através da alteração do parâmetro level. Este componente é
semelhante ao já existente VolumeRaycaster, entretanto utiliza os valores de curvatura
normalizados que foram calculados no módulo de curvatura e interpolação tricúbica como
parâmetro para função de transferência enquanto o primeiro utiliza os valores de intensidades da função em cada um dos pontos. A imagem 5.4 traz o processador em destaque.
Os dados de entrada do processador são texturas contendo o volume a ser visualizado além
das texturas contendo os pontos de entrada e saı́da mencionadas no capı́tulo anterior.

Figura 5.4: CurvatureRaycaster
A figura 5.5 mostra um resultado obtido com o processador CurvatureRaycaster,
além disso traz algumas de suas propriedades principais como aquelas que permitem editar
a função de transferência (figura 5.6) e os que permitem escolher o nı́vel a ser visualizado,
o tipo de cálculo de derivadas (forward differences, central differences, trilinear e tricubic)
e avaliação de função. Uma outra propriedade que pode ser alterada é o tipo de curvatura
visualizada (Gaussiana, Média, k1 e k2 ).

5.5. MÓDULO DE RENDERIZAÇÃO DE VOLUMES

44

Figura 5.5: Resultado obtido com o CurvatureRaycaster e suas propriedades

Figura 5.6: Editor para funções de transferência

5.5

Módulo de Renderização de Volumes

O módulo de renderização de volumes segue o algoritmo definido no capı́tulo 2. Ele está
implementado no processador CurvatureVolumeRaycaster e é semelhante ao processador
descrito na seção anterior, recebe como dados de entrada o volume e os pontos de entrada
e saı́da. A saı́da desse processador é uma textura com a imagem do volume em diferentes valores de opacidade que dependem da superfı́cie de nı́vel. As cores, como no caso
anterior, são determinadas pelos valores de curvatura normalizados. A figura 5.7 mostra
o processador em destaque já a figura 5.8 traz o resultado obtido mostrando o efeito de
transparência mencionado.

5.5. MÓDULO DE RENDERIZAÇÃO DE VOLUMES

45

Figura 5.7: CurvatureRaycaster

Figura 5.8: Resultado obtido com o CurvatureVolumeRaycaster mostrando diversas
porções do volume com efeito de transparência e levando em consideração os valores de
curvatura para o cálculo das cores.

Capı́tulo 6
Resultados
Este capı́tulo é dedicado à descrição dos resultados obtidos a partir dos módulos implementados e descritos no capı́tulo anterior.

(a) Diferenças Pra Frente

(b) Diferenças Centrais

(c) Filtrado

(d) B-Splines

Figura 6.1: Resultados obtidos utilizando o mapa de cores JET para visualização de
valores de curvatura. Vermelho indica curvaturas positivas, azul curvaturas negativas,
verde valores próximos de zero, e em preto está delimitada a curva de nı́vel zero.
Na imagem 6.1 podemos observar que há uma variação na simetria da linha preta
que indica os valores onde a curvatura é nula sobre o toro. Na imagem 6.1(a), que
utiliza diferenças finitas para frente para o cálculo das derivadas, aparece um padrão
46

47
triangular, em 6.1(b), que foi gerada utilizando diferenças centrais para o cálculo das
derivadas, aparece um padrão em forma de pentágono. Na figura 6.1(c) foi utilizada
uma filtragem que calcula a derivada através de uma ponderação dos valores nos voxels
vizinhos e na 6.1(d) foi utilizada interpolação tricúbica para o cálculo das derivadas,
em ambas o resultado é mais próximo da realidade, entretanto a interpolação tricúbica
oferece melhores resultados. Nas três primeiras imagens a função foi avaliada utilizando
interpolação trilinear e na última foi utilizada a interpolação tricúbica.

(a) Diferenças Pra Frente

(b) Diferenças Centrais

(c) Filtrado

(d) B-Splines

Figura 6.2: Resultados obtidos utilizando o volume smile.
A imagem 6.2 mostra um resultado semelhante ao primeiro aplicado ao volume smile.
A imagem 6.3 mostra o resultado do ray casting de curvaturas aplicado ao volume cow
que sofreu uma suavização em um passo de pré-processamento. A figura 6.3(a) mostra
os valores da curvatura Gaussiana sobre a superfı́cie enquanto a figura 6.3(b) exibe os
valores de curvatura média. As cores das figuras abaixo das duas primeiras representam
os valores das curvaturas principais. Assim como nos resultados do inı́cio do capı́tulo, a cor
azul representa curvatura negativa, verde curvatura nula e vermelho curvatura positiva.
Observe que há uma predominância da cor verde ou de valores muito próximos de zero,
isso deve-se ao passo de suavização do volume. Algumas inconsistências também ocorrem
devido aos ruı́dos provenientes do volume original.
A figura 6.4 exemplifica os problemas que surgem quando se utiliza diferentes formas

48

(a) K

(b) H

(c) k1

(d) k2

Figura 6.3: Resultados obtidos utilizando o volume cow.
de representação numérica para visualização da curvatura Gaussiana. A figura 6.4(a) foi
obtida utilizando inteiro de 8 bits para representação enquanto em 6.4(b) foram utlizados
16 bits para representação. Os melhores resultados ocorrem quando utilizamos 32 bits
para representação de dados (figura 6.4(c)).

(a) int8

(b) int16

(c) float

Figura 6.4: Erros que ocorrem utilizando diferentes formas de representação de valores.
A figura 6.5 mostra o raycasting de diferentes volumes utilizando a curvatura Gaussi-

49
ana para definir as cores. O verde indica curvatura Gaussiana nula.

(a) teardrop

(b) hyperbolic paraboloid

(c) cylinder

(d) cyclide

Figura 6.5: Visualização da curvatura Gaussiana em volumes diferentes.

(a) nı́vel 0.011

(b) nı́vel 0.015

(c) nı́vel 0.720

Figura 6.6: Diferentes nı́veis do volume cubo.
Para finalizar, a figura 6.6 representa a visualização de diferentes superfı́cies de nı́vel
do volume cubo onde as cores representam a curvatura Gaussiana. Observe a coerência
do resultado da figura 6.6(a) onde apenas valores positivos de curvatura aparecem mesmo
que muito próximos de zero.

Capı́tulo 7
Considerações Finais
Este trabalho fez um estudo dos problemas envolvendo a visualização de volumes utilizando informações de curvatura. Para tanto foi implementado um método de interpolação
tricúbica que utiliza funções B-splines para avaliar funções amostradas garantido que elas
sejam contı́nuas e duas vezes diferenciáveis. O cálculo de curvatura utilizou esses resultados para obter derivadas de segunda ordem e foi baseado no algoritmo descrito por
Kindlmann et. al. [9].
Um framework de visualização volumétrica chamado Voreen foi utilizado na implentação dos algoritmos. Para este sistema foram criados processadores aos quais demos os
nomes de CubicSplineProcessor, CurvatureRayCaster, CurvatureVolumeRaycaster e
FunctionSelector. O primeiro realiza a conversão dos valores de intensidade contidos
no volume de dados em coeficientes B-spline utilizados no processo de interpolação tricúbica. O processador CurvatureRayCaster possibilita a visualização do volume através
algoritmo de traçado de raios utilizando os valores de curvatura sobre a superfı́cie para a
definição de cores. O processador CurvatureVolumeRaycaster é semelhante ao anterior
exceto pelo fato de que possibilita a visualização diversas porções do volume através da
definição de valores de opacidade para as superfı́cies de nı́vel. Por fim, o processador FunctionSelector gera volumes a partir de funções predefinidas eliminando a necessidade de
armazenamento desses volumes em arquivos.
Os algoritmos de renderização foram implementados utilizando a linguagem GLSL e
fazendo uso do alto poder de processamento das placas gráficas para proporcionar bons
resultados em tempo real.
Como trabalhos futuros destacamos:
• A implementação dos cálculos de curvatura e a avaliação de funções utilizando a
linguagem CUDA que oferece ferramentas para programação paralela em GPU.
• A criação de módulos para o processamento das imagens obtidas como resultado da
visualização.
50

51
• A busca por boas aplicações que utilizem as informações obtidas com o cálculo de
curvatura.

Referências Bibliográficas
[1] Blinn, J. (1977), ‘Models of light reflection for computer synthesized pictures’, Computer Graphics 11(2), 192–198.
[2] Carmo, M. P. D. (2005), Geometria Diferencial de Curvas e Superfı́cies, SBM.
[3] Drebin, R. A., Carpenter, L. & Hanranhan, P. (1988), ‘Volume rendering’, Computer
Graphics 22(4), 65–74.
[4] FSF (1985), ‘Free software foundation’. URL http://www.fsf.org, último acesso
em fevereiro de 2010.
[5] Glassner, A. S. (2002), An Introduction to Ray Tracing, Morgan Kaufmann Publishers, Inc.
[6] GPGPU (2002), ‘Gpgpu: General-purpose computation on graphics hardware’. URL
http://gpgpu.org, último acesso em fevereiro de 2010.
[7] Hadwiger, M., Ljung, P., Salama, C. R. & Ropinski, T. (2009), Advanced Illumination
Techniques for GPU Volume Raycasting, Communications of the ACM.
[8] Hladuvka, J., Konig, A. & Gröller, E. (2000), ‘Curvature-based transfer functions for
direct volume rendering’, Spring Conference on Computer Graphics 2000.
[9] Kindlmann, G., Whitaker, R., Tasdizen, T. & Möller, T. (2003), ‘Curvature-based
transfer function for direct volume rendering: Method and aplications’, IEEE Visualization pp. 513–520.
[10] Levoy, M. (1990), ‘Efficient ray tracing of volume data’, ACM Transactions on
Graphics 9(3), 205–261.
[11] Lima, E. L. (2008), Curso de Análise, Vol. 2, Associação Instituto Nacional de Matemática Pura e Aplicada.
[12] Lorensen, W. E. & Cline, H. E. (1987), ‘Marching cubes: A high resolution 3d surface
construction algorithm’, ACM SIGGRAPH Computer Graphics 21(4), 163–169.

52

REFERÊNCIAS BIBLIOGRÁFICAS

53

[13] Marroquin, R. & Maximo, A. (2009), ‘Introdution to gpu programming with glsl’,
Sessão de Tutoriais do SIBGRAPI 2009.
[14] Mello, V. M. (1999), Representação multiescala de objetos implı́citos, Master’s thesis,
Universidade Federal da Bahia.
[15] Meyer-Spradow, J., Ropinski, T., Mensmann, J. & Hinrichs, K. (2009), ‘Voreen: A
rapid-prototyping environment for ray-casting-based volume visualizations’, IEEE
Computer Graphics and Applications 29(6), 6–13.
[16] Möller, T., Mueller, K., Machiraju, R. & Yagel, R. (1998), ‘Design of accurate and
smooth filters for function and derivative reconstruction’, Proceedings of the 1998
IEEE symposium on Volume visualization pp. 143 – 151.
[17] NVIDIA (2009), Cuda programming guide, Technical report, Nvidia Corporation.
URL http://www.nvidia.com/cuda.
[18] NVIDIA (2010), ‘Nvidia’. URL http://www.nvidia.com.br, último acesso em fevereiro de 2010.
[19] Paiva, A. C., Seixas, R. B. & Gattas, M. (1999), Introdução à Visualização Volumétrica, Pontifı́cia Universidade Católica do Rio de Janeiro.
[20] Phong, B. T. (1975), ‘Illumination for computer generated pictures’, Communications
of the ACM 18(6), 311–317.
[21] Rost, R. J. (2006), OpenGL Shading Language, Addison Wesley Professional.
[22] Ruijters, D., Romeny, B. & Suetens, P. (n.d.), ‘Acurracy of gpu-based b-spline evaluation’.
[23] Shreiner, D., Woo, M., Neider, J. & Davis, T. (2006), OpenGL Programming Guide,
Addison Wesley Professional.
[24] Sigg, C. & Hadwiger, M. (2005), ‘Fast third-order texture filtering’, GPU Gems 2
pp. 313–317.
[25] Stegmaier, S., Strengert, M., Klein, T. & Ertl, T. (2005), ‘A simple and flexible volume rendering framework for graphics-hardware-based raycaster’, Volume Graphics
pp. 187–241.
[26] Voreen (2009), ‘Vooren - volume rendering engine’. URL http://www.voreen.org,
último acesso em fevereiro de 2010.