Автоматика и телемеханика, № 7, 2023

Управление в технических системах

(АО «Концерн Радиоэлектронные технологии», Москва),

А.А. МОЗОЛЬ, канд. техн. наук (amozol@bk.ru)

(АО «ВНИИ «Градиент», Ростов-на-Дону)

ПЕРИОДО-ВРЕМЕННОЙ МЕТОД ПАРАМЕТРИЧЕСКОЙ

ИДЕНТИФИКАЦИИ ДЛЯ РЕШЕНИЯ

ЛОКАЦИОННЫХ И НАВИГАЦИОННЫХ ЗАДАЧ

Применительно к локационным и навигационным задачам для однопо-

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

тификации параметров полиномиальной модели движения объекта с уче-

том эволюции невязки между периодическим излученным и принятым

квазипериодическим сигналом. Рассматривается прохождение сигнала в

произвольной физической среде, при этом не требуются знание периода

излученного сигнала и традиционное оценивание текущей частоты До-

плера, вызывающей невязку. Метод основан на подсчете числа периодов

принимаемого сигнала в заданном интервале наблюдения. Рассмотрены

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

метода и его точностными характеристиками. Даны полезные практиче-

ские рекомендации и иллюстративный пример.

Ключевые слова: излучающая цель, периодический сигнал, квазиперио-

дический сигнал, однопозиционный пассивный наблюдатель, беспеленго-

вый метод, временная невязка, периодо-временной метод, полиномиаль-

ное движение, параметрическая идентификация, наблюдаемость метода,

корреляционная матрица ошибок оценивания, адаптация.

DOI: 10.31857/S0005231023070048, EDN: FCTVPA

1. Введение

Методы пассивной локации и навигации излучающей цели на базе одно-

позиционного пассивного наблюдателя широко отражены в известной лите-

ратуре (см., например, [1-20]). Среди них весьма популярны доплеровско-

временные беспеленговые методы, оперирующие с периодическими сигнала-

ми и ориентированные на возможность измерения непрерывного смещения

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

ленного движением цели (для задач локации) или движением наблюдате-

ля (для задач навигации); в [6] на стр. 169-173 дан исчерпывающий список

литературы по данному вопросу, доступной в открытой печати. При этом

измерения могут осуществляться на любой характерной частоте из спектра

излучаемого периодического сигнала (например, на центральной) либо моду-

лирующей функции, а также путем сопоставления моментов прихода фрон-

тов последовательных импульсов с учетом известного периода. Указанные

методы основаны на идее «синтеза базы», что приводит в конечном итоге

к формированию нескольких точек наблюдения на траектории движения и

возможности использования известных методов многопозиционной локации

и навигации (например, триангуляционного, разностно-дальномерного, три-

латерационного и их комбинаций [21, 22]). При этом, как правило, рассмат-

риваются такие траектории, которые на участке наблюдения либо извест-

ны (например, орбитальные с известными параметрами движения), либо с

достаточной для практики точностью аппроксимируются моделью прямоли-

нейного равномерного движения (как с известными, так и неизвестными па-

раметрами движения). При этом принципиальным моментом является учет

априорной информации о величине скорости цели или наблюдателя, что для

практики зачастую является неприемлемым.

В [20] развивается периодо-временной метод (ПВМ), который снимает

ограничение, связанное с получением априорной информации о величине

скорости, а также рассматривается вопрос параметрической идентификации

применительно к модели криволинейного движения, учитывающего возмож-

ный маневр цели или наблюдателя. При этом не требуется предварительная

текущая оценка доплеровской частоты, что эквивалентно нахождению произ-

водной от временной невязки между периодами излучаемого (периодическо-

го) и принимаемого (квазипериодического) сигналов. Однако результаты, по-

лученные в [20], распространяются лишь на радиосигналы (распространяю-

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

риодом, а также не исследована зависимость возникающей временной невяз-

ки от параметров движения цели. Настоящая статья является дальнейшим

развитием известного ПВМ в части устранения указанных недостатков при-

менительно к сигналам, распространяющимся в произвольных физических

средах. Не снижая общности рассуждений, а также с целью компактного опи-

сания предлагаемого метода ограничимся решением задачи пассивной лока-

ции излучающей цели (ИЦ) на базе стационарного однопозиционного пассив-

ного наблюдателя (СОПН). Применение полученных результатов к решению

задачи навигации не вызовет принципиальных затруднений.

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

Пусть движущаяся ИЦ формирует в текущем времени t периодический

сигнал S₀(t) (период сигнала T_S = const может быть неизвестным), распро-

страняющийся в заданной физической среде в виде волны со скоростью v_S

(речь может идти о различных волнах, например электромагнитных или аку-

стических). В точке наблюдения, связанной с СОПН, на интервале наблюде-

ния [0, T ] принимается квазидетерминированный сигнал S(t) с переменным

периодом.

Согласно ПВМ отрезок наблюдения представляется в виде

⋃

(2.1)

[0, T ] =

[tn-1, t_n] , t_n > tn-1, t₀ = 0, t_N

≤T,

n=1

где t₀ = 0 — фиксированный момент времени, соответствующий началу при-

нимаемого сигнала (например, поступлению первого импульса), t_n — фикси-

∑_n

рованный момент времени поступления M_n =

ΔM_p периодов принимае-

p=1

мого квазипериодического сигнала (ΔM_p - количество периодов, подсчиты-

ваемых на отрезке [tp-1, t_p]), при этом в момент времени t_n число M_n периодов

целиком укладывается в отрезок [0, t_n].

Теоретические и практические вопросы, связанные с подсчетом указанных

периодов, решаются с использованием электронных цифровых частотомеров

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

[23, стр. 148-161]).

