УДК 004.921

Визуализация больших сегментированных объемов в web-приложениях

Беляев Сергей Юрьевич – кандидат физико-математических наук, доцент кафедры прикладной математики Санкт-Петербургского политехнического университета Петра Великого

Кацман Надежда Игоревна – магистрант Санкт-Петербургского политехнического университета Петра Великого

Аннотация: Подлежащие визуализации данные могут занимать сотни и даже тысячи гигабайт. В каждый момент времени web приложение выбирает требуемые для визуализации данные, которые пересылаются по сети и визуализируются используя WebGL. Проблема этой схемы при визуализации объемов состоит в том, что для формирования одного изображения требуется переслать несколько гигабайт данных. В случае сегментированных данных эта проблема может быть решена за счет упаковки данных. Однако на алгоритм упаковки накладывается ограничение – упакованный формат должен позволять выполнять визуализацию на WebGL. В работе предложен такой формат в предположении, что вектор наблюдения направлен вдоль одной из осей объема. Кроме того, разработан алгоритм, реализующий визуализацию данных в этом формате на WebGL в режиме реального времени.

Ключевые слова: визуализация объемов, большие данные, упаковка данных, web приложения.

1. Введение

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

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

В случае объемной визуализации для формирования требуемого изображения требуется на порядок больше данных, чем, например, для карты. Такой объем данных передать невозможно за требуемое время по крайней мере по сетям 4G. Для сегментированных данных задача упрощается, так как такие данные поддаются сжатию, поскольку сегменты полностью определяются своими границами. Однако, на формат компрессированных данных накладывается ограничение – он должен позволять выполнять визуализацию на GPU. Предлагаемый в данной работе метод решает эту проблему в предположении, что объем наблюдался в заранее определенном направлении и для визуализации используется параллельное проецирование. В этом случае пользователю доступны только операции масштабирования и сдвига. Это предположение часто является естественным для протяженных объемов (см. рисунок 1).

Рисунок 1. Направление просмотра множества данных.

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

2. Связанные работы

Существуют два подхода к  упаковки объемных сегментированных данных. Первый подход основан на триангуляции поверхностей сегментов [1, 2]. Так как сегменты являются прозрачными, визуализация в этом случае выполняется с помощью алгоритма трассировки лучей. Трассировка лучей на GPU реализуется с помощью вычислительных шейдеров [3], но скорость визуализации оказывается приемлемой только для очень небольшого числа сегментов.

Другой подход к упаковке основан на сохранении только тех вокселей объема, которые расположены на границах сегментов. Сжатие этих границ с помощью цепного кода позволяет получить очень высокую степень компрессии [4-6]. К сожалению, формат полученных данных не может быть использован для визуализации на GPU, поэтому он в основном используется для хранения и передачи больших сегментированных данных.

Различные реализации алгоритма визуализации сегментированных объемных данных разработаны в работах [7-9]. Эти реализации в качестве исходных данных используют исходную трехмерную матрицу, содержащую номера сегментов, но могут применяться только к относительно небольшим объемам данных. В случае больших данных реализуются различные механизмы динамической подгрузки данных [10], что существенно сказывается на производительности.

3. Метод

3.1 Краткое описание

Пусть исходные данные представлены трехмерной матрицей V, элементы которой содержат индексы сегментов сегментированного объема. Эти индексы определяют цвет и прозрачность сегментов через палитру, которая определяется пользователем. Обозначим через X, Y и Z размерности матрицы V по осям x, y, z, a через С – число сегментов различного типа содержащихся в V. Предлагаемый метод не накладывает ограничений на значения X и Y, но значения  Z должны удовлетворять условию:

              (1)

В случае Z=1 задача сводится к визуализации изображения большого размера. Для решения этой задачи в google map и других подобных приложениях по заданному изображению, строится пирамида [11] изображений меньшего размера. Когда пользователь выбирает масштаб просмотра, используется изображение подходящего масштаба. Для того, чтобы не пересылать по сети изображения целиком, они делятся на блоки. В зависимости от заданной пользователем области просмотра, пересылаются только те блоки, которые попадают в требуемую область.

