Особенности применения некоторых методов дискретной фильтрации для задач сопровождения баллистических объектов, совершающих спиральный маневр Текст научной статьи по специальности «Физика»

_Доклады БГУИР_

2015 № 1 (87)

УДК 621.396.96

ОСОБЕННОСТИ ПРИМЕНЕНИЯ МЕТОДОВ ДИСКРЕТНОЙ ФИЛЬТРАЦИИ ДЛЯ ЗАДАЧ СОПРОВОЖДЕНИЯ БАЛЛИСТИЧЕСКИХ ОБЪЕКТОВ, СОВЕРШАЮЩИХ СПИРАЛЬНЫЙ МАНЕВР

АС. СОЛОНАР, ПА. ХМАРСКИЙ

Военная академия Республики Беларусь Независимости, 220, Минск, 220057, Беларусь

Поступила в редакцию 28 июля 2014

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

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

Введение

Алгоритмы сопровождения баллистических объектов (БО) реализуются с большей вычислительной сложностью, чем алгоритмы сопровождения аэродинамических летательных аппаратов [1-4]. Движение БО за пределами атмосферы хорошо предсказуемо, так как там объект не может маневрировать, поскольку единственной силой, действующей на него, является гравитация [1, 5-8]. Наиболее сложным участком сопровождения является атмосферный, где поведение БО становится менее предсказуемым из-за асимметрии в конструкции БО, которая может вызвать движение по спирали [1, 5, 6]. Кроме того, возникающие замедления могут восприниматься как маневр [1, 6]. Обеспечение более высокой точности измерения траекторных параметров БО достигается совершенствованием алгоритмов сопровождения.

Из открытых публикаций зарубежных и отечественных авторов [1-10] известно, что для сопровождения БО наиболее распространено использование модификаций фильтра Калмана при косвенных измерениях (ФК) [1, 8, 11-13]. Такие фильтры успешно применяются как для сопровождения аэродинамических летательных аппаратов, так и БО. Однако эти фильтры могут иметь большие ошибки сопровождения и несостоятельную оценку при оценивании траекторных параметров БО, совершающих сложный маневр [1, 4, 7]. Учет особенностей динамики БО возможен при использовании ансцентного фильтра Калмана (АФК), реализующего метод нелинейной фильтрации [14, 15].

Целью статьи является проведение сопоставительного анализа результатов фильтрации АФК и модификаций ФК для случая сопровождения БО, совершающего спиральный маневр при возвращении в атмосферу.

Постановка задачи

Начальные условия. Будем считать, что наблюдается БО, входящий в атмосферу; интервал обновления данных равен Т; вектор наблюдаемых параметров БО 8 содержит разовые

оценки дальности г, азимута в и угла места 8, формируемые трехкоординатной радиолокационной станцией (РЛС). Ошибки наблюдения распределены по гауссовскому закону и определяются корреляционной матрицей ошибок измерения:

(1)

где аг, Ор, ое - СКО гауссовских ошибок наблюдения г, в и 8 .

Задачи: 1) рассмотреть модель движения БО, совершающего спиральный маневр при возвращении в атмосферу; 2) рассмотреть особенности построения модификаций фильтров Калмана при косвенных измерениях и ансцентного фильтра Калмана при сопровождении БО; 3) методами математического моделирования провести сопоставление суммарных ошибок измерения местоположения БО для АФК и модификаций ФК.

2 0 0

*е = 0 2 0

0 0 о2

Модель движения БО, совершающего спиральный маневр

При синтезе фильтров сопровождения БО чаще всего рассматриваются следующие модели движения:

1) модель движения со случайным маневром в виде гауссовского белого шума [2, 8, 12, 13];

2) модель движения со случайным маневром в виде винеровского процесса [2];

3) модель движения БО, входящего в атмосферу [2, 6];

4) модель движения БО, совершающего спиральный маневр при вхождении в атмосферу [4].

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