В точке наблюдения (где расположен СОПН) с учетом движения ИЦ сиг-

нал становится квазипериодическим, поскольку возникает временная невяз-

ка δ(t) между периодами излученного и принятого сигналов

(2.2)

δ(t) = v-1SΔ_R(t) = v-1S [R(t) - R₀

t ∈ [0,T],

где R(t) — текущая дальность до ИЦ, R₀ = R(0) — начальная дальность.

В прямоугольной декартовой системе координат XY Z (в центре которой

находится СОПН) движение ИЦ описывается полиномиальной моделью (для

упрощения выкладок и наглядности метода вместо обобщенного конечного

полинома с произвольными базисными функциями ограничимся степенным

полиномом второй степени с начальным условием r₀ = r (0), ∥r₀∥ = R₀)

(2.3)

r(t) = r₀ + v₀t + 2-1a₀t²

t ∈ [0,T],

где r(t) = r = [x, y, z]^T — вектор положения (∥r(t)∥ = R(t)),

v₀ = [vx0,vy0,vz0]^T — вектор начальной скорости (v₀ = ∥v₀∥ — величина

скорости),

a₀ = [ax0,ay0,az0]^T — вектор ускорения (a₀ = ∥a₀∥ — величина ускоре-

ния), при этом векторы r₀, v₀ и a₀ априорно неизвестны.

Если в качестве измеряемого параметра принять величину t_n, то можно

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

(2.4)

h=t+ξ=^t

+ δ + ξ,

где

[

]_T

[

]_T

h_n,n = 1,N

t_n,n = 1,N

t₌

^[t_n,n = 1,N]T ,

[

]_T

[

]_T

δ=

δ_n,n = 1,N

ξ=

ξ_n,n = 1,N

h_n = h(t_n), ξ_n = ξ(t_n) .

[

]_T

В (2.4) под ξ =

ξ_n,n = 1,N

понимается гауссовская погрешность из-

мерений с нулевым математическим ожиданием и корреляционной матри-

цей K_ξ, а измеряемый параметр t_n связан с числом подсчитываемых периодов

соотношением

(2.5)

t_n = M_nT_S + δ_n = t_n + δ_n = t_n + v-1S [R_n - R₀

где δ_n = δ(t_n) — неизвестная временная невязка, t_n = M_nT_S , R_n = R (t_n),

t₀ = 0.

Формулу (2.5) можно прокомментировать так (см., например, [6, стр. 154]):

за время t_n = M_nT_S дальность между ИЦ и СОПН изменится на величину

ΔR_n = R_n - R₀, что соответствует временной невязке δ_n = v-1SΔR_n между

периодами излученного и принятого сигналов. Если бы цель была стацио-

нарной или двигалась по окружности (в центре которой находится СОПН),

то приращение дальности отсутствовало бы и δ_n = 0 для всех n. Именно про-

хождение волной дополнительного участка пути длиной ΔR_n со скоростью

v_S является причиной возникновения невязки δ_n.

Напомним, что при известном периоде T_S в качестве измеряемого парамет-

ра можно было принять величину δ_n = t_n - M_nT_S (именно так формирова-

лось уравнение наблюдения в [6, 20]), при неизвестном периоде T_S доступны

для измерения только величины t_n и M_n.

Если расстояние между ИЦ и СОПН уменьшается, то δ_n < 0, в против-

ном случае δ_n > 0. Появление временной невязки δ_n = δ(t_n) обусловлено эф-

фектом сжатия или растяжения исходного периодического сигнала в точке

наблюдения за счет движения ИЦ.

Требуется с учетом (2.1)-(2.5) разработать метод параметрической иден-

тификации ИЦ с криволинейным (полиномиальным) движением на базе

периодо-временного СОПН, не требующий знания периода T_S излучаемого

сигнала и вычисления текущей частоты Доплера. Метод должен включать в

себя решение следующих вопросов:

— получение зависимостей, позволяющих оценить характер эволюции пе-

риода принимаемого сигнала (обусловленную движением ИЦ), что является

принципиальным для данного метода;

— формирование алгоритма идентификации наклонной дальности и ряда

характерных параметров движения ИЦ на точных данных (принимая ξ_n = 0,

n = 1,N);

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

(т.е. определение условий наблюдаемости метода);

— учет случайных погрешностей измерений;

— решение задачи идентификации на избыточных данных (h) с учетом

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

ратов (МНК)) и получение соотношения для расчета корреляционной матри-

цы ошибок идентификации;

— проведение вычислительного эксперимента с целью демонстрации воз-

можностей метода.

3. Исследование эволюции периода сигнала

Набег δ(t) описывается выражением (для случая прямолинейного равно-

мерного движения)

}