В данной работе используется тот же подход, но каждый блок вместо одного изображения содержит Z двухмерных матриц. Это приводит к существенному увеличению передаваемого по сети объема данных. Например, если Z = 2¹² и C = 256 , то для формирования изображения размером 1024*1024 потребуется передать по сети и загрузить в видеопамять 4Gb данных. Для уменьшения этого объема на этапе подготовки данных выполняется двухуровневая упаковка данных.  Упаковка первого уровня использует когерентность данных в направлении оси z, обусловленную тем, что данные сегментированы. Данные этого уровня компрессии предназначены для загрузки в видео память. Формат этих данных позволяет шейдерам выполнять визуализацию без распаковки. Второй уровень сжатия использует комбинацию форматов jpeg8 и jpeg12 и предназначен для сокращения объема передаваемых по сети данных. Полученные в этом формате данные распаковываются. Результатом распаковки являются данные упакованные до первого уровня, которые, как сказано выше, загружаются в видео память. Описанная схема показана на рисунке 2.  

Рисунок 2. Двух уровневая компрессия данных.

• Подготовка данных

По заданному объему V строится пирамида объемов V^K. Объем V^K получается из объема V^K-1 путем уменьшения его размера в два раза в плоскости x, y. Каждый объем V^K делится на блоки V^K_IJ фиксированного размера равного r*r*Z (см. рисунок 3). Выбор величины r обсуждается в пункте “анализ метода”. 

Рисунок 3. Пирамида объемов.

Объемы V^K_IJ  упаковываются независимо друг от друга с помощью следующего алгоритма.

Пусть V^* – трехмерная матрица размерности r*r*Z, элементы V^*[i,j,k] которой содержат индексы сегментов. Обозначим через M максимальное число сегментов, содержащихся в столбцах V^*[i,j], а через  h_ijm длину отрезка, соответствующего m-му сегменту в столбце  V*[i,j] (см. Рис. 4)

Рисунок 4. Длины сегментов в столбце (i, j) объема V^*.

Для всех значений i и j строятся списки величин h_ijm , которые cохраняются в матице S размерности r*r*M (см. рисунок 5). Дополнительно в матрице P той же размерности сохраняются индексы сегментов.

 

Рисунок 5. Трехмерная матрица S содержит длины сегментов в направлении k.

Особенностью матрицы S является то, что с ростом номера m возрастает число элементов S[i,j,m] равных нулю. Эту особенность использует алгоритм визуализации (см. пункт 3.4). Для работы этого алгоритма требуется знание среднего значения длины списков, содержащихся в S. Далее это значение обозначается через n. Вычисление n является частью алгоритма упаковки. Данный алгоритм применяется ко всем объемам V^K_IJ  и результатом его применения являются трехмерные матрицы S^K_IJ  и P^K_IJ , a также значения n^k_ij. Эти данные представляют первый уровень компрессии.

Второй уровень компрессии создается следующем образом. К двухмерным матрицам S[m], составляющим трехмерную матрицу S, применяется jpeg компрессия. В зависимости от числа бит, требуемых для представления элементов матрицы S[m], используется формат jpeg8 или jpeg12. Формат jpeg позволяет выполнять сжатие с потерями, что увеличивает степень сжатия. Степень сжатия зависит также от когерентности данных. Данные, располагающиеся в верхних слоях пирамиды объемов, имеют меньшую когерентность. С другой стороны, при их визуализации требуется меньшая точность. Это позволяет компрессировать их с большими потерями, добиваясь примерно одинаковой степени сжатия.

В отличие от матрицы S, компрессия матрицы P должна выполняться без потерь.

• Распаковка данных

По заданному пользователем масштабу и смещению определяются и скачиваются по сети необходимые для визуализации блоки. Каждый блок после выполнения декомпресии содержит две трехмерные матрицы – S и P. Дополнительно по сети передается число n. Матрицы S и P объединяются в одну матрицу T, элементы которой вычисляются по формуле: 

                                                                                               (2)

Для уменьшения объема данных передаваемых в видеопамять, каждая матрица T разбивается на две – Т₁ и T₂. Матрица T₁ имеет размерность r*r*n и совпадает с T при m меньших и равных n. Из этой матрицы создается 3D текстура, которая используется во фрагментном шейдере во время визуализации.

Матрица T₂ содержит элементы матрицы T при k>n. Ненулевые столбцы этой матрицы представляются списками L[i,j], структура которых показана в таблице 1.