Уравнения движения БО. Уравнения движения объекта приводятся в системе координат Север-Восток-Низ ('ЫоПк-Еаз^Вомт), изображенной на рис.1, а. Земля считается сферической и не вращающейся.

Ускорения БО ах, а2, а^ представляются как сумма гравитационной и аэродинамической составляющих: ах = Ох + Д, аг = Ог + Д, ау = С + Ау где Ох у г и Ах г у - гравитационные и аэродинамические составляющие ускорения х, г, у.

Гравитационная составляющая ускорения БО учитывается следующим образом [1,2]:

х

СУ

М

( х2 + г2 + (у - ЯТ )2 )3/2

У - Х

(2)

где м = 3,986 • 1014 м3/с2 - гравитационная постоянная Земли, Х « 6371 км - радиус Земли.

Аэродинамические силы, воздействующие на БО, могут быть разделены на две основные составляющие: силу сопротивления вдоль вектора скорости и нормальную силу, воздействующую на БО перпендикулярно вектору скорости. Спиральное движение может быть аппроксимировано нормальным вектором ускорений, вращающимся с некоторой скоростью в перпендикулярной плоскости. Аэродинамические составляющие ускорения представляются при помощи следующих уравнений [1]:

Ax cos у cos x - sinx sin у cosx

A = cos у sin x cosx sin у sinx

A sin у 0 cos у

р V Г S

(

где у = arctan

\

ylV+V

■2 7

, х = arctan

2m

гуЛ

v v у

-Cd

CN cos ф Cn sin ф

горизонтальным и вертикальным углы наклона

вектора скорости (см. рис.1); ф - угол вращения БО, зависящий от времени так, что скорость

его изменения равна круговой частоте вращения со: ф=со; Ут - полный вектор скорости БО; и

т - площадь сечения миделя и масса БО; р - атмосферная плотность воздуха; Св, Сы -

коэффициенты сопротивления и подъемной силы. Псевдокод алгоритма расчета траектории БО приведен на рис. 1, б.

Для различных типов БО круговая частота ю и коэффициенты С0, См зависят от формы объекта, размеров и массы. Для наблюдателя эти параметры являются неизвестными.

-Vy s V

> ^ .У VT

Y Vx „

/ X / / / /

1

l|T

[at+,] = trajectoryBT ]

T

Вход: И» = IX Vx Zk V7k Ук Vyk Фк||

'■k ||xk у xk V zk

Выход: <L, = IV^s z

l|T

к+1 = ||А(к+1) кх(к+1) ^(к+1) Кг(к+1) У(к+1) "у (к+1) •

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

I I V, 1 |тг |2

Y = arctan I

= arctan I v7k I Vt

\2

V Tk I

,/(Vxk )2 + Vzk )2 )

2. Определение атмосферной плотности воздуха для высоты yk : р(у )=c0el

3. Расчет угла вращения: фА+1=фА+юТ .

4. Вычисление составляющих ускорения баллистического объекта. cos у cos х - sinx sin у cos х

cos у sinx cosx sin у sinx

y -R siny 0 cosy

(xl + z2+ (yk - R )2)

Pfa) iV-k S

2m

CN c0s ф2+1 CN sin фм

5. Определение прямоугольных координат и скоростей их изменения:

1 Т Т2/2

yCk+l

V

0 l Т 0 0 l

6. Вернуть значения ¡.

а б

Рис. 1. Система координат БО (а) и псевдокод алгоритма расчета траектории БО (б)

Модель движения БО для алгоритма фильтрации

Модель движения, приведенная на рис.1, б, учитывает траекторные параметры БО (местоположение и скорость) и априорно неизвестные параметры: аэродинамические коэффициенты Св, Сы, массу т, площадь сечения миделя S, угол вращения ф и круговую

частоту вращения ю. Использование радиолокационных измерений не позволяет однозначно оценить все эти параметры. Для снижения количества неизвестных параметров вводят модифицированные коэффициенты сопротивления и подъемной силы [1]:

СОт = С / т СЬт = / т .

Более того, вводятся две новые составляющие вектора состояния: ^, Z2, называемые гармоническими параметрами состояния, которые определяются следующим образом -Zl =СЬт со^сЫ + ф(1), Z2 =СЬт зт(со/ + ф0) . Эти параметры описываются соответствующим

уравнением динамики:

Зависимость атмосферной плотности воздуха от высоты у определяется экспоненциальной функцией р(у)=е0 еСУ, где с0 = 1,293 кг/м3, с = 1,694 хШ^м-1.

С учетом введенных коэффициентов расширенный вектор состояния будет иметь вид:

0 -go Zi

¿2 ш 0 Z2

а =

X V

Zl Z2 CDm ®

x 7 V У V Ускорение БО будет вычисляться согласно выражению

(4)

z

T

a.

(x2 + z2 + (y — RT )2 )3/2

x cos у cos X —sinx sin у cosx

z + cos у sin X cosx sin у sinx

y—rt sin у 0 cos у

2