{[

]1/2

(3.1)

δ(t) = v-1S R²⁰ + 2tR₀v₀ cos γ₀ + t²v²⁰

-R₀

t ≥ 0, δ(0) = 0,

где γ₀ — угол между векторами r₀ и v₀.

При 0 < γ₀ ≤ π/2 функция δ(t) является неотрицательной, гладкой и стро-

го выпуклой, δ⁽¹⁾(t) = dδ(t)/dt = 0 в точке t = 0. При π/2 < γ₀ < π функ-

ция δ(t) является гладкой и строго выпуклой, имеет два корня (t = 0 и

t = -2R0 cosγ₀/v₀), в точке t = -R0 cosγ₀/v₀ достигает минимального значе-

ния (v-1SR₀ (sin γ₀ - 1)). При γ₀ = 0 имеем δ(t) = (v₀/v_S ) t, т.е. набег является

линейной неотрицательной функцией, не зависящей от R₀. При γ₀ = π имеем

δ(t) = - (v₀/v_S ) t для 0 ≤ t ≤ R₀/v₀, т.е. δ(t) является линейной функцией и

достигает своего минимума (-R₀/v_S ) в точке t = R₀/v₀. Поскольку при γ₀ = 0

и γ₀ = π набег δ(t) не зависит от R₀, то для этих некорректных случаев, свя-

занных с движением ИЦ вдоль линии визирования, определение дальности с

учетом эволюции периода сигнала в точке наблюдения невозможно.

Для более детального исследования δ(t) найдем несколько первых произ-

водных по времени (в точке t = 0):

⎧

⎪

δ⁽¹⁾⁰ = v-1Sv_R,

⎨

(2)

(3.2)

= (v_S R₀)-1 v2τ,

⎪

⎩

(

)

δ⁽³⁾⁰ = -3

v_SR²⁰

−1 v2τv_R,

где v_R = R₀ cos γ₀ и v_τ = v₀ sin γ₀ — соответственно величины радиальной и

тангенциальной скорости.

В итоге можно воспользоваться разложением на основе ряда Тейлора

(

)

v2τt

v2τv_Rt²

δ(t) = v-1St v_R +

+...

2R₀

2R²

(3.3)

[

(

)

]

v2τt

v_Rt

=v-1St v_R +

+... ,

2R₀

R₀

из которого следует, что спектральный состав функции δ(t) существенно за-

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

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

бенно для больших интервалов наблюдения и малых дальностей.

Формулы (3.1)-(3.3) весьма полезны при обосновании возможности прак-

тической реализации развиваемого ПВМ в каждом конкретном случае с уче-

том принятых исходных данных.

4. Построение алгоритма параметрической идентификации

на точных данных

С учетом (2.3) можно воспользоваться следующей зависимостью

(

)

(4.1)

R²(t) - R²⁰ = 2t〈r₀,v₀〉 + t²

v²⁰ + 〈r₀,a₀〉

+ t3 〈v₀,a₀〉 + 4-1t⁴a²⁰,

где 〈·, ·〉 — символ скалярного произведения двух векторов, ∥·∥ — символ

нормы вектора.

Формула (4.1) представляет собой первое базовое соотношение развивае-

мого ПВМ.

Второе базовое соотношение следует непосредственно из формулы (2.2):

(4.2)

R²(t) - R²⁰ = 2v_SR₀δ(t) + v2Sδ²

(t).

Приравнивая выражения (4.1) и (4.2), после несложных преобразований

получаем уравнение

(4.3)

-2v_S δ(t)χ₁ + 2tχ₂ + t²χ₃ + t³χ₄ + 4-1t⁴χ₅ = v2S δ²

(t),

где

⎧

⎪

χ₁ = R₀,

⎪

χ₂ = 〈r₀,v₀〉 ,

⎨

(

)

(4.4)

χ₃ =

v²⁰ + 〈r₀,a₀〉

⎪

χ₄ = 〈v₀,a₀〉,

⎪

⎩ χ₅ =a²⁰

— неизвестные коэффициенты, имеющие понятный физический смысл и под-

лежащие идентификации.

Поскольку величины δ_n неизвестны, то с учетом (2.5) для дискретного вре-

мени запишем уравнение относительно неизвестных величин T_S и χ_i, i = 1, 5:

-2v_S (t_n - M_nT_S ) χ₁ + 2t_nχ₂ + t2nχ₃ + t3nχ₄ + 4-1t4nχ₅ =

(4.5)

= v2S [(t_n - M_nT_S)]2 .

После несложных, но громоздких преобразований формулу (4.5) можно

представить в виде нового уравнения (относительно коэффициентов A_i)

∑

(4.6)

B_inA_i = D_n,

i=1

где

⎧

⎪

A₁ = (v_Sχ₁ - χ₂)v-2ST-1S, A₂ = -χ₁v-1S,

⎪

(

)(

)-1

⎪

A₃ =

v2S - χ₃

2v2S T_S

A₄ = 2-1T_S,

⎨

(

A₅ = -χ₄

2T_S v2S

)-1 , A₆ = χ₅ (8T_S v2S)-1 ,

(4.7)

⎪

B1n = t_n, B2n = M_n, B3n = t2n,

⎪

B4n = M2n, B5n = t3n, B6n = t4n,

⎩

D_n = M_nt_n.

Соотношения (4.6) и (4.7) являются основой для идентификации пара-

метров криволинейного движения ИЦ при неизвестном периоде излучаемого

сигнала. В (4.6) неизвестными являются коэффициенты A_i, i = 1, 6, кото-

рые однозначно связаны с искомыми параметрами движения ИЦ и перио-

дом излучаемого сигнала. Если уравнение (4.6) записать для всех n = 1, N,

где N ≥ 6, то получим систему линейных алгебраических уравнений (СЛАУ)

(с прямоугольной матрицей B)

(4.8)

BA = D,

[

]

[

]_T

[

]_T

где B =

B_in,n = 1,N,i = 1,6

, A=

a_i,i = 1,6

,D=

D_n,n = 1,N

Данная СЛАУ позволяет решать задачу оценивания указанных коэффи-

циентов и параметров, а также периода сигнала при избыточных измерениях.

При N > 6 речь идет о задаче сглаживания на основе МНК с использованием

ортогонально-сингулярного разложения [24].

Рассмотрим частный случай, когда ИЦ движется прямолинейно и равно-

мерно, а период сигнала неизвестен. Теперь вместо (4.3) имеем уравнение

(4.9)

-2v_S δ(t)χ₁ + 2tχ₂ + t²χ₃ = v2Sδ²

(t),

где

⎧

⎨ χ1 = R0,

(4.10)

χ₂ = 〈r₀,v₀〉,

⎩

χ₃ = v²⁰.

При этом вместо (4.6) в этом случае имеем

(4.11)

t_nA₁ + M_nA₂ + t2nA₃ + M2nA₄ = M_nt_n.

Если предположить, что период сигнала известен, т.е. известны величи-

ны δ_n, то с учетом (4.9) для нахождения параметров прямолинейного равно-

мерного движения ИЦ достаточно решить СЛАУ (относительно χ_i, i = 1, 3)

(4.12)

-2v_S δ_nχ₁ + 2t_nχ₂ + t²χ₃ = v2S δ2n,

n = 1,N.

При этом находим дальность R₀, величину скорости v₀ = ∥v₀∥ и угол γ₀ меж-

ду векторами r₀ и v₀ с учетом очевидных соотношений:

⎧

⎪

R₀ = χ₁,

⎨

v₀ =

^√χ₃,

(4.13)

⎪

[

]

^⎩ γ₀ = arccos χ2 (R₀v₀)-1

В случае прямолинейного равноускоренного движения ИЦ (когда векторы

v₀ и a₀ являются коллинеарными) необходимо решить СЛАУ (относитель-

но χ_i, i = 1, 5)

(4.14)

-2v_S δ_nχ₁ + 2t_nχ₂ + t2nχ₃ + t3nχ₄ + 4-1t4nχ₅ = v2S δ2n.

Теперь имеем:

⎧

⎪

χ₁ = R₀,

⎪

χ₂ = R₀v₀ cos γ₀,

⎨

(

)

(4.15)

χ₃ =

v²⁰ + R₀a₀ cos γ₀

⎪

χ₄ = v₀a₀,

⎪

⎩ χ₅ =a²⁰.

По найденным значениям χ₁, . . . , χ₅ вычисляем следующие параметры

движения ИЦ:

⎧

⎪

R₀ = χ₁,

⎪

⎨

a₀ =

^√χ₅,

(4.16)

⎪

v₀ = χ₄a-10,

⎪

[

]

^⎩ γ₀ = arccos χ2 (R₀v₀)-1

Выражения (4.1)-(4.16) составляют математическую основу развиваемого

ПВМ.

В следующем разделе проанализируем условия наблюдаемости развивае-

мого метода, т.е. выявим ситуации, при которых он становится некорректным

с вычислительной точки зрения.

5. Анализ наблюдаемости метода

Развиваемый ПВМ можно реализовать на любом наборе узлов из мно-

жества {t₁, . . . , t_N }, что позволяет не только уменьшить объем вычислений,

но и в ряде случаев повысить надежность формируемых оценок (особенно

при отсутствии достоверной априорной информации о весовых коэффициен-

тах, необходимых для реализации МНК). Для этого введем векторы времен-

[

]_T

ных узлов t[l] =

t[l]p,p = 1,P[l]

, где l = 1, L, t[l]p ∈ {t₁, . . . , t_N }, t[l]p+1 > t[l]p.

Здесь L — число наборов, P[l] — количество узлов в l-м наборе, t[l]p - узел с

номером [l] p (это натуральное число, принадлежащее множеству {1, . . . , N}).

На основе (4.12) сформируем следующую СЛАУ:

(5.1)

C[l]χ[l] = Y[l],

[

]_T

[

где Y[l] = δ2[l]p, p = 1, P[l]

, χ[l] =

χi[l],i = 1,5]T , а матрица C[l] (разме-

(

)

ром P[l] × 5) образована строками v-2S -2vS δ[l]p, 2t[l]p, t2[l]p, t3[l]p, 4-1t4

,p=

[l]p

= 1, P[l].

Введение t[l] позволяет с учетом геометрии наблюдения, характеристик

ИЦ и СОПН находить такие наборы узлов, на которых вопрос идентифика-

ции решается наиболее качественно (это относится к известной задаче пла-

нирования эксперимента [25]).

Не снижая общности рассуждений, ограничимся плоским случаем (при-

нимая z = 0) и сигналом с известным периодом, а также зададимся P[l] = 5,

что соответствует квадратной матрице C[l]. Очевидно, что для корректного

применения развиваемого метода, связанного с решением СЛАУ (5.1), необхо-

димо и достаточно выполнение условия det C[l] = 0, что приводит к искомому

результату χ[l] = C-1[l]Y[l]. Для выявления случаев, при которых это условие

нарушается, запишем столбцы матрицы C[l] в виде векторов:

[

]_T

[

]_T

C[l]1 =

-2v_S δ[l]p, p = 1, 5

C[l]2 =

2t[l]p, p = 1, 5

[

]_T

[

]_T

[

]_T

C[l]3 = t2[l]p,p = 1,5

C[l]4 = t3[l]p,p = 1,5

C[l]5 =

4-¹t4[l]p, p = 1, 5

Несложно заметить, что столбцы C[l]2, C[l]3 и C[l]4 линейно независимы, сле-

довательно, для проверки условия det C[l] = 0 достаточно показать, что стол-

бец C[l]1 не может быть представлен в виде линейной комбинации этих столб-

цов.

[

]-2 (

(

)

Так как R[l]p = x2[l]p + y2[l]p

где R[l]p = R

t[l]p

, x2[l]p = x₀ +vx0t[l]p+

)₂

(

)₂ )

