способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени

Формула изобретения

1. Способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени, включающий фрагментацию всего объема скалярных трехмерных данных воксельного массива на множество суб-объемов, определение подмножества суб-объемов, расположенных вдоль луча наблюдения, определение цвета пикселя как максимальное значение интенсивности из подмножества суб-объемов, принадлежащих лучу наблюдения, отличающийся тем, что суб-объемом является воксель, воксельный массив приводят к виду, при котором длина ребра куба вокселя равна длине стороны квадрата пикселя, из совокупности которых в дальнейшем формируют изображение, возможный диапазон значений интенсивностей вокселей воксельного массива берут равным количеству элементов палитры цветов, используемых для изображения пиксельного массива, длину каждого ребра воксельного массива приводят к размеру 2ⁿ в вокселях, вначале визуализируют границы воксельного массива в виде ребер прямоугольного параллелепипеда толщиной в 1 пиксель по 8 наборам координат, используя алгоритм Брезенхэма, при этом цвет пикселя ребра прямоугольного параллелепипеда определяют как максимально возможное значение интенсивности вокселя, умноженное на весовой коэффициент, который обратно пропорционален расстоянию от вокселя до плоскости построения изображения вдоль луча наблюдения, и имеет значение от 0 до 1, затем, используя массив значений трехмерных координат вокселей, принадлежащих ребрам прямоугольного параллелепипеда, определяют лучи наблюдения, пересекающие воксельный массив, далее, используя алгоритм Брезенхэма в плоскостях, которым принадлежат лучи наблюдения, определяют координаты вокселей, являющихся точками пересечения соответствующего луча наблюдения и стороны прямоугольного параллелепипеда, ограничивающего воксельный массив, затем посредством обратного преобразования угловых координат для каждого из этих вокселей получают координаты входа и выхода луча наблюдения в воксельный массив, далее вычисляют приращение каждой из координат, соответствующее «единичному шагу» вдоль луча наблюдения, после чего, последовательно сравнивая значения интенсивностей вокселей, умноженные на весовой коэффициент, вдоль луча наблюдения, определяют цвет пикселя как максимальное из этих значений, однако выборку этих значений производят по направлению от наблюдателя «в глубину» воксельного массива, при этом сравнивают текущее значение пикселя с максимально возможным значением вокселя, умноженным на весовой коэффициент, и если текущее значение пикселя больше или равно вычисленному значению, то выборку для данного луча завершают, к тому же воксельный массив визуализируют с возможностью произвольного изменения ракурса наблюдения путем изменения значений углов Эйлера - угла чистого вращения и угла нутации, при этом при нулевых значениях этих углов ось вращения проходит параллельно плоскости построения изображения и вертикально относительно наблюдателя через центр вращения - центр воксельного массива, к тому же визуализируют воксельный массив с возможностью получения произвольного среза.

2. Способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени по п.1, отличающийся тем, что приведение длины каждого ребра воксельного массива к размеру 2ⁿ в вокселях осуществляют методом поиска максимального значения объединяемых данных.

3. Способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени по п.1, отличающийся тем, что приведение длины каждого ребра воксельного массива к размеру 2ⁿ в вокселях осуществляют методом усреднения.

4. Способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени по п.1, отличающийся тем, что приведение длины каждого ребра воксельного массива к размеру 2ⁿ в вокселях осуществляют методом адаптивной фильтрации.

5. Способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени по п.1, отличающийся тем, что приведение длины каждого ребра воксельного массива к размеру 2ⁿ в вокселях осуществляют методом сплайн-интерполяции.

6. Способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени по п.1, отличающийся тем, что приведение длины каждого ребра воксельного массива к размеру 2ⁿ в вокселях осуществляют методом пространственной спектральной передискритизации.

7. Способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени по п.1, отличающийся тем, что строят «параллельную» проекцию, при которой лучи наблюдения ортогональны к плоскости построения изображения.

8. Способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени по п.1, отличающийся тем, что строят «перспективную» проекцию, при которой лучи наблюдения сходятся перед плоскостью построения изображения в предполагаемой точке расположения наблюдателя.

9. Способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени по п.1, отличающийся тем, что произвольный срез воксельного массива визуализируют в дополнительном окне в виде плоского среза в B-режиме.

10. Способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени по п.1, отличающийся тем, что формируют стереоизображение путем расчета не одного, а двух изображений, составляющих стереопару.

Описание изобретения к патенту

Изобретение относится к обработке данных изображения, а именно к визуализации трехмерного массива данных, и может быть использовано для решения задачи визуализации внутренней структуры «мутных» сред.

Несмотря на огромный прогресс в развитии алгоритмов и аппаратной поддержки компьютерной визуализации и анимации, рассматриваемая задача визуализации трехмерных данных является сложной для большинства современных средств, как программных, так и аппаратных. Это связано с тем, что «классическое» трехмерное моделирование и визуализация, в основном, рассматривают задачу визуализации поверхностей произвольной формы с различным «цветовым» заполнением и с различными отражающими способностями. Это удовлетворяет потребностям большинства прикладных задач. Для решения этих задач визуализируемые данные представляются в виде различных «векторных» объектов, таких как элементы поверхности, источники освещения, камера наблюдения и так далее. Разбиение визуализируемых данных на ограниченное количество объектов с простой геометрией сильно сокращает требуемое количество вычислений. Кроме того, развитие аппаратной поддержки именно этих вычислений и соответствующая оптимизация архитектуры современных графических процессоров (GPU) позволяет успешно решать эти задачи. Однако данные, получаемые при спектральной оптической когерентной томографии, используемые, в частности, для решения задачи визуализации внутренней структуры "мутных" сред, как правило, плохо подходят для обработки «классическими» методами векторной графики. Аппаратной поддержки визуализации воксельных массивов пока тоже нет. Эта задача приобрела в последнее время значительную актуальность в оптической когерентной томографии и в оптической когерентной голографии.

Известен способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных по статье «А fast maximum-intensity projection algorithm for generating magnetic resonance angiograms», IEEE TRANSACTIONS ON MEDICAL IMAGING, Vol.12, No. 1, pages 50-57, March 1993, авторов Steven Schreiner, Robert L. Galloway.

Данный способ используют, в частности, для построения магнитно-резонансных (MP) ангиографических изображений. Кровоток в интересующем объеме представляется в виде распределения яркости MP данных. Применяя алгоритмы MIP (метода проекции максимальной интенсивности), производят поиск максимальной интенсивности вдоль параллельных лучей, пересекающих объем MP данных. Проецируемое изображение формируется на основе этих максимальных значений интенсивностей. Таким образом, поток в сосудистой системе виден в проецируемой плоскости. Данный способ направлен на уменьшение времени, необходимого для создания проекций, что позволяет сделать технологию более доступной для пользователей систем малой вычислительной мощности. Каждый срез изображения предварительно сортируется в виде совокупностей интенсивностей. Благодаря соотнесению интенсивности пикселей с их положением на плоскости проекции уменьшается общее число пикселей, необходимое для построения плоскости проекции, и, таким образом, экономится время вычислений. Только наиболее яркие интенсивности относятся к кровотоку и используются в проекции. Дополнительная экономия времени расчетов получается в результате предварительного расчета "шаблонов" проекций и в то же время в результате заполнения сразу нескольких проекционных плоскостей. Алгоритм написан на языке Си на базе систем 80386. Этот алгоритм демонстрирует более чем шестикратное увеличение скорости вычислений относительно других подобных алгоритмов.

Известен способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных по статье «А fast hierarchical algorithm of maximum intensity projection», Image Processing, 2000. Proceedings. 2000 International Conference on. Vol.1. Pages 745-748. Sep 2000, авторов Keun Ho Kim, Min Jeong Kwon and Hyun Wook Park.