Таблица 1. Структура списка L[i,j].

Число элементов

i

j

T[i,j,n+1]

T[i,j,n+2]

…

Эти списки сохраняются в массивах фиксированной длинны подходящего размера. Используются массивы B_t с размерами элементов 4*2^t байт. Список L[i,j] помещается в массив B_t, если его длина находится промежутке (4^t-1, 4^t]. Структура массива B_t показана в таблице 2.

Таблица 2. Структура массива B_t.

	4*2t байт	4*2t байт
Число элементов	L[i1,j1]	L[i2,j2]	…

Таким образом,  матрица T₂ представляется множеством массивов B_t. Из этих массивов создаются вершинные буферы, которые  пересылаются в видеопамять для обработки в вершинных шейдерах.

Aппаратным ограничениям GPU устанавливают максимальное значение t равным 4. Как следствие, слишком длинные списки L[i,j] могут помещаться в массивах B_t. В этом случае значение n увеличивается уменьшая число элементов в L[i,j] за счет увеличения размера трехмерной матрицы T₁.

4. 4. Визуализация

Визуализация выполняется с помощью алгоритма ray-casting во front-to-back порядке [10]. Пусть на пути луча встречается М сегментов. Обозначим через C_m и A_m цвета и прозрачности участков пересекаемых сегментов. Согласно ray-casting алгоритму цвет пикселя вычисляется по рекуррентной формуле

                                                                                (3)

                                                                               (4)

где C^*_m и A^*_m накопленные цвета и прозрачности. Как показано в работе [10], если материал имеет прозрачность , то прозрачность отрезка этого материала длинной h определяется соотношением 

                                                                                                  (5)

Следовательно, уравнение (4) может быть заменено на следующее:

                                                          (6)

Расчеты по этим формулам выполняются до достижения величины m значения M или величины A^*_m значения близкого к 1.

В нашем случае первые n шагов выполняются во фрагментном шейдере, который выбирает данные из 3D текстуры без фильтрации. Результат сохраняется в RGBA текстуре, которая используется как начальные данные для следующих шагов. Эти шаги выполняются в вершинных шейдерах, использующих данные из вершинных буферов, которые были созданы из массивов B_t. В вершинных шейдерах дополнительно вычисляются координаты пикселей. Для этого используются индексы i, j, содержащиеся в буферах вершин, и координаты блоков, которые передаются в вершинные шейдеры через буфер констант. Общая схема визуализации показана на рисунке 6.

Рисунок 6. Схема визуализации.

Описанный алгоритм имеет сложность O(n log₂(Z)) в среднем. Наличие логорифма в этом выражении объясняется ростом требуемой памяти для хранения величин h_m входящих формулу (6), что увеличивает число обращений к памяти.  Коэффициент ускорения по сравнению со стандартным методом имеет порядок

                                                                                                 (7)

в среднем. Это следует из того, что формулы (5) и  (6) имеют примерно одинаковую вычислительную сложность, но число итераций разное. Достигнутое ускорение позволяет выполнять визуализацию в реальном времени.

5. Анализ метода

Подлежащие пересылки по сети блоки представляют собой трехмерные матрицы размерности r²Z. Предположим, что C < 256. Тогда блок занимает r²Z байт. После упаковки первого уровня блок заменяются двумя трехмерными матрицами S и P размерности r²M. Элементы матрицы P при C < 256 являются однобайтовыми, так что P занимает r²M байт. Элементы матрицы S занимают log₂g бит, где g – максимальная длина отрезков h_ijm в блоке. Так как h_ijm < Z, то S занимает не более r²М log₂(Z/256) байт. Учитывая, что Z удовлетворяет условию (1), это значение ограничено величиной 4 r²М байт. Таким образом, коэффициент сжатия первого уровня  определяется соотношениями

                                                                                    (8)

На втором уровне сжатия к двухмерным матрицам S[m] и P[m] применяется jpeg компрессия. Так как матрицы S[m] компрессируются с потерями, а P[m] без потерь, то основной вклад в итоговый результат вносят матрицы P[m]. С ростом m число нулевых элементов в P[m] возрастает, поэтому коэффициент сжатия P и второго уровня в целом определяется формулой

                                                                                                             (9)

где d – коэффициент jpeg сжатия двухмерных матриц P[m] без потерь. Следовательно, общий коэффициент сжатия может быть вычислен по формуле

                                                                                                            (10)

Пусть, например, r = 1024, Z= 2¹². Тогда исходный размер блока занимает 4Gb памяти, а упакованный – n/d Mb.  Если n = 5, а d = 10, то упакованный блок будет занимать 0.5 Mb. Такой объем может быть передан по  сети 4G за время меньшее 0.1 сек обеспечивая визуализацию в интерактивном режиме в окне визуализации 1024*1024. Если же n=20 , а d = 5, то упакованный блок будет занимать 4 Mb исключая возможность визуализации в интерактивном режиме. Эта проблема может быть решена, если выбрать r=128. При таком выборе для отрисовки в окне того же размера потребуется 64 блока объема  4 Mb /64 = 0.0625 Mb. Это не изменит время пересылки данных для первого кадра визуализации, но при сдвиге точки наблюдения потребуется передать по сети не более 14 блоков (см. Рис. 7), то есть всего 14 * 0.0625 = 0.875 Mb.  

Рисунок 7. Сдвиг окна просмотра. Синим цветом показаны загружаемые по сети блоки.

6. Заключение

Разработан алгоритм упаковки объемных сегментированных данных большого размера, который позволяет создать web приложения для визуализации таких данных. Это приложение может работать в реальном времени при небольших изменениях масштаба и/или области наблюдения и в интерактивном режиме в противном случае. Пользователь может управлять цветом и прозрачностью сегментов. Алгоритм упаковки предполагает, что направление просмотра данных задано.  При необходимости просмотра данных под другим углом, необходимо подготовить новый набор упакованных данных. Еще одним ограничением метода является отсутствие возможности использования центрального проецирования, то есть данные могут просматриваться только в параллельной проекции.

Список литературы

1. Kaufman and K. Mueller, “Overview of volume rendering,” The Visualization Handbook. 2005. vol. 7. P. 127–174.
2. Lorensen, and H. Cline, “Marching cubes: a high resolution 3D surface construction algorithm,” Proc. of SIGGRAPH’87. 1987. P. 163-169.
3. Marrs, A., Shirley, P. and Wald, I., "Ray tracing Gems II: next generation real-time rendering with DXR, Vulkan, and OptiX." Springer Nature, Berlin. 2021. P. 281–299.
4. E Aldemir, OAT Dueñas, AE Kavur. Chain code strategy for lossless storage and transfer of segmented binary medical data. Elsevier. 2023. Vol. 219.
5. Kavur, A., Emre, N. S., Gezer, M. B., Aslan, S., Conze, P. H., Groza, V., … Alper Selver, M. CHAOS challenge - Combined (CT-MR) healthy abdominal organ segmentation. Medical Image Analysis. 2021. P. 69-75.
6. Selver, M. A., & Emre Kavur, A. Implementation and use of 3D pairwise geodesic distance fields for seeding abdominal aortic vessels. International Journal of Computer Assisted Radiology and Surgery. 2016. N. 11(5). P. 803–816.
7. Helwig Hauser, Lukas Mroz, Gian-Italo Bischi, and Eduard Gr¨oller. “Two-level volume rendering-fusing MIP and DVR.” In Proceedings of IEEE Visualization. 2000. P. 211–218.
8. Helwig Hauser, Lukas Mroz, Gian Italo Bischi, and Eduard Gr¨oller. “Two-Level Volume Rendering.” IEEE Transactions on Visualization and Computer Graphics. 2001. N. 7. P. 242–252.
9. Cs´ebfalvi, L. Mroz, H. Hauser, A. K¨onig, and E. Gr¨oller. “Fast Visualization of Object Contours by Non-Photorealistic Volume Rendering.” In Proceedings of Eurographics. 2001. P. 452–460.
10. Engel K, Hadwiger M, Kniss J, Rezk-Salama C, Weiskopf D. Real-time volume graphics. AK Peters. 2006. P. 444-457.
11. Fisher, R., Breckon, T.P., Dawson-Howe, K., Fitzgibbon, A., Robertson, C., Trucco, E., Williams, K.I. Dictionary of computer vision and image processing. Second edition. New Jersey: John Wiley & Sons. 2013. P. 382- 391.