Известия РАН. Теория и системы управления, 2023, № 1, стр. 148-163

МЕТОД СОГЛАСОВАНИЯ ПРИБОРНЫХ И СВЯЗАННОЙ СИСТЕМ КООРДИНАТ НА ЛЕТАТЕЛЬНОМ АППАРАТЕ

В. М. Лисицын a*, К. В. Обросов a, Г. Г. Себряков a

a ФАУ “ГосНИИАС”
Москва, Россия

* E-mail: lvm@gosniias.ru

Поступила в редакцию 01.08.2022
После доработки 12.08.2022
Принята к публикации 26.09.2022

Полный текст (PDF)

Аннотация

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

Введение. Современный многофункциональный летательный аппарат (ЛА) (самолет или вертолет) не может обходиться без обзорно-прицельных оптико-электронных систем (ОЭС). Их преимуществом является возможность работы в пассивном режиме для решения задач обеспечения ситуационной осведомленности экипажа. ОЭС часто используются для решения задач мониторинга подстилающей поверхности (ПП) и обнаружения наземных объектов. ОЭС может состоять из нескольких разнесенных по корпусу и крыльям ЛА постов (подсистем). Все посты ОЭС должны работать в единой системе координат (СК), согласованной со связанной СК, которая задается инерциальной системой ЛА. Выполнение этого требования необходимо для обеспечения адекватного обмена информацией об угловых координатах объектов информационными системами между собой и комплексом бортового оборудования. Отсюда очевидна актуальность задачи согласования приборных СК разных постов ОЭС между собой и со связанной СК.

В настоящее время согласование СК осуществляется на основе так называемой статической юстировки, которая производится на заводах-изготовителях при первоначальной установке всех информационных систем на ЛА, а также в процессе эксплуатации при замене оборудования и после регламентных работ. Статическая юстировка заключается в следующем: вывешивая ЛА на домкратах, устанавливают его в строго горизонтальное положение и определяют положение строительной оси ЛА относительно направления на север. После этого производится выставка бесплатформенной инерциальной навигационной системы (БИНС) и согласование осей БИНС со связанной СК ЛА. Для юстировки ОЭС применяются юстировочные мишени. С помощью теодолита, устанавливаемого на специальных узлах крепления ЛА, выставляется центральная точка мишени, соответствующая направлению строительной оси. Положение реперных точек, соответствующих направлению горизонтальных осей приборных СК, определяется по заданным в чертежах ЛА смещениям ОЭС относительно связанной СК. Для юстировки тепловизионных систем в реперных точках устанавливаются излучатели. Юстировка производится путем получения изображения реперных точек и определения их угловых координат в приборной СК (угловых поправок), которые фиксируются в памяти вычислителя и применяются для корректировки измеренных угловых координат объектов в полете при переходе от приборной СК к связанной и от связанной СК к приборной. С помощью статической юстировки взаимное согласование СК ОЭС и БИНС производится не непосредственно между ними, а относительно строительных осей ЛА.

Такие методы хорошо отработаны на практике и используются как основной способ первоначальной юстировки ОЭС и БИНС. Существенным недостатком методов статической юстировки является их высокая трудоемкость, необходимость применения специального оборудования квалифицированным персоналом и трудности проведения в условиях аэродромов при эксплуатации. Кроме того, такие методы не позволяют компенсировать ошибки, которые возникают из-за различия изгибных деформаций корпуса ЛА в полете и при размещении его на домкратах.

Поиск способов, позволяющих отказаться от предполетной юстировки по наземной мишени, привел к разработке метода так называемой динамической юстировки СК на основе векторного согласования. Под ним понимается измерение трех ортогональных проекций как минимум двух неколлинеарных векторов в обеих СК и определение их взаимной ориентации по разности одноименных проекций [1]. При движении ЛА в качестве векторов могут быть использованы данные об угловой скорости движения ЛА, определяемой с помощью датчиков угловых скоростей (ДУС), и измерение суммарного значения линейных ускорений ЛА и вектора силы тяжести Земли с помощью измерителей линейных ускорений (ИЛУ). Выбор векторов, по которым будет производиться согласование СК, может быть различным в зависимости от множества факторов, включая возможности ЛА, особенности оборудования и пр. В частности, это могут быть векторы ускорений и угловых скоростей; векторы скоростей, проинтегрированных скоростей; вектора, задающего угловое положение БИНС, а также варианты, включающие комбинации указанных векторов. При этом необходимо на каждом посту ОЭС устанавливать приборные инерциально-измерительные блоки (ИИБ), оборудованные ДУС и ИЛУ. Маневры носителя и нежесткость конструкции приводят к тому, что приборный ИИБ ОЭС движется относительно БИНС, что требует коррекции при сопоставлении измерений. Результаты измерений, как правило, подвергаются фильтрации Калмана [2]. Выбор используемых измерений в математической модели фильтрации Калмана во многом определяет и маневры ЛА, на которых будет происходить динамическая юстировка СК.

В отношении согласования СК ОЭС и БИНС очевидно, что динамическая юстировка требует оснащения всех постов ОЭС ИИБ, что приводит к увеличению массогабаритных характеристик и стоимости, а в некоторых случаях это невозможно в силу ограничений по габаритам. Поэтому в настоящее время поиск путей решения проблемы выставки СК не потерял актуальность [3]. Авторами предлагается новый подход к согласованию приборной СК ОЭС с СК БИНС, который не имеет указанных недостатков.

1. Постановка задачи. Определим системы координат:

начало координат связанной СК находится в центре масс ЛА;

ось OX направлена вдоль строительной оси;

ось OZ перпендикулярна плоскости симметрии ЛА и направлена в его правую часть;

ось OY перпендикулярна плоскости ХОZ и дополняет СК до правой тройки.

Оси приборной СК каждого поста ОЭС имеют направление аналогично связанной СК. В идеале оси приборных СК коллинеарны осям связанной СК. Однако в силу объективных причин, как правило, имеется некоторое рассогласование между направлениями осей. Введем также СК кадра, связанную с датчиком изображения. Центр изображения, соответствующий главной оптической оси датчика и направлению OX приборной СК, принимается за начало координат. Ось OZ направлена вправо, ось OY – вверх.

Предлагаемый метод решения проблем выставки СК заключается в том, что векторами, подлежащими согласованию, являются векторы скорости ЛА: линейной и угловой, проекции которых на оси юстируемых ОЭС измеряются в результате независимой обработки изображений, формируемых этими системами, и сопоставления результатов с информацией БИНС. Измерение ориентации указанных векторов производится не одновременно: измерение ориентации вектора линейной скорости производится при отсутствии угловой, а измерение ориентации вектора угловой скорости – при отсутствии линейного перемещения ЛА.

Рассмотрим случай, когда ЛА совершает прямолинейный горизонтальный полет на относительно небольшой высоте и в поле зрения ОЭС с тепловизионным или телевизионным каналом находится ПП. Если на изображении ПП имеются характерные точки, которые могут быть выделены в результате обработки, то такие точки при стабилизированном в инерциальном пространстве поле зрения ОЭС в каждый момент времени двигаются на изображении по лучам, расходящимся из одной точки (точки схождения). Эта точка является пересечением вектора скорости ЛА с фокальными плоскостями датчиков изображений независимо от их ориентации в пространстве. Другими словами, положение точки схождения на изображении соответствует направлению вектора скорости ЛА в приборной СК (связанной с полем зрения ОЭС). Для определения координат точки схождения необходимо сопровождать хотя бы две характерные точки. Если предположить, что вектор линейной скорости совпадает с направлением строительной оси ЛА (углы атаки и скольжения равны нулю), то положение рассматриваемой точки схождения относительно центра изображения определяет рассогласование систем координат относительно вертикальной и поперечной осей связанной СК (т.е. угловые поправки по курсу и тангажу) при условии отсутствия крена. Отметим, что рассмотренная процедура может быть выполнена на любом типе ЛА. Для определения параметров таких траекторий и положения точки схождения в [4, 5] был использован метод двух засечек на вложенных интервалах времени с усреднением. В данной работе предложен более адекватный подход, основанный на методах оптимального оценивания.