Авторами реализован иерархический метод проекции максимальной интенсивности (MIP) и проанализирована сложность вычислений с помощью данного метода. Трудно реализовать метод MIP при получении трехмерных медицинских изображений в реальном времени, так как время обработки зависит от количества вокселей трехмерных данных. Для того чтобы преодолеть трудности, связанные с большим количеством вычислений для MIP, авторами был разработан иерархический метод MIP, который позволяет на первом этапе выполнить процедуру MIP с субдискретизированным изображением. На втором этапе такое MIP изображение низкого разрешения может использоваться для экономии времени вычислений MIP изображения с полным разрешением.

Недостатком описанных выше способов является то, что рассматриваются варианты реализации MIP без использования эффекта «затемнения изображения в глубину» (немодифицированный MIP). К тому же рассматривается вариант реализации MIP только с прямой параллельной проекцией. Тогда как для адекватного восприятия больших изображений желательна не параллельная, а перспективная проекция. Кроме того, применение прямой проекции приемлемо только в случае, когда размер вокселей исходных данных много меньше размера пикселей на плоскости проекции. В противном случае при некоторых ракурсах возможны артефакты, которые принято называть артефактами типа «шахматная доска». Чтобы избежать этого, необходимо применять не прямую, а обратную проекцию.

К тому же описанные выше аналоги рассматривают варианты реализации MIP, при которых рассчитываются только ракурсы вращения визуализируемых данных вокруг центральной вертикальной оси и не предполагается обеспечения интерактивного просмотра пользователем проекций других произвольных ракурсов.

Ближайшим аналогом разработанного способа является способ проекции максимальной интенсивности для визуализации скалярных трехмерных данных, известный по патенту US 8149238, МПК G06T 15/50 от 03.04.2012 г.

Сначала осуществляют фрагментацию всего объема трехмерных данных на множество суб-объемов, определяют подмножества суб-объемов, расположенных вдоль луча наблюдения. Цвет пикселя определяют как максимальное значение интенсивности из подмножества суб-объемов, принадлежащих лучу наблюдения.

Недостатком данного способа является то, что при построении изображения не применяют «эффект затемнения изображения в глубину», т.е. не учитывают весовые коэффициенты, характеризующие расстояние от вокселя до плоскости проекции, что не дает хорошее зрительное восприятие изображения. К тому же данный способ проекции максимальной интенсивности для визуализации трехмерных данных для обеспечения интерактивности просмотра пользователем получаемых изображений требует большой вычислительной мощности, доступной, в настоящее время, только в компьютерных кластерных системах.

Задачей, на решение которой направлено данное изобретение, является разработка способа оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени, который позволяет улучшить восприятие объемности получаемого изображения, определить ракурс визуализируемого изображения, а также повысить скорость вычислений за счет уменьшения объема оперативной памяти, затрачиваемой на построение изображения, что обеспечит интерактивность работы с получаемым изображением, т.е. даст возможность изменения ракурса изображения пользователем в реальном времени.

Указанный технический результат достигается благодаря тому, что разработанный способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени так же, как и способ, который является ближайшим аналогом, включает фрагментацию всего объема скалярных трехмерных данных воксельного массива на множество суб-объемов, определение подмножества суб-объемов, расположенных вдоль луча наблюдения, определение цвета пикселя как максимальное значение интенсивности из подмножества суб-объемов, принадлежащих лучу наблюдения.

