CC BY
2
УДК 004.932.2
DOI: 10.24412/2071-6168-2023-7-562-563
АЛГОРИТМ ВЫЧИСЛЕНИЯ ОБЪЁМА ПРИМЕСЕЙ В ЖИДКОСТИ ПО ОПТИЧЕСКИМ ДАННЫМ
Л.О.А. Салех, С.В. Хлопин, Л.В. Черненькая
Работа посвящена разработке методов и их программной реализации для обработки бинарных изображений, вычисления суммарного объема содержания примесей в жидкости. Для вычисления производится поиск контуров на бинарном изображении, выделение замкнутых областей, анализ иерархии и вычисление площади и объема содержания примесей. Предложен авторский метод поиска характерных геометрических размеров контуров примесей и применения этих размеров для вычисления объема примесей.
Ключевые слова: компьютерное зрение, поиск контуров, анализ иерархии контуров.
Проблема детектирования, классификации и анализа объектов посредством обработки графической информации требует совершенствования методов исследования изображений. Применение полученных методов может быть применено, например, для поиска включений, содержащих нефтепродукты и подсчета суммарной составляющей примесей.
Постановка задачи. Для проведения вычисления объема содержания примесей нефтепродукта в жидкости требуется провести обработку и анализ фотоснимков потока жидкости. Алгоритм следует разделить на несколько этапов:
1. Предварительная обработка изображения, фильтрация (сглаживание);
2. Коррекция цвета и выделение границ объектов;
3. Выделение объектов в определенном цветовом диапазоне;
4. Фильтрация высокочастотных шумов на изображении;
5. Определение иерархии замкнутых контуров границ объектов;
6. Вычисление площади детектированных объектов;
7. Определение характерных геометрических размеров объектов;
8. Вычисление объема примесей.
Пункты 1-4 алгоритма описаны в статье [1].
В качестве исходных данных для задачи будем использовать изображения жидкости с наличием включений - пузырьков нефтепродуктов. Инородные включения можно определить оптическим образом.
Рис. 1. Примеры исходных изображений
Алгоритм определения иерархии замкнутых контуров и извлечения массива точек контуров бинарных изображений. После применения предварительной обработки, выделения объектов с спектром цвета, находящимся в определенном цветовом диапазоне и фильтрации шумов алгоритмом определения концентраций примесей [1] для рис. 1 мы получаем следующие бинарные изображения: (рис. 2).
На рис. 2 представлены бинаризованные изображение жидкости с примесью тяжелых нефтепродуктов при естественном освещении (2. а) и при облучении коротковолновым ультрафиолетом (2.б) размером M х N.
Для определения иерархии замкнутых контуров границ будем частично использовать алгоритм suzuki85 [2]. Исходное бинарное изображение jBmary ^х,y^, х N) содержит значения цветности пикселей 0 для обозначения фона (background) и 1 для объекта (foreground) соответственно [3]. Введем матрицу jReglons ^х, y^, х N) , содержащую информацию об областях объектов изображения:
562
б
а
0 background, IRegions (x, y) = J 1 foreground, (1)
2,3,... region label.
TT jRegions j-Binary
На этапе инициализации значения матрицы I равны значениям матрицы I для всех элементов.
а б
Рис. 2. Результат предварительной обработки изображения
Определение иерархии замкнутых контуров границ объектов. Иерархия замкнутых контуров описывается структурой Contours:
Contours
{
Id - уникальный номер контура;
Type - идентификатор принадлежности типа контура к внешним границам объектов или контурам отверстий;
Points[] - массив пикселей принадлежащих текущему контуру;
Parent - индекс внешнего по отношению к текущему контуру (граница изображения имеет нулевой индекс номера контура);
}
Для определения принадлежности точек к иерархии замкнутых контуров, необходимо определить к какому типу контуров относится рассматриваемая точка. Точки объектов lRegions (x, y) = 1 могут принадлежать к контурам объектов внешней границы (OB - outer border) или контурам отверстий внутренней границы (HB - hole border).
Для хранения номера контура в иерархии контуров применяются две переменные: номер текущего контура NBD (newest border number), по умолчанию NBD имеет значение 1 и номер предыдущего найденного алгоритмом контура LNBD (last newest border number). По умолчанию LNBD также имеет значение 1.
Предварительным шагом алгоритма suzuki85 проводится инициализация границы исходного изображения нулем -увеличения размера изображения IRegions путем добавления по одной нулевой строки и столбца слева, справа, сверху и снизу. Новое изображение будет иметь размеры ((M + 2) х (N + 2)) .
Шаг 1. производится сканирования изображения построчно слева направо со второй до N — 1 строки со второго до M — 1 столбца (сканирование в размерах исходного изображения ^M х N) ) для поиска принадлежности точки к контуру объекта (foreground) iReglcns (x, y) = 1 и поиск к какому Id из контуров структуры Contours принадлежит текущая точка. Если ^^ют- (x, y) = о или
^IRegicns(x — 1, y) ф о и iRegicns(x, y) ф 0 и IRegicns(x + 1,y) Ф 0;) , то текущая точка не принадлежит контурам,
а является внутренней частью фона или объекта. Необходимо переместиться к следующей координате точки, увеличивая индекс столбца или строки, или завершить алгоритм при достижения правого нижнего угла изображения.
Определение принадлежности к контуру производится вычислением соотношений:
F d „ Job, iRegicns(x,y)—iReglcns(x — 1,y) = 1
Found type = < (2)
" [HB, IRegions (x, y) — IRegions (x +1, y) > 1
x x+1 x-1 x
У I 0 | 1 | y | >1 | 0
Контуры объектов Контуры отверстий
В случае обнаружения при сканировании входного изображения принадлежности текущей точки к внешней или внутренней границе контура производится увеличение значения переменной NBD на единицу и создается новый объект структуры Contours . При этом элементы структуры получают значения:
Contours {
Type = Found _ type; Id = NBD;
Point[] = Point (x, y);
Parent = 0;}
На начальном этапе происходит детектирование внешней границы OB . Определение на начальном этапе внутренней границе невозможен, так как на первом шаге мы увеличили внешние границы бинарного изображения на 1 элемент и заполнили внешние строки и столбцы нулями. Появление отличного от нуля значения элемента массива
IRegi°m проявится сначала в элементе IRegions (x + 1, y) , элемент IRegions (x, y) при этом будет иметь значение ноль. Дальнейшее смещение индекса x вправо приведет к тому, что IRegions (x, y) примет значение единицы и выполнится первое условие (2). Будет найдена внешняя граница. Contours Type первого элемента структуры будет иметь
значение OB. Contours.Id первого элемента структуры будет иметь значение 2. Contours.Parent первого элемента структуры примет значение 0.
В случае обнаружения принадлежности контура в границе объекта OB, необходимо запомнить значение точки начала обхода контура PointTraversalStart = (x -1, y)
В случае обнаружения принадлежности контура в границе отверстия HB, необходимо запомнить значение точки начала обхода контура PointTraversalStart = (x +1, y) и заменить значения переменной LNBD, которое имело значение 1, значением max ^LNBD, IRegions (x, y)^ .
В алгоритме suzuki85 иерархия строится только в случае обязательного расположения контура объекта в контуре отверстий (или внешней рамки). Описание выбора родительского контура приведено в [2, табл 1]. В случае если рамка объекта по ширине равна одному элементу (точке), внутри такого объекта не будет отверстия и любые объекты, размещенные в области внутри такого объекта, получат значение Contours .Parent = 0. Алгоритм не детектирует внутренних контуров, у которых верхняя грань представляет из себя горизонтальную прямую линию толщиной в 1 элемент.
В случае не обнаружения контуров переходим к шагу 3.
Шаг 2. При создании очередного объекта структуры Contours необходимо произвести вычисление принадлежности созданного объекта к объектам из существующей иерархии контуров. Обозначим центральную точку обхода PointTraversalCenter = (x, y).
2.1. Производится обход по часовой стрелке вокруг точки PointTraversalCenter начиная с точки PointTraversalStart и поиск ненулевого значения. Координаты найденной точки сохраняются PointTraversalBack . В случае, если не было найдено ненулевого значения, мы определили, что нашли единичный пиксель (контур, состоящий из одной точки), присваиваем IRegions (PointTraversalCenter) = -NBD, переходим к шагу 3.
2.2. Так как мы определили в предыдущем шаге что был найден контур, а не единичный пиксель, продолжаем обход контура в обратном направлении для поиска всех точек, принадлежавших ему. Производится обход против часовой стрелке вокруг точки PointTraversalCenter начиная с точки PointTraversalBack и поиск первого ненулевого значения. Координата найденной точки сохраняется PointTraversalForward и добавляется в массив точек контура структуры Contours .Point[] = PointTraversalForward.
2.3. Если мы достигли правую границу объекта IRegions (PointTraversalCenter .x +1, PointTraversalCenter. y) = 0 то присваиваем IRegions (PointTraversalCenter) = -NBD.
2.4. Если IRegions(PointTraversalCenter.x +1,PointTraversalCenter.y) ф 0, а текущее значение
IRegions (PointTraversalCenter) = 1, то присваиваем IRegions (PointTraversalCenter) = NBD.
2.5. Если после выполнения шага 2.2 координата точек конца обхода контура PointTraversalForward совпала с точкой (x, y) и PointTraversalCenter = PointTraversalBack то это означает, что мы обошли весь контур против часовой стрелке. Переходим к шагу 3. Иначе переопределяем координаты точек PointTraversalStart = PointTraversalCenter и PointTraversalCenter = PointTraversalForward. Переходим к шагу 2.2.
Шаг 3. Если значение IRegions (x, y) ф 1 то присваиваем значение LNBD = |IRegions (x, y)|. Переходим к
следующей координате точки увеличивая индекс столбца или строки и переходим к шагу 1. Алгоритм завершается по достижению правого нижнего угла изображения.
В результате работы алгоритма поиска контуров мы получаем иерархическую структуру, содержащую тип контура (контур объекта или контур отверстия) и массив точек, принадлежащий контуру. Для анализа объема примесей нам понадобятся внешние контуры объектов с индексом родителем 0.
Вычисление площади детектированных объектов. Для вычисления площади контура нам необходимо используя массив точек контура [4] вычислить площадь фигуры, ограниченную точками.
0 1 2 3 4 5 6 7 Я 9 10 и 12 13
1
2 ■ ■ а а
3 ■ ■ ■
4 ■ в
5 ■ в
6 ■ - • ■
7 ■ ■
К
10
Ро1п1 [х,у]— {
{5,2}, {4,3} ,{3,4},
{3,5},{4,6},{5,6},
{6,6М7,7М8,7},
{9,6},{10,5},{10,4}, {10,3}, {9,3}, {8,2}, {7,2}, {6,2}, {5,2}
Рис. 3. Соответствие контура объекта массиву вершин
Для вычисления площади требуется рассмотреть задачу определения площади фигуры [5] - зная координаты вершин, необходимо определить площадь фигуры Р1-Р4. Для вычисления площади необходимо вычислить площади четырех треугольников ОР1Р2, ОР2 Р3, ОР3 Р4, ОР4Р1 и произвести арифметические действия.
О[0,0]
Р|[Х],У1]
Р.1[х3,у3]
Р4|Х4,У4]
Рис. 4. Вычисление площади фигуры
Вычисление площади треугольников производится по правилу вычисления площади параллелограмма векторному произведению двух его сторон [6]:
' 1 ] к
= I
= \а х =
X У 0
Х2 У2 0
У1 0 х1 0 + к х1 У1 (4)
0 - у х2 0 =
У 2 Х2 У 2
= 1 • у1 • 0 -1 • у2 • 0 - у • х • 0 + у • х2 • 0 + к • х • У2 - к • х2 • у1 = = к (х1 у 2 - х2 У1)
где 1 — единичный вектор оси абсцисс, у — единичный вектор оси ординат, к — единичный вектор оси аппликат.
Рис. 5. Площадь параллелограмма
565
Так как наш контур Contours.Point[] находится в двухмерном вещественном плоскости □ . Единичный вектор k можно исключить и тогда площадь параллелограммы SParalldogram будет вычисляться следующим образом:
SParall'bgram =| xy2 - x2y1 | (5)
Так как диагональ параллелограмма делит его на два равных треугольника, то площадь треугольника STriangle, который нам и требуется вычислить, равна половине площади параллелограмма SParalldogram :
1 ,
ST,
= -I У2 - *2У1
(6)
Для вычисления площади фигуры, заданной/и точками контура Contours.Point[] расположенных массиве V = (v1,V2,...,vm) [7], воспользуемся формулой Гаусса [8] для определения площади многоугольника:
SCO"'alr = 1 * У 2 + Х2Уэ + • - + Ут + ХтУ1 - Х2 У - Х3 У 2 - • - - ХтУт-1 - X Ут[
£ ХУ+1 + ХтУ1 -£ ^+1 Ук + X1Уm
k=1 k=1
(7)
где m - количество точек контура Contours.Point[], SContour - площадь контура Contours.Point[].
Определение характерных геометрических размеров объектов. Следующим шагом будет вычисление объема фигуры, образованной контуром. Будем предполагать, что фигура примесей является либо эллипсоидом, либо шаром. Для вычисления объема эллипсоида и шара нам необходимо знать геометрические размеры (радиусы).
С
а'\
—г"
__\У J-
а'
у Where _ _
3
Kdipsu"J =-л-а-Ь-с = 3
4
_ ^clipgoid
Рис. 6. Объём фигуры
Для вычисления объемов фигур нам необходимо знать радиусы. В случае со сферой, все радиусы равны. В случае с эллипсоидом радиус вращения равен минимальному радиусу образующего эллипса - половине малой полуоси, большая полуось равна максимальному диаметру эллипса. Из массива точек контура Contours.Point[] необходимо определить точки такие, что расстояние [9] между этими двумя точками будет максимальным. Это будет большая полуось.
a = — max (R{Contours.Point[i],Contours.Point[у']))
1 2 '7--2-2 (8)
alk = — max y(Roint[i].x-Point[y'].x) + (Point[i].y-Point[y'].y) ,
где l, k - индексы точек начала и конца большой полуоси. Применение метод перебора всех точек
Contours.Point[] мы найдем точки большой полуоси за время O2.
После нахождения точек начала и конца большой полуоси мы можем найти точку центра контура
ОСопьоиг (она будет находиться в центре большой полуоси), используя следующие соотношения:
г^сопитг | РоШ[1 ].х + РоМ[к].х РоМ[1 ].у + РоМ[к].у О ~
Теперь мы можем найти радиус малой полуоси эллипса:
b = min (r (GContour, Contours.Point[i]j) b = min ^(Roint[i].x - GContour.x)2 + (Point[i].y - GContour.y)2
(9)
(10)
Вычисление большой полуоси, центра и малой полуоси приведёнными методами (8) - (10) требует полного перебора попарно всех точек массива контура. Пользуясь свойством, что Contours.Point[] заполнен последовательно точками контура в порядке их обхода, можно предложить метод нахождения большой и малой полуоси со сложностью порядка O1. Будем перебирать не все точки массива, а только находящиеся на расстоянии M/2. Найдем индексы l, k точек находящихся на максимальном удалении, и далее, для точек l, k проверим изменение
566
расстояния для соседних точек массива I +1, к +1. В случае получения полученного расстояния большего предыдущего переобозначим I, к новыми значениями и итерационно проверим новый максимум расстояния.
Точку центра контура G
будем определять по формуле (9), а затем для нахождения радиуса малой
полуоси будем проводить вычисление расстояния не для всех точек, а только для тех, что расположены в окрестностях точек s =i+, t=к+. В случае выхода индексов s, t точек контура за границу массива необходимо
скорректировать индексы, взяв остаток от деления на размер массива Contours.Point[]: s = s modM, t = t modM . Находим два радиуса до полуосей от центра контура до точек с индексами s, t и итеративно проверяем изменение расстояния от центра контура до соседних с s, t точками (s +1, t +1). В результате
получаем два радиуса, из которых следует выбрать наибольший.
В результате приведенных вычислений мы определили радиус
С _ b _ max( R (GContour, Contours. Point[s, t ])) вращения эллипсоида (он равен радиусу малой полуоси эллипса).
Вычисление объема примесей. Зная геометрические параметры контуров, мы можем вычислить объемы фигур (эллипсоидов):
Velipsoid _ gC
(11)
Bi oelipsoid
формуле площадь эллипса S , вычисляемая через геометрические размеры эллипса, заменена на
о Contour
численное вычисление площади контура объекта S .
Для вычисления объёмной концентрации примесей необходимо знать суммарный объем всех примесей и реальный объем кадра.
Зная геометрические размеры установки w х h х d (ширина, высота, глубина) и соответствующую им матрицу изображения W х H необходимо вычислить физический размер одной ячейки изображения.
Для соответствия размера кадра (установке) изображения реальному физическому размеру мы сначала должны вычислить DPM (Dots per meter) [10] цифрового изображения (матрицы), используя следующую формулу:
DPM _{wW (12)
Зная величину DPM , мы можем преобразовывать найденные нами цифровые величины длин, площадей и объемов к реальным размерам. Объем примесей Velipsoid в метрических размеров можно вычислить по формуле:
V^elipsoid _ veliPsoid , DPM3 (13)
Чтобы вычислить объёмную концентрацию нефтепродуктов в жидкости на основе снимков потока жидкости и соответствия размера кадра изображения реальному физическому размеру необходимо:
ZV eliP
real
psoid
•100%,
(14)
w • h • d
где (- объёмная процентная концентрация нефтепродуктов в жидкости.
Программная реализация и анализ результатов. Для изображения 1 .а определим иерархию замкнутых контуров границ объектов. На рис. 7.а. показаны все иерархии контуров. Выберем среди иерархии внешние контура (родитель Parent = 0) рис. 7.б.
Рис. 7. Контуры объектов
567
б
а
Среди найденных внешних контуров проведем вычисление их площади, а также характерных размеров
(рис. 8).
>156 J1>1S6
21,25 69) ГХ?в21,2570)
а б
Рис. 8. Площадь контура, радиус малой и большой полуоси эллипса
На левом рисунке вычисления производились полным перебором всех точек (8)-(10), на правом рисунке вычисления производились авторскими оптимизированными методами.
Используя метод полного перебора всех точек контура, а также авторский методы проведем вычисление объемной концентрации примесей. Результаты вычислений приведены на рис.9.
image Size: 1870x3062 Total piséis: 5725940 All Q$je¡Et detected = 63
ImfSurjjy Concentration = 0.31869064 59047308 % VOL
4* J á7 9 -SU2
4-3
Image Size: 1870x3062 Total pixels: 5725940 All aBjefít detected - 63
ImiSusj^ Concentration = 0.37676292029722847 % VOL
if „
é3
45
а б
Рис. 9. Вычисление объемной концентрации примесей
Различия в вычисленных концентрациях примесей для полного перебора всех точек (8)-(10) и авторским методом составляют 15%. Предложенный метод требует дополнительного уточнения данных вычислительного алгоритма с экспериментальными данными.
Заключение. Предложенный алгоритм позволяет определять контуры таким образом, чтобы они содержали только внешние границы объектов, вычислять геометрические размеры объектов, переводить размеры объектов в метрическую систему и вычислять объемную концентрацию примесей. В дальнейшем обработка серии фотоизображений (видеорядов) позволит создавать системы динамического контроля уровня загрязненности жидкости. По результатам работы возможно создание датчиков для определения количественного состава примесей.
Работа выполнена при поддержке Правительства Российской Федерации (гос. задание 0784-2020-0021)
Список литературы
1. Салех Л.О.А., Хлопин С.В., Черненькая Л.В., Тарасевский Ф.Г., Царев М. М., Алгоритм определения концентрации примесей в жидкости по оптическим данным. Известия Тульского государственного университета. Технические науки. 2023. № 1. С. 247-256.
2. Satoshi S. Computer Vision, Graphics, and Image processing N30, 32-46 (1985) Topological Structural Analysis of Digitized Binary Images by Border Following.
3. Burger W., Burge M. J., Digital Image Processing: An Algorithmic Introduction Using Java. Springer; 2nd ed. 2016. P. 210 - 211.
4. Kenny A. H. The Art of Image Processing with Java. A K Peters, 2016. P.311 - 317.
5. Геометрические олимпиады им. И.Ф. Шарыгина / Сост. А. А. Заславский, В. Ю. Протасов, Д. И. Шары-гин. М.: МЦНМО, 2007. С. 15 - 19.
6. Axler, Sheldon Jay. Linear Algebra Done Right (3rd ed.). Springer. 2015. P. 309 - 328.
568
7. Wilhelm Burger, Mark J. Burge. Principles of Digital Image Processing. Core Algorithms. Springer. 2009. P.
17 - 48.
8. Bart Braden. The Surveyor's Area Formula. The College Mathematics Journal, Volume 17, Number 4. September 1986. P. 326 - 337.
9. Козлов В. Н. Функциональный анализ (с приложениями в энергетике). СПб.: Издательство Политехнического университета. 2011. C. 36 - 43.
10. Sara Campos, Understanding DPI and Pixel Dimensions, iPrintfromHome.com, 2007 [Электронный ресурс]. URL: https://www.iprintfromhome.com/mso/understandingdpi.pdf (дата обращения: 18.07.2023).
Салех Лабиб Омар Абдулла, аспирант, labib_saleh_92@yahoo.com, Россия, Санкт-Петербург, Санкт-Петербургский политехнический университет Петра Великого,
Хлопин Сергей Владимирович, канд. техн. наук, доцент, x@spbstu.ru, Россия, Санкт-Петербург, Санкт-Петербургский политехнический университет Петра Великого,
Черненькая Людмила Васильевна, д-р техн. наук, профессор, ludmila@qmd.spbstu.ru. Россия, Санкт-Петербург, Санкт-Петербургский политехнический университет Петра Великого
ALGORITHM FOR CALCULATION OF THE VOLUME OF IMPURITIES IN A LIQUID FROM OPTICAL DATA
L.O.A. Saleh, S.V. Khlopin, L.V. Chernenkaya
The work is devoted to the development of methods and their software implementation for processing binary images, calculating the total volume of impurities in a liquid. For the calculation, the search for contours in a binary image, the selection of closed areas, the analysis of the hierarchy, and the calculation of the area and volume of the content of impurities are performed. An author's method is proposed for finding the characteristic geometric dimensions of impurity contours and using these dimensions to calculate the volume of impurities.
Key words: computer vision, finding contours, contour hierarchy analysis.
Saleh Labib Omar Abdullah, postgraduate, labib_saleh_92@yahoo.com. Russia, St. Petersburg, Peter the Great St. Petersburg Polytechnic University,
Khlopin Sergey Vladimirovich, candidate of technical sciences, docent, x@spbstu.ru, Russia, St. Petersburg, Peter the Great St. Petersburg Polytechnic University,
Chernenkaya Lyudmila Vasilievna, doctor of technical sciences, professor, ludmila@qmd.spbstu.ru. Russia, St. Petersburg, Peter the Great St. Petersburg Polytechnic University
УДК 004.42
DOI: 10.24412/2071-6168-2023-7-569-70
АЛГОРИТМ ОПЕРАТИВНОЙ ИДЕНТИФИКАЦИИ ДИНАМИЧЕСКИХ СИСТЕМ С ПЕРИОДИЧЕСКИМИ ВХОДНЫМИ СИГНАЛАМИ
Я.Н. Гусеница, М.И. Колоколов, А.А. Донченко
В работе предложен алгоритм оперативной идентификации динамических систем с периодическими входными сигналами. Алгоритм базируется на гипердельтной аппроксимации аналогов автокорреляционной и взаимно-корреляционной функций входного и выходного сигналов, которая предполагает вычисление начальных моментов для их аналогов плотностей распределения. Новый алгоритм позволяет вычислять начальные моменты только для одного периода и повысить оперативность расчетов в количество раз, равное периодам сигналов.
Ключевые слова: динамическая система, периодический сигнал, непараметрическая идентификация, автокорреляционная функция.
Идентификация представляет собой получение или уточнение математической модели исследуемой системы по ее измеряемым входным и выходным сигналам [1].
Важное место в идентификации отводится решению уравнения
да
Кх(0 = K(t-T)h(T)dx, 0)
0
где т - момент подачи входного x(t) сигнала относительно начала отсчета; Kyx(t) - взаимно-корреляционная функция входного x(t) и выходного y(t) сигналов; Kxx(t) - автокорреляционная функция входного x(t) сигнала; h(t) - неизвестная весовая функция, позволяющая при произвольной x(t) определять y(t).
569