Для определения рассогласования СК вертолета относительно продольной оси связанной СК (т.е. по крену) предлагается использовать специальный маневр – вращение вертолета, в результате которого измеряется ориентация вектора угловой скорости вращения: параметры углового вращения формируется в БИНС, а оценка направления вектора относительно датчиков изображений производится в результате обработки видеопотока в каждом посту ОЭС. При этом, как уже отмечалось, предполагается отсутствие линейного перемещения аппарата, что при наличии РЛС на борту вертолета обеспечивается с высокой точностью путем идентификации и оценки дальностей до элементов аэродромной инфраструктуры при вращении вертолета [68].

Рассмотрим ситуацию, когда вертолет после взлета зависает на небольшой высоте и совершает вращение вокруг вертикальной оси OY связанной СК. При этом аналогично рассмотренному случаю в поле зрения ОЭС находится ПП и элементы инфраструктуры аэродрома. Будем полагать, что поле зрения ОЭС стабилизировано и в течение всего маневра положение центра масс вертолета и направление вектора угловой скорости не изменяется. На изображении ПП и инфраструктуры аэродрома выделенные характерные точки при вращении вертолета будут перемещаться по экрану ОЭС по разным траекториям. Требуется по результатам идентификации траекторий характерных точек изображений ПП определить направление вектора угловой скорости вертолета в приборной СК ОЭС. Методы детектирования и сопровождения характерных точек будут рассмотрены далее.

2. Метод измерения ориентации вектора линейной скорости ЛА. При прямолинейном полете и компенсации стабилизационных колебаний ЛА системой стабилизации ОЭС, как уже отмечалось, характерные точки ПП на изображении будут двигаться по лучам, исходящим из точки схождения, положение которой соответствует направлению вектора скорости. При этом не имеет значения, в какие промежутки времени формируются траектории точек. Другими словами, траектории точек могут регистрироваться и одновременно, и последовательно, а положение точки схождения должно определяться после набора заданного количества таких точек, что повышает точность определения ее координат. Хотя формально для нахождения положения точки схождения достаточно сопровождать лишь две точки. В процессе сопровождения характерных точек могут наблюдаться ошибки определения текущего положения точек вследствие стохастического характера работы алгоритмов в условиях изменения локальных контрастов и масштабов фрагментов изображения в процессе сближения ЛА сопровождаемыми участками ПП. Таким образом, полученные при сопровождении какой-либо точки координаты, вообще говоря, не лежат на одной прямой. Для оценки ее параметров необходимо аппроксимировать полученные координаты прямой линией. Считая, что ошибки измерения текущих координат точки имеют нулевое математическое ожидание, а сами измерения некоррелированные, аппроксимацию можно осуществить с помощью метода наименьших квадратов [9]. Учитывая, что траектории движения сопровождаемых точек могут иметь всевозможные направления в экранной СК, в том числе и вертикальное, уравнение прямой необходимо задать в общем виде, которое может быть сведено к следующему:

$Ny + Mz = 1.$

Пусть за время движения точки по полю зрения ОЭС произведено n измерений ее координат (y1, z1), (y2, z2), …, (yn, zn). Тогда может быть составлена система линейных уравнений относительно коэффициентов уравнения прямой:

$\left\{ \begin{gathered} N{{y}_{1}} + M{{z}_{1}} = 1, \hfill \\ N{{y}_{2}} + M{{z}_{2}} = 1, \hfill \\ ...................... \hfill \\ N{{y}_{n}} + M{{z}_{n}} = 1. \hfill \\ \end{gathered} \right.$

Метод наименьших квадратов состоит в минимизации следующей суммы:

$\sum\limits_{i = 1}^n {(N{{y}_{i}}} + M{{z}_{i}} - 1{{)}^{2}} \to \mathop {\min }\limits_{N,\;M} .$

Взяв производные по параметрам N и M и приравняв их к нулю, получим систему

$\left\{ \begin{gathered} N\sum\limits_{i = 1}^n {y_{i}^{2}} + M\sum\limits_{i = 1}^n {z_{i}^{{}}{{y}_{i}}} = \sum\limits_{i = 1}^n {{{y}_{i}}} , \hfill \\ N\sum\limits_{i = 1}^n {y_{i}^{{}}z_{i}^{{}} + M\sum\limits_{i = 1}^n {z_{i}^{2} = \sum\limits_{i = 1}^n {{{z}_{i}}.} } } \hfill \\ \end{gathered} \right.$

Отсюда находим значения коэффициентов, задающих прямую:

(2.1)
$N = \frac{{\sum\limits_{i = 1}^n {{{y}_{i}}\sum\limits_{i = 1}^n {z_{i}^{2} - \sum\limits_{i = 1}^n {{{z}_{i}}\sum\limits_{i = 1}^n {{{z}_{i}}{{y}_{i}}} } } } }}{{\sum\limits_{i = 1}^n {y_{i}^{2}\sum\limits_{i = 1}^n {z_{i}^{2} - {{{\left( {\sum\limits_{i = 1}^n {{{z}_{i}}{{y}_{i}}} } \right)}}^{2}}} } }},\quad M = \frac{{\sum\limits_{i = 1}^n {{{z}_{i}}\sum\limits_{i = 1}^n {y_{i}^{2} - \sum\limits_{i = 1}^n {{{y}_{i}}\sum\limits_{i = 1}^n {{{z}_{i}}{{y}_{i}}} } } } }}{{\sum\limits_{i = 1}^n {y_{i}^{2}\sum\limits_{i = 1}^n {z_{i}^{2} - {{{\left( {\sum\limits_{i = 1}^n {{{z}_{i}}{{y}_{i}}} } \right)}}^{2}}} } }}.$

Если взять произвольную пару прямых с номерами  j и k, то точка их пересечения, т.е. оценка положения точки схождения, определяется на основании решения системы

$\left\{ \begin{gathered} {{N}_{j}}y + {{M}_{j}}z = 1, \hfill \\ {{N}_{k}}y + {{M}_{k}}z = 1. \hfill \\ \end{gathered} \right.$

Координаты точки схождения Z, Y, найденные с использованием  j-й и k-й траекторий, удобно обозначить соответственно zjk, yjk с индексным указанием рассматриваемых прямых. Очевидно,

(2.2)
${{z}_{{jk}}} = \frac{{{{N}_{j}} - {{N}_{k}}}}{{{{M}_{j}}{{N}_{k}} - {{M}_{k}}{{N}_{j}}}},\quad {{y}_{{jk}}} = \frac{{{{M}_{j}} - {{M}_{k}}}}{{{{M}_{j}}{{N}_{k}} - {{M}_{k}}{{N}_{j}}}}.$

Для уменьшения погрешности определения положения точки схождения (Y, Z) следует произвести усреднение полученных координат по всем парам прямых. Если таких прямых m, то количество пар составляет $Q = C_{m}^{2} = m(m - 1){\text{/}}2$. Однако, как будет показано ниже, чем меньше острый угол γ между прямыми, тем больше ошибка определения положения точки схождения. Поэтому целесообразно проводить весовое суммирование и использовать в качестве весовых коэффициентов величину $\sqrt {1 - \cos \gamma } $ или, что почти эквивалентно, значение модуля синуса угла γ. Формула для определения усредненных значений координат $\overline Z $, $\overline Y $ в этом случае принимает вид

(2.3)
$\overline Z = \frac{{\sum\limits_{i = 1}^{m - 1} {\sum\limits_{j = i + 1}^m {{{z}_{{ij}}}\left| {\sin {{\gamma }_{{ij}}}} \right|} } }}{{\sum\limits_{i = 1}^{m - 1} {\sum\limits_{j = i + 1}^m {\left| {\sin {{\gamma }_{{ij}}}} \right|} } }},\quad \overline Y = \frac{{\sum\limits_{i = 1}^{m - 1} {\sum\limits_{j = i + 1}^m {{{y}_{{ij}}}\left| {\sin {{\gamma }_{{ij}}}} \right|} } }}{{\sum\limits_{i = 1}^{m - 1} {\sum\limits_{j = i + 1}^m {\left| {\sin {{\gamma }_{{ij}}}} \right|} } }},$
где ${{\gamma }_{{ij}}}$ – угол между прямыми с номером i и j, который может быть найден [10] как

${{\gamma }_{{ij}}} = \arccos \frac{{{{N}_{i}}{{N}_{j}} + {{M}_{i}}{{M}_{j}}}}{{\sqrt {N_{i}^{2} + M_{i}^{2}} \sqrt {N_{j}^{2} + M_{j}^{2}} }}.$

Таким образом найдены оценки рассогласований связанной и приборной СК относительно осей OY и OZ, т.е. обусловленных поворотом СК вокруг этих осей (обозначим их Δφ(Y) и Δφ(Z)). При отсутствии рассогласований осей относительно оси OX, т.е. при Δφ(X) = 0, полученные величины полностью характеризуют рассогласование СК. Однако при Δφ(X) ≠ 0 это не так.

3. Метод измерения ориентации вектора угловой скорости вертолета. Метод предполагает стабильное положение начал всех приборных систем координат вращающегося вертолета относительно местности. Прецизионное выполнение этого требования гарантируется путем радиолокационного дальнометрирования элементов инфраструктуры и управления положения вертолета [68]. При вращении вертолета анализ возможных траекторий движения изображения точек ПП по экрану показывает, что точки, находящиеся в плоскости, перпендикулярной вектору угловой скорости, и проходящей через главную точку объектива, будут двигаться по прямой. Траектории остальных точек будут представлять собой кривые. Чтобы определить тип этих кривых, рассмотрим, как отображаются элементы ПП на экране ОЭС при вращении вертолета. Воспользуемся принципом обращения движения, т.е. представим, что не вертолет вращается относительно оси OY, а ПП вращается вокруг вертолета. Тогда любая прямая, соединяющая точку на ПП и вектор угловой скорости, формирует относительно оси вращения коническую поверхность. Для прямых, пересекающих матричное фотоприемное устройство (МФПУ), фокальная плоскость ОЭС по сути является фрагментом секущей плоскости конуса, а угол между осью конуса, т.е. между вектором угловой скорости и секущей плоскостью, – углом отклонения приборной СК от связанной СК по тангажу. Поскольку этот угол, как правило, мал, то характерные точки ПП, за исключением точек, лежащих в упомянутой выше плоскости, будут двигаться на изображении по гиперболам. Точки могут двигаться по разным гиперболам, но все гиперболы имеют общий центр. Анализ траекторий характерных точек на изображении ПП показывает, что при наличии рассогласования приборной и связанной СК по тангажу центр будет соответственно перемещаться по вертикали. Наличие бокового смещения ОЭС относительно строительной горизонтали фюзеляжа будет приводить к смещению центра гипербол в боковом направлении. Если есть смещение по крену, то гиперболы будут повернуты относительно приборной СК. Таким образом, в общем случае в уравнениях гипербол будут присутствовать все члены, и любое уравнение может быть сведено к следующему виду:

(3.1)
$A{{z}^{2}} + 2Bzy + C{{y}^{2}} + 2Dz + 2Ey = 1.$

Чтобы определить отклонение приборной СК от связанной по тангажу и крену, достаточно вычислить коэффициенты уравнения (3.1) и найти такую СК (O1Z1Y1), в которой уравнение гиперболы имеет канонический вид. При этом приводить уравнение к каноническому виду нет необходимости: положение искомой СК по вертикальной оси приборной СК определяет отклонение по тангажу, а угол поворота соответствует отклонению по крену. Положение начала координат новой СК (z0, y0) может быть найдено с помощью решения системы линейных уравнений [10]:

(3.2)
$\left\{ \begin{gathered} A{{z}_{0}} + B{{y}_{0}} + D = 0, \hfill \\ B{{z}_{0}} + C{{y}_{0}} + E = 0. \hfill \\ \end{gathered} \right.$

Угол поворота – с помощью выражения

(3.3)
$\alpha = \frac{1}{2}\operatorname{arctg} \frac{{2B}}{{A - C}}\,.$

Для определения коэффициентов уравнения (3.1) необходимо измерить координаты сопровождаемой характерной точки не менее чем 5 раз за время ее движения по полю зрения ОЭС. Если количество измерений более пяти, то коэффициенты могут быть найдены с помощью метода наименьших квадратов [9]. Пусть за время движения точки по полю зрения ОЭС произведено n ≥ 6 измерений ее координат (z1, y1), (z2, y2), …, (zn, yn). Тогда может быть составлена система линейных уравнений относительно коэффициентов уравнения гиперболы:

(3.4)
$\left\{ \begin{gathered} Az_{1}^{2} + 2B{{z}_{1}}{{y}_{1}} + Cy_{1}^{2} + 2D{{z}_{1}} + 2E{{y}_{1}} = 1, \hfill \\ Az_{2}^{2} + 2B{{z}_{2}}{{y}_{2}} + Cy_{2}^{2} + 2D{{z}_{2}} + 2E{{y}_{2}} = 1, \hfill \\ ............................................................ \hfill \\ Az_{n}^{2} + 2B{{z}_{n}}{{y}_{n}} + Cy_{n}^{2} + 2D{{z}_{n}} + 2E{{y}_{n}} = 1. \hfill \\ \end{gathered} \right.$

Обозначим через M матрицу, составленную из соответствующих измерений:

${\mathbf{M}} = \left( \begin{gathered} \begin{array}{*{20}{c}} {z_{1}^{2}}&{2{{z}_{1}}{{y}_{1}}}&{y_{1}^{2}}&{2{{z}_{1}}}&{2{{y}_{1}}} \\ {z_{2}^{2}}&{2{{z}_{2}}{{y}_{2}}}&{y_{2}^{2}}&{2{{z}_{2}}}&{2{{y}_{2}}} \end{array} \hfill \\ ............................. \hfill \\ \begin{array}{*{20}{c}} {z_{n}^{2}}&{2{{z}_{n}}{{y}_{n}}}&{y_{n}^{2}}&{2{{z}_{n}}}&{2{{y}_{n}}} \end{array} \hfill \\ \end{gathered} \right),$
а через f – вектор, компонентами которого являются коэффициенты уравнения гиперболы:

${{{\mathbf{f}}}^{{\text{Т}}}} = \left( {A,B,C,D,E} \right).$

Тогда система (3.4) примет вид Mf = I, где I – вектор размерности n, все элементы которого равны единице.

Можно считать все измерения равноточными и некоррелированными. Необходимо найти минимум суммы квадратов разностей левой и правой частей этого уравнения, т.е. ${{({\mathbf{Mf}}\, - \,{\mathbf{I}})}^{{\text{Т}}}}({\mathbf{Mf}}\, - \,{\mathbf{I}})\, \to \,\min .$ Дифференцируя эту функцию по вектору параметров и приравняв производную к нулю, получаем следующее матричное уравнение:

(3.5)
${{{\mathbf{M}}}^{{\text{Т}}}}{\mathbf{Mf}} = {{{\mathbf{M}}}^{{\text{Т}}}}{\mathbf{I}}.$

Выполнив умножение матриц, получаем (3.5) в виде

$\left( {\begin{array}{*{20}{c}} {\sum {z_{i}^{4}} }&{2\sum {z_{i}^{3}{{y}_{i}}} }&{\sum {z_{i}^{2}y_{i}^{2}} }&{2\sum {z_{i}^{3}} }&{2\sum {z_{i}^{2}{{y}_{i}}} } \\ {2\sum {z_{i}^{3}{{y}_{i}}} }&{4\sum {z_{i}^{2}y_{i}^{2}} }&{2\sum {{{z}_{i}}y_{i}^{3}} }&{4\sum {z_{i}^{2}{{y}_{i}}} }&{4\sum {{{z}_{i}}y_{i}^{2}} } \\ {\sum {z_{i}^{2}y_{i}^{2}} }&{2\sum {{{z}_{i}}y_{i}^{3}} }&{\sum {y_{i}^{4}} }&{2\sum {{{z}_{i}}y_{i}^{2}} }&{2\sum {y_{i}^{3}} } \\ {2\sum {z_{i}^{3}} }&{4\sum {z_{i}^{2}{{y}_{i}}} }&{2\sum {{{z}_{i}}y_{i}^{2}} }&{4\sum {z_{i}^{2}} }&{4\sum {{{z}_{i}}{{y}_{i}}} } \\ {2\sum {z_{i}^{2}{{y}_{i}}} }&{4\sum {{{z}_{i}}y_{i}^{2}} }&{2\sum {y_{i}^{3}} }&{4\sum {{{z}_{i}}{{y}_{i}}} }&{4\sum {y_{i}^{2}} } \end{array}} \right)\left( \begin{gathered} A \hfill \\ B \hfill \\ C \hfill \\ D \hfill \\ E \hfill \\ \end{gathered} \right) = \left( {\begin{array}{*{20}{c}} {\sum {z_{i}^{2}} } \\ {2\sum {{{z}_{i}}{{y}_{i}}} } \\ {\sum {y_{i}^{2}} } \\ {2\sum {{{z}_{i}}} } \\ {2\sum {{{y}_{i}}} } \end{array}} \right).$

В скалярном виде уравнение (3.5) имеет вид системы линейных уравнений:

(3.6)
$\left\{ \begin{gathered} A\sum {z_{i}^{4} + 2B\sum {z_{i}^{3}{{y}_{i}} + C\sum {z_{i}^{2}y_{i}^{2} + D2\sum {z_{i}^{3} + 2E\sum {z_{i}^{2}{{y}_{i}} = \sum {z_{i}^{2},} } } } } } \hfill \\ 2A\sum {z_{i}^{3}{{y}_{i}} + 4B\sum {z_{i}^{2}y_{i}^{2} + 2C\sum {{{z}_{i}}y_{i}^{3} + 4D\sum {z_{i}^{2}{{y}_{i}} + 4E\sum {{{z}_{i}}y_{i}^{2}} = 2\sum {{{z}_{i}}{{y}_{i}}} ,} } } } \hfill \\ A\sum {z_{i}^{2}y_{i}^{2} + 2B\sum {{{z}_{i}}y_{i}^{3} + C\sum {y_{i}^{4}} + 2D\sum {{{z}_{i}}y_{i}^{2}} + 2E\sum {y_{i}^{3}} = \sum {y_{i}^{2},} } } \hfill \\ 2A\sum {z_{i}^{3} + 4B\sum {z_{i}^{2}{{y}_{i}} + 2C\sum {{{z}_{i}}y_{i}^{2} + 4D\sum {z_{i}^{2}} + 4E\sum {{{z}_{i}}} y = 2\sum {{{z}_{i}},} } } } \hfill \\ 2A\sum {z_{i}^{2}{{y}_{i}} + 4B\sum {{{z}_{i}}y_{i}^{2} + 2C\sum {y_{i}^{3}} + 4D\sum {{{z}_{i}}} {{y}_{i}} + 4E\sum {y_{i}^{2}} = 2\sum {{{y}_{i}}.} } } \hfill \\ \end{gathered} \right.$

Здесь все суммирования производятся для i = $\overline {1,n} $. Решение системы (3.6) может быть получено, например, с помощью метода Жордана. Подставив найденные значения компонентов вектора f в уравнения (3.2) и (3.3), находим положение начала координат новой СК (z0, y0) и угол поворота α в приборной СК. Угол α соответствует рассогласованию приборной и связанной СК по крену. Угол β рассогласования СК по тангажу определяется как

$\beta = {{y}_{0}}\delta ,$
где δ – угловое разрешение соответствующего канала ОЭС.

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

Необходимо отметить, что область определения гипербол рассматриваемого типа расположена от –∞ до ∞, а наблюдаемая часть траектории сопровождаемой точки представляет собой лишь малый фрагмент. Кроме того, наблюдаемые гиперболы имеют очень малую кривизну, особенно для ОЭС с узким полем зрения, и на непрерывную кривую наложен шум дискретизации, связанный с матричной структурой МФПУ, т.е. измеряемые координаты траектории сопровождаемой точки округлены до целых чисел. Перечисленные факторы могут приводить к тому, что искомые коэффициенты в уравнении гипербол будут определяться с недостаточной точностью. Это в свою очередь может привести к возникновению недопустимо больших ошибок при определении рассогласования СК относительно оси OX. В связи с этим был предложен альтернативный подход к нахождению упомянутого рассогласования.

4. Альтернативный подход к определению рассогласования систем координат относительно продольной оси связанной СК. Подход заключается в следующем. Производится аппроксимация наблюдаемых траекторий сопровождаемых точек прямыми линиями. Если ОЭС не имеет бокового смещения относительно строительной оси вертолета, то гиперболы на экране будут симметричны относительно центральной вертикали. В этом случае угол наклона получаемых прямых будет соответствовать углу рассогласования СК относительно оси OX. Если есть боковое смещение, то, как указывалось выше, это приведет к смещению центра гипербол в боковом направлении и аппроксимирующие прямые будут иметь наклон, который должен быть скомпенсирован.

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

(4.1)
$y = az + b.$

Тогда, если имеется n наблюдений сопровождаемой точки, т.е. зафиксированы координаты ${{\left\{ {{{z}_{i}},{{y}_{i}}} \right\}}_{{i = }}}_{{\overline {1,n} }}$, то метод наименьших квадратов для нахождения параметров прямой состоит в минимизации выражения:

$\sum\limits_{i = 1}^n {({{y}_{i}}} - a{{z}_{i}} - b{{)}^{2}} \to \mathop {\min }\limits_{a,b} .$

Нетрудно показать, что параметры a и b принимают следующие значения:

(4.2)
$a = \frac{{\sum\limits_{i = 1}^n {{{z}_{i}}{{y}_{i}} - n{{m}_{z}}{{m}_{y}}} }}{{\sum\limits_{i = 1}^n {z_{i}^{2} - nm_{z}^{2}} }};{\kern 1pt} \quad b = {{m}_{y}} - a{{m}_{z}},\quad где\quad {{m}_{z}} = \frac{1}{n}\sum\limits_{i = 1}^n {{{z}_{i}}} ;\quad {{m}_{y}} = \frac{1}{n}\sum\limits_{i = 1}^n {{{y}_{i}}} .$

Угол наклона прямой в экранной СК определяется выражением

$\alpha = \operatorname{arctg} a.$

Если есть боковое смещение L датчика изображения ОЭС относительно строительной оси вертолета или сам вертолет в процессе вращения имел боковую составляющую вектора линейной скорости, то центр гипербол сместится в боковом направлении на соответствующую величину ΔZ  :

$\Delta Z = \frac{{\operatorname{arctg} (L{\text{/}}R)}}{{\Delta \varphi }},$
где R – расстояние от датчика ОЭС до сопровождаемой точки, Δφ – угловое разрешение датчика изображения. В этом случае во избежание ошибки определения рассогласования относительно оси OX необходимо брать для аппроксимации прямой часть гиперболы, симметричную относительно вычисленного положения ее центра. При этом появляется необходимость измерения дальности до сопровождаемой точки, что обеспечивается в результате обработки информации РЛС [6].

5. Анализ точности метода определения ориентации вектора линейной скорости ЛА. В данном разделе решается задача получения аналитических выражений для СКО определения углов ориентации вектора линейной скорости ЛА изложенным выше методом (разд. 2).

Для получения аналитических выражений будем считать измерения координат сопровождаемой характерной точки равноточными с СКО σ на всей траектории ее наблюдения, которая аппроксимируется прямой, описываемой уравнением $y = az + b$. Необходимо найти дисперсии оценок вектора G = (b, a)T, задающего положение аппроксимирующей прямой в СК кадра. В соответствии с [9] корреляционная матрица ошибок оценок вектора G

(5.1)
${{{\mathbf{R}}}_{{\mathbf{G}}}} = \left( {\begin{array}{*{20}{c}} {\sigma _{b}^{2}}&{r_{{ab}}^{{}}} \\ {r_{{ab}}^{{}}}&{\sigma _{a}^{2}} \end{array}} \right) = {{({\mathbf{B}})}^{{ - 1}}} = {{({{{\mathbf{T}}}^{{\text{T}}}}{{{\mathbf{M}}}^{{ - 1}}}{\mathbf{T}})}^{{ - 1}}},$
где T = (I, Z) – матрица размерности n × 2, n – количество измерений, I – единичный вектор, Z – вектор горизонтальных координат измерений, M – корреляционная матрица ошибок измерения. Поскольку измерения координат сопровождаемой точки в каждый момент времени независимы, то корреляционная матрица (5.1) диагональна. Тогда

(5.2)
${\mathbf{B}} = \left( {\begin{array}{*{20}{c}} {\frac{n}{{{{\sigma }^{2}}}}}&{\frac{1}{{{{\sigma }^{2}}}}\sum\limits_{i = 1}^n {{{z}_{i}}} } \\ {\frac{1}{{{{\sigma }^{2}}}}\sum\limits_{i = 1}^n {{{z}_{i}}} }&{\frac{1}{{{{\sigma }^{2}}}}\sum\limits_{i = 1}^n {z_{i}^{2}} } \end{array}} \right),\quad {{{\mathbf{R}}}_{{\mathbf{G}}}} = \frac{1}{{\frac{n}{{{{\sigma }^{4}}}}\sum\limits_{i = 1}^n {z_{i}^{2}} - {{{\left( {\frac{1}{{{{\sigma }^{2}}}}\sum\limits_{i = 1}^n {{{z}_{i}}} } \right)}}^{2}}}}\left( {\begin{array}{*{20}{c}} {\frac{1}{{{{\sigma }^{2}}}}\sum\limits_{i = 1}^n {z_{i}^{2}} }&{ - \frac{1}{{{{\sigma }^{2}}}}\sum\limits_{i = 1}^n {z_{i}^{{}}} } \\ { - \frac{1}{{{{\sigma }^{2}}}}\sum\limits_{i = 1}^n {z_{i}^{{}}} }&{\frac{n}{{{{\sigma }^{2}}}}} \end{array}} \right).$

Найдем входящие в (5.2) выражения для $\sigma _{a}^{2},\;\sigma _{b}^{2}$, rab:

$\sigma _{a}^{2} = \frac{{n{\text{/}}{{\sigma }^{2}}}}{{\frac{n}{{{{\sigma }^{4}}}}\sum\limits_{i = 1}^n {z_{i}^{2}} - {{{\left( {\frac{1}{{{{\sigma }^{2}}}}\sum\limits_{i = 1}^n {{{z}_{i}}} } \right)}}^{2}}}} = \frac{{{{\sigma }^{2}}}}{{\sum\limits_{i = 1}^n {z_{i}^{2}} - \frac{1}{n}{{{\left( {\sum\limits_{i = 1}^n {{{z}_{i}}} } \right)}}^{2}}}} = \frac{{{{\sigma }^{2}}}}{{n\left( {\frac{{\sum\limits_{i = 1}^n {z_{i}^{2}} }}{n}} \right) - {{{\left( {\frac{{\sum\limits_{i = 1}^n {{{z}_{i}}} }}{n}} \right)}}^{2}}}}.$

Для получения аналитических выражений будем считать, что промежутки между измеряемыми горизонтальными координатами сопровождаемой точки одинаковы, т.е. zizi – 1 = Δz для ∀i. Такая последовательность может быть сформирована с помощью некоторого прореживания измерений. Поскольку ошибки определения угла ориентации прямой не зависят от положения начального измерения по z, можно сдвинуть начало координат экранной СК так, чтобы z1 = Δz. Тогда

$\sigma _{a}^{2} = \frac{1}{n}\frac{{{{\sigma }^{2}}}}{{\left( {\frac{{\sum\limits_{i = 1}^n {{{{(i\Delta z)}}^{2}}} }}{n}} \right) - {{{\left( {\frac{{\sum\limits_{i = 1}^n {i\Delta z} }}{n}} \right)}}^{2}}}} = \frac{1}{n}\frac{{{{\sigma }^{2}}}}{{\left( {\frac{{\Delta {{z}^{2}}}}{n}\sum\limits_{i = 1}^n {{{i}^{2}}} } \right) - {{{\left( {\frac{{n(n + 1)\Delta z}}{{2n}}} \right)}}^{2}}}} = \frac{1}{n}\frac{{{{\sigma }^{2}}}}{{\Delta {{z}^{2}}\left[ {\frac{1}{n}\sum\limits_{i = 1}^{{{n}_{j}}} {{{i}^{2}}} - {{{\left( {\frac{{n + 1}}{2}} \right)}}^{2}}} \right]}}.$

Используя выражение

$\sum\limits_{i = 1}^n {{{i}^{2}}} = \frac{{n(n + 1)(2n + 1)}}{6},$
находим

(5.3)
$\sigma _{a}^{2} = \frac{{{{\sigma }^{2}}}}{{\Delta {{z}^{2}}n\left[ {\frac{{n(n + 1)(2n + 1)}}{{3n}} - {{{\left( {\frac{{n + 1}}{2}} \right)}}^{2}}} \right]}} = \frac{{12{{\sigma }^{2}}}}{{\Delta {{z}^{2}}n({{n}^{2}} - 1)}}.$

Согласно предыдущим рассуждениям, можно получить

(5.4)
$\sigma _{b}^{2} = \frac{1}{{{{n}_{j}}}}\frac{{\sigma _{1}^{2}\sum\limits_{i = 1}^{{{n}_{j}}} {x_{i}^{2}} }}{{\sum\limits_{i = 1}^{{{n}_{j}}} {x_{i}^{2}} - {{{\frac{{\left( {\sum\limits_{i = 1}^{{{n}_{j}}} {{{x}_{i}}} } \right)}}{{{{n}_{j}}}}}}^{2}}}}.$

Относительно координаты z0 = z1 – Δz запишем (5.4) в виде

(5.5)
$\begin{gathered} \sigma _{b}^{2} = \frac{1}{n}\frac{{\sigma _{{}}^{2}\sum\limits_{i = 1}^n {{{{({{z}_{0}} + i\Delta z)}}^{2}}} }}{{\sum\limits_{i = 1}^n {{{{({{z}_{0}} + i\Delta z)}}^{2}}} - \frac{{{{{\left( {\sum\limits_{i = 1}^n {({{z}_{0}} + i\Delta z)} } \right)}}^{2}}}}{n}}} = \frac{{12\sigma _{{}}^{2}}}{n}\frac{{z_{0}^{2} + (n + 1){{z}_{0}}\Delta z + \frac{{\Delta {{z}^{2}}}}{{12}}(4{{n}^{2}} + 6n + 2)}}{{\Delta {{z}^{2}}({{n}^{2}} - 1)}} = \\ \, = \frac{{\sigma _{{}}^{2}}}{n}\left[ {12\frac{{{{{\left( {{{z}_{0}} + \frac{{n + 1}}{2}\Delta z} \right)}}^{2}}}}{{\Delta {{z}^{2}}({{n}^{2}} - 1)}} + 1} \right]. \\ \end{gathered} $

Аналогично можно получить выражение для коэффициента корреляции:

(5.6)
${{r}_{{ab}}} = - \frac{1}{n}\frac{{\sigma _{{}}^{2}\sum\limits_{i = 1}^n {{{z}_{i}}} }}{{\sum\limits_{i = 1}^n {z_{i}^{2}} - {{{\frac{{\left( {\sum\limits_{i = 1}^n {{{z}_{i}}} } \right)}}{n}}}^{2}}}}.$

Сместим начало координат так, чтобы абсцисса z = 0 находилась в середине траектории сопровождаемой точки. В этом случае в выражениях (5.5), (5.6) числитель дроби равен нулю, и таким образом коэффициент корреляции rab = 0, а $\sigma _{b}^{2} = {{\sigma }^{2}}{\text{/}}n$. Определим приведенную дисперсию $\sigma _{{{\text{прив}}}}^{2}$ отклонения прямой, аппроксимирующей траекторию точки, в районе точки схождения. Учитывая (5.2), получаем

(5.7)
$\sigma _{{{\text{прив}}}}^{2} = \sigma _{a}^{2}\Delta {{L}^{2}} + {{\sigma }^{2}}{\text{/}}n = {{\sigma }^{2}}\left( {\frac{{12\Delta {{L}^{2}}}}{{\Delta {{z}^{2}}n({{n}^{2}} - 1)}} + \frac{1}{n}} \right),$
где ΔL – расстояние от середины траектории до ожидаемого положения точки схождения. По сути, $\sigma _{{{\text{прив}}}}^{2}$ определяет СКО точки, лежащей на аппроксимированной прямой в районе точки схождения в направлении, перпендикулярном этой прямой. Если две прямых, имеющих равные $\sigma _{{{\text{прив}}}}^{2}$, пересекаются под острым углом γ, то большая ось эллипса рассеяния точек пересечения прямых совпадает с биссектрисой угла γ (в силу симметрии задачи), а максимальное СКО σmax точки пересечения по этому направлению определяется выражением

(5.8)
${{\sigma }_{{\max }}} \simeq \frac{{{{\sigma }_{{{\text{прив}}}}}}}{{\cos (90 - \gamma + \gamma {\text{/}}2)}} = \frac{{{{\sigma }_{{{\text{прив}}}}}}}{{\sin (\gamma {\text{/}}2)}} = \frac{{\sqrt 2 {{\sigma }_{{{\text{прив}}}}}}}{{\sqrt {1 - \cos \gamma } }}.$

Если при нахождении положения точки схождения было выделено не менее 10 характерных точек, то ошибки определения ее координат будут на порядок меньше ошибок выделения характерных точек, т.е. меньше углового разрешения ОЭС.

6. Алгоритмы выделения характерных точек. Для реализации метода измерения ориентации вектора угловой скорости вертолета необходимо разработать алгоритмы выделения характерных точек на ПП и инфраструктуре аэродромов. При этом следует учитывать следующие особенности, которые могут возникнуть при обработке изображений ОЭС. На территории аэродромов постоянно присутствуют движущиеся объекты. Это могут быть ЛА во время руления, топливозаправщики, прочая техника, а также люди, которые на тепловизионных изображениях имеют высокий контраст. При попадании в поле зрения ОЭС фрагменты их образов могут быть выделены как характерные точки, что приведет к некорректному результату работы предлагаемого метода. Поэтому алгоритмы выделения характерных точек должны обнаруживать в первую очередь верхние (наиболее высоко расположенные) фрагменты стационарных наземных объектов.

Используемые детекторы характерных точек (Моравека, Харриса, Форстнера, Ши-Томаси-Кэнэда и др.) ориентированы на решение других задач (создание панорамных изображений, привязки изображений и т.п.) [1113]. Эти детекторы, формируя так называемую “меру отклика угла”, находят не собственно углы, а любые информативно насыщенные участки изображения. Поэтому их использование при обработке изображений, полученных в передней полусфере ЛА, приводит при уменьшении порога для “меры отклика угла” к увеличению количества обнаружений в окрестности нижних фрагментов наземных объектов, где, как правило, наблюдается достаточное разнообразие градиентов яркостного поля. В результате на изображениях в первую очередь выделяются характерные точки, которые с большой вероятностью могут принадлежать движущимся объектам, селекция которых может быть проведена доплеровскими методами [6], и идентификацией радиолокационных (РЛ) и оптико-электронных образов. Поэтому было необходимо создать детектор характерных точек, который имеет свои отличия:

выделяет в основном верхние углы объектов, реагируя на их геометрические свойства;

инвариантен к величинам контрастов, формирующих эти углы;

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

Разработанный детектор анализирует интенсивности pij в скользящем окне 9 × 9 с вертикальной (направленной вниз) координатой y и горизонтальной координатой z (i – номер строки,  j – номер столбца). Далее индексы приведены в системе координат скользящего окна. Для удобства введем следующие обозначения: M[…] – среднеарифметическое из чисел в квадратных скобках,

$\Delta \mu _{j}^{Н} = {\text{M}}[{{p}_{{7,j + 6}}},{{p}_{{8,j + 6}}},{{p}_{{9,j + 6}}}] - {\text{M}}[{{p}_{{7,j}}},{{p}_{{8,j}}},{{p}_{{9,j}}}],$
$\Delta \mu _{j}^{В} = {\text{M}}[{{p}_{{1,j + 6}}},{{p}_{{2,j + 6}}},{{p}_{{3,j + 6}}}] - {\text{M}}[{{p}_{{1,j}}},{{p}_{{2,j}}},{{p}_{{3,j}}}].$

Алгоритм детектора.

1. Находим нижний максимальный перепад (НМП):

${\text{НМП}} = \max \left\{ {\left| {\Delta \mu _{1}^{Н}} \right|;\left| {\Delta \mu _{2}^{Н}} \right|;\left| {\Delta \mu _{3}^{Н}} \right|} \right\}.$

2. Вычисляем верхний средний перепад (ВСП):

${\text{ВСП}} = \left| {\frac{{\Delta \mu _{1}^{В} + \Delta \mu _{2}^{В} + \Delta \mu _{3}^{В}}}{3}} \right| = \left| {M[\Delta \mu _{1}^{В},\Delta \mu _{2}^{В},\Delta \mu _{3}^{В}]{\kern 1pt} } \right|.$

3. Оцениваем модуль производных по координате z по первой разности в центральных точках верхних трех (В1, В2, В3) и нижних трех (Н1, Н2, Н3) строк скользящего окна:

$\begin{gathered} {\text{В}}1 = \left| {\frac{{{{p}_{{16}}} - {{p}_{{14}}}}}{2}} \right|;\quad {\text{В}}2 = \left| {\frac{{{{p}_{{26}}} - {{p}_{{24}}}}}{2}} \right|;\quad {\text{В}}3 = \left| {\frac{{{{p}_{{36}}} - {{p}_{{34}}}}}{2}} \right|, \\ {\text{H}}1 = \left| {\frac{{{{p}_{{76}}} - {{p}_{{74}}}}}{2}} \right|;\quad {\text{H}}2 = \left| {\frac{{{{p}_{{86}}} - {{p}_{{84}}}}}{2}} \right|;\quad {\text{H}}3 = \left| {\frac{{{{p}_{{96}}} - {{p}_{{94}}}}}{2}} \right|. \\ \end{gathered} $

4. Если

$\begin{gathered} \left( {{\text{H}}1 > \eta \,{\text{НМП}}} \right) \wedge \left( {{\text{H}}2 > \eta \,{\text{НМП}}} \right) \wedge \left( {{\text{H}}3 > \eta \,{\text{НМП}}} \right) \wedge \\ \, \wedge \left( {{\text{В}}1 < \eta \,{\text{ВСП}}} \right) \wedge \left( {{\text{В}}2 < \eta \,{\text{ВСП}}} \right) \wedge \left( {{\text{В}}3 < \eta \,{\text{ВСП}}} \right), \\ \end{gathered} $
то выполняется п. 5. Иначе переход к новому положению “скользящего окна” по координате z. Здесь $ \wedge $ – логическое умножение, η – коэффициент, зависящий от аберрационных свойств объектива. В нашем случае экспериментально определено η = 0.2.

5. Определяем наклон k некоторой пеленгационной характеристики:

$k = \frac{1}{3}{\text{M}}[\Delta \mu _{1}^{Н},\Delta \mu _{2}^{Н},\Delta \mu _{3}^{Н}].$

Находим смещение по нижним строкам:

$\begin{gathered} \Delta {{z}_{7}} = \frac{1}{k}\left( {{\text{M}}[{{p}_{{74}}},{{p}_{{75}}},{{p}_{{76}}}] - {\text{M}}[{{p}_{{71}}},{{p}_{{72}}},{{p}_{{73}}}]} \right), \\ \Delta {{z}_{8}} = \frac{1}{k}\left( {{\text{M}}[{{p}_{{84}}},{{p}_{{85}}},{{p}_{{86}}}] - {\text{M}}[{{p}_{{81}}},{{p}_{{82}}},{{p}_{{83}}}]} \right), \\ \Delta {{z}_{9}} = \frac{1}{k}\left( {{\text{M}}[{{p}_{{94}}},{{p}_{{95}}},{{p}_{{96}}}] - {\text{M}}[{{p}_{{91}}},{{p}_{{92}}},{{p}_{{93}}}]} \right). \\ \end{gathered} $

Среднее смещение Δz от центрального пиксела окна по координате z равно

$\Delta z = {\text{M}}[\Delta {{z}_{7}},\Delta {{z}_{8}},\Delta {{z}_{9}}].$

6. Производим операции 1–5, рассматривая вместо строк столбцы и вместо нижней границы левую границу или правую границу скользящего окна. Это позволяет позиционировать угол объекта в окне и определять величины Δy7, Δy8, Δy9 и Δy как в случае левого, так и правого верхнего угла здания.

7. Близость (±1 пиксел) результатов косвенных измерений Δz7, Δz8, Δz9 и Δy7, Δy8, Δy9 используется как еще один признак правильности выделения стороны угла.

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

При отсутствии на анализируемом изображении выраженных характерных точек для формирования траекторий могут применяться другие алгоритмы, в частности, корреляционного типа. Это очевидно, поскольку за время сопровождения фрагмента изображения не меняются условия наблюдения. Для обеспечения эффективного сопровождения характерных точек (недопущения перепутывания точек) или фрагментов изображения необходимо использовать информацию, поступающую с БИНС и РЛ-системы. Построение области допоиска характерной точки, обнаруженной на предыдущем изображении, осуществляется на основе данных о величине текущей угловой скорости, промежутке времени между моментами формирования изображений и возможных углах наклона траекторий в экранной СК.

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

7. Результаты моделирования. Моделирование проводилось для ОЭС с полем зрения 40 × 30° и МФПУ 640 × 480 пикселов. При этом предполагалось, что дисторсия объектива, если она присутствует, программно скомпенсирована. Сначала было проведено исследование возможности использования метода, в котором определяются все коэффициенты уравнения гиперболы. Для этого синтезировалась траектория сопровождаемой точки в нижней части экрана, т.е. траектории с возможно большей кривизной с заданными величинами угла поворота и отклонения центра гиперболы по вертикали от центра экрана. Для значительного количества наблюдений (несколько сот) координаты фиксировались с высокой точностью. Затем зафиксированные координаты загружались в программу определения коэффициентов гиперболы, ее центра и угла поворота в соответствии с (3.6), которые сопоставлялись с заданными при формировании гиперболы. Действия производились для одной и той же траектории, но при этом координаты округлялись до все меньшего числа знаков после запятой, вплоть до округления до целого числа. В качестве примера в табл. 1 приведены задаваемые и определяемые коэффициенты уравнения гиперболы A, B, C, D, E и оценки угла поворота гиперболы α (в градусах) от количества знаков после запятой координат сопровождаемой точки для случая, когда задаваемый угол поворота гиперболы равен нулю.

Таблица 1
Количество знаков после запятой в определяемых параметрах A B C D E α, град
Округление до целого –3 × 10–6 –1.3 × 10–5 –0.00021 0.001575 0.014466 –3.50003
1 знак 0.000317 –1.3 × 10–5 –0.00026 –0.1009 0.206549 –1.27107
2 знака 0.024669 –9.8 × 10–6 –0.00379 –7.89361 14.81396 –0.01968
3 знака 0.096909 –3.8 × 10–7 –0.01425 –31.0109 58.146 –1.9 × 10–4
4 знака 0.10011 –5 × 10–9 –0.01472 –32.0353 60.06611 –2.5 × 10–6
5 знаков 0.099999 –2.1 × 10–10 –0.0147 –31.9997 59.99946 –1.1 × 10–7
6 знаков 0.100001 –7.8 × 10–11 –0.0147 –32.0004 60.00037 –3.9 × 10–8
Заданные параметры 0.1 0 –0.0147 –32 60 0

Этот простой пример показывает, что для успешного применения рассмотренного метода необходимо увеличить точность определения координат сопровождаемой точки, по крайней мере, на два порядка, что в настоящее время невозможно. Аналогичное моделирование было проведено для альтернативного метода определения рассогласования СК по крену на основе аппроксимации наблюдаемых траекторий сопровождаемых точек прямыми линиями с помощью (4.1), (4.2). В табл. 2 приведена оценка угла поворота гиперболы с помощью аппроксимации от количества знаков после запятой. Была задана гипербола со следующими параметрами: A = 0.01, B = = 0.001, C = –0.025, D = –3.44, E = 59.68, z0 = 320,  y0 = 240, α = 13.2222'.

Таблица 2
Количество знаков после запятой, учитываемых при аппроксимации 4 3 2 1 Округление до целого
Оценка угла поворота, угл. мин 14.0442 14.0442 14.0461 14.07 14.1095

На рис. 1 представлена эта гипербола при округлении вертикальной координаты до целых значений и ее аппроксимация прямой линией (за начало координат в данном примере выбран левый нижний угол изображения). Для наглядности масштаб по оси OY  во много раз больше, чем по оси OZ. Из табл. 2 и рис. 1 очевидно, что данный метод намного более устойчив к уменьшению точности измерения координат сопровождаемой точки. Дальнейшее моделирование проводилось для округления измеряемых координат до целых значений.

Рис. 1.

Аппроксимация гиперболы прямой линией при округлении координат до целых чисел

На рис. 2 приведена зависимость ошибки измерения рассогласования по крену от величины этого рассогласования, т.е. ошибка измерения угла поворота рассмотренной выше гиперболы от этого угла. При этом относительная ошибка в диапазоне измерений от 10 до 130' не превышала 4–6%.

Рис. 2.

Зависимость ошибки измерения поправки относительно продольной оси БИНС от величины этой поправки

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

(7.1)
$k = \mathop {\max }\limits_i \frac{{{{y}_{i}} - {{y}_{1}} - {{x}_{i}}({{y}_{n}} - {{y}_{1}}){\text{/}}({{x}_{n}} - {{x}_{1}})}}{{({{x}_{n}} - {{x}_{1}})}}.$

На рис. 3 представлена указанная зависимость для случая, когда рассогласование СК равно нулю.

Рис. 3.

Зависимость ошибки измерения рассогласования относительно продольной оси БИНС от степени кривизны траектории сопровождаемой точки

Полученная зависимость носит линейный характер. Отсюда очевидно, что при наборе заданного количества траекторий необходимо усреднение результатов измерения за время всего маневра. При этом углам рассогласования, найденным по траекториям с меньшей кривизной, надо придавать больший вес. Если было проведено сопровождение m точек на ПП, то для определения среднего угла наклона траекторий $\overline \alpha $ может использоваться следующее выражение:

(7.2)
$\bar {\alpha } = \operatorname{arctg} \frac{{\sum\limits_{j = 1}^m {{{a}_{j}}(1 - q{{k}_{j}})} }}{{n - q\sum\limits_{j = 1}^m {{{k}_{j}}} }},$
где q – коэффициент, величина которого может быть выбрана в результате моделирования с реальной аппаратурой.

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

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

Рис. 4.

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

Статистический анализ результатов обработки изображений промышленных и городских сцен показал, что средние квадратические ошибки оценок горизонтальных и вертикальных координат характерных точек равны 0.36–0.71 пиксела в зависимости от условий проведения регистрации изображений (от простых до ограниченно сложных метеоусловий).

8. Комплексная обработка видеоинформации и навигационных параметров. Рассматривая влияние возмущающих факторов на точность согласования СК, необходимо отметить, что идеально выполнить маневр (полет по прямой или вращение) невозможно независимо от того, выполняется маневр автопилотом или летчиком вручную. При выполнении маневра на вертолет будут действовать возмущающие воздействия, приводящие к появлению стабилизационных колебаний, которые будут проявляться в виде отклонений вектора угловой или линейной скорости от среднего положения. При отсутствии стабилизации поля зрения ОЭС в экранной СК это будет проявляться в виде отклонений сопровождаемых точек от идеальных траекторий. Такие отклонения могут быть скомпенсированы с помощью информации, получаемой из БИНС. Для этого необходимо вычислять приращения углов ориентации вертолета по крену, тангажу и курсу за время, прошедшее между моментами формирования предыдущего и текущего изображений, и пересчитывать координаты характерной точки в экранной СК, компенсируя возмущения. Однако реальные временные задержки, связанные с формированием и передачей информации, а также неопределенность моментов обновления параметров могут приводить к недопустимо большим ошибкам при построении траекторий точек. Поэтому для компенсации возмущающих воздействий необходимо использовать возможности стабилизации поля зрения ОЭС. В этом случае погрешность стабилизации изображений существенно меньше углового разрешения датчиков ОЭС. В такой постановке компенсация стабилизационных колебаний сводится к необходимости определения основных требований к выполнению маневров.

Для определения направления линейной скорости и оценок поправок относительно боковой и вертикальной осей БИНС ЛА должен:

двигаться прямолинейно под управлением автопилота;

продолжительность полета по прямой должна быть не менее нескольких периодов стабилизационных колебаний (зависит от характеристик ЛА);

поле зрения ОЭС стабилизировано;

должно быть детектировано и построено траекторий характерных точек на изображении ПП не менее заданного количества (порядка 8–10).

Для оценки поправки относительно продольной оси БИНС (по крену) вертолету необходимо:

совершить вращение на месте под управлением автопилота;

продолжительность вращения должна быть не менее нескольких периодов стабилизационных колебаний (зависит от характеристик вертолета);

поле зрения ОЭС стабилизировано по крену и тангажу;

должно быть детектировано и построено траекторий характерных точек на ПП не менее заданного количества (порядка 8–10).

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

При определении рассогласования приборных СК ОЭС со связанной СК нужно учесть, что БИНС выдает параметры крен, тангаж, курс с заданной частотой, а для определения направления движения по экрану ОЭС сопровождаемых характерных точек изображения ПП необходимо некоторое время. При выполнении маневра на вертолет будут действовать возмущающие силы, приводящие к появлению упомянутых стабилизационных колебаний. Поэтому для определения искомой угловой поправки следует проводить усреднение выдаваемых БИНС значений крена за время выполнения маневра. Этот усредненный крен нужно сопоставить со средним углом наклона прямых (7.2) с учетом кривизны (7.1), аппроксимирующих траектории точек, сопровождаемых на экране ОЭС за все время выполнения маневра. Разность указанных величин определяет поправку приборной СК относительно продольной оси БИНС Δγ:

(8.1)
$\Delta \gamma = \frac{{\sum\limits_{i = 1}^n {{{\gamma }_{i}}\Delta t} }}{T} - \bar {\alpha },$
где γi – текущее значение крена, получаемое из БИНС, Δt – период выдачи информации из БИНС, T – время выполнения маневра, n – количество периодов за время T, m – количество построенных траекторий характерных точек на ПП, $\bar {\alpha }$ – средний угол наклона прямых.

Комплексная обработка сигналов БИНС и видеоинформации ОЭС для векторного согласования СК относительно боковой и вертикальной осей связанной СК заключается в следующем. Положение точки схождения в приборной СК $(\overline Y ,\overline Z )$, найденное изложенным выше способом (2.1)–(2.3), корректируется на поправку Δγ относительно продольной оси БИНС, которая определена на предыдущем этапе согласования СК (8.1):

(8.2)
$\left( \begin{gathered} {{Y}_{{{\text{кор}}}}} \hfill \\ {{Z}_{{{\text{кор}}}}} \hfill \\ \end{gathered} \right) = \left( {\begin{array}{*{20}{c}} {\cos \Delta \gamma }&{\sin ( - \Delta \gamma )} \\ { - \sin ( - \Delta \gamma )}&{\cos \Delta \gamma } \end{array}} \right)\left( \begin{gathered} \overline Y \hfill \\ \overline Z \hfill \\ \end{gathered} \right).$

Навигационные параметры выдаются БИНС в нормальной подвижной СК (НПСК). Поэтому с момента начала и до окончания прослеживания траекторий характерных точек информация, выдаваемая БИНС, преобразуется следующим образом. Определяются проекции вектора путевой скорости, задаваемого текущими значениями вертикальной VH, северной VN и восточной VE составляющих, на оси связанной СК (ССК) с помощью выражений [14]:

${{\left( \begin{gathered} {{R}_{x}} \hfill \\ {{R}_{y}} \hfill \\ {{R}_{z}} \hfill \\ \end{gathered} \right)}_{{{\text{ССК}}}}} = \left( {\begin{array}{*{20}{c}} {\cos \vartheta \,\cos \psi }&{\sin \vartheta }&{ - \cos \vartheta \,\sin \psi } \\ { - \cos \gamma \,\sin \vartheta \,\cos \psi + }&{\cos \gamma \,\cos \vartheta }&{\cos \gamma \,\sin \vartheta \,\sin \psi + } \\ { + \sin \gamma \,\sin \psi }&{}&{ + \sin \gamma \,\cos \psi } \\ {\sin \gamma \,\sin \vartheta \,\cos \psi + }&{ - \sin \gamma \,\cos \vartheta }&{ - \sin \gamma \,\sin \vartheta \,\sin \psi + } \\ { + \cos \gamma \,\sin \psi }&{}&{ + \cos \gamma \,\cos \psi } \end{array}} \right){{\left( \begin{gathered} {{V}_{N}} \hfill \\ {{V}_{H}} \hfill \\ {{V}_{E}} \hfill \\ \end{gathered} \right)}_{{{\text{НПСК}}}}},$
где Ψ – курс, $\vartheta $ – тангаж, γ – крен. Полученные значения проекций вектора путевой скорости в связанной СК преобразуются в вектор направления R (у-пеленг phiY и z-пеленг phiZ) с помощью выражений:

$\operatorname{phi} Y = - \operatorname{arctg} \left( {{{R}_{Z}}{\text{/}}{{R}_{X}}} \right);\quad \operatorname{phi} Z = - \operatorname{arctg} ({{R}_{Y}}\cos (\operatorname{phi} Y){\text{/}}{{R}_{X}}).$

Значения у-пеленга и z-пеленга в связанной СК усредняются на время сопровождения характерных точек. Соответствующие разности между проекциями скорректированного положения точки схождения на вертикальную и боковую оси приборной СК, найденные с помощью (8.2), и усредненными значениями пеленгов определяют поправки приборной СК относительно боковой Δϑ и вертикальной осей ΔΨ БИНС:

$\begin{gathered} \Delta \vartheta = \overline {phiZ} - {{Y}_{{{\text{кор}}}}}, \\ \Delta \Psi = \overline {phiY} - {{Z}_{{{\text{кор}}}}}, \\ \end{gathered} $
где $\overline {\operatorname{phi} Y} $ и $\overline {\operatorname{phi} Z} $ – усредненные значения у-пеленга и z-пеленга вектора направления R.

Заключение. Предложен не имеющий аналогов метод согласования приборных СК разных постов ОЭС вертолета между собой и с инерциальной навигационной системой вертолета. Метод не требует предполетной подготовки и включения в состав каждого поста ИИБ. Метод основан на обработке потоков видеоинформации, формируемых тепловизионными и телевизионными каналами ОЭС, и использовании информации инерциальной навигационной системы. За время совершения специальных маневров, представляющих собой вращение вертолета на небольшой высоте и горизонтальный полет по прямой, выполнение которых можно автоматизировать, на изображениях формируются траектории характерных точек ПП. С помощью анализа этих траекторий определяются направления векторов угловой и линейной скоростей вертолета в приборных системах координат, которые сопоставляются с информацией инерциальной навигационной системы.

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

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

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

  1. Lipton A.H. Alignment of Inertial Systems on a Moving Base. Cambridge: Electronics Research Center, 1967.

  2. Schneider A. M. Kalman Filter Formulation for Transfer Alignment of Strapdown Inertial Units // J. of the Institute of Navigation. 1983. V. 30. № 1.

  3. Бельский А.Б. Основные задачи и требования к бортовым ОЭС для современных и перспективных вертолетов // Тр. XXV Междунар. научно-техн. конф. и школы по фотоэлектронике и приборам ночного видения. Т. 1. М.: ОФСЕТ МОСКВА, 2018. С. 2–23.

  4. Мужичек С.М., Обросов К.В., Ким В.Я. и др. Определение направления полета по сигналам оптико-электронной системы переднего обзора // Вестн. компьютерных и информационных технологий. 2013. № 5. С. 8–14.

  5. Мужичек С.М., Обросов К.В., Лисицын В.М. и др. Способ измерения курса летательного аппарата. Пат. 2556286 Российская Федерация, МПК G01C 21/12 C1. № 2014115385/28; заявл. 17.04.2014; опубл. 10.07.2015, Бюл. № 19. 13 с.

  6. Бобин А.В., Лисицын В.М., Обросов К.В., Сикачева М.И. Доплеровская селекция наземных объектов, движущихся со случайными изменениями ориентации вектора скорости // Изв. РАН. ТиСУ. 2021. № 5. С. 143–151.

  7. Клочко В.К. 3D-радиовидение на базе бортовой доплеровской РЛС // Ракетно-космическое приборостроение и информационные системы. 2015. Т. 2. Вып. 2. С. 53–57.

  8. Денисов П.В., Зайцев С.Э., Костюк Е.А. и др. Вопросы дешифрирования радиолокационных снимков при радиовидении // Радиотехника. 2014. № 7. С. 7–14.

  9. Тихонов В.И. Статистическая радиотехника. 2-е изд. перераб. и доп. М.: Радио и связь, 1982. 624 с.

  10. Беклемишев Д.В. Курс аналитической геометрии и линейной алгебры. М.: Наука, 1971. 328 с.

  11. Harris C.G., Stephens M.J. Combined Corner and Edge Detector // Proc. Fourth Alvey Vision Conf. Manchester, 1988. P. 147–151.

  12. Foerstner W. A Feature Based Correspondence Algorithm for Image Matching // ISPRS. Commission III Sympos. Rovaniemi. Finland, 1986. V. 26-3/3. P. 150–166.

  13. Сергунов А.А. Перспективы применения детекторов характерных точек для обнаружения движущихся малоразмерных объектов на сложном фоне // Научная сессия ГУАП: сб. докладов: В3. Ч. II. Технические науки // СПб.: ГУАП, 2010. С. 60–62.

  14. ГОСТ 20058–80. Динамика летательных аппаратов в атмосфере.

Дополнительные материалы отсутствуют.