Новым в разработанном способе оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени, является то, что суб-объемом является воксель, воксельный массив приводят к виду, при котором длина ребра куба вокселя равна длине стороны квадрата пикселя, из совокупности которых в дальнейшем формируют изображение, возможный диапазон значений интенсивностей вокселей воксельного массива берут равным количеству элементов палитры цветов, используемых для изображения пиксельного массива, длину каждого ребра воксельного массива приводят к размеру 2ⁿ в вокселях, вначале визуализируют границы воксельного массива в виде ребер прямоугольного параллелепипеда толщиной в 1 пиксель по 8 наборам координат, используя алгоритм Брезенхэма, при этом цвет пикселя ребра прямоугольного параллелепипеда определяют как максимально возможное значение интенсивности вокселя, умноженное на весовой коэффициент, который обратно пропорциональней расстоянию от вокселя до плоскости построения изображения вдоль луча наблюдения, и имеет значение от 0 до 1, затем, используя массив значений трехмерных координат вокселей, принадлежащих ребрам прямоугольного параллелепипеда, определяют лучи наблюдения, пересекающие воксельный массив, далее, используя алгоритм Брезенхэма в плоскостях, которым принадлежат лучи наблюдения, определяют координаты вокселей, являющихся точками пересечения соответствующего луча наблюдения и стороны прямоугольного параллелепипеда, ограничивающего воксельный массив, затем посредством обратного преобразования угловых координат для каждого из этих вокселей получают координаты входа и выхода луча наблюдения в воксельный массив, далее вычисляют приращение каждой из координат, соответствующее «единичному шагу» вдоль луча наблюдения, после чего, последовательно сравнивая значения интенсивностей вокселей, умноженные на весовой коэффициент, вдоль луча наблюдения, определяют цвет пикселя как максимальное из этих значений, однако выборку этих значений производят по направлению от наблюдателя «в глубину» воксельного массива, при этом сравнивают текущее значение пикселя с максимально возможным значением вокселя, умноженным на весовой коэффициент, и если текущее значение пикселя больше или равно вычисленному значению, то выборку для данного луча завершают, к тому же воксельный массив визуализируют с возможностью произвольного изменения ракурса наблюдения путем изменения значений углов Эйлера - угла чистого вращения и угла нутации, при этом при нулевых значениях этих углов ось вращения проходит параллельно плоскости построения изображения и вертикально относительно наблюдателя через центр вращения - центр воксельного массива, к тому же визуализируют воксельный массив с возможностью получения произвольного среза.

В первом частном случае реализации разработанного способа оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени, приведение длины каждого ребра воксельного массива к размеру 2ⁿ в вокселях осуществляют методом поиска максимального значения объединяемых данных.

Во втором частном случае реализации разработанного способа, приведение длины каждого ребра воксельного массива к размеру 2ⁿ в вокселях осуществляют методом усреднения.

В третьем частном случае реализации разработанного способа, приведение длины каждого ребра воксельного массива к размеру 2ⁿ в вокселях осуществляют методом адаптивной фильтрации.

В четвертом частном случае реализации разработанного способа, приведение длины каждого ребра воксельного массива к размеру 2ⁿ в вокселях осуществляют методом сплайн-интерполяции.

В пятом частном случае реализации разработанного способа, приведение длины каждого ребра воксельного массива к размеру 2ⁿ в вокселях осуществляют методом пространственной спектральной передискритизации.

В шестом частном случае реализации разработанного способа, строят «параллельную» проекцию, при которой лучи наблюдения ортогональны к плоскости построения изображения.

В седьмом частном случае реализации разработанного способа, строят «перспективную» проекцию, при которой лучи наблюдения сходятся перед плоскостью построения изображения в предполагаемой точке расположения наблюдателя.

В восьмом частном случае реализации разработанного способа, произвольный срез воксельного массива визуализируют в дополнительном окне в виде плоского среза в B-режиме.

В девятом частном случае реализации разработанного способа формируют стереоизображение путем расчета не одного, а двух изображений, составляющих стереопару.

На фиг.1 представлено положение визуализируемого воксельного массива (вид сбоку) относительно плоскости экрана при повороте на угол .

На фиг.2 представлено положение оси вращения относительно границ визуализируемого воксельного массива и направление вращения на угол .

На фиг.3 представлено трехмерное ОКТ-изображение (томограмма участка кожи пальца), как пример применения разработанного способа оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени.

