УДК 681.5.08
ОЦЕНКА ПАРАМЕТРОВ ЭЛЛИПТИЧЕСКОГО КОНТУРА
© 2015 Р.Р. Диязитдинов
Поволжский государственный университет телекоммуникаций и информатики, г.Самара
Поступила в редакцию 30.07.2015
В данной статье представлен вывод алгоритма оценки параметров эллиптического контура: значения полуосей, координаты центра и угол поворота. Данная задача возникает при контроле диаметра трубы и оценки ее отклонения от номинала. Исходными данными являлись уравнение эллипса и совокупность точек на плоскости. Алгоритм базируется на связи уравнений эллипса и координатами точек, которая описывается через матрицу преобразования (смещения и поворота). Решение сводится к формированию системы нелинейных уравнений, решаемой итерационным способом. Ключевые слова: оптический способ измерения, профильный сканер, эллиптический контур, линеаризация, итерационный алгоритм, нелинейное уравнение, метод Ньютона, метод простой итерации
ВВЕДЕНИЕ
Бесконтактный оптический способ измерения получил широкое распространение в различных технологических процессах и измерительных системах. Оптические устройства применяются для измерения расстояния - либо в виде дальномеров или в виде профильных сканеров, которые позволяют строить 2Б контур (профиль).
При производстве труб возникает вопросы контроля диаметра трубы и оценки ее отклонения от номинала. В современном производстве для этих целей используют профильные системы, которые позволяет измерять контур объекта, представляя его в виде совокупности точек (см. рис. 1 а). Измеренный контур трубы, что естественно, будет представлять собой эллипс (см. рис. 1 б).
Задача оценивания параметров эллипса интересна не только с практической, но с научной точки зрения. С точки зрения науки интерес представляет вывод алгоритма оценивания совокупности его параметров (значению полуосей, координат центра эллипса и угла поворота). А практическая значимость работы заключается в применении разработанного способа для автоматизации сбора параметров по контролю качества труб и конических поверхностей.
1. ПОСТАНОВКА ПРОБЛЕМЫ
Исследования в этой области, а именно оценивание параметров эллипса, проводилось в основном зарубежными авторами [1]-[7]. Основополагающей работой является [1], в которой автором предложен следующий подход: уравнение эллипса записывается в форме:
Диязитдинов Ринат Радмирович, кандидат технических наук, доцент кафедры системы связи. E-mail: [email protected]
Г (а, х) = а • х = ах2 + Ьху + + су2 + йх + еу + / = 0 где а = [а Ь с й е / ],
х = [2 ху у2 х у 1].
Проводится минимизация выражения
N
^(а) = ]ТГ (х г )2
1=1 ,
при котором выполняется условие а2 + 0,5Ь2 + с2 = 1 (алгебраическая инвариантная константа, верная для уравнения эллипса). После нахождения коэффициентов [а Ь с й е /], определяются канонические параметры эллипса - величина большой и малой полуоси, угол поворота и координаты центра.
В работах [2]-[7] алгоритмы строятся аналогичным образом, с той лишь разницей, что выбираются другие инварианты.
В данной статье предлагается новый способ к решению указанной проблемы. Он отличается более общим методологическим подходом при выводе алгоритма оценивания параметров, не предполагает наличие инвариантов и использует только аналитические уравнения, описывающие контур.
Основная идея заключалась в оценке параметров эллипса с помощью итерационных методов решения нелинейных уравнений, содержащие функции синуса и косинуса. Эта идея была подсказана устройством фазовой автоподстройки частоты (ФАПЧ). В частности, похожий прием был использован в работе [8] для оценивания параметров компенсатора искажений квадратурного амплитудного модулированного сигнала (КАМ-сигнала). Представленная работа является продолжением работы [9] и [10], где были ис-
а б
Рис. 1. К задаче оценивания параметров эллипса
пользованы похожие приемы для оценивания параметров контура на плоскости.
Эллиптические контуры представляют собой линии второго порядка - окружности и эллипсы, которые описываются уравнением вида:
т + — 2 b
■ = 1.
(1)
В ходе измерений производится оценка параметров контура, по которым определяется качество его изготовления - отклонение от номинальных параметров и дефекты изготовления.
Преобразование, связывающее между собой координаты эллипса (x, z) (см. формулу 1) и координаты контура, получаемого в профильной системе (xs, zs), имеет вид:
Гx = x • cos(/)-z • sin (/) + X, {z = x • sin (/)+ z • cos(/)+ Z. (2)
Задача оценивания состоит в оценке параметров: значению полуосей, координат центра эллипса и угла поворота по измеренным координатам (xs, zs).
2. ОПИСАНИЕ АЛГОРИТМА
Разработанный алгоритм состоит из двух последовательных этапов.
На первом этапе производится грубая оценка параметров a, b, X, Z, /, где a - большая полуось, b - малая полуось, (X, Z) - координаты центра эллипса, 3 - угол поворота.
На втором этапе производится точная оценка параметров a, b, X, Z, /. Этап точной оценки параметров позволяет почти на порядок уменьшить погрешность оценивания параметров.
Для вывода алгоритма в данной статье были использованы два метода:
- решение нелинейного уравнения методом Ньютона, которое сводится к решению уравнения f (x)= 0 итерационным способом
х - * -А) ■
^ г'{*„
- решение системы нелинейных уравнение методом простой итерации, которая сводится к следующему: для системы уравнений Г х -1 (х, z )-0 [ z - g (х, z )-0
Хи+1 - 1 (хп , zn ), zn+1 - g (п , zn ) .
, итерационное решение -
3. АЛГОРИТМ ОЦЕНИВАНИЯ КООРДИНАТ ЦЕНТРА ЭЛЛИПСА
Для вывода алгоритма оценки координат его центра будет использовано уравнение окружности с центром в точке (X, Z):
(xs - X)2 + (zs - Z )2 = R2, (3)
где xs, zs - координаты точек контура, R - радиус.
В уравнении (3) присутствует три неизвестных - R, X, Z, которые необходимо оценить.
Для решения использовался метод наименьших квадратов, и для этого была использована функция, которую необходимо минимизировать:
f (X, Z, R) = X К - X )2 + (zs, - Z )2 - R2 f ^ min, (4)
г=1
где xs, zs - координаты точек контура, N - количество точек контура.
Наиболее очевидный способ поиска минимума функции - это нахождения частных производных по искомым переменным и сведения его к системе нелинейных уравнений. Однако это способ не дал практических результатов - итерационный способ решения уравнения приводит к некорректным результатам (не выполняется условие сходимости).
Вместо этого использовался следующий прием.
При фиксированном значении R=r можно оценить параметры X, Z:
2
2
f (X, Z) = £ [ - X )2 + (zs, - Z)2 - г2 f min, (5)
/=1
Решение определяется системой нелинейных уравнений, найденных путем дифференцирования выражение (5) по параметрам (X, Z) :
А2 - X3 + B2 - X2 + С2 • X + В2 = 0; Аъ - Z3 + B3 • Z2 + C3 • Z + Д = 0.
(6)
где А2 = 4N, B2 =-12
xs
C2 = X {12 -
xs2 - 4 - r2 +
4(Z - zs,, )2}
/=i
D2 = £ {-4-xs, ( -r2 + (Z-zs,)2)}
A3 = 4N B3 =-12-£
zs.
C3 = £ {12 - zs,2 - 4 - r2 + 4( -zs,)2} ,=1
D3 = X {- 4 - zs, (xs,2 - r2 + (x - xs, )2)}.
,=1
N - количество точек в контуре (профиле). Итерационное решение имеет вид:
x А2 - x" + B2 - x2 + D2
^ n+1 C
C2
Z = A3 - z" + B3 - z2 + D3
Zn+1 C
(7)
При фиксированном значении X = x0, Z = z0 можно оценить параметр R:
f (R) = X [(xs - x0)2 + (zs, - z0)2 - R2 f min .(8)
/=1
Решать уравнение (9) будем искать в следующем виде:
q(R)= f (R)- f (R + S)* 0, (9)
где S - малая величина (в расчетах было использовано S = 106). Тогда
q(Rn)
R, = R - * \ = R -
f (R )-f (R + S)
"+1 "" q'(R„) "" f'(R)-f(R+*)' (10)
f'(R) = -4R -X {(xs,. - x0) + (zs, - z0)2 - R2},(И)
Объединяем описанные выражения в единую процедуру оптимизации. Алгоритм имеет следующий вид:
1. Выбираем первоначальные значения Xо,Zо,R0 и п=0.
2. Определяем Хй+1, Zn+1 при г = Rn по формуле (8).
3. Вычислить /R ) /'(Rn ) при x0 = Хи+1, z0 = Zn+1 по формулам (8) и (11).
4. Вычислить /((n +j), /'((п +S) при x0 = Xn+1, z0 = Zn+1 по формулам (8) и (11).
5. Вычислить Rn+1 по формуле (10).
6. n = n + 1 и перейти к шагу №2.
Переходов от шага №6 к шагу №2 определяет
количество итераций алгоритма - n.
На выходе работы алгоритма имеется грубая оценка параметров, который обозначим как
(X' , Z').
Зная оценку координат центра (X', Z'), можно оценить параметры - a, b, р. Для этого рассчитывается расстояние от (X', Z') до точек контура (xs., zs.) и угол, который формируют прямые, проходящие через (X', Z') и (xs., zs.):
Г = - X)2 +(zsi - Z )2,
a, = atan! ((zs, - Z), (( - X)) .(12)
Максимальное значение max(r.) будет являться грубой оценкой параметра af, минимальное значение min(r¿) - Ъ. В направлении, где Г. достигает своего максимума, будет определять угол поворота р. Для подтверждения этого вывода на рис. 2 представлена геометрическая иллюстрация.
4. ТОЧНАЯ ОЦЕНКА ПАРАМЕТРОВ
Последний этап - это точное оценивание параметров - a, b, X, Z, р. На этом этапе одновременно используется поиск по сетке (фиксируется значения параметров a, b), также итерационная процедура - по параметрам X, Z, р.
-60 -40 -20 0 20 40 60 80 100
Рис. 2. К оценке параметров величины полуосей и угла поворота
,=1
,=1
,=1
,=1
Преобразуем координаты (xs.,zsi) в ((1., zs\i) согласно найденным J, (XZ'):
( xs1t} zs1i
V1
(cos(j) - sin(JJ') 0^
V
sin J) 0
cos(jJ') 0
(1 0 - X 4 (xs. ^
0 1 - Z'
V00 1 J l1 У
zs,.
(17)
A - AJ7 + B -A J6 +Cj- AJ5 + Dj • AJ4 + + £1 • AJ3 + F-AJ + Gj-AJ + H = 0, A2 - AX3 + B2 - AX2 + C2 - AX1 + D2 = 0, A3 -AZ3 + B3 -AZ2 + C3 -AZ1 + D3 = 0. Итерационное решение имеет вид:
(20)
(13)
AJ
A - A J + Bj - AJ + Cj - A J + Dj - A J + £ - A Pi + F - AJj + Hl
Обозначим как (AX, AZ ) и A/ - поправки к центру эллипса и углу разворота. Координаты точек (xs1¿, zs1;) связаны с координатами точек «идеального» контура эллипса (, zi) через преобразование:
Г xs = x • cos(A/)- z • sin (A/)+AX, [zs = x• sin(A/) + z • cos(A/)+AZ. (14) Тогда:
Гx = (xs - X0) cos(A/)+(zs - AZ)• sin (A/), [z = -(xs - X0)sin(A/)+(zs - AZ)-cos(A/). (15) Подставим (15) в выражение (1), получим:
((xs - X 0) cos(/)+(zs - Z 0) sin (3))2 b2 + (- (xs - X0) sin(/)+(zs - Z0) cos(/))2a2 - (16)
- a2b2 = 0
Полагая, что A/ мал (- 5 < A/ < 5 градусов), заменим функции косинуса и синуса на их линеаризованные выражения:
sin (A/)^A/, cos(A/)~ 1 - 0,5A/2.
G
= A2 -AXi + B2 -AXj + D2
AX1+1 =--с-:
C 2
= AA3Z + B3 - Z + D3
AZ"+1 = C3 •
(21)
(22)
Выражения для коэффициентов A1, B1 и т.д. в формулах (21) - (22) представляют собой полином с большим числом переменных. Для их вычисления в аналитическом виде целесообразно использовать математические пакеты (например, Matlab, MathCad и т.д.).
Так как в данном алгоритме используется поиск по сетке параметров a, b, то каждому паре этих значений ставится в соответствие вычисленные параметры (АХ, AZ ), A3 и метрика:
metric = ]Г( 1-АХ)cos(A3) + ( 1-AZ)sin(A3))2-b2 +
¡=1
(-( 1i -AX )• sin(A3) + ( 1i-AZ )• cos(A3))2- a2 - a 2b2 }2
(23)
Таким образом, отыскивая минимум метрики (28), будут оценены параметры a, b. А используя найденные значения (AX, AZ ) , A3 , определятся параметров X, Z, 3, которые связаны через следующую матрицу:
Подставляя (19) в (18) получаем: ((1 - ДХ) • (1 - 0,5 АЗ2 )+(1 -М) • ДЗ)2 62 + (- (1 - АХ) • З + (1 - А2) • ^(1 - 0,5АЗ2 ))2 а2 -
- а2Ь2 - 0 .(18)
При фиксированных значениях a, Ь, можно найти параметры (АХ, Д2 ) и ДЗ , для этого используем метод наименьших квадратов, минимизирую функцию:
1 (АХ, Д2, З) - |ж1, - АХ) ■ (1 - 0,5 ДЗ2)+ (\-Д2) ■ Дз)2 62 +
1-1
(- (ж1, - АХ)-ДЗ+(.■¡1,-62)■ с^- 0,5 ДЗ2)) а2 - а2Ь2 ^ .(19)
Решение определяется системой нелинейных уравнений, найденных путем дифференцирования выражение (19) по искомым параметрам
(ДХ, м, З).
После приведения выражений к единому виду, можем записать:
f cos(J) - sin(J) X^ sin(J) cos(J) Z 0 0 1
cos(AJ) - sin(AJ) AXy1 sin(AJ) cos(AJ) AZ 0 0 1
V
f cos(J') - sin(J') 0 y1
sin(J')
0
cos(J') 0
(
1 0 - X 0 1 - Z' 0 0 1
Л
(24)
ЗАКЛЮЧЕНИЕ
В статье представлен алгоритм оценивания параметров эллиптического контура: значения полуосей, координат центра и угла поворота. Данная задача была сведена к системе нелинейных уравнений, решаемой итерационным способом. Результаты работы планируется использовать для автоматизации сбора параметров по контролю качества труб и конических поверхностей, конту-
X
X
X
-1
X
х
ры которых измеряются с помощью профильных систем. Дальнейшие исследования вопросов совмещения и оценки параметров будут связаны со сложными контурами, отдельные части которых описываются различными уравнениями.
СПИСОК ЛИТЕРАТУРЫ
1. Bookstein F.L. Fitting Conic Sections to Scattered Data // Computer Graphics and Image Processing. 1979, № 9. P. 56-71.
2. Ellipse Detection and Matching With Uncertainty / Т. Ellis, A. Abbood, B. Brillault // Image and Vision Computing.1992, № 2. P. 271-276.
3. Least-Square Fitting of Circles and Ellipses / W. Gander, G.H. Golub, R. Strebel// BIT. 1994, № 43. P. 558-578.
4. Rosin P.L. A Note on the Least Squares Fitting of Ellipses // Pattern Recognition Letters. 1994. № 14. P. 799-808.
5. Rosin P.L. West G.A. Nonparametric Segmentation of Curves Into Various Representations // IEEE Trans.
Pattern Analysis and Machine Intelligence. 1995, № 12. P. 1140-1153.
6. Sampson P.D. Fitting Conic Sections to Very Scattered Data: An Iterative Refinement of the Bookstein Algorithm // Computer Graphics and Image Processing. 1982, 18. P. 97-108.
7. Taubin G. Estimation of Planar Curves, Surfaces and NonPlanar Space Curves Defined by Implicit Equations, With Applications to Edge and Range Image Segmentation // IEEE Trans. Pattern Analysis and Machine Intelligence. 1991, № 11. P. 1115-1138.
8. Поборчая Н.Е. Анализ работы компенсатора искажений КАМ-сигнала, наблюдаемого на фоне аддитивного шума // Электросвязь. 2014, № 5.C. 20-25.
9. Диязитдинов Р.Р. Оценивание параметров положения контура кривой в профильной системе // Инфокоммуникационные технологии. 2014, № 2. C. 70-73.
10. Васин Н.Н., Куринский В.Ю. Расширение функциональных возможностей систем видеонаблюдения // Информационные технологии. 2013, № 6. С. 63-66.
PARAMETERS ELLIPSE'S COUNTER RATING
© 2015 R.R. Diyazitdinov
Povolzhskiy State University of Telecommunications and Informatics, Samara
In this article is described developing algorithm for parameters ellipse counter's rating: the value of semimajor axis and the semi-minor axis, coordinate's centre, rotate. This problem turns out for control tube's diameter and rating deviation from the nominal value. Initial data are ellipse's equation and sequence points on a plane. Algorithm bases on connection ellipse's equation and sequence points, which is showed by transform matrix (offset and rotate). The solution is found as nonlinear equation's system, which is solved iteration method.
Keywords: optical method of measurement, profile scanner, ellipse's counter, linearization, iteration algorithm, nonlinear equation, Newton's method, simple iteration's method.
Rinat Diyazitdinov, Candidate of Technics, Associate Professor at the of Communications System Department. E-mail: [email protected]