—C

Dm

Z

(5)

Законы изменений во времени оценок модифицированного коэффициента сопротивления СВт и круговой частоты ю считаются винеровскими. Псевдокод расчета модели движения БО для алгоритма фильтрации приведен на рис. 2.

[а к+1 ] = езФ-ВТ [а к ]

вход: ак = ||х У1к 2к Ул ук Уук 21к 21 к Сшк ю1||Т.

Выход'. ак+1 = ||Хк+1 ^Х(к+1) гк+1 ^г (к+1) ук+1 ^у (к+1) ^1(к+1) ^2(к+1) ССм(к+1) Ю(к+1)|| 1. Вычисление горизонтального и вертикального углов наклона вектора скорости:

Y = arctanl

ЩШ? J" X = arctan(, N* = ^)2 + V)2 + (Vyk)2

2. Модифицированный коэффициент сопротивления считается постоянным: Сс

3. Нахождение гармонических параметров состояния:

4. Вычисление составляющих ускорения баллистического объекта:

0U

(x2 + + (y - Rt )2 )3/

xk cos у cos X —sinX slnуcosx

zk + cosуsin X cosX sin уsinx

yt — RT sin у 0 cos у

5. Определение прямоугольных координат и скоростей их изменения:

x(k+1) xk Z(k+1) zk y(k+1) yk 1 T T2/2

Vx(k+1) = B Kk Vz(k+1) = B Kk Vy (k+1) = B Vyk , B = 0 1 T

ax(k+1) axk az(k+1) azk ay(k+1) ayk 0 0 1

6. Вернуть значения акл