На фиг.4 представлены трехмерное ОКТ-изображение с произвольным плоским срезом (a) и произвольный плоский срез в B-режиме (б).

На фиг.5 представлено трехмерное ОКТ-изображение как способ визуализации внутренних структур массива.

На фиг 6 представлено анаглифное стереоизображение (для просмотра необходимо использовать анаглифные очки).

Разработанный способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени предназначен для визуализации воксельного массива. При этом результатом вычислений изображения произвольного ракурса наблюдения является пиксельный массив. Для оптимизации вычислений вводятся следующие ограничения на исходные данные. Визуализируемый воксельный массив и получаемый пиксельный массив имеют одинаковый и изотропный пространственный масштаб. То есть каждый воксель имеет кубическую форму, а каждый пиксель - квадратную форму. При этом длина ребра куба вокселя равна длине стороны квадрата пикселя. Кроме того, диапазон значений визуализируемых данных равен количеству элементов палитры цветов, используемых для изображения пиксельного массива. Например, если визуализируемые данные предполагается «рисовать» 256-ю градациями серого цвета, то значения каждого элемента исходных данных (вокселя) должны быть целочисленными и находиться в пределах от 0 до 255.

В общем случае визуализируемый объем будет иметь форму прямоугольного параллелепипеда. Для дальнейшей оптимизации вводится следующее ограничение: длина каждого ребра визуализируемого прямоугольного параллелепипеда в вокселях должна быть 2ⁿ. Если размеры визуализируемого массива будут равны степени двойки, то преобразования трехмерных координат в линейный индекс исходного массива, требующие две операции умножения и две операции сложения, можно легко заменить комбинацией сдвиговых операций и операций побитового сложения. А это на 30-40% ускоряет процедуру этого преобразования и, таким образом, существенно влияет на общую скорость вычислений. Если исходные данные не соответствуют этим ограничениям, то производятся предварительные вычисления, которые привели бы исходные данные к требуемому виду. Так как для визуализации используется метод проекции максимальной интенсивности, в этих расчетах проще использовать метод поиска максимального значения объединяемых данных. Однако, если этого требует специфика задачи, можно применить методы усреднения, фильтрации, сплайн-интерполяции или пространственной спектральной передискритизации.

Для предоставления пользователю возможности произвольного изменения ракурса наблюдения визуализируемых данных достаточно, чтобы пользователь имел возможность изменять две угловые координаты. В терминах углов Эйлера одна из этих угловых координат будет соответствовать углу чистого вращения, а другая угловая координата будет соответствовать углу нутации. При этом считаем, что при нулевых значениях этих углов ось вращения проходит параллельно плоскости экрана визуализации и вертикально относительно наблюдателя через центр визуализируемого объема, а центр вращения находится в центре визуализируемого объема. При изменении угла нутации ось вращения наклоняется в сторону наблюдателя. Таким образом, пользователь имеет возможность при помощи клавиатуры, мыши или других средств управления произвольно или последовательно изменять значения этих углов. Назовем текущее значение этих углов текущим ракурсом наблюдения.

Для того чтобы пользователь ясно представлял текущее положение ракурса наблюдения, а так же для того, чтобы пользователь понимал, как и в какую сторону происходит текущее последовательное изменение ракурса наблюдения, предлагается в обязательном порядке визуализировать границы воксельного массива. Эти границы предлагается изображать в виде непрерывных отрезков, которые являются ребрами соответствующего прямоугольного параллелепипеда. Кроме того, эти отрезки должны быть толщиной в один пиксель и должны отображаться теми же цветами, что и визуализируемый массив. То есть согласно методам модифицированной MIP предполагается, что точки этих отрезков изначально имеют максимальную интенсивность и изменяют свой цвет в соответствии с установленным весовым коэффициентом «затемнения в глубину». Точки этих отрезков при визуализации согласно методам MIP должны «скрываться» более яркими вокселями визуализируемого массива и «скрывать» менее яркие. Наличие такой визуализации границ визуализируемого массива даже при наблюдении «статичной» картинки дает пользователю лучшее представление о пространственной суперпозиции элементов визуализируемого объема.