+2-1ax0t2[l]p

и y2[l]p = y₀ +vy0t[l]p +2-1ay0t2[l]p

, то нарушение усло-

[

]_T

вия det C[l] = 0 эквивалентно тому, что векторы μ[l] = x2[l]p, p = 1, 5

[

]_T

η[l] = y2[l]p,p = 1,5

не связаны условием коллинеарности: μ[l] = kη[l], где

k — коэффициент пропорциональности. В противном случае имеем

[

]-2

[

]-2



_y

(5.2)

R[l]p = x2[l]p + y2[l]p

= k²y2[l]p + y2[l]p

[l]p

,

[

]

[



]

(5.3)

−2v_S δ[l]p = -2

R[l]p - R₀

= -2

^y[l]p

- R₀

(

)-2

где q =

k² + 1

Из (5.2) и (5.3) следует, что координаты вектора C[l]1 можно представить

линейной комбинацией из координат векторов C[l]2, C[l]3 и C[l]4. Физический

смысл условия μ[l] = kη[l] (условие вычислительной некорректности метода)

состоит в том, что ИЦ движется прямолинейно вдоль линии визирования

СОПН.

Таким образом, для корректности метода необходимо исключить случаи,

когда ИЦ движется вдоль указанной линии или в ее окрестности. Это накла-

дывает определенные ограничения на условия наблюдения ИЦ, что необхо-

димо предусмотреть на практике.

Если ограничиться моделью(прямолинейного равномерного движения и

сигналом с известным периодом в этом случае в (5.1) надо положить p = 1, 3

)

[

иt[l] =

t[l]1,t[l]2,t[l]3

, то решение СЛАУ (5.1) при корректном применении

метода позволяет определить искомые параметры движения ИЦ

⎧

(

)

δ2[l]1Δt[l]23 - δ2[l]2Δt[l]13 + δ2[l]3Δt[l]12

⎪

R0[l] = 2-1v_S

⎪

-δ[l]1Δt[l]23 + δ[l]2Δt[l]13 - δ[l]3Δ^t

⎪

[l]12

⎪

(

)

⎪

t2[l]1Δδ[l]23 - t2[l]2Δδ[l]13 + t2[l]3Δδ[l]12

⎪

〈r₀,v₀〉[l] = 2-1v²

⎨

-δ[l]1Δt[l]23 + δ[l]2Δt[l]13 - δ[l]3Δ^t

[l]12

(5.4)

[

]1/2

⎪

t[l]3Δδ[l]12 - t[l]2Δδ[l]13 + t[l]1Δδ[l]23

⎪

v0[l] =

⎪

δ[l]1Δt[l]23 - δ[l]2Δt[l]13 + δ[l]3Δ^t

⎪

[l]12

⎪

[

]

⎪

〈r₀,v₀〉[l]

⎩ γ0[l] = arccos

R0[l]v0[l]

(

)

(

)

где Δt[l]12 = t[l]1t[l]2

t[l]1 - t[l]2

, Δδ[l]12 =δ[l]1δ[l]2

δ[l]1 - δ[l]2

и, если не учиты-

вать ошибки измерений и вычислений, R0[l] = R₀, v0[l] = v₀, 〈r₀, v₀〉[l] = 〈r₀, v₀〉,

γ0[l] = γ₀.

Следователь

(

[

])

ния R₀, v₀ и γ₀ где R₀ = χ₁, v₀ =

^√χ₃, γ₀ = arccos χ₂ (R₀v₀)-1

, не прибе-

гая к численному решению СЛАУ, что является несомненным достоинством

развиваемого ПВМ.

6. Учет случайных погрешностей измерений

Полагая период сигнала известным, для оценки влияния случайных по-

грешностей измерений ξ на точностные характеристики метода восполь-

зуемся традиционной процедурой расчета элементов корреляционной мат-

рицы Kχ[l] ошибок оценивания координат вектора χ в линейном при-

ближении [26]. Для этого с учетом СЛАУ (5.1) (полагая для простоты

матрицу C[l] квадратной размером 5 × 5) воспользуемся представлением

[

(

)

]_T

[

]_T

χ[l] = C-1[l]Y[l] =

χ_k

δ[l]

,k = 1,5

(где δ[l] =

δ[l]p,p = 1,5

) и частными

(

)

производными следующего вида: ∂χk[l]

δ[l]

/∂δ[l]p. Корреляционная матри-

ца находится по правилу

(6.1)

Kχ[l] = Fχ[l]K_ξFTχ[l],

[

(

)

]

где Fχ[l] =

∂χk[l]

δ[l]

/∂δ[l]p, k = 1, 5, p = 1, 5

Выражение (6.1) позволяет априорно на математических ожиданиях изме-

ряемых параметров оценить потенциальные возможности развиваемого ПВМ

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

при конкретных условиях наблюдения ИЦ, а также обоснованно подходить

к выбору основных параметров метода (длины интервала наблюдения (T ),

количества узлов (N) и временных наборов (t[l]) ). Так, номер l^∗ ∈ {1, . . . , L}

оптимального набора δ[l∗_], обеспечивающего минимизацию ошибки оценива-

ния, находится по следующему адаптивному правилу:



(6.2)

l^∗ = arg min

Kχ[l],



где

Kχ[l]

— любая из норм матрицы Kχ[l], применяемая в задачах оценива-

ния.

При практической реализации развиваемого ПВМ следует учитывать тот

фактор, что при больших значениях v_S (например, когда v_S = c, где c —

скорость света) решение квадратной СЛАУ (4.8) при наличии случайных по-

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

ним этот факт для случая N = 6 на примере вычисления скорости v₀. По-

скольку v₀ = c√1 - 2TSA3, то ошибка Δ3 =Aˆ3 - A3 (где Aˆ3 - вычисленное

значение коэффициента A₃ путем решения СЛАУ (4.8) с учетом ошибок из-

√

мерений) приводит к следующей оценке скорости: v₀ =

v²⁰ + 2c²T_SΔ₃. Тo

есть корректная оценка скорости возможна лишь при выполнении условия

(

)-1

Δ₃ > -v²⁰

2c²T_S

, что накладывает очень жесткое ограничение на величи-

ну погрешности Δ₃. Данный эффект относится также ко всем коэффициен-

там СЛАУ (4.8), кроме A₂ иA₄.

Для преодоления указанной некорректности (при больших скоростях v_S)

рекомендуется двухэтапный подход к идентификации. На первом этапе ре-

шается СЛАУ (4.8), из которой потребуются лишь оценка

A₄ для A₄. Это

позволяет сформировать искомую оценк

T_S =

A₄ для периода T_S, а на ее

основе и оценки для невязокδn = t_n - M_n

T_S. Все оценки параметров движе-

ния ИЦ находятся на базе СЛАУ (5.1), в которую вместо δ_n подставляется

величинаδn.

7. Учет избыточных измерений

Теперь рассмотрим случай избыточных измерений, когда матрица C[l]

и вектор Y[l] в (5.1) имеют произвольное число строк P[l] ≤ N, которое,

как правило, значительно превышает количество оцениваемых параметров.

Для упрощения расчетов будем рассматривать в СЛАУ (5.1) составляющую

[

]_T

Y[l] =

v2Sδ2[l]p,p = 1,P[l]

в качестве вектора вторичных измеряемых пара-

метров, а первичные измерения h[l]1, . . . , h

полагаем некоррелированны-

[l]P[l]

ми. С учетом этого по аналогии с (5.1) корреляционную матрицу ошибок

измерений координат вектора Y[l] можно представить так

(7.1)

KY[l] = Fδ[l]K_ξFTδ[l].

В предположении, что матрица K_ξ является диагональной, имеем_[

]

KY[l] = diag 4δ2[l]1, 4δ2[l]2, . . . , 4δ2

. При условии достаточно малых ошибок

[l]P[l]

измерений для построения сглаженной оценки χ∗[l] вектора χ можно восполь-

зоваться методом наименьших квадратов [25]

(

)-1

(7.2)

χ∗[l] = CT[l]K-1Y[l]C[l]

CT[l]K-1Y[l]hY[l],

[

]_T

где hY[l] =

hY[l]p,p = 1,P[l]

— вектор вторичных измерений.

Корреляционную матрицу ошибок оценивания находим так:

(

)-1

(7.3)

K_χ∗

= CT[l]K-1Y[l]C[l]

[l]

Для выбора оптимального набора с номером l^∗ ∈ {1, . . . , L} применяем

адаптивный алгоритм типа (6.2).

Следует отметить, что подход (7.1)-(7.3) не является строго оптимальным,

поскольку элементы матрицы C[l] зависят от результатов наблюдений. Но

при определенных ограничениях на погрешности измерений он дает вполне

приемлемый результат.

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

ры нелинейного оптимального оценивания, которые приводят на практике к

трудоемким рекуррентным вычислительным алгоритмам, предполагающим

задание достаточно качественного начального условия.

Другой наиболее простой и достаточно надежный способ построения сгла-

женной оценки χ∗[l] состоит в предварительном сглаживании первичных изме-

рений h[l]1, . . . , h[l]P[l] соответствующим полиномом δ∗[l](t) и применении полу-

ченных результатов к решению СЛАУ (5.1). Кроме того, можно найти сгла-

женную оценку дальности для любого t ∈ [0, T ], а именно,

(7.4)

R[l](t) = R∗0[l] + cδ∗[l]

(t).

Здесь в качестве оптимального принимаем набор с номером l = l^∗ ∈

∈ {1,... ,L}.

8. Некоторые обобщения и практические рекомендации

Выше рассматривался случай оценивания начальной дальности R₀ = R (0)

для момента времени t = 0. Однако если ряд Тейлора, используемый для

описания криволинейного движения ИЦ, записать относительно не началь-

ного, а любого произвольного t = t_∗ ∈ [0, T ], то по аналогии с вышеизложен-

ным можно решить задачу идентификации именно для момента времени t^∗,

в частности найти дальность R_∗ = R (t_∗).

Развитый метод несложно реализовать в виде следующих алгоритмов: по

выборке нарастающего объема, на «скользящей сетке» или в виде филь-

тра [25]. При этом движение ИЦ на интервале наблюдения можно рассматри-

вать как кусочно-полиномиальное (в [20] оно рассматривалось как кусочно-

линейное).

При практической реализации метода возникают вопросы (например, вы-

бор степени полинома, описывающего движение ИЦ или количества подсчи-

тываемых импульсов), связанные с организацией измерительного экспери-

мента. В [25] даются практические рекомендации для решения этих вопросов

в полном объеме. Очевидно, что развиваемый метод наиболее эффективен в

случае, когда речь идет о больших пройденных расстояниях (т.е. «синтези-

руется база» достаточного размера), а это задает определенные ограничения

на тип ИЦ (в частности на его скорость, на возможности маневра и т.д.), на

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

и на технические характеристики СОПН.

Для случаев, связанных с движением ИЦ по линии визирования, можно

предложить гибридный вариант использования развитого и известного энер-

гетического метода [27]. Доказано, что данный метод, оперирующий с отно-

сительным уровнем принимаемого сигнала, реализует при движении ИЦ по

линии визирования свои потенциальные возможности. В некотором смысле

развитый и энергетический методы «ортогональны» друг другу в плане точ-

ности. Следовательно, комбинируя эти методы, можно выровнять рабочую

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

для различных условий наблюдения ИЦ.

Для более эффективного применения энергетического метода следует ис-

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

ции и отсева ненадежных измерений.

9. Иллюстративный пример

Предположим, что ИЦ осуществляет плоскостное движение: x(t) = x₀ +

+vx0t, y(t) = y₀ + vy0t, где x₀ = y₀ = 11 · 10³, vx0 = -5 · 10², vy0 = 6 · 10²,

γ₀ = 85. Здесь и далее время и погрешности измерений временных интер-

валов задаются в секундах (с), координаты и дальность — в метрах (м), ско-

рость — в м/с, ускорение — в м/с², частота — в герцах (Гц), угол — в градусах,

относительная погрешность — в процентах.

ИЦ формирует импульсный радиосигнал

∑

[

]

S₀(t) =

rect

(t - kT_S)τ-1

cos (2πf₀t),

k=1

0,208

0,206

0,204

0,202

0,200

0,198

0,196

0,194

0,192

Время, с

Относительная погрешность оценивания дальности.

где T_S = 10-2, τ = 10-5, f₀ = 10¹⁰. Параметры работы СОПН: T = 18,

v_S = c = 3 · 10⁸, L = 1 (т.е. использован один единственный набор узлов),

[

]

P_[1] = 4 (размер набора), ΔM_p = ΔM = 10, K_ξ = diag

σ²,... ,σ²

, при этом

погрешности измерений временного положения фронтов импульсов полага-

лись некоррелированными и задавались по нормальному закону распределе-

ния с нулевым математическим ожиданием и значением среднеквадратиче-

ского отклонения σ = 10-9.

Реализация метода осуществлялась в два этапа с применением датчи-

ка случайных чисел и усреднением по тысяче экспериментов. На первом

этапе решалась СЛАУ (4.8) с квадратной матрицей B размером 4 × 4 (по-

лагалось A₅ = A₆ = 0, так как рассматривается ИЦ с нулевым ускорени-

ем), при этом для расчета элементов матрицы B и столбца D исполь-

[

]_T

зован вектор t_[1] = t_[1] + δ_[1] =

t[1]p,p = 1,4

с номерами узлов: [1]1 = 12,

]_T

[

]_T

[1]2 = 65, [1]3 = 118, [1]4 = 171, t_[1] =

[t_[

[1]p · 10-1, p = 1, 4

1]p, p = 1, 4

Из всех четырех оценок неизвестных коэффициентов выбирается только

оценка периода сигнала

TS[1] = 9,999999731646 · 10-3 (получена на основе

набора t_[1] =t[1] + δ_[1]), которой соответствует относительная погрешность

δTS[1] = 2,683540941544882 · 10-6.

На втором этапе с учетомδ[1]p = t[1]p - M[1]p

TS[1] = t[1]p - [1] pΔ

TS[1] ре-

шалась СЛАУ (5.1) с квадратной матрицей C[l] (размером 3 × 3), при этом

χ₄ = χ₅ = 0 и использовался наборt[1] =

^[t

]^T = [1,1; 9,1; 17,1]^T.

[1]p, p = 1, 3

(

)

Сама матрица образована строками c-2

-2cδ[1]p, 2t[1]p, t²

, p = 1,3. В ито-

[1]p

ге истинной дальности R₀ = 1,555634918 · 10⁴ соответствуют оценка

0[1] =

= 1,559672203 · 10⁴ и погрешность ΔR_0[1] = 0,259526489, истинной скорос-

ти v₀ = 7,810249675 · 102 — оценка v0[1] = 7,821417156 · 102 и погрешность

Δv_0[1] = 0,142984942, истинному углу γ₀ = 84,805571092 — оценка γ0[1] =

= 84,761511501 и погрешность Δγ_0[1] = 0,051953650.

На рисунке представлен график зависимости относительной погрешности

оценивания дальности, полученный с учетом (7.4).

Для более эффективного использования разработанного в статье метода

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

а также их согласования с динамикой движения ИЦ и величиной погрешно-

стей измерений следует решать в оптимизационной постановке. При решении

СЛАУ следует привлекать известные методы регуляризации. Результаты чис-

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

используемой временной сетки, тем меньшее влияние оказывают случайные

погрешности измерений на результирующую точность оценивания. Это рас-

стояние должно быть согласовано с динамикой ИЦ, а именно: чем меньше

скорость движения ИЦ, тем больше должны быть шаг этой сетки и длитель-

ность интервала наблюдения.

10. Заключение

Разработанный ПВМ позволяет идентифицировать модель криволиней-

ного полиномиального движения ИЦ по результатам регистрации времен-

ной невязки между периодами излученного сигнала и этими же периодами,

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

нала и предварительной оценки текущей частоты Доплера, а также знания

каких-либо априорных данных о параметрах принятой модели движения ИЦ.

Исследованы наблюдаемость и основные ограничения метода, условия его

наиболее эффективного применения. Получены аналитические соотношения,

позволяющие оценить эволюцию временной невязки с учетом характеристик

ИЦ и СОПН, а также точностные характеристики метода для различных

условий наблюдения.

Метод можно реализовать в различных вариантах: по фиксированной вы-

борке измерений, по выборке измерений нарастающего объема, в виде алго-

ритма динамической фильтрации (линейной, квазилинейной или нелинейной)

и др.

Метод может быть реализован как самостоятельно, так и в составе гибрид-

ного метода, объединяющего другие известные подходы пассивной однопози-

ционной и многопозиционной локации и навигации ИЦ. Поскольку разви-

тый метод позволяет определять дальность, то он может быть использован в

дальномерно-дальномерных системах многопозиционной локации при реше-

нии известной задачи трилатерации [21, 22].

При наличии в периодо-временных измерениях не только флуктуацион-

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

подвергнуть эти измерения процедуре обобщенного инвариантно-несмещен-

ного оценивания [28], обеспечивающей компенсацию этих ошибок, достиже-

ние эффекта сглаживания и оптимальное оценивание различных числовых

характеристик (линейных функционалов, например производных, интегра-

лов, спектральных коэффициентов и т.д.), полезных не только для повыше-

ния вычислительной устойчивости метода, но и оценки его эффективности.

Для решения СЛАУ с использованием процедуры регуляризации можно при-

менять известный подход [29].

СПИСОК ЛИТЕРАТУРЫ

Основы маневрирования кораблей / Под ред. М. Скворцова. М.: Воениздат, 1966.

Шебшаевич В.С. Введение в теорию космической навигации. М.: Сов. радио,

1971.

Громов Г.Н. Дифференциально-геометрический метод навигации. М.: Радио и

связь, 1986.

Хвощ В.А. Тактика подводных лодок. М.: Воениздат, 1989.

Соловьев Ю.А. Спутниковая навигация и еe приложения. М.: Экотрендз, 2003.

Мельников Ю.П., Попов С.В. Радиотехническая разведка. М.: Радиотехника,

2008.

Ярлыков М.С. Статистическая теория радионавигации. М.: Радио и связь, 1985.

Сосулин Ю.Г., Костров В.В., Паршин Ю.Н. Оценочно-корреляционная обра-

ботка сигналов и компенсация помех. М.: Радиотехника, 2014.

Булычев Ю.Г., Манин А.П. Математические аспекты определения движения

летательных аппаратов. М.: Машиностроение, 2000.

10.

Булычев Ю.Г., Васильев В.В., Джуган Р.В. и др. Информационно-измеритель-

ное обеспечение натурных испытаний сложных технических комплексов. М.: Ма-

шиностроение - Полет, 2016.

11.

Гельцер А.А. Однопозиционный метод определения местоположения источника

радиоизлучения с использованием отражений сигналов от множества элементов

рельефа и местных предметов // Автореф. дисс. Том. гос. ун-т систем управле-

ния и радиоэлектроники. 2012.

12.

Сиренко И.Л., Донец И.В., Рейзенкинд Я.А. и др. Однопозиционное определе-

ние координат и вектора скорости радиоизлучающих объектов // Радиотехника.

2019. № 10 (16). С. 28-32.

13.

Булычев Ю.Г., Булычев В.Ю., Ивакина С.С., Насенков И.Г. Пассивная локация

группы движущихся целей одним стационарным пеленгатором с учетом апри-

орной информации // АиТ. 2017. № 1. С. 152-166.

Bulychev Yu.G., Bulychev V.Yu., Ivakina S.S., Nasenkov I.G. Passive Location of

a Group of Moving Targets with One Stationary Bearing with Prior Information //

Autom. Remote Control. 2017. V. 78. No. 1. Р. 125-137.

14.

Булычев Ю.Г., Булычев В.Ю., Ивакина С.С., Николас П.И. Оценка наклонной

дальности до цели с полиномиальным законом движения // Вестн. Казан. гос.

ун-та. 2013. № 1. С. 67-74.

15.

Булычев Ю.Г. Некоторые аспекты идентификации динамических объектов при

некорректных условиях наблюдения // АиТ. 2020. № 6. С. 131-152.

Bulychev Yu.G. Some Aspects of Identification of Dynamic Objects under Incorrect

Observation Conditions // Autom. Remote Control. 2020. V. 81. No. 6. Р. 1073-1090.

16.

Дятлов А.П., Дятлов П.А. Доплеровские обнаружители подвижных объектов

с использованием «постороннего» источника излучения // Спец. техника. 2010.

№ 5. С. 16-22.

17.

Aidala V.J., Nardone S.C. Biased Estimation Properties of the Pseudolinear

Tracking Filter // IEEE Transact. Aerospas. Electron. Syst. 1982. V. 18. No. 4.

P. 432-441.

18.

Amelin K.S., Miller A.B. An Algorithm for Refinement of the Position of a Light

UAV on the Basis of Kalman Filtering of Bearing Measurements // J. Commun.

Techn. Electron. 2014. V. 59. No. 6. Р. 622-631.

19.

Miller A.B. Development of the Motion Control on the Basis of Kalman Filtering

of Bearing-Only Measurements // Autom. Remote Control. 2015. V. 76. No. 6.

Р. 1018-1035.

20.

Булычев Ю.Г., Мозоль А.А. Однопозиционная пассивная локация источника

излучения с криволинейным движением и учетом эволюции периода сигнала в

точке приема // РЭ. 2023. Т. 68. № 2. С. 131-137.

21.

Кондратьев В.С., Котов А.Ф., Марков Л.Н. Многопозиционные радиотехниче-

ские системы. М.: Радио и связь, 1986.

22.

Черняк В.С. Многопозиционная радиолокация. М.: Радио и связь, 1993.

23.

Нефедов В.И., Сигов А.С., Битюков В.К., Самохина Е.В. Электрорадиоизме-

рения. М.: Форум: Инфра-М. 2018.

24.

Лоусон Ч., Хенсон Р. Численное решение задач метода наименьших квадратов.

М.: Наука, 1986.

25.

Жданюк Б.Ф. Основы статистической обработки траекторных измерений. М.:

Сов. радио, 1978.

26.

Вентцель Е.С. Теория вероятностей. М.: Высш. шк., 1999.

27.

Булычев Ю.Г., Ивакина С.С., Насенков И.Г. Метод пассивно-энергетической

локации и навигации в стационарной и нестационарной постановках // Радио-

техника. 2015. № 6. С. 107-115.

28.

Булычев Ю.Г., Елисеев А.В. Вычислительная схема инвариантно несмещенно-

го оценивания значений линейных операторов заданного класса // ЖВМиМФ.

2008. Т. 48. № 4. С. 580-592.

29.

Булычев Ю.Г., Бурлай И.В. Метод параметрической идентификации систем

управления при неточном задании входных данных // АиТ. 1997. № 11. С. 56-65.

Статья представлена к публикации членом редколлегии О.Н. Граничиным.

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

После доработки 04.05.2023

Принята к публикации 21.05.2023