зависимости для прицельных расстояний 2 и 0 и относительной скорости 0,1 в численном эксперименте получились сходного вида, но дали вероятности образования иона 83,7 и 84,28% соответственно. Таким образом, с уменьшением прицельного расстояния растет вероятность образования молекулярного иона. Видно также, что в области межчастичных расстояний, соответствующих молекулярному иону, форма функции распределения межчастичных расстояний сохраняется. В частности, растояние, соответствующее минимуму функции распределения, во всех случаях одно и то же и соответствует минимуму потенциальной энергии С/г-
СПИСОК ЛИТЕРАТУРЫ
1. О z hi gov Yu. I. Selection of quantum states in the simulation of quantum mechanics // Quantum Computers and Computing. 2007. 7. N 1. P. 88-110.
2. Ozhigov A., Arakelov K., Ozhigov Yu. Principles of the numerical simulation of many body quantum dynamics // Quantum Computers and Computing. 2006. 6. N 1. P. 137-148.
3. Блохинцев Д. В. Принципиальные вопросы квантовой теории // М.: Наука, 1977.
4. Ландау J1. Д. Ассоциация двухатомных молекул // Советский физический журнал. 1932. 2. С. 46-52.
5. Zener С. Molecular terms and chemical reactions // Proc. R. Soc. Ser. A. 1932. 137. P. 515-546.
6. Никитин E. E., Смирнов Б.М. Атомно-молекулярные процессы в задачах с решениями. М.: Наука, 1989.
Поступила в редакцию 12.02.08
УДК 519.688
А.Ю. Дижевский
ВИЗУАЛИЗАЦИЯ ТРЕХМЕРНЫХ ОБЪЕКТОВ, ВЛОЖЕННЫХ В ПОЛУПРОЗРАЧНЫЙ ОБЪЕМ
(кафедра вычислительной математики механико-математического факультета,
e-mail: [email protected])
В данной работе описаны основные алгоритмы объемной визуализации (volume rendering) и методы поверхностной визуализации (surface rendering) трехмерных объектов. Представлены методы гибридной визуализации:
• поверхности, заданной триангуляцией, вложенной в полупрозрачный объем;
• неявно заданной поверхности, вложенной в полупрозрачный объем.
Введение. Объемная и поверхностная визуализация применяются в медицине (визуализация снимков, полученных с помощью компьютерной томографии, ультразвука или магнетического резонанса), вычислительной динамике жидкости (аэродинамические симуляторы), метеорологии (прогнозы), астрофизике и других областях.
При поверхностной визуализации изображаемая модель представляется в виде поверхности (набор полигонов). При объемной визуализации модель задана поточечно в некоторой области. Алгоритмы поверхностной визуализации вычислительно более просты, чем алгоритмы объемной визуализации. В то же время объемное изображение хранит больше информации, поскольку включает информацию о всех деталях объекта.
Визуализация поверхностей, вложенных в полупрозрачный объем, играет важную роль как в медицине (хирургическое планирование, обучающие системы), так и во многих других областях исследований (например, при создании различных симуляторов).
Целью данной работы являются описание основных алгоритмов гибридной визуализации и разработка алгоритма, способного изображать поверхности, вложенные в полупрозрачный объем. Изображаемые в этом алгоритме поверхности задаются неявными функциями, т. е. уравнениями вида f(x,y,z) — к = 0. Здесь / — некоторая функция, принимающая вещественные значения; к — пороговое значение; физически f(x,y,z) обычно представляет собой плотность изучаемой модели в точке (.x,y,z). Изображаемая поверхность представляет собой поверхность уровня функции /.
Как правило, триангуляции в гибридных алгоритмах синтезированы или получены с помощью порогового значения. Поверхность, заданная пороговым значением, может быть отображена с помощью данного алгоритма без дополнительных операций построения ее триангуляции. Многие синтезированные данные также могут быть представлены неявными функциями.
1. Алгоритм объемной визуализации. В медицине для представления полной картины с учетом всех органов используется объемная визуализация. Изображения, созданные в результате объемной визуализация цифровых данных, похожи на классические рентгеновские снимки и поэтому привычны для докторов. При таком способе визуализации модель задается поточечно в некоторой области. Существует четыре типа алгоритмов объемной визуализации: испускание лучей [1-4], разложение на слои со сдвигающей деформацией [5], проектирование вокселей [6], отображение текстуры [7]. В данной работе дается подробное описание алгоритма испускания лучей.
Процесс построения изображения можно представить следующим образом: из некоторого источника испускаются лучи света с определенной яркостью параллельно заданному вектору. Лучи при пересечении с исследуемой областью (обычно куб или прямоугольный параллелепипед) теряют яркость. Интенсивность луча, падающего в заданную точку плоскости изображения, определяет яркость этой точки. Исследуемый объект представляется набором срезов (например, томографических снимков), а значения плотностей между срезами получаются с помощью линейной интерполяции по значениям на соседних срезах (возможна также кубическая интерполяция).
На рис. 1 представлена картина объемной визуализации при помощи данного алгоритма. Картина отображается на плоскости ii, которая перпендикулярна оси Z. Центр куба находится в начале координат. Точка наблюдения находится слева, туда же помещен и источник излучения (в точку с положительной координатой Z и нулевыми X, Y). Чтобы найти яркость какой-нибудь точки (жо,уо) на плоскости, проводится луч из этой точки в направлении вектора с координатами (0,0, —1) (по направлению взгляда обозревателя). При пересечении луча с кубом получается отрезок. В первой точке этого отрезка определяется плотность. Затем вычисляются плотности остальных точек отрезка, взятых с определенным шагом. Все плотности учитываются при вычислении яркости точки (жо,уо) с помощью определенной процедуры. Процедура может вычислять среднее, максимальное (MIP — Maximal Intensity Projection) или минимальное значение плотности на отрезке.
Когда точка отрезка попадает между срезами, определяется объемлющий куб с известными плотностями на его вершинах (A, i?, С, D, Ai, i?i, Ci, D\).
Точка находится между слоями I и I + 1. Допустим, расстояние между слоями равно точка расположена на расстоянии р от слоя следовательно, она на расстоянии к — р от слоя /+ 1. Сторона квадрата равна /. Таким образом, плотность точек будет вычислена по следующей формуле:
density0 = (b ■ А + (/ - Ь) ■ В) ■ а + (Ъ ■ D + (/ - Ь) ■ С) ■ (/ - а), density х = (b-A1 + (f-b)-B1)-a+(b-D1 + (f-b)- С г) ■ (/ - а),
где density0, density-l — плотности точек на слоях I + 1, полученных проекциями искомой точки. Плотность искомой точки получается по формуле density = density-l + (к—р) -density0, где density — плотность искомой точки.
2. Алгоритм поверхностной визуализации. Триангуляция трехмерных поверхностей является одним из удобных методов визуализации. Заданная область пространства (обычно прямоугольный параллелепипед), в которой требуется построить триангуляцию, разбивается на маленькие ячейки (кубические, тетраэдрические и др.). При пересечении поверхности и ячейки образуются треугольники. Более подробно про триангуляцию кубическими ячейками можно посмотреть в [8, 9], про триангуляцию тетраэдрическими ячейками — в [10-12]. Полученный набор треугольников для всех ячеек представляет собой требуемую триангуляцию.
Поверхность задается функцией /(ж, у, ¿) — к = 0 (/ — функция плотности, к — пороговая величина). Значения функции вычисляются только в вершинах ячейки. Точки пересечения вычисляются с помощью линейной интерполяции. Вершины треугольников хранятся в одном массиве, а треугольники — в другом как ссылки на эти вершины. Треугольник представляется набором из трех индексов в массиве вершин. Для визуализации триангулированной поверхности вычисляется расстояние до каждого треугольника. Рисование происходит с помощью Я-буфера для непрозрачных поверхностей и начиная с самого дальнего к точке наблюдения для полупрозрачных поверхностей. При рисовании также учитывается вектор нормали в каждой вершине треугольника.
3. Визуализация триангуляции, вложенной в полупрозрачный объем. Обычно вложение триангуляции или нескольких триангуляций в полупрозрачный объем возникает в задачах визуализации медицинских данных. Триангуляции синтезированы или получены с помощью порогового значения. Поэтому можно предположить, что в каждой вершине треугольника триангуляции известна нормаль (вычисляемая как градиент функции при использовании порогового значения).
Можно предложить два способа визуализации набора треугольников, вложенных в полупрозрачный объем:
1) вокселизация треугольников и визуализация одним из существующих алгоритмов объемной визуализации;
2) непосредственное вычисление точки пересечения луча и треугольника в процессе построения объемного изображения методом испускания лучей. При этом для уменьшения количества рассматриваемых треугольников применяется алгоритм двумерного "хеширования" на матрице изображения.
3.1. Вокселизация триангуляции. Алгоритм объемной визуализации рисует изображение с учетом всех внутренних вокселей (трехмерных пикселей) в объекте. Вокселизация треугольника из триангуляции подразумевает поиск вокселей, покрывающих треугольник, с применением одного из критериев покрытия [13]. Существует несколько критериев покрытости вокселя: воксель покрыт, если он покрыт шестью (по общим граням), восемнадцатью (по граням и ребрам) или двадцатью шестью (по грани, ребрам и вершинам) вокселями.
Треугольник покрывается вокселями следующим образом: для вершин треугольника строим ограничивающие сферы радиуса Д, ребра покрываем цилиндрами длиной, равной длинам ребер, и радиуса Д, для внутренностей треугольника строим две ограничивающие плоскости на расстоянии Д от треугольника и три плоскости, перпендикулярные треугольнику и пересекающие стороны треугольника. Для критерия "6 покрывающих вокселей" Д = 1./2. где Ь — диагональ вокселя. Для критерия "26 покрывающих вокселей Д = тД-Ь/2 (с м. [13]). Все воксели, полностью покрытые сферами и цилиндрами либо ограниченные пятью плоскостями, и будут покрывающими вокселями для треугольника.
Для визуализации с учетом отражающих свойств поверхности необходим поиск нормалей в векселях поверхности. При известных значениях градиента поверхности на вершинах треугольника вектор нормали каждого вокселя треугольника можно вычислять как среднее арифметическое градиентов в вершинах треугольника. Хранение нормалей можно организовать в массиве.
Данный процесс строится для каждого треугольника триангуляции. С целью уменьшения количества вычислительных операций для треугольников, имеющих общее ребро или вершину, строим общую покрывающую сферу или цилиндр для смежной вершины или ребра. Для этого за линейное время надо обойти все треугольники и с каждой вершиной сопоставить список треугольников, у которых одна из вершин совпадает с ней. После этого для ребер достаточно обойти треугольники, у которых две вершины совпадают с вершинами ребер, т. е. сопоставляются списки треугольников для обеих вершин ребра. Совпадающие треугольники имеют общее ребро.
После окончания вокселизации поверхности надо сделать сопоставление или наложение значений в вокселях на поверхности и в объеме, представленном вокселями. Как правило, значения в векселях объема ограничены, например, для томографических снимков 5000. Наложение значений можно проводить следующим образом: для всех общих вокселей поверхности и объема к значениям объема прибавляем пороговое т = 5000. При визуализации идет проверка каждого вокселя на превышение значения т. Если значение выше т, то в вокселе хранится поверхностное и объемное значение, и для получения чистого значения объемного вокселя производится вычитание значения т.
3.2. Алгоритм двумерного хеширования. Прямоугольная плоская область //. на которой строится изображение, разбивается на квадратные ячейки заданного размера. Каждой ячейке сопоставляется список треугольников, проекции которых могут пересекать ячейку. Рассматриваем последовательно все треугольники. Каждый треугольник проецируем на плоскость и для его проекции строим объемлющий ее прямоугольник со сторонами, параллельными сетке. Треугольник добавляется в списки тех ячеек, которые пересекаются с построенным прямоугольником. Матрица списков номеров треугольников, индексами в которой являются координаты ячеек сетки, представляет собой двумерную хеш-таблицу.
Списки треугольников могут быть отсортированы по расстоянию от Н до центров треугольников, если используется полупрозрачный режим рисования. Данный режим используется в следующих работах [14-17].
Алгоритм нахождения оптимального размера ячейки можно представить следующим образом. Каждый спроецированный треугольник окружается прямоугольником наименьшей диагонали. Среднее арифметическое диагоналей всех прямоугольников и будет искомой диагональю ячейки.
Испускается луч из области Н в направлении к объему. На луч попадают воксели объема и треугольники из ячейки, содержащей данный луч. Для вокселей, лежащих на луче, вычисляется яркость как максимальное или среднее значение плотностей всех вокселей. При прохождении первого треугольника процесс останавливается. Итоговая плотность пикселя в Н состоит из двух компонент — плотности объема и плотности треугольника. Процесс рисования вокселя аналогичен процессу рисования неявной поверхности, вложенной в полупрозрачный объем. Нормаль для треугольника будет равна среднему арифметическому нормалей его вершин.
Для определения точки пересечения луча с поверхностью, представленной триангуляцией, рассматриваются все треугольники из списка, соответствующего квадратной ячейке, в которую попадает луч. Для каждого такого треугольника вычисляется точка пересечения с лучом. В случае непрозрачной поверхности находится ближайшая к плоскости изображения точка пересечения, и луч на этом обрывается. В случае полупрозрачной поверхности луч ослабевает при переходе каждой точки пересечения пропорционально степени непрозрачности поверхности. Основным отличием представленного алгоритма от алгоритмов гибридного рендеринга из работ [14, 15, 17] является применение хеш-таблицы для треугольников. В [14, 16, 17] треугольники группируются в плоские "плиты" для каждого шага по лучам, перпендикулярным плоскости наблюдателя.
4. Алгоритм визуализации неявной поверхности, вложенной в полупрозрачный объем.
В данной работе предложен алгоритм, визуализирующий неявную поверхность, вложенную в полупрозрачную среду. Такой тип задач возникает в задачах хирургического планирования и в обучающих медицинских системах, когда поверхности заданы как неявные функции. Обычно полупрозрачные среды визуализируются с помощью объемного рендеринга. Поверхность обладает свойством отражения света, в свою очередь отраженный свет, проходя через среду, теряет интенсивность пропорционально ее плотности. Модель можно представить как рассеивающую свет поверхность, вложеную в полупрозрачный объем. При этом считается, что свет, падающий на поверхность, не теряет яркости при прохождении объема (т. е. проходит сквозь среду, как будто она абсолютно прозрачна). Среда частично поглощает лишь свет, отраженный поверхностью.
Поскольку поглощение среды влияет на яркость рассеиваемого поверхностью света, то при ее визуализации одним и тем же цветом может быть потеряна информация либо об объеме, либо о поверхности. Эту проблему можно решить путем смешивания независимых цветов с учетом неизменности яркости до и после помещения поверхности в объем. Например, поверхность отражает красный свет, а среда имеет зеленый цвет.
Перейдем к описанию алгоритма. Как и в алгоритме объемной визуализации, есть плоскость H, на которой изображается вся картина.
Идет поточечный процесс на каждом испущенном луче, как и в обычном алгоритме построения объемного изображения. Единственное отличие состоит в том, что процесс останавливается, как только точка на луче достигает поверхности. Предполагается, что поверхность непрозрачная. Создаются две матрицы изображения, содержащие яркость точек в диапазоне от 0 до 1: первая для объема, вторая для поверхности. Заполнение первой матрицы происходит так же, как и в алгоритме объемной визуализации.
Для заполнения второй матрицы вычисляются значения функции, задающей поверхность, для двух соседних точек на луче. Если произведение этих значений неположительно, то точка пересечения луча с поверхностью находится между этими точками. Точка пересечения вычисляется по следующей формуле: point = pointy + (point2 — point-J • (а/(a + ft)), где point — искомая точка; pointy point2 — соседние точки на луче.
Значение яркости в точке point на поверхности вычисляется следующим образом: value = coeffx х (grad(point), г>), где coeff— заданный коэффициент, grad(point) — градиент поверхности в данной точке, V — единичный вектор направления луча, (ж, у) — евклидово скалярное произведение. Если луч не пересекает поверхность, то значение в точке матрицы изображения будет равно нулю. Чтобы яркость точки на плоскости H была неизменной для поверхности и после изъятия из среды (обнуление первой матрицы изображения), значения во второй матрице изображения вычисляются по формуле
bitmapSurfacefi] = value — transparency • bitmapVolumefi],
где transparency — коэффициент прозрачности объема, bitmapSurface и bitmap Volume — матрицы изображения для поверхности и объема соответственно, value получено выше. Такое вычисление яркости цвета поверхности соответствует физической модели. Полупрозрачный объем зеленого цвета будет затенять поверхность красного цвета. Может случиться и так, что поверхность не будет видна, например при низком коэффициенте прозрачности объема.
5. Результаты программной реализации
5.1. Общее понятие о входных данных. Данные, полученные с томографа, можно представить набором матриц плотностей каждого среза. Матрицы имеют размерность 512x512. Используемые программой плотности пробегают значения от —1024 для воздуха до 4000-5000 наиболее твердых костей и зубов. По полученным матрицам требуется восстановить объемное изображение. Срезы получены в параллелепипеде, в общем случае наклоненном к плоскости изображения.
5.2. Метод построения объемного изображения. Программно реализованы два способа вычисления плотности. Первый способ вычисляет максимальное значение плотности при пересечении объема лучами (MIP). Ниже представлены некоторые изображения, полученные этим способом. Рисование поточечное без поддержки OpenGL (рис. 2).
Рис. 2. Объемный рендеринг с помощью М1Р: а — изображение грудной клетки пациента; б — изображение легких, полученное с помощью удаления пикселов вне контуров легких
а б
Рис. 3. Поверхностная визуализация: а — поверхность непрозрачна; б — поверхность полупрозрачна
5.3. Визуализация поверхности с помощью триангуляции. На рис. 3 справа изображена поверхность для порога 1150 со вставленными цилиндрическими штифтами. Рисование происходит с помощью библиотеки OpenGL. Поверхность полупрозрачна и рисуется по треугольникам, начиная с самого дальнего. Полупрозрачность получена функциями
glEnable(GL_BLEND),
glBlendFunc(GLJ3RC .ALPHA, GL_ONE_MINUS_SRC_ALPHA).
На рис. 3 слева изображена эта же поверхность без полупрозрачности. Используется Z-буфер, поддерживаемый функциями OpenGL:
glEnable(GL_DEPTH_TEST), glDepthFunc(GL_LEQUAL).
Рис. 4. Визуализация поверхности, вложенной в объем с помощью
MIP
5.4. Метод визуализации объекта в полупрозрачной среде. На рис. 4 показаны изображения поверхности, вложенной в полупрозрачный объем. Значения точек матрицы изображения для среды вычисляются с помощью алгоритма испускания лучей, где берется максимальное значение по лучу до поверхности объекта. Поверхностью является эллипсоид. Среда показана зеленым цветом, поверхность — красным.
6. Заключение. В данной работе описаны и программно реализованы алгоритмы построения объемного изображения для максимального и среднего значения по лучам. Описан алгоритм визуализации поверхности с помощью триангуляции. Рассмотрены алгоритмы визуализации триангуляции, вложенной в полупрозрачный объем. Предложен способ ускорения визуализации посредством использования двумерного "хеширования". Разработан алгоритм наложения полупрозрачного объема
и поверхности, заданной уравнением f(x,y,z) = 0. Данный алгоритм ранее не был представлен в
литературе.
СПИСОК ЛИТЕРАТУРЫ
1. Levoy M. Display of surfaces from volume data // IEEE Comp. Graph, and Appl. 1988. 8. N 3. P. 29-37.
2. Levoy M. Efficient ray tracing of volume data // ACM Trans. Comp. Graph. 1990. 9. N 3. P. 245-261.
3. Hoehne K.H., Pfiesser В., Pommert A., Riemer M., Schiemann T., Schubert R., Tiede U. A virtual body model for surgical education and rehearsal // IEEE Comp. 1996. 29. N 1. P. 25-31.
4. Bhalerao A., Pfister H., Halle M., Kikinis R. Fast re-rendering of volume and surface graphics by depth, color, and opacity buffering // Med. Image Analysis. 1999. 4. N 3. P. 235-251.
5. Lacroute P., Levoy M. Fast volume rendering using a shear-warp factorization of the viewing transformation // Proc. SIGGRAPH. Orlando, USA: ACM, 1994. P. 451-458.
6. Westover L. Footprint evaluation for volume rendering // Proc. of the 17th Annual Conf. on Comp. Graph, and Inter. Tech. Dallas, USA: ACM, 1990. P. 367-376.
7. Cabrai В., Cam N., Foran J. Accelerated volume rendering and tomographic reconstruction using texture mapping hardware // Proc. ACM/IEEE Symp. on Vol. Vis. Tysons. Corner, USA: ACM, 1994. P. 91-98.
8. Lorensen W. E., С line H. Marching cubes: A high resolution 3D surface construction algorithm // Proc. of the 14th ACM Siggraph Annual Conf. on Comp. Graph, and Inter. Tech. N.Y., USA: ACM, 21. N 4. 1987. P. 163-169.
9. Nat araj an В. K. On generating topologically consistent isosurfaces from uniform samples // The Vis. Сотр.: Intern. J. of Сотр. Graph. 1994. N 11. P. 52-62.
10. G uez iec A., H ummel R. Exploiting triangulated surface extraction using tetrahedral decomposition // IEEE Trans, on Vis. and Сотр. Graph. 1995. 1. N 4. P. 328-342.
11. Skala V. Precision of iso-surface extraction from volume data and visualization // Algoritmy'2000 Int. Conf. on Sei. Сотр. Slovakia: AMUC, 2000. P. 368-378.
12. Carneiro B.P., Kaufman A. E. Tetra-cubes: an algorithm to generate 3D isosurfaces based upon tetrahe-dra // Proc. of SIGGRAPH. New Orleans, USA: ACM, 1996. P. 205-210.
13. Huang J., Yagel R., Filippov V., Kurzion Y. An accurate method for voxelizing polygon meshes // Proc. IEEE Sym. on Vol. Vis. Research Triangle Park, USA: ACM, 1998. P. 119-126.
14. Kreeger K., Kaufman A. Mixing translucent polygons with volumes // Proc. of the Conf. on Vis. California, USA: IEEE Сотр. Soc. Press, 1999. P. 191-198.
15. Tietjen C., Isenberg T., Preim В. Combining silhouettes, surface, and volume rendering for surgery education and planning // EUROGRAPHICS — IEEE VGTC Sym. on Vis. Aire-la-Ville, Switzerland: Eurographics Association, 2005. P. 303-310.
16. Levoy M. A hybrid ray tracer for rendering polygon and volume data // IEEE Сотр. Graph, and Appl. 1990. 10. N 3. P. 33-40.
17. Kaufman A., Yagel R., Cohen D. Intermixing surface and volume rendering // 3D-Imaging in Medicine: Algorithms, Systems, Applications. Berlin, Germany: Springer-Verlag, 1990. P. 217-227.
Поступила в редакцию 06.12.07