Для реализации метода MIP в нашем случае необходимо использовать только обратную проекцию, так как пространственный масштаб вокселей и пикселей одинаковый и при применении прямой проекции для некоторых ракурсов могут возникать артефакты типа «шахматная доска». Для описания метода введем понятие лучей наблюдения. Это прямые, которые пересекают экран, на котором мы строим изображение, в центре каждого пикселя этого экрана. В воображаемом полупространстве за экраном они пересекают визуализируемый массив. В рамках метода MIP можно использовать как «параллельную», так и «перспективную» проекцию. Если используется «параллельная» проекция, то лучи наблюдения ортогональны к плоскости экрана, а если «перспективная», то в полупространстве перед экраном они сходятся в предполагаемой точке расположения наблюдателя. Таким образом, значение пикселя на экране определяется значениями вокселей, которые пересекаются соответствующим лучом наблюдения.

Необходимость использования «перспективной» или «параллельной» проекции определяется субъективными угловыми размерами получаемого изображения. Например, если визуализируемый массив отображается на стандартном компьютерном мониторе с разрешением 96 dpi и имеет кубическую форму с длиной ребра в 128 точек, то использование «перспективной» проекции необязательно. Если длина составляет 256 точек, то желательно использовать «перспективную» проекцию. Если же эта длина 512 точек и более, то использование «перспективной» проекции обязательно для правильного восприятия пользователем «объемности» картинки.

Для реализации метода MIP предстоит рассчитать, какие лучи наблюдения пересекают визуализируемый объем, и какие воксели пересекает каждый из них. Сначала произведем расчет для визуализации границ воксельного массива и воспользуемся промежуточными результатами этого расчета для дальнейших вычислений. Для каждой из восьми точек, являющихся вершинами прямоугольного параллелепипеда, произведем расчет координат в воображаемом трехмерном пространстве.

На фиг.1 представлено положение визуализируемого воксельного массива (вид сбоку) относительно плоскости экрана при повороте на угол . На фиг.2 представлено положение оси вращения относительно границ визуализируемого воксельного массива и направление вращения на угол . Пусть текущему ракурсу наблюдения соответствуют углы и , где - угол наклона оси вращения 1, а - угол вращения.

Тогда расчет координат сводится к преобразованию поворота координат в виде:

где x, y, z - исходные координаты визуализируемой точки,

x', y', z' - координаты визуализируемой точки в результате поворота.

Если мы используем «перспективную» проекцию, то нам необходимо дополнительное преобразование:

где x', y', z' - координаты визуализируемой точки в результате преобразований (1), (2), (3),

x'', y'', z'' - координаты визуализируемой точки в результате поворота при построении «перспективной» проекции,

а - размер наибольшего ребра прямоугольного параллелепипеда, ограничивающего визуализируемый объем,

r - расстояние от начала отсчета координат (центра 2 прямоугольного параллелепипеда) до плоскости экрана 3.

и - угловые координаты соответствующего пикселя относительно точки наблюдения.

Следует отметить, что все необходимые тригонометрические функции и их произведения при расчете текущего ракурса считаются один раз и на протяжении всего расчета ракурса фактически являются константами. При использовании «перспективной» проекции тангенсы углов и рассчитываются один раз для каждого пикселя и используются во время расчетов в виде таблицы.

Полученные восемь наборов координат позволяют построить контур, который будет отображать границы визуализируемого объема. То есть необходимо попарно соединить соответствующие точки двенадцатью отрезками. Для этого построения предлагается воспользоваться алгоритмом Брезенхэма. При этом координату z, которая будет определять цвет соответствующего пикселя отрезка, следует считать в виде действительного числа путем последовательного приращения координаты z на величину z/N, где z - разность координат z начальной и конечной точки отрезка, а N - число пикселей, из которых состоит отрезок. Если за начало отсчета координаты z считать плоскость экрана, то цвет пикселя отрезка будет соответствовать интенсивности I(z):

