CC BY
2
Субпиксельное накопление и обнаружение смазанных изображений звёзд, полученных астроинерциальным датчиком ориентации на фоне дневного неба
Н.Н. Василюк1
1 ООО «НПК Электрооптика», 107076, г. Москва, ул. Стромынка, д.18, к.1
Аннотация
Для определения ориентации по изображению звёздного поля, полученному цифровой камерой, необходимо обнаружить изображения отдельных звёзд и определить их координаты с субпиксельным разрешением. Эта задача относительно легко решается, если изображение получено неподвижной камерой, наблюдающей ясное ночное небо. Затруднение возникает при наблюдении звёзд поворачивающейся камерой на фоне дневного неба, когда смазанные изображения звёзд практически невозможно обнаружить в фоновом шуме в одном кадре. Из этого затруднения можно выйти за счёт накопления некоторого количества последовательных кадров, в которых предварительно выполнена согласованная фильтрация смаза. В работе рассматривается алгоритм накопления последовательности кадров, реализованный в виде многоканального накопителя, состоящего из нескольких однотипных накопительных каналов, работающих параллельно. Плоскость накопленного изображения в каждом канале совпадает с плоскостью первого кадра накапливаемой последовательности. Субпиксельное разрешение в определении координат обнаруженных изображений звёзд достигается за счёт введения внутрипиксельной дискретизации в первом кадре. Ориентация текущего кадра относительно первого кадра вычисляется из измерений гироскопов астрои-нерциального датчика. При синтезе импульсных реакций фильтров, согласованных со сма-зом в очередном кадре, учитываются ориентация этого кадра и расположение узлов внутри-пиксельной дискретизации первого кадра. Приводятся расчётные формулы для определения числа накапливаемых кадров в зависимости от угловой скорости вращения камеры, интенсивности фона и яркости наблюдаемых звёзд. Результаты моделирования показывают возможность обнаружения смазанных изображений относительно тусклых звёзд на фоне дневного неба и определения их координат с субпиксельным разрешением.
Ключевые слова: астродатчик, астрокорректор, астроинерциальная навигационная система, коррекция смаза, накопление изображений.
Цитирование: Василюк, Н.Н. Субпиксельное накопление и обнаружение смазанных изображений звёзд, полученных астроинерциальным датчиком ориентации на фоне дневного неба / Н.Н. Василюк // Компьютерная оптика. - 2024. - Т. 48, № 2. - С. 303-311. - DOI: I0.18287/2412-6179-C0-I309.
Citation: Vasilyuk NN. Subpixel stacking and detection of blurred star images observed by an astroinertial attitude sensor against the background of the daytime sky. Computer Optics 2024; 48(2): 303-311. DOI: 10.18287/2412-6179-C0-1309.
Введение
Данная работа продолжает работы [1, 2], посвященные алгоритмам обработки изображений звёзд в бесплатформенных астроинерциальных датчиках ориентации (астродатчиках). Рассматриваемый астродатчик представляет собой жёсткую конструкцию, состоящую из цифровой камеры, вычислительного устройства и инерциального измерительного модуля (ИИМ).
Если астродатчик неподвижен относительно звёзд, то объектив цифровой камеры строит изображение звезды на поверхности матричного фотоприёмника (МФП) в виде рассеянного пятна. При вращении астродатчика цифровое изображение звезды смазывается и из пятна превращается в произвольно ориентированную полоску. Форма полоски зависит от угловой
скорости камеры, что вносит непредсказуемую систематическую ошибку в определение координат изображения звезды [3, 4]. Отношение сигнал/шум (С / Ш) в смазанном изображении ухудшается, причём это ухудшение не зависит от оптических условий наблюдения, а определяется исключительно параметрами вращения камеры [5, 6].
В [1, 2] рассмотрен подход к улучшению отношения С / Ш в смазанном изображении, полученному в одном кадре. Улучшение возникает за счёт использования при обработке кадров дополнительной информации, получаемой от ИИМ. Измерения ИИМ позволяют синтезировать цифровой фильтр, согласованный со смазом, непосредственно из параметров движения астродатчика [7]. При таком синтезе фильтра снимается неопределённость в координатах смазанного изображения: максимум профильтрованного
изображения возникает в месте расположения пятна в момент начала экспозиции фильтруемого кадра.
Согласованная фильтрация внутри одного кадра позволяет обнаруживать на фоне дневного неба только яркие звёзды [2]. Таких звёзд очень мало, и они редко используются для определения ориентации астродатчика. Вероятность попадания звезды в поле зрения камеры растёт с увеличением звёздной величины, т.е. с уменьшением яркости звезды. На практике для определения ориентации астродатчика используются звёзды 5-й - 7-й величины, для которых согласованная фильтрация не даёт улучшения С / Ш, достаточного для их обнаружения на фоне дневного неба в одном кадре.
Дальнейшее улучшение С / Ш в изображениях тусклых звёзд возможно за счёт цифрового сложения (накопления) нескольких последовательных кадров. При накоплении кадров уровень полезного сигнала от звёзды растёт пропорционально числу накопленных кадров, а СКО фонового шума и прочих шумов, добавляемых к каждому кадру при его считывании из МФП [8], растёт медленнее, как квадратный корень из числа кадров. Повороты кадров, полученных в различных положениях вращающейся камеры, определяются из измерений гироскопов ИИМ [9,10].
При накоплении смазанных изображений возникает дополнительная проблема. Смазы в изображениях одной и той же звезды в разных кадрах, вообще говоря, ориентированы по-разному. Смазанные изображения, даже после компенсации поворотов кадров, не полностью накладываются друг на друга. То есть вклад отдельных изображений в общую сумму оказывается неполным, а положение максимума этой суммы - неопределённым.
Согласованная фильтрация смазов решает эту проблему - корреляционные максимумы в профильтрованном кадре возникают в тех пикселях, в которых находились изображения звезд в момент начала экспозиции. Если кадры фильтруются независимо, без сквозного учёта траектории движения изображения звезды по всем кадрам накапливаемой последовательности, то положение корреляционного максимума определяется с точностью до пикселя. Максимум возникает в пикселе независимо от того, в каком месте фотоприёмной ячейки находилось изображение звезды в момент начала экспозиции. Поворот независимо профильтрованных кадров ограничивает разрешение в определении координат звезды по накопленному изображению до пиксельного уровня.
В данной работе рассматривается процедура накопления профильтрованных кадров, позволяющая определять координаты накопленного изображения звезды с субпиксельным разрешением. Для этого предлагается организовать накопление кадров в виде многоканального накопителя, состоящего из нескольких однотипных накопительных каналов, работающих параллельно [11, 12]. Каждый накопительный
канал вычисляет своё собственное накопленное изображение, характеризующееся априорным предположением о попадании изображения звезды в начальную точку, одинаково расположенную внутри всех фотоприёмных ячеек в момент начала экспозиции первого кадра накапливаемой последовательности. Для каждого последующего кадра из траектории движения камеры вычисляются новые точки, в которые перешли начальные точки из первого кадра, и для новых точек синтезируется неоднородная импульсная реакция согласованного фильтра. Очередной кадр обрабатывается своим согласованным фильтром и попиксельно складывается с изображением в накопительном канале.
Смазанное изображение звезды считается обнаруженным, если значение пикселя в изображении, накопленном в каком-то из каналов, превысило пороговое значение. Если превышение порога обнаруживается в одном и том же пикселе в нескольких накопленных изображениях, выбирается канал, в котором получено максимальное значение пикселя. За координаты обнаруженного изображения звезды принимаются координаты начальной точки, соответствующей выбранному каналу, а за эпоху наблюдения -момент начала экспозиции первого кадра. Количество накопительных каналов равно количеству дискретных начальных точек внутри фотоприёмных ячеек.
1. Основные понятия, условные обозначения и системы координат
Термины, определения и системы координат, используемые в работе, введены и обоснованы в [1, 2]. Для ясности изложения некоторые определения здесь воспроизводятся, но без подробного объяснения.
МФП цифровой камеры представляет собой прямоугольную матрицу, составленную из фотоприёмных ячеек квадратной формы с длиной стороны a. Размер матрицы Н^, где H, W - число строк и столбцов матрицы соответственно. Положение фотоприёмной ячейки внутри матрицы задаётся парой целочисленных индексов (И, V), где И = 0...Н—1 - номер строки и V = 0. W—1 - номер столбца, в котором расположена ячейка. Значение индексов (0, 0) соответствует левой верхней ячейке матрицы.
В плоскости МФП задана растровая система координат с началом в левом верхнем углу ячейки (0, 0). Положение точки на поверхности МФП задаётся парой действительных чисел [И w]T, где 0 ^ И < Н, 0 ^ w < W. Точка с растровыми координатами [И w]T расположена внутри ячейки с целочисленными индексами
И =Аоог(И), V = Аоог^), (1)
где йоог(...) - функция округления в направлении к —да.
Положение точек на плоскости МФП относительно главной точки снимка [Ио wo]T задаётся векторными координатами [х у]т. Растровые [И w]T и векторные [х у] координаты одной и той же точки связаны:
х 1 И - Иа ' И' = а х + Иа
1_ у _ а W - М!а w 1_ у _ _ Wа _
С конструкцией камеры жёстко связана правая ортогональная система координат камеры (СКК). Точке с векторными координатами [х у]т в плоскости МФП соответствует единичный направляющий вектор в пространстве предметов с координатами в СКК
8 =•
■у Р]т
(2)
•у/х2 + у2 + Р2 '
где Р - фокусное расстояние объектива.
2. Многоканальное накопление кадров
2.1. Определение ориентации кадров
В накоплении участвуют К последовательных кадров с порядковыми номерами к =0...К-1. Экспозиция очередного кадра с номером к начинается в момент времени ¡эк = ¡э0 + кТк начала экспозиции и завершается в момент времени ¡эк + Тэ окончания экспозиции, где ¡э0 -момент начала экспозиции кадра с номером к = 0; ТК -постоянный период получения кадров, Тэ - время экспозиции, ТК > Тэ. Накопление выполняется в плоскости начального кадра с к = 0, то есть все последующие кадры перепроектируются на эту плоскость. Координаты звезды определяются на момент времени /Э0.
Ориентация кадра к задаётся при помощи инерци-альной системы координат СККк, которая в точности совпадает с мгновенным положением СКК поворачивающейся камеры в момент времени ¡эк: СККк = СКК(/эк). Преобразование координат из СКК0 в СККк для кадров с номерами 0 и к задаётся ортогональной матрицей 8°:
8СКК0 (. ) =
ДСКК Мэк ) —
С21 Лк
с22
с23 к
8 0 = т
8 0 - т 3
где s'kJ, /, ] = 1 ... 3 - скалярные компоненты матрицы 80; 13 - единичная матрица 3x3. Верхний и нижний индексы в обозначении матрицы 80 нужно понимать так: если один и тот же вектор ? задан своими координатами 80 в СКК0 и 8к в СККк, то 8к = 8080 .
Оси чувствительности трёхосного гироскопа в составе ИИМ сонаправлены с одноимёнными осями СКК. Гироскоп измеряет интегралы
9т = [9хт 9ут 9хт Г = | ®(т) <Т
от проекций ю(0 = юсск(0 вектора ю(/) угловой скорости камеры на СКК между эквидистантными моментами времени ¡Гт = ¡Г0 + ТГт, где /го - начальный момент времени измерения угловой скорости, ТГ -постоянный период измерения гироскопа, Тг < Тэ / 2 < ТК. За ¡Г0 принимается момент измерения
гироскопа, непосредственно предшествующий моменту ¡э0: ¡Г0 ^ ¡э0 < ¡Г1 (рис. 1). Количество М измерений гироскопа, необходимое для интегрирования, определяется по моменту начала экспозиции последнего кадра: ¡гм-2 ^ ¡эк-1 < ¡гм-1. Каждый момент ¡эк располагается между двумя последовательными измерениями гироскопа с номерами шк-1 и тк: ¡Гтк-
2 ^ ¡эк < ¡Гтк.
¿эО *Э0 + ^э ^эк Ьк + Тэ
Рис.1. Взаимное расположение моментов времени измерения гироскопа и начала экспозиции очередного кадра.
Стрелки направлены из начальной точки поворота, заданного соответствующим кватернионом, в конечную
Матрица 80 вычисляется из решения дифференциального уравнения для нормированного кватерниона ц = [^0 Ц1 Ц2 Цз]\ ||я|| = 1, задающего поворот камеры от СКК в момент времени ¡Г0 к СКК в момент ¡ [13, 14]
Я = 2Я ою(/), ц(/Г0) = Ц0,
где ^0 - скалярная, д\, ^2, qъ - векторные части кватерниона; Я0 - начальное значение кватерниона в момент /Г0; о - символ кватернионного умножения. Распространённым методом численного решения кватерни-онных уравнений ориентации на коротких промежутках времени является метод последовательных приближений Пикара 3-го порядка [13] с линейной интерполяцией угловой скорости внутри шага интегрирования [¡Гш-1, ¡Гш]:
ю(т) = ^ + А9„
ТГ
Т 2
V ТГ
1
2Тт
Г
где т = ^ - ¡Гш-1, те[0, Тг]; А9ш = 9ш - 9ш-1. Линейная интерполяция позволяет аналитически проинтегрировать отдельные слагаемые в приближениях Пикара и рассчитать ненормированное приращение
АЯт (т) = [А$0т (т) А#1т (т) Д#2т (т) А^т (т)]Т
кватерниона ориентации камеры 1
АЦ_1т (Т) А^2т (Т) А%т (Т)
0„
+ 9т хА9т | А9и
24
2 т2
Т 2 V ТГ
Тт
Г У
Ад0т (т) = 1 —т^Т —
¡т-1
При х = ТГ это приращение соответствует приращению кватерниона
АС1 т = [Д?0т Д?1т Д?2т Д?3т Г = А^т (Тг )
ориентации камеры за весь промежуток времени [/Гт -1, /Гт] между измерениями ИИМ:
A?0m = 1
Aqim A?2m
A^3m
=i íi - ai)+Ss^
2 l 24 J 24
Начальное значение qo = я(/го) кватерниона ориентации камеры выбирается из условия 80 = 13, то есть q(tэо) = [1 0 0 0]т. При таком выборе кватернион ql = q(/п) рассчитывается в виде
qi = Aqi(tgo - tro) ° Aqi, qi =
qi
llqi
где знак ~ обозначает комплексное сопряжение кватерниона, q = [^0 - Ц1 - - ?з]т. Кватернион ориентации камеры в момент времени /Гт, qm = q(íГm), т ^ 2, вычисляется после получения очередного измерения от гироскопа:
qт = qт-1 ° Дqт , qm =
Кватернионы {qm }"=-0 вычисляются независимо от состояния экспозиции камеры. Кватернион fk ориентации камеры в момент /эк, = q(tэk), вычисляется из кватерниона ориентации камеры qmk-l = q(tгmk-l) и приращения кватерниона за промежуток времени [ктк-1, /эк]:
?к = qтк-1 ° Дqmk (/эк - т-1),
fk = [ f0k fik f2k f3k ]T = 77"S_
Из компонентов кватерниона fk рассчитываются элементы искомой матрицы S0 ориентации кадра:
sii = i - 2( f22k + fik), s'2 = 2(fokf3k + ff k), Sl3 = 2( fikf3k — f0kf2k ), Sk = 2( fikf2k — f0kf3k ),
s22 = i - 2(fi2 + f32k), s23 = 2(f0kfik + f20kf3k ),
sii = 2( f0kf2k + fikf3k ), = 2( f3kf2k - f0kfik ),
S? = i - 2(fi2 + f2k ).
2.2. Внутрипиксельная дискретизация кадров
Для определения координат изображения звезды с субпиксельным разрешением внутри чувствительной области каждой фотоприёмной ячейки МФП выделяются Ыд ^ 1 несовпадающих точек - узлов внутрипиксель-ной дискретизации. Положение узла с порядковым номером п = 0...Ыд - 1 внутри ячейки задаётся вектором смещений [Ах(п) Ау(п)]т относительно её левого верхнего угла. Расположение узлов минимизирует матема-
тическое ожидание {й(А, $)) расстояния между истинным положением звезды и узлом внутри ячейки (А, $):
а а
{й(И, 1)) = | йх| йу й(х, у)р(х, у | И, 14"),
0 0
й(х,у) = шш {((х- Дх(п))2 + (у - Ду(п))2)1/2},
п=0...Ыд-1
где й(х,у), 0 ^ х, у < а - расстояние между геометрическим изображением звезды и узлом с номером п; р(х, у | А, - плотность распределения изображений звезды внутри ячейки (А, $). Координаты узлов внут-рипиксельной дискретизации для равномерной плотности р(х,у | А, = а 2 и Ыд = 3 перечислены в табл. 1. Для дальнейших расчётов будет использоваться именно эта конфигурация узлов.
Табл. 1. Координаты узлов внутрипиксельной дискретизации в единицах длины стороны пикселя а
Номер Ах(и) A.V(«)
i 0,2 i 0,50
2 0,69 0,76
3 0,69 0,23
Каждый кадр в накапливаемой последовательности подвергается внутрипиксельной дискретизации по отдельности. Узел с номером nk = 0 ... N - 1 внутри пикселя (hk, Wk), hk = 0...H-1, Wk = 0...W-1, расположенный в кадре k = 0.. .K-1, далее будет обозначаться (nk, hk, wk). Векторные координаты
П(Пк, hk, l) = [Ц (nk, hk, Wk) Цу (nk, hk, l )]T этого узла на полноразмерной матрице:
, hk, w k) = -a
hk - h0 + Ax(«k)
_ Hk - w0 _ Ay(«k) _
(3)
В момент времени /Э0 узлу по кадра 0 соответствует некоторое постоянное направление в пространстве предметов. Координаты
S0 («0, h?0,1^0) =
S^0(«0, h), И0)
Sy0 («0 , hb 1!0)
Sz 0(^0, 4,1^0)
единичного направляющего вектора ?(п0, /0, И0) этого направления в СКК0 получаются из (2):
S0 («0, К, ÍV0) =-
i
^f («0, ha, ÍV0)
-п(«0, К, ÍV0) F
(4)
где пКп, Ак, ш) = (|п(пк, Ак, ^к)|2 2)1/2. В СККк у этого вектора получается новый столбец координат
Sk («0, h¡0, ÍV0) =
Sxk («0, h0, ÍV0) Syk («0, К, w0) Szk («0, h¡0, TÍ-0)
= S° s0(«0,К,H0). (5)
k
Из компонент 8к (пк, кк, $к) находятся векторные координаты
4к («0, 4,1^0) = [§ хк («0, 4,1^0) | ук (П0, 4, 1^0 )]т
точки в плоскости кадра к, в которую проецируется направление ?(п0, И0,Н0) в момент времени ¡эк
(«о, Ид, Но) = --
F
„W («о, Ио, Но)
„Хк («о, Ио, Wo) («о, Иго, Но)
(6)
У вектора 4к («к, к к, $к), в отличие от (3), появился подстрочный индекс к номера кадра. Выражение (3) даёт векторные координаты дискретного узла (пк, кк, $к) в плоскости «своего» кадра к, т.е. множество значений вектора ц(пк, кк, $к) действительно, но дискретно. Выражение (6) позволяет найти векторные координаты точки в плоскости кадра к, в которую перешёл дискретный узел («, кк, $к) из кадра 0. Множество значений вектора 4к (пк, кк, $к) действительно и непрерывно. После подстановки (4 - 5) в (6) эти координаты записываются в форме известных проективных преобразований одной плоскости в другую
(«о, Ио, Но) = -F
Векторные координаты п(«0, Ьо, $0) узла (п0, Й0, $0) в плоскости «своего» кадра 0, в правой части этого преобразования, рассчитанные согласно (3).
Точка 4к («к, кк, $к) попадает в кадре к в пиксель с растровыми индексами (1):
Г Vй к „12 к „13 к П(«о, Ио, Но)
_ „21 „ 22 к с 23 к _ - F
[„31 „32 к „33] П(«о, Ио, Но)
-F
Ик («о, Ио, Но)
W к («о, Ио, Но)
(
= floor
\ к («о, Ио, Но)
Ио
и заменяется ближайшим узлом внутрипиксельной дискретизации с номером
шк = arg min {| dk(no,h0,Wo | m)|},
m= 0...Лд-1
где dk (n0, ho, w0 | ш) - вектор, отложенный из точки tk (nk, hk, w£) в узел ш = 0...Лд—1 внутри пикселя
(hk(По,4,Wo), Wk(По,ho,WWo)),
л , i -I ч ^к («о, Ио, Но) , Ак(«о,Ио, Н |m) = --+
Ио
Wo
Ик («о,Ио, Но) Н к («о, Ио, Но)
Ax(m) Ay(m)
Таким образом, узел (п0, И0, Н0) из кадра 0 переходит в узел (тк, Ик, Нк | п0, И0, Н0) в кадре к. Далее для сокращения обозначений растровые индексы начального пикселя (4, Н0) и текущего пикселя
(Ик, ) будут обозначаться одинаковыми строчными буквами. То есть запись (тк, Ик, Нк | п0) следует читать как (тк, /гк, Нк | п0, И0, Н0).
2.3. Обработка последовательных кадров в многоканальном накопителе
Каждый канал многоканального накопителя настраивается на определённый номер п0 = 0 ... Ыд - 1 узла внутри всех пикселей кадра с номером 0. Таким образом, в накопителе содержится Ыд независимых, параллельно работающих накопительных каналов.
Новый узел (тк,Ик, Нк | п0) рассматривается в накопительном канале п0 как начальная точка траектории смаза для возможного изображения звезды в кадре к. Вдоль этой траектории в канале п0 синтезируется ядро смаза Итк (р,ц;Ик, Нк | п0) [1] с использованием измерений ИИМ, накопленных за промежуток времени [¡эк, э + Тэ ]. Здесь индексы (Ик, Нк) и номер тк являются параметрами ядра смаза, задающими пиксель и узел внутри этого пикселя, в котором начинается траектория смаза. Индексы (р, ц) являются аргументами дискретной параметрической функции Итк (р, ц; Ик, Нк | п0). Ядро смаза неоднородно в плоскости изображения (при переходе от одного пикселя (Ик, Нк) к другому структура ядра изменяется), поэтому с ядром Итк (р,ц;Ик, Нк | п0) рассчитывается единственное значение корреляции, определённое только в том пикселе (Ик, Нк), для которого синтезировано ядро:
Л (4, Нк | п,) =
И-1 (7)
= XX Итк(p, Ик, | п0) ь(p,
р=0 =0
где {1к (р, д)}И,д=!Г-1 - матрица цифрового изображения, полученная в кадре к. Если возможное изображение звезды реализуется в рассматриваемом узле, то значение 1к(Ик, Нк | п0) содержит корреляционный максимум с наилучшим отношением сигнал / шум. В противном случае в этом пикселе реализуется смесь ослабленного сигнала со сглаженным фоновым шумом, или просто сглаженное значение фонового шума.
Очевидно, что для каждого узла (па к0, $0) существует единственный узел (тк,Ик, Нк | п0) (обратное неверно). Отсюда, при вычислении корреляции (7) в накопительном канале п можно сразу учесть поворот этого кадра относительно кадра 0. Для этого результат вычисления Л (кк, $к | нужно присвоить не элементу (кк, $к), как это следует из определения корреляции, а элементу (Й0, $0) в канале щ
Jk(ко,Но | «о) = 1к(Ик, Нк | «о);
(8)
где Jk (Ь0, $0 | п0) - значение пикселя промежуточного изображения ук (4, Н0 | п0)}И0—Н0,101.
Для накопления изображений в каждом накопительном канале выделен отдельный аккумулятор
{А(/0,Н0 | Яо)}^1. Все аккумуляторы инициализируются результатами фильтрации (8) начального кадра 0
АоФо, 14-0 | По) = Jо(/?о, 14-0 | Яо),
где А)(/г0,Н01 п0) - содержимое аккумулятора А(/0, Н01 п0) после обработки кадра 0. После обработки очередного кадра к > 0 к аккумуляторам прибавляются значения из промежуточных буферов:
Ак (4, Но I По) = Ак^фо, Но I По) + Jk (К, Но | По),
где Ак (4, Н01 я0) - содержимое аккумулятора А(/, Н01 я0) после обработки кадра к. После получения последнего кадра с номером К-1 в накапливаемой последовательности все аккумуляторы нормируются для получения Ыд накопленных изображений:
и (/о, Но | По) = -1 Ак _!(/,, Но | По) .
к
3. Моделирование многоканального накопления
Моделирование многоканального накопления обобщает моделирование однокадровой согласованной фильтрации в [2] на случай, когда яркости звезды недостаточно для её обнаружения в одном кадре. Модельные звёзды наблюдаются камерой, параметры которой перечислены в [2, табл. 1].
Модельное накопление выполняется многоканальным накопителем с координатами узлов, перечисленными в табл.1. В момент tэо начала экспозиции кадра к=0 геометрическое изображение модельной звезды находится в узле № 1 внутри фотоприёмной ячейки с растровыми индексами
(Ац, ^ц) = (3240, 2430), в которой расположена главная точка изображения. В этом описываемая модель многокадрового накопления отличается от модели [2] однокадровой фильтрации, в которой модельная звезда располагалась точно в центре ячейки (Ац, в момент начала экспозиции. Моделирование проводится для угловой скорости ю = [0,63 1,90 0]т °/с, |ю| = 2 °/с при частоте кадров 1/ Тк = 20 Гц (межкадровый период ТК = 50 мск).
В качестве опорной звезды для моделирования выбрана Вега со звёздной величиной т = 0,03. Спектр излучения Веги был линейно отмасштабирован для получения модельной звёзды со звёздной величиной т = 3,39. Моделирование излучения фона выполняется для зенитного угла оптической оси цифровой камеры в 0° (камера смотрит в зенит) и зенитного угла Солнца в 40° [15]. Начальной высотой наблюдения выбирается /нбл = 8 км, т.к. на этой высоте спектральные плотности освещённости от фона и от Веги в полосе пропускания объектива очень близки (рис.2). Время экспозиции для этого случая, рассчитанное согласно методике [2] и использованное для моделирования однокадровой фильтрации, ТЭ = 5 мс.
400 500 600 700 800 900 1000
Длина волны, нм
Рис. 2. Графики спектральной плотности освещённости входного зрачка объектива, создаваемые атмосферным фоном и модельными звёздами. Полоса пропускания объектива выделена светлым прямоугольником
Эффект от накопления последовательности кадров демонстрируется для К = 30, что эквивалентно времени измерения одного значения ориентации (К-1)ТК = 1,45 с. Численные значения энергетических параметров наблюдаемой сцены приведены в табл. 2 с использованием обозначений [2]. Максимальное число электронов в одном пикселе кадра со смазанным изображением звезды обозначено
/ гшах \ vепикс/ ■
Табл. 2. Энергетические параметры сигналов от модельных звёзд. Параметры фона приведены для /нбл = 8 км
т 0,03 (Вега) 3,39
Еезв >[е_ /(см2 Х с)] 232489 10529
< -и Х[е" ] 10144 459
< -егмкс >,[е" ] 647 29
Еефон >[е_ /(см2 Х с)] 202126 202126
< 4юН Х[е_] 8819 8819
°ефон >[е_ ] 94 94
Смазанные изображения звезды, полученные в одном кадре в различных условиях наблюдения, показаны в [2, рис. 3, 4]. На рис. 3 показаны результаты накопления изображения одной и той же звезды, полученной при |ю| = 2 / с, в трёх накопительных каналах при различных фоновых шумах. Черным прямоугольником выделен ярчайший пиксель изображения, линия - траектория смаза, жирная точка - начало траектории. Истинное положение изображения звезды в момент времени ^о - узел 0 в центральном пикселе, т. е. звезда должна быть обнаружена в центральном пикселе изображения, накопленного в канале 0.
На рис. 3а - в показаны результаты накопления кадров без фонового шума для звезды с т = 3,39 в трёх накопительных каналах. Фон в этих кадрах присутствует только в качестве постоянной «фоновой подставки», равной математическому ожиданию фона на высоте /НБЛ = 8 км. Уровни сигналов в ярчайших пикселях изображений, накопленных в отдельных каналах, перечислены в табл. 3. Из этой таблицы видно, что в идеальных условиях, вследствие плохих корре-
ляционных свойств ядра смаза, абсолютный корреляционный максимум, соответствующий истинному положению звезды в канале 0, мало отличается от максимумов в других накопительных каналах.
На рис. 3г - е показаны те же самые кадры, накопленные с фоновым шумом на высоте /НБЛ = 8 км. Здесь звезда обнаруживается во всех трёх каналах, абсолютный максимум проявляется в канале 1 (табл. 3). При этом ни один максимум в каналах обнаружения не проявляется в пикселе с истинным положением звезды. Т.е. при К = 30 в таких условиях наблюдения координаты изображения звезды будут определены с погрешностью порядка размера пикселя.
На рис. 3ж - и показаны накопленные изображения той же самой звезды, но наблюдаемой на высоте /НБЛ = 25 км. Здесь абсолютный максимум проявляется в канале 0, что соответствует истинному положению звезды (табл. 3).
Табл. 3. Значения ярчайших пикселей (в единицах младшего разряда) накопленных изображений для т = 3,39
Условие Без шума НнБЛ = 8 км Ннбл = 25 км
Канал 0 3512,49 3513,15 76,72
Канал 1 3512,36 3513,22 76,53
Канал 2 3512,35 3513,20 76,64
Заключение
Многоканальное накопление последовательных кадров, рассмотренное в работе, позволяет обнаруживать смазанные изображения тусклых звёзд в сильных фоновых шумах и измерять координаты этих изображений с субпиксельным разрешением. Обязательным этапом обработки каждого кадра, предшествующим накоплению, является цифровая фильтрация, согласованная со смазом. Предлагаемый алгоритм обобщает известные корреляционные алгоритмы многоканального обнаружения точечной цели на случай, когда нельзя пренебречь движением изображения цели на интервале экспозиции, а направление этого движения на межкадровом интервале может произвольно изменяться.
Движение изображения учитывается за счёт измерения вектора угловой скорости камеры при помощи гироскопов. Из этих измерений рассчитывается ориентация плоскости очередного кадра относительно плоскости накапливаемого изображения. Ориентация позволяет связать каждый пиксель накапливаемого изображения с пикселем в очередном кадре и указать на узел внутрипиксельной дискретизации, в котором нужно синтезировать ядро смаза и выполнить согласованную фильтрацию.
Если в каком-либо пикселе очередного кадра начинается смазанное изображение звезды, то после согласованной фильтрации в этом пикселе проявляется корреляционный максимум. Если сигнал от звезды в пикселе отсутствует, то его профильтрованное значение представляет собой корреляцию между ядром смаза и реализацией изображения фона.
Ядро смаза обладает плохими корреляционными свойствами: значения пикселя, содержащего максимум автокорреляционной функции ядра, мало отличается от значений других пикселей из ближайшей окрестности этого максимума. Поэтому число накапливаемых кадров, требуемое для измерения координат смазанного изображения звезды с субпиксельным разрешением, заметно превышает число накапливаемых кадров, достаточное для простого обнаружения этого изображения.
Плохие корреляционные свойства ядра смаза не позволяют разделить смазанные изображения двух близких звёзд. Если расстояние между начальными пикселями смазанных изображений меньше половины длины смаза, то два накопленных изображения сольются в одно общее пятно. Максимум в этом пятне, скорее всего, удастся обнаружить. Но невозможно указать, к какому из двух слившихся изображений он относится.
Максимальное количество накапливаемых кадров ограничено временем нахождения звезды в поворачивающемся поле зрения цифровой камеры. Кроме указанного геометрического ограничения, максимальное время накопления может ограничиваться нарастающими погрешностями определения ориентации очередного кадра относительно начального. Эти погрешности обусловлены как погрешностями измерения компонент вектора угловой скорости, так и методическими погрешностями алгоритма численного интегрирования этих измерений. Поэтому параметры алгоритма накопления и условия эксплуатации аст-роинерциального датчика нужно выбирать так, чтобы минимизировать количество накапливаемых кадров.
References
[1] Vasilyuk NN. Synthesis of the rotational blur kernel in a digital image using measurements of a triaxial gyroscope. Computer Optics 2022; 46(5): 763-773. DOI: I0.18287/2412-6179-C0-I081.
[2] Vasilyuk NN. Correction of rotational blur in images of stars observed by an astroinertial attitude sensor against the background of the daytime sky. Computer Optics 2023; 47(1): 79-91. DOI: 10.18287/2412-6179-C0-1141.
[3] Liebe CC, Gromov K, Meller DM. Toward a stellar gyroscope for spacecraft attitude determination. J Guid Control Dyn 2004; 26(1): 91-99. DOI: 10.2514/1.9289.
[4] Samaan MA, Mortari D, Pollock TC, Junkins JL Predictive centroiding for single and multiple FOVs star trackers. J Astronaut Sci 2002; 50(1): 113-123. DOI: 10.1007/BF03546333.
[5] Aretskin-Hariton E, Swank AJ. Star tracker performance estimate with IMU. AIAA Guidance, Navigation and Control Conf 2015: AIAA2015-0334. DOI: 10.2514/6.20160334.
[6] Wang Z, Jiang J, Zhang G. Global field of view imaging model and parameter optimization for high dynamic star tracker. Opt Express 2018; 26(25): 33314-33332. DOI: 10.1364ЮЕ.26.033314.
[7] Wang S, Zhang S, Ning M, Zhou B. Motion blurred star image restoration based on MEMS gyroscope aid and blur
kernel correction. Sensors 2018; 18(8): 2662. DOI: 10.3390/s18082662.
[8] Loke SC. Astronomical image acquisition using an improved track and accumulate method. IEEE Access 2017; 5: 9691-9698. DOI: 10.1109/ACCESS.2017.2700162.
[9] Ni Y, Dai D, Tan W, Wang X, Qin S. Attitude-correlated frames adding approach to improve signal-to-noise ratio of star image for star tracker. Opt Express 2019; 27(11): 15548-15564. DOI: 10.1364/0E.27.015548.
[10] Ma L, Bernelli-Zazzera F, Qin S, Wang X, Ma L. Performance analysis of the attitude-correlated frames approach for star sensors. IEEE Metrology for Aerospace 2016: 81-86. DOI: 10.1109/MetroAeroSpace.2016.7573190.
[11] Kirichuk VS, Kosykh VP, Kurmanbek UT. Algorithm of detection of moving small-scale objects in a sequence of images. Optoelectronics, Instrumentation and Data Processing 2009; 45(1): 8-13. DOI 10.1134/S8756699009010026.
[12] Kirichuk VS, Kosykh VP. Construction of a multichannel filter for detecting point targets in the image formed by a matrix photodetector. Optoelectronics, Instrumentation and Data Processing 2012; 48(5): 497-505. DOI 10.3103/S875669901205010X.
[13] Branets VN, Shmyglevski IP, Applications of quaternions in problems on the motion of rigid bodies [In Russian], Moscow: "Nauka" Publisher; 1973.
[14] Molodenkov AV, Sapunkov YG, Perelyaev SE, Moloden-kova TV, Analytical solutions in the Darboux problem, the Bortz equation and the approach to orientation algorithm of sins based on them [In Russian]. Applied Mathematics and Mechanics 2019, 83(4): 586-596. DOI 10.1134/S003282351904009X.
[15] Bessonov RV, Zhukov BS, Karavaeva ES, et al. The basic principles of design of astrocorrector for endoatmospheric vehicles [In Russian]. Current Problems in Remote Sensing of the Earth from Space 2018; 15(6): 21-30. DOI: 10.21046/2070-7401-2018-15-6-21-30.
ж)
з)
и)
Рис.3. Результаты накопления изображений одной звезды с т = 3,39 в трёх накопительных каналах по =0,1,2 при \т\ = 2°/с. (а-в) - накопление без фонового шума при Инбл = 8 км; (г-е) - накопление с фоновым шумом при Инбл = 8 км;
(ж-и) - накопление с фоновым шумом при Инбл = 25 км
Сведения об авторе
Василюк Николай Николаевич, 1976 года рождения, кандидат физико-математических наук, заместитель главного конструктора ООО «НПК Электрооптика». Область научных интересов: оптико-электронные следящие системы, системы навигации и ориентации, методы комплексирования в навигационных системах. E-mail: nik-vasilyuk@yandex. ru ORCID: https://orcid.org/0000-0003-2317-8066.
ГРНТИ: 28.23.15
Поступила в редакцию 28 марта 2023 г. Окончательный вариант - 15 августа 2023 г.
Subpixel stacking and detection of blurred star images observed by an astroinertial attitude sensor against the background of the daytime sky
N.N. Vasilyuk1 1 Electrooptika, LLC, 107076, Moscow, Russia, Stromynka, d.18, k.1
Abstract
To determine the attitude from the image of a star field obtained by a digital camera, it is necessary to detect images of individual stars and determine their coordinates with subpixel resolution. This problem is relatively easy to solve if the image of the star field is obtained by a stationary camera observing a clear night sky. The difficulty arises when the stars are observed by the rotating camera against the background of the daytime sky. The blurred images of stars are almost impossible to detect in the background noise in one frame. This difficulty can be overcome by accumulating a certain number of consecutive frames after the matched blur filtering has been performed. This paper considers an algorithm for accumulating the sequence of frames, implemented as a multichannel accumulator, consisting of several accumulative channels of the same type, and operating in parallel. The plane of the accumulated image in each channel coincides with the plane of the first frame of the accumulated sequence. Subpixel resolution in determining the coordinates of detected star images is achieved by introducing intrapixel discretization in the first frame. The orientation of the current frame relative to the first one is calculated from the measurements of the gyroscopes of the astroinertial sensor. When synthesizing the impulse responses of the blur-matched filters for the current frame, the attitude of this frame and the location of the intrapixel discretization nodes of the first frame are taken into account. Calculation formulas are given for determining the number of accumulated frames depending on the angular velocity of the camera rotation, the background intensity and the brightness of the observed stars. The simulation results show the possibility of detecting blurred images of fairly dim stars against the background of the daytime sky and determining their coordinates with subpixel resolution.
Keywords: star tracker, astrocorrector, astro inertial navigation system, blur correction, images stacking.
Citation: Vasilyuk NN. Subpixel stacking and detection of blurred star images observed by an astroinertial attitude sensor against the background of the daytime sky. Computer Optics 2024; 48(2): 303-311. DOI: I0.18287/2412-6179-C0-I309.
Author's information
Nikolay Nikolaevich Vasilyuk (b. 1976), PhD in Physics, deputy chief design officer at Electrooptika LLC. Research interests: optoelectronic tracking and detection systems, navigation and orientation systems, navigation sensors fusion. E-mail: [email protected] ORCID: https://orcid.org/0000-0003-2317-8066
Received 28 March, 2023. The final version - 15 August, 2023.