Рис. 2. Псевдокод алгоритма расчета модели движения БО для алгоритма (

жльтрации

Особенности построения некоторых модификаций фильтров Калмана при

сопровождении БО

1. Фильтр Калмана при наблюдении прямоугольных координат с коррелированными ошибками наблюдения [11]. В этой модификации оценка вектора состояния

а =

x К ax z К az y Vy ay

формируется по результатам наблюдения разовых

оценок прямоугольных координат 0х2у=||х^ у. . Для пересчета оценки вектора

наблюдения из полярной системы координат в прямоугольную применяется вектор -

функция Ь0 хг(0) :

К = х/

zf y А

(6)

h0xz(0) = r sin (в) sin (s) r cos (в) sin (e) rcos ( s )

При этом необходимо пересчитать корреляционную матрицу ошибок текущих измерений R0 в прямоугольную систему координат при помощи выражения [8, 11-13]:

os (в) sin (s) —r sin (в) sin (s) r cos (в) cos (s)

in (в) sin (s) r cos (в) sin (s) r sin (в) sin (s) cos (s) 0 —rsin ( 8 )

R11 R12 R

R f = HxzyReHTzy = R12 R 22 R

R13 R23 R

H

~xzy

где Rn = o^2 (cos2 (в) cos2 (s) + sin2 (в)sin2 (s)) + o2 cos2 (в)sin2 (s), R

= O^ 2sin

in2(s)-

o^os2 (s),

R12 =

o^2 sin (2в) cos (2в) o2 sin (2в) sin (2S) cos (в), 2 . --V +-0 sin2 (s), Rl3 =-o-(Q — o2),

R22 = or (cos2 (в)sin2 (s) + sin2 (в)cos2 (s)) + o2 sin2 (в)sin2 (s), R23 = sin(2s)sin(в) (o2 — o2). Неслучайная динамическая матрица пересчета приращения вектора состояния:

a

x

a

z

2

T

Z1(k+1) = Z2k®T + Z1k; Z2(k+1) = Z1(k+1)®T + Z

c, ecyk K

Z

a

zk

2

a

T

b 03x3 03x3 0 0 0 1 T T2 /2

B = 03x3 b 03x3 , где 03x3 = 0 0 0 , b = 0 1 T

03x3 03x3 b 0 0 0 0 0 1

S =

Матрица пересчета случайных приращений вектора состояния:

Г3/6 Г2/2 Г 0 0 0 0 0 0 0 0 0 Г3/6 Т2 / 2 Т 0 0 0 0 0 0 0 0 0 Г3/6 Г2/2 Г

100000000

Матрица статистического пересчета: н =

000100000 000000100

маневра

характеризует случайный маневр по трем

Дисперсия случайного прямоугольным координатам.

2. Фильтр Калмана при наблюдении прямоугольных координат с независимыми ошибками наблюдения [11]. В данной модификации оценки полярных координат также пересчитываются в прямоугольные с помощью (6), но при этом принято, что ошибки наблюдения x^, z^, y^ между собой не коррелированы. Подобное допущение позволяет

фильтровать прямоугольные координаты x, y, z раздельно при помощи трех идентичных ФК. В этих фильтрах:

- корреляционная матрица флуктуационных ошибок текущих измерений будет представлена матрицами единичного размера: для ФК (канал x) R. = R; для ФК (канал z)

R/z = R22; для ФК (канал y) R/y = R33;

- матрица статического пересчета: H = |

- неслучайная динамическая матрица пересчета приращения вектора состояния соответствует выражению (7);

ф 0 о||;

T 3/6 0 0

T2/2 0 0

T 0 0

- матрица пересчета случайных приращений вектора состояния S =

- дисперсия случайного маневра оМ будет одинаковой для х, у, z.

3. Ансцентный фильтр Калмана. Реализация ансцентного фильтра Калмана подробно описана в работах [14, 15]. Особенности реализации ансцентного фильтра Калмана при сопровождении БО определяют: вектор состояния (4); корреляционную матрицу ошибок измерений ^ ; функцию экстраполяции (рис. 2); функцию пересчета вектора состояния в

систему координат вектора наблюдения: 0 =

r -y/x2 + z2 + y2

ß = arctan( z / x) + (n/2) [2 - sgn( x) - sgn(x) sgn( z)]

£ arctan(Vx2 + z2 / y) + n [1 - sgn( y)]

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

M = S0QX + SiQzST + S2QmST, где so =

T2/2 T 0 00 0 0000 0 0 T2 / 2 T 0 0 0000 0 0 0 00 T2/2 T 000

Si =

Qa

0 0 0 0 0 0

0 0 0 0 0 0

oa 0 0

0 oa 0 , Qz

0 0 oa

T

-fflT2 /2

roT2 /2

T

S2 =

0 0 T2/2 T 0 0 0 0 0 0 0 0

T 0

ö,

0

O,

'Cdm

2

0

ü

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

модифицированному коэффициенту сопротивления соответственно [1].

С

От

и круговой частоте ю,

Результаты сопоставительного анализа

Сопоставительное моделирование результатов фильтрации АФК и модификаций ФК проводилось при помощи комплекса математического моделирования. Комплекс формировал входное воздействие в виде суммы задающего и возмущающего воздействия. Задающее воздействие представляло собой атмосферный участок траектории БО, совершающего спиральный маневр. Возмущающее воздействие на каждом интервале обновления данных задавалось в виде случайных ошибок измерений полярных координат, с нулевым математическим ожиданием и корреляционной матрицей ошибок измерения (1). Условия проведения моделирования соответствуют работе [3] для случая с постоянной круговой частотой (модельный эксперимент № 1) и круговой частотой, изменяющейся по линейному закону (модельный эксперимент № 2).

Показатели качества. Показателем качества АФК и модификаций ФК являлось СКО суммарной ошибки измерения местоположения А (к), и ее среднее значение ДЕ(к), которые вычислялись

согласно выражениям: АЕ (к) =

7Т £ ([ ^ (к)" ХФ-< (к )]2 ^ (к)" ^ (к)]2 ^ (к)" (к)]2)

оп 7=1

Д = к=1

ДЕ =

£ А, (к)

N

где N - количество проведенных опытов (под одним опытом понимается

результат фильтрации траектории в течение 350 интервалов обновления данных); к - номер шага наблюдения; хТ ■ (к), - (к), у7] (к) - требуемые значения фазовых координат;

Хф7 (к), (к), уф. (к) - фильтрованные значения фазовых координат; N - количество

измерений в одном опыте. Значения (к) нормировались к СКО текущих оценок для 5000

проведенных опытов. На рис. 3, а-б представлены значения суммарной ошибки измерения местоположения для модельных экспериментов № 1 и № 2.

28 20"

Фильтр Калмана при наблюдении прямоугольных координат с коррелированными ошибками наблюдения

Ансцентный фильтр Калмана

Фильтр Калмана при наблюдении прямоугольных координат с независимыми ошибками наблюдения

/ «V

Фильтр Калмана при наблюдении

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

Ансцентный фильтр Калмана

Фильтр Калмана при наблюдении

прямоугольных координат с независимыми ошибками наблюдения

номер шага обновления данных

номер шага обновления данных

а б

Рис. 3. Результаты сопоставительного моделирования: а - модельный эксперимент № 1; б - модельный

эксперимент № 2

Результаты моделирования. Анализ рис. 3, а, б и данных таблицы показывает, что все предложенные алгоритмы фильтрации позволяют уменьшить ошибку измерения. Ансцентный фильтр Калмана имеет лучшие показатели качества, чем модификации ФК (в случае с постоянной частотой вращения: до 15 % относительно фильтра Калмана при наблюдении прямоугольных координат с коррелированными ошибками наблюдения и до 10 % относительно фильтра Калмана при наблюдении прямоугольных координат с независимыми ошибками наблюдения; в случае с переменной частотой вращения: до 21 % относительно

д. %

д. %

80

60

О 60

г 40

О 40

20

300 к

к

фильтра Калмана при наблюдении прямоугольных координат с коррелированными ошибками наблюдения и до 18 % относительно фильтра Калмана при наблюдении прямоугольных координат с независимыми ошибками наблюдения). Это объясняется использованием в АФК алгоритма экстраполяции, учитывающего спиральное вращение БО при возвращении в плотные слои атмосферы.

Результаты математического моделирования

Среднее значение суммарной ошибки Д£, м Разовая ошибка Ансцентный фильтр Калмана Фильтр Калмана при наблюдении прямоугольных координат с независимыми ошибками наблюдения Фильтр Калмана при наблюдении прямоугольных координат с коррелированными ошибками наблюдения

Модельный эксперимент №1. Постоянная круговая частота 62,366 9,74 17,186 15,035

Модельный эксперимент №2. Переменная круговая частота 61,681 10,473 27,02 25,349

Заключение

В ходе проведенных исследований была рассмотрена модель движения БО, совершающего спиральный маневр при вхождении в атмосферу; представлены особенности построения модификаций ФК и АФК при сопровождении БО. Модель БО рассмотрена в первой части, особенности построений фильтров перечислены во второй части. В третьей части проанализировано качество измерений координат БО в различных модификациях фильтров. Результаты математического моделирования демонстрируют преимущество АФК в сравнении с рассматриваемыми модификациями ФК при наблюдении маневрирующего БО.

FILTERING FEATURES FOR TRACKING OF SPIRALING REENTRY VEHICLES

A.S. SOLONAR, PA. KHMARSKI Abstract

The target dynamics model for tracking of spiraling reentry vehicles is considered. The features of Extended Kalman filter modification and Unscented Kalman filter are listed. The estimation performance of Extended Kalman filter modification and Unscented Kalman filter is compared.

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

1. Kim J., Vaddi S., Menon P., Ohlmeyer E. // IEEE transaction on A&ES. 2012. Vol. 48. №1. P. 313-328.

2. Ohlmeyer E., Menon P., Kim J. // AIAA Guidance, Navigation and Control Conference. 2010.

3. Kim J., Vaddi S., Menon P., Ohlmeyer E. // AIAA Guidance, Navigation and Control Conference. 2010.

4. Vaddi S., Menon P., Ohlmeyer E. // AIAA Guidance, Navigation and Control Conference, 2007.

5. Moon K.R., Kim T.H., Song T.L. // JASS. 2012. Vol. 29. № 4. P. 363-374.

6. Zarchan P. // IEEE Position, Location, and Navigation Symposium, San Diego, CA, 2000.

7. Farina A., Ristic B., Benvenuti D. // IEEE transaction on A&ES. 2002. Vol. 38. № 3. P. 854-867.

8. Кузьмин С.З. Основы теории цифровой обработки радиолокационной информации. М., 1974.

9. Солонар А.С., Хмарский П.А. // Матер. МНТК, приуроченной к 50-летию МРТИ-БГУИР. Минск, 1819 марта 2014. С. 78-79.

10. Солонар А.С.,Хмарский П.А. // Тез. докл. 30 НТК «Агат-Системы управления», 2013.

11. Солонар А.С., Хмарский П.А. // Докл. БГУИР. 2012. № 7. С. 47-53.

12. Ширман Я.Д. Радиоэлектронные системы: Основы построения и теория. Справочник. М., 2007.

13. Blackman S. Design and analysis of modern tracking systems. Boston, London, 1999.

14. Julier S., Uhlmann J. // IEEE Trans. on Automatic Control. 2000. Vol. 45. № 3. P. 477-482.

15. Солонар А.С.,Хмарский П.А. // Докл. БГУИР. 2013. № 2. С. 79-86.