где I_max - значение цвета используемой палитры, соответствующее максимальной интенсивности,

d(z) - весовой коэффициент «затемнения в глубину».

где A - длина диагонали прямоугольного параллелепипеда, ограничивающего визуализируемый объем.

Трехмерные координаты всех точек, из которых состоят полученные отрезки, можно использовать для определения лучей наблюдения которые пересекают визуализируемый объем. Это можно произвести также при помощи алгоритма Брезенхэма, но на плоскости, на которой лежат соответствующие лучи наблюдения. Соответствующие обратные преобразования координат дадут координаты, которые будут соответствовать индексам элементов визуализируемого массива (вокселям), через которые проходит начало и конец луча наблюдения. Зная эти координаты, можно вычислить, какому приращению каждой из координат будет соответствовать «единичный шаг» вдоль соответствующего луча наблюдения. Далее, используя эти значения для последовательного приращения соответствующих координат, мы будем последовательно определять значение интенсивности каждого вокселя, умножать его на величину d(z) и сравнивать с предыдущим максимальным значением для данного луча. Если оно больше, то это и будет текущее значение соответствующего пикселя. Для оптимизации вычислений следует производить каждую последующую выборку значения вокселя вдоль соответствующего луча наблюдения по направлению от наблюдателя «в глубину» визуализируемого объема. При этом текущее значение пикселя следует сравнивать с величиной I(z) (7). Если текущее значение пикселя больше или равно I(z), то перебор для данного луча следует прекратить и текущее значение этого пикселя будет окончательным для данного ракурса.

Введенные выше ограничения на формат исходных данных позволяют значительно ускорить вышеописанные вычисления. Чем меньше объем исходных данных, тем выше скорость обращения к ним. Это связано с использованием в современных CPU и GPU кэш памяти. Кроме того, для обращения к буферу памяти, содержащему исходные данные, в любом случае требуется пересчет трехмерных координат в линейные. Если все линейные размеры визуализируемого массива являются степенью двойки, то эти вычисления сводятся к побитовым логическим и сдвиговым операциям. А это значительно ускоряет процедуру вычислений. Вычисления для каждого луча наблюдения не зависят друг от друга. Поэтому могут выполняться параллельно. Это позволяет эффективно использовать многоядерную и многопоточную архитектуру современных CPU и GPU. В результате авторами были разработаны программы, которые в реальном времени могут качественно отобразить данные ОКТ, MPT, КТ, УЗИ, акусто-тепловидения, конфокальной микроскопии, дефектоскопии и других технических и научных исследований и расчетов. На фиг.3 представлено трехмерное ОКТ изображение - томограмма участка кожи пальца.

Возможность «мгновенной» визуализации любого ракурса наблюдения в зависимости от изменения пользователем положения манипулятора "мышь", создает у пользователя иллюзию того, что он сам вращает в нужные ему стороны и рассматривает некий трехмерный объект. Дополнительные возможности программы по визуализации в отдельном дополнительном окне любого плоского среза в B-режиме, позволяют разглядеть все содержимое визуализируемого массива (фиг.4а и фиг.4б). Возможность «срезать» любые внешние слои и углы, позволяет пользователю получить дополнительную информацию о структуре и суперпозиции визуализируемых данных (фиг.5).

