+380 57 755 34 05 team@fulcrum.software

Автор: сотрудник компании Fulcrum Software.

Одной из реальных рабочих задач, порученных мне компанией Fulcrum Software, было предоставление докторам (пользователям) информации о различных метриках (размерах, площадях, объемах) органов пациента по сходным медицинским изображениям (рентген, узи, прочее). В ходе выполнения этой задачи мне пришлось познакомиться с различными алгоритмами из области обработки многоугольников. В этой статье я хотел бы описать интересный опыт выполнения поставленной задачи. А начал я, как принято, со знакомства с теорией.

Приложения, связанные с алгоритмами обработки многоугольников

Алгоритмы обработки многоугольников находят применения в различных практических приложениях [“a fast trapezoidation …”]. К ним относятся приложения компьютерной графики (задача видимости), геоинформационные системы (задача закраски областей), робототехника (планирование траектории движений) и машинное обучение (распознавание образов).

1. Задача видимости (visibility problem)

В общей постановке задача видимости [Preparata, Fujimura] звучит следующим образом. Пусть на плоскости заданы набор многоугольников и точка наблюдения O. Требуется определить, какие участки многоугольников видимы из точки наблюдения O. Точка Pвидима из точки O в случае, если отрезок прямой OP не пересекает ни один из заданных многоугольников.

2. Задача о картинной галерее (Art Gallery Problem)

Понятие видимости используется в одной из классических задач, которая называется «задача о картинной галерее» (Art Gallery Problem) [“Vis in Plane”]. Эта задача звучит следующим образом. Требуется определить, сколько нужно установить охранников для того, чтобы каждая точка интерьера картинной галереи находилась под наблюдением. То есть, каждая точка интерьера должна быть видимой для одного из охранников так, чтобы все картины в галерее находились под охраной.

3. Закраска областей (region filling)

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

4. Планирование (траектории) движения (path planning)

Задача планирования движений является одной из важнейших задач, встречающихся в робототехнике. В самом простом виде эта задача звучит следующим образом.

Пусть на плоскости заданы:

  • движущийся объект A, изначально расположенный исходной позиции s;
  • множество неподвижных объектов Bi, расположенных в заданных позициях и имеющих форму многоугольников;
  • допустимая целевая позиция g.

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

4. Распознавание образов

Одним из приложений, связанных с разбиением многоугольников, является задача распознавания образов [D.T. Lee]. В подобных задачах фигура вначале разбивается на более простые компоненты, называемые примитивами, а затем сравнивается с имеющимися образцами по какому-либо критерию подобия. Такие примитивы, как правило, имеют выпуклую форму.

Проведение измерений на медицинском изобзажении

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

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

Измерительные инструменты

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

Среди базовых/основных типов измерительных инструментов можно назвать следующие:

Кронциркуль (caliper). Может быть двухточечным, либо состоять из набора точек. Используется для вычисления расстояний на 2D изображении;

Контур. Контуры используются для измерения площади и объема на 2D изображении.

Измерение произвольной замкнутой области

Мне предстояло реализовать инструмент, позволяющие измерять площадь замкнутой области на 2D-изображении.

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

Во время выделения контура пользователь может сбиться и некорректно отметить фрагмент области. Для подобных случаев предусмотрена возможность удаления фрагментов контура. Пользователь просто ведет трекбол в обратном направлении вдоль отмеченного контура, и этот фрагмент контура удаляется.

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

Требования к реализации вычисления площади замкнутой области

Выделение области на 2D-снимке реализовано таким образом, что контур области аппроксимирован многоугольником. Каждая сторона многоугольника имеет определенную длину, достаточно малую (несколько пикселей), чтобы контур на экране выглядел как гладкая кривая линия.

Расчет периметра области реализовать просто. Достаточно просуммировать длины всех сторон многоугольника.

Рассчитать площадь многоугольника значительно сложнее. В данном случае алгоритм расчета должен учитывать следующие особенности:

  • пересчет по мере выделения (рисования) области; когда пользователь выделяет область на снимке, он должен видеть текущее значение  площади, поэтому алгоритм должен пересчитывать площадь по мере добавления вершин в многоугольник;
  • пересчет площади при удалении части контура; пользователь может удалить часть контура, после чего значение площади должно быть быстро пересчитано;
  • большое количество вершин многоугольника при достаточно малой длине сторон; данная особенность может вызвать проблемы с точностью расчетов.

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

Существующие подходы к вычислению площади многоугольника

Разбиение на треугольники. Алгоритм отсечения ушей

Среди алгоритмов триангуляции многоугольника одним из самих простых для реализации является алгоритм отсечения ушей (Ear cutting) [“comp. geometry in C”], или алгоритм Ван-Гога [Скиена]. Суть данного алгоритма состоит в последовательном отсекании «ушей» многоугольника.

Для определения понятия уха многоугольника требуется ввести понятие диагонали многоугольника.

Диагональю многоугольника P называется отрезок s, соединяющий две его вершины viи vj (ij) такой, что s не пересекает никаких граней многоугольника и лежит внутри многоугольника P.

Ухом многоугольника P называется треугольник, образованный последовательными вершинами a,b,c такой, что отрезок ac является диагональю многоугольника. Вершина bназывается вершиной уха.

Например, ниже на рисунке 1 треугольник, образованный вершинами A,B и H, образуют ухо, так как BH является диагональю многоугольника. Треугольник D,C,E не может быть ухом, т.к. его основание EC пересекает отрезок FG и поэтому не является диагональю многоугольника.

Рисунок 1. Отсечение «уха» у многоугольника.

В работе [“comp. geometry in C”] приводится и доказывается следующая теорема.

Теорема о двух ушах. Любой многоугольник, содержащий n≥4 вершин, содержит, по крайней мере, 2 непересекающихся уха.

Из этой теоремы можем получить следующий алгоритм (листинг 1), позволяющий определить площадь многоугольника.

Листинг 1. Алгоритм отсечения ушей.

Шаг 1. Обнулить начальное значение площади многоугольника.Шаг 2. Выполнять шаги 3-6, пока количество вершин многоугольника больше трех.Шаг 3.     Найти вершину многоугольника, определяющую вершину уха.Шаг 4.     Образовать ухо и вычислить его площадь (треугольник).Шаг 5.     Увеличить текущее значение площади многоугольника.Шаг 6.     Удалить вершину из многоугольника.Шаг 7. Вычислить площадь оставшегося многоугольника, состоящего из трех вершин.Шаг 8. Увеличить текущее значение площади многоугольника.

Разбиение на треугольники. Простейший метод

Существует довольно простой метод расчета площади многоугольника [Порублев]. Он основан на теореме, которая приводится и доказывается в работе [Лопшиц].

Теорема. Пусть дана упорядоченная совокупность точек a1, a2, …, an, произвольно расположенных на плоскости и определяющих некоторый ориентированный n-угольник. Для произвольной точки c составим следующую сумму:

Эта сумма определяет сумму площадей ориентированных треугольников, для которых общей вершиной является произвольно выбранная точка c, а основаниями – направленные стороны заданного ориентированного n-угольника. Тогда сумма  не зависит от положения точки с.

Пусть задан простой многоугольник (выпуклый или невыпуклый), состоящий из n вершин (и n сторон). Из вышеприведенной  теоремы следует, что если из любой вершины aiэтого многоугольника провести отрезки во все остальные вершины, то получим n-2 треугольника aia0a1, aia1a2, …, an-1an-2ai. Сумма ориентированных площадей этих треугольников равна ориентированной площади всего многоугольника.

Рисунок 2. Площадь многоугольника как сумма площадей треугольников

Рассмотрим рисунок 2. Область Ω в результирующей площади многоугольника не учитывается, так как в результате вычисления ориентированных площадей треугольников a0a3a4 иa0a2a3 эта область сначала учитывается с отрицательным знаком, а затем – с положительным. Область Ψ учитывается один раз в результирующей площади многоугольника. Эта область два раза учитывается с положительным знаком  (как часть треугольников a0a1a2 и a0a3a4) и один раз – с отрицательным (как часть треугольника a0a2a3).

Разбиение на трапеции

Рассмотрим многоугольник, приведенный на рисунке 3 ниже.

Рисунок 3. Разбиение многоугольника на трапеции.

Возьмем две любые смежные вершины многоугольника, например, C  и D. Опустим из этих точек проекции на ось Ox, получим точки C ‘ и  D ‘. Точки C, D, C ‘, D ‘  образуют трапецию с основаниями  и  и боковыми сторонами (С, D) и (C ‘, D ‘ ) . Ориентированная площадь данной трапеции вычисляется по формуле:

Здесь Сy и Dy – это длины оснований трапеции, которые представляют собой у‑координаты точек C и D соответственно. Разность (Dx— Cx) представляет собой высоту трапеции (разностьy-координат оснований трапеции). Эта разность может быть положительной либо отрицательной (что и требуется, как будет показано далее).

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

Простой алгоритм, основанный на значениях координат вершин

В работе [D. Chandler] приводится очень простой алгоритм, позволяющий получить площадь произвольного многоугольника. В качестве исходных данных для алгоритма требуются значения координат всех вершин многоугольника.

Допустим, имеется многоугольник с вершинами (x1, y1), (x2, y2), …, (xn, yn). Перечислим эти вершины в столбик в порядке их обхода на многоугольнике. Первая вершина должна быть и последней.

x1y1
x2y2
xnyn
x1y1

Теперь, перемещаясь вниз от первой вершины до n-й, перемножим каждую x‑координату текущей вершины с y-координатой следующей вершины, и запишем результат справа от y-координаты. И точно так же, перемещаясь вниз,  перемножим каждую y-координату текущей вершины с x-координатой следующей вершины, и запишем результат слева от x-координаты. В результате к перечисленному выше списку вершин добавятся столбцы слева и справа, в которых находятся перемноженные значения:

 x1y1 
y1*x2x2y2x1*y2
yn-1*xnxnynxn-1*yn
yn*x1x1y1xn*y1

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

 x1y1 
y1*x2x2y2x1*y2
yn-1*xnxnynxn-1*yn
yn*x1x1y1xn*y1
Σ1 Σ2

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

Рассмотрим, как работает данный алгоритм на примере многоугольника, приведенного на рисунке 4 ниже.

Рисунок 4. Пример многоугольника для расчета площади.

Проделав шаги, описанные выше, получим следующее:

 14 
4 * 3 = 12371 * 7 = 7
7 * 11 = 771153 * 5 = 15
5 * 9 = 459111 * 1 = 11
1 * 4 = 4429 * 2 = 18
2 * 1 = 2144 * 4 = 16
Σ1 = 140 Σ2 = 67

Тогда площадь рассматриваемого многоугольника равна:

Этот алгоритм и был выбран мной для реализации мной поставленной задачи.