Стоит упомянуть о возможности использования стереоскопичности нашего зрения. Наблюдаемая картинка приобретает «истинную» объемность для нашего восприятия, если для левого и правого глаза сформировать и показать отдельные соответственно преобразованные изображения. Самой простой и дешевой возможностью показа разных изображений для левого и правого глаза является использование анаглифных очков и соответствующее формирование анаглифного изображения. Но этот способ имеет удовлетворительный результат только при использовании черно-белой палитры. В противном случае цвета сильно искажаются. В последнее время бурно развиваются и широкодоступны другие, более совершенные аппаратные средства для стереовизуализации. Все эти средства легко сочетаются с заявленным способом визуализации скалярных трехмерных данных. Достаточно рассчитать не один ракурс наблюдаемой картинки, а два. Причем при расчете "левого" и "правого" ракурса достаточно изменить соответствующим образом массив углов у для лучей наблюдения. То есть углы для "левого" и "правого" ракурсов должны отличаться на постоянную величину, которая зависит от предполагаемого расстояния между наблюдателем и экраном. Сформировав таким образом "левую" и "правую" картинку, мы имеем возможность показать их в виде стереопары любым из доступных способов (фиг.6). В итоге мы получаем довольно впечатляющий результат.

Таким образом, как правило, результатом метода MIP является фильм, состоящий из последовательности предварительно рассчитанных кадров, которые отображают «вращение» визуализируемого объекта вокруг некоторой оси. Это обусловлено тем, что, с одной стороны, недостатки метода MIP не дают пользователю ясно представить пространственную суперпозицию элементов визуализируемого объекта. С другой стороны, при использовании метода MIP для визуализации данных MPT, КТ или УЗИ приходится обрабатывать большой объем исходной информации. Поэтому, несмотря на простоту вычислений, расчет изображения каждого ракурса, без должной оптимизации, занимает значительное время.

Предлагаемая методика позволяет, с одной стороны, улучшить восприятие пользователем результатов визуализации даже в «статическом» режиме. С другой стороны предлагаемая оптимизация вычислений позволяет достаточно быстро рассчитать изображение произвольного ракурса визуализируемого объема благодаря повышению скорости вычислений за счет уменьшения объема памяти, затрачиваемой на построение изображения, что делает возможным обеспечение «интерактивности» просмотра данных пользователем даже при использовании мощности современного настольного или переносного персонального компьютера. То есть пользователь сможет «поворачивать» визуализируемые данные и рассматривать их с любой стороны в реальном времени.

В процессе оптимизации вышеизложенных алгоритмов для современных CPU и GPU с использованием методов параллельных многопоточных вычислений достигнуты следующие результаты. Массив размером 256×256×256 рассчитываем в среднем за 60-70 миллисекунд на CPU, аналогичном Intel Core i3 (~2.8 GHz 4-threads). Пропорциональное увеличение тактовой частоты и количества потоков у используемого процессора дает примерно пропорциональное уменьшение среднего времени расчета. Но при увеличении объема исходных данных скорость расчета падает непропорционально. Это связано с уменьшением скорости доступа к памяти из-за увеличения количества кэш-промахов. Использование параллельных вычислений на GPU значительно расширяет наши возможности и улучшает результаты. Опыты по использованию GPU Nvidia с интерфейсом CUDA показали следующие тенденции. При размерах исходных данных 512×512×512 ракурс рассчитывается в среднем за 60-70 миллисекунд на GPU Nvidia GeForce GTX 260. На GPU Nvidia GeForce GTX 470 за 40-50 миллисекунд. И так далее. Существенным в этом случае является следующее требование к объему видеопамяти. Особенности разработанного способа оптимизации требуют, чтобы количество доступной видеопамяти было в три раза больше, чем объем визуализируемого массива. Имеются широкие возможности использования вышеизложенных алгоритмов для других аппаратных средств.

Таким образом, разработанный способ оптимизации метода проекции максимальной интенсивности для визуализации скалярных трехмерных данных в статическом режиме, в интерактивном режиме и в реальном времени позволяет улучшить восприятие объемности получаемого изображения, определить ракурс визуализируемого изображения, а также повысить скорость вычислений за счет уменьшения объема оперативной памяти, затрачиваемой на построение изображения, что обеспечит интерактивность работы с получаемым изображением, т.е. даст возможность изменения ракурса изображения пользователем в реальном времени.