Виконано дослідження середньоквадратичних помилок (СКО) спільних оцінок широти, довготи і частоти наземного джерела вимірювання (ІІ). Оцінки формуються на виході квазіоптимального фільтра, вхідними сигналами якого є вимірювання поточної несучої частоти високочастотних сигналів ІІ, прийнятих бортовими приймачами одного (або кількох) КА, що мають високоеліптичного орбіти. Величини СКО квазіоптимальних оцінок порівнюються з СКО оптимальних Байесови оцінок при квадратичної функції втрат.

Анотація наукової статті з фізики, автор наукової роботи - Савін Олександр Олександрович, Тисленко Володимир Ілліч


Область наук:

  • фізика

  • Рік видавництва: 2006


    Журнал: Доповіді Томського державного університету систем управління і радіоелектроніки


    Наукова стаття на тему 'Квазіоптимальний оцінка координат наземного джерела випромінювання в космічній системі з вимірами частоти радіосигналів'

    Текст наукової роботи на тему «Квазіоптимальний оцінка координат наземного джерела випромінювання в космічній системі з вимірами частоти радіосигналів»

    ?УДК 621.396

    А.А. Савін, В.І. Тисленко

    Квазіоптимальний оцінка координат наземного джерела випромінювання в космічній системі з вимірами частоти радіосигналів

    Виконано дослідження середньоквадратичних помилок (СКО) спільних оцінок широти, довготи і частоти наземного джерела випромінювання (ІІ). Оцінки формуються на виході квазіоптимального фільтра, вхідними сигналами якого є вимірювання поточної несучої частоти високочастотних сигналів ІІ, прийнятих бортовими приймачами одного (або кількох) КА, що мають високоеліптичного орбіти. Величини СКО квазіоптимальних оцінок порівнюються з СКО оптимальних Байесови оцінок при квадратичної функції втрат.

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

    Геометрія задачі визначення місця розташування показана на рис. 1, де ОХУ - декартова геоцентрична система координат (ГСК). Її центр розташований в центрі Землі, вісь ОХ знаходиться в площині екватора і спрямована в точку весняного рівнодення; вісь 02 паралельна осі обертання Землі і орієнтована в напрямку півночі; напрямок осі ОУ відповідає правій системі координат [1].

    Мал. 1 - Геоцентрична система координат

    Поточне положення на орбіті? -Го КА (спостерігача) в момент часу t визначено завданням радіус-вектора roi (t), положення ІІ визначає радіус-вектор rg (i). Таким чином, вектора

    Л1; (О = ro? (*) "Г. (*) І Av? (T) = roi {t) ~ r3 {t) = voi (t) -vs (t) (1)

    визначають в момент часу t відповідно поточне взаємне положення і швидкість спрямованої лінії візування «джерело - спостерігач» в ГСК.

    Нехай ІІ має частоту радіосигналу f3 і розташований в точці на поверхні Землі, широта і довгота якої дорівнюють відповідно <ра і Xs. У грінвічській системі координат [1], що обертається разом з Землею, широту і довготу точки визначає постійний в часі вектор грінвічського координат джерела

    Tsg = [^ Е 'cos Ф8 "COSК&Е | C0SФа 's * Nк ^ Е "S" 1 Ф8Г' (2)

    де Re = 6378136 м - радіус Землі. Таким чином, в грінвічській системі координат вектор швидкості ІІ rsg = v $ g = 0. Вектора положення і швидкості в ГСК і грінвічській системі координат пов'язані відомими [1] співвідношеннями

    (О = [у. (О (0 z * (0] Г = і (0 • Гв,; Ve (t) = Q? • [- ", (t) * 8 (t) 0] r, (3 )

    де матриця Х1 (<) Враховує добове обертання Землі; &Е = 0,729211-1 (Г4 з-1 - середня

    кутова швидкість обертання Землі.

    Введемо тривимірний вектор змінних стану

    = [* 1 * 2 * з] = [/. До Ф.], (4)

    координати якого постійні в часі і є невідомими випадковими величинами. Очевидно, що для нерухомого ІІ в дискретному часі для вектора х справедлива система різницевих рівнянь

    х (А) = х (А-1), І = 1,2, .... (5)

    Випадкові початкові умови для рівняння (5) визначені завданням апріорної щільності розподілу ймовірностей (ПРВ) Ж (х (0)). Апріорне математичне очікування

    = V,

    о, має діагональ-

    М [х (0)] = ш0. Коваріаційна матриця М (х (0) - т0) (х (0) - т0) Т

    (2 + 2 2 \

    о01 о02 ст031. Далі по тексту М - оператор статистичного усереднення.

    Припустимо, що рух 1-го спостерігача відбувається по кеплерова орбіті. Вектор фазових координат, що визначає поточний стан 1-го спостерігача, визначимо у вигляді

    х<м (г) = [Гог (Про УТОГ • Відомо, що в цьому випадку хи (<) Є рішенням системи нелінійних диференціальних рівнянь. У векторній формі вони мають вигляд [1]

    (0. «Про (*) _

    (6)

    де ас (г) = -ц-Гог (г) / | го; (г) | 3 - гравітаційне прискорення і ц = 0,3986013 Ю15 кг м / с2 - добуток маси Землі на гравітаційну постійну. Початкові умови хо? (0) для

    системи (6) визначають орбіту 1-го спостерігача.

    Похибка вимірювання несучої частоти радіосигналу на борту КА обумовлена ​​впливом кількох чинників. Це, зокрема, середовище поширення радіохвиль (РРВ), власний шум приймальних пристроїв, неточність знання навігаційних параметрів КА, інструментальна і методична помилки, а також помилка, обумовлена ​​нестабільністю частоти излученного сигналу.

    Домінуючий внесок в похибка вимірювання частоти вносять канал РРВ, власний шум приймача і нестабільність частоти ІІ. Тропосфера і іоносфера є джерелами систематичних і випадкових похибок вимірювання частоти. Вплив систематичної помилки може бути значно зменшено за рахунок обліку регулярних (рефракційних) явищ.

    В роботі досліджується похибка визначення координат ІІ, пов'язана з наявністю тільки випадкових помилок. Ця компонента похибки має інтервал кореляції не більше одиниць секунд. При рознесенні пунктів прийому на відстані близько декількох сотень кілометрів і більше випадкові похибки, зумовлені впливом каналу РРВ, некорре-ліровать між собою.

    При темпі надходження даних з виходу цифрового частотоміра порядку секунди можна вважати, що послідовності випадкових помилок вимірювання частоти радіосигналу на борту кількох КА некорреліровани в часі в кожному з пунктів прийому. Взаємна кореляція помилок вимірювання залежить, очевидно, від співвідношення дисперсії ст ^ СА помилки, зумовленої впливом каналу РРВ, і дисперсії про ^ помилки, пов'язаної з нестабільністю частоти ІІ. Очевидно, що ці помилки можна вважати статистично незалежними між собою.

    Викладені міркування і співвідношення (1-5) дозволяють уявити в ГСК сукупність т спостережуваних сигналів у вигляді вектора

    >(*) = [Г1 (А) ... 2г (А) ... 2т (й)] Т = Ь [х; х0 (А), А] + п (й), (7)

    де

    1-4

    з

    г

    | Дгг (*) |

    с - швидкість світла; (•, |) - скалярний добуток двох векторів, зазначених в круглих

    дужках, визначає взаємну радіальну швидкість? -го КА і ІІ; п (й) - / га-мірний вектор випадкових дискретних послідовностей помилок вимірювань частоти.

    Вектор помилок п (&) Має гауссову ПРВ з нульовим середнім значенням і ковариационной

    матрицею В ". Елементи головної діагоналі цієї матриці є А2 {= а2сЛ + О23, де - СКО

    про

    вимірювання несучої частоти на борту КА. Елементи поза головною діагоналі рівні. Взаємна ковариация М (0) - т0) • п (к) Т = 0, оскільки при о01 «/ 5 і ст02 = ст03« | ДГ? | можна знехтувати фізичними передумовами наявності статистичної залежності між координатами ІІ в апріорної зоні і помилками вимірювання частоти.

    Алгоритм оцінки координат

    Завдання Байесови фільтрації станів динамічної системи в постановці марковської теорії нелінійної фільтрації [2, 3] при заданої функції втрат визначена завданням лінійних рівнянь (5) і нелінійних рівнянь (7). Для квадратичної функції втрат оптимальна Байесови оцінка х (&) Поточного стану х (&) Реалізується у вигляді оператора апостеріорного середнього по розподілу ймовірностей Ту (х (а) | 7, к), тобто.

    оо оо

    х (А) = ...] "х (Л) | IV (х (А) | г *) <<х (Л), (8)

    -оо -оо

    де 2! * = | г (1), г (2 послідовність спостережень. В силу нелінійності рівнянь (7) апостериорная ПРВ в даній задачі не є гаусом. Оператор (8) не може бути строго представлений в замкнутій формі, тобто . суворо оптимальна оцінка не може бути реалізована на виході дискретного лінійного фільтра Калмана (ФК).

    У цьому завданню можливе отримання квазіоптимальний поточної оцінки стану. Зокрема, оцінка х (А) може бути реалізована за допомогою алгоритму узагальненого (розширеного) фільтра Калмана (РФК). Відомо [2, 3], що РФК виходить на основі алгоритму лінійного ФК при лінійної апроксимації нелінійних функцій спостережень г (А) = Ь [х; до \ поруч Тейлора в точці відомої поточної оцінки стану. Рівняння оцінки та її ковариационной матриці мають вигляд

    х (*) =? "(*) + К (*) [« (*) - Г (А)]; (9)

    УЖ (Л) = [Е-К (*) НЖ (А)] \ ^ (А), (10)

    де х ~ (к) - екстрапольована оцінка стану; V- (А) - ковариационная матриця помилок екстраполяції; К (к) - матричний коефіцієнт посилення фільтра; Нх (к) = =

    = | Е / ь ( «) /'хл __ > - матриця Якобі розмірністю ТПХ 3. екстраполювати оцінка х (А)

    [ 'Ж-X (*)]

    спостережень в рівнянні (9) має вигляд

    Г (А) = ь [х-(А); А, ж0 (А)]. (І)

    Відзначимо дві обставини, пов'язаних з реалізацією алгоритму (9-11). По-перше, обчислення матриці К (А) і оцінки? ~ (А) виконуються наближено на основі лінійної апроксимації, що веде до збільшення помилок оцінок стану при зростанні похибки вимірювання частоти. І, по-друге, складність виразів, що визначають функції у виразі (7), призводить до істотного збільшення обсягу обчислень при розрахунку мат-

    Ріди Hx (fe). Бажання використовувати наближення 2-го порядку, як відомо, вимагає розрахунку матриці Гесса, що призводить до додаткових громіздким розрахунками.

    У зв'язку з цим доцільно застосувати, розвинений в роботах [4, 5], метод наближеного обчислення середнього і ковариаций величин при багатовимірних нелінійних перетвореннях, який істотно скорочує обсяг обчислень. Застосування цього методу до задачі фільтрації призводить до алгоритму UKF (unscented Kaiman Filter). У дослівному перекладі - не має відношення ( «Не пахне») до фільтру Калмана алгоритм.

    Суть методу в тому, що для обчислення (оцінки) моментів нелінійної функції z = h [x] використовується кінцеве безліч точок ( «сигма-точки») Х "= {хот |, i = 0, ..., 2пх; Хст е Кх, де Мх - простір станів і пх - його розмірність. Вибірка з «сигма-точок» використовується для розрахунку всіх необхідних моментів вектор-функції z = h [x] по її вибірковими значеннями на безлічі Х ". Розташування точок в задачі фільтрації на (А - 1) -му кроці визначається виразом

    (Xai (A-l)} = {x (A-l), х (А-1) ± уа; (* - 1)}, (12)

    де ог (k -1) - i-й стовпець матриці <JVX (k -1); у = <Jnx + ^! ^ - складовою параметр масштабу. Знак плюс в рівнянні (12) приймаємо для i = 1, ..., пх і знак мінус - для I = (пх + 1),...

    ..., 2пх, Квадратний корінь з матриці ^ Ух (k -1) (в даному випадку симетричної) випливає з її подання у вигляді

    1 = 1

    Для кожної «сигма-точки» (в даному завданні на етапі прогнозу спостереження) обчислюється відповідна точка в просторі спостережень М2, тобто.

    z ^ (A) = h [x ^ (ft)]. (14)

    Підсумкова екстрапольована оцінка спостережень формується за всіма «сигма-точкам» (14) у вигляді вагового середнього

    z- (fe) = 2X • (15)

    ?= 0

    де сог0 = Х / (пх + X) - вага центральної (нульовий) точки хо0 (ft -1) = х (k -1); (0zi = l / [2 (nx + Л)]; i = 1, ..., 2nx - вагові коефіцієнти нецентральних точок.

    Для розрахунку матриць К (а) і V "(a) в цьому випадку використовується інша (еквівалентна) форма, яка пов'язана з обчисленням ковариационной матриці Vs (й) оновлюючого процесу z (&) =? Z (a) - г (a)] і ковариационной матриці Vxi. Зазначені матриці також обчислюються з використанням «сигма-точок». Зокрема:

    2 пх у

    va (*) = Е С0уг • [ь (*) - Z- (*)] [z ^ (*) - Г (А)], (16)

    ?= 0

    де зі у 0 = Xj [nx + X) + (l - А2 +? 2) - вага центральної (нульовий) точки; avi = cozi для всіх i * 0. Масштабний параметр X = А2 (пх + l) -nx і параметри а,?, I забезпечують збіг моментів ПРВ нелінійної вектор-функції z = h [x] з моментами апроксимуючої її гаусом ПРВ. Детальний аналіз з цього питання виконаний в [4, 5].

    Результати дослідження імовірнісних характеристик оцінок координат

    Моделювання алгоритму обробки було виконано для випадку трьох КА, що мають високоеліптичного орбіти. Спостережувані сигнали (7) формувалися після рішення рівнянь (6) і розрахунку xQi (k) для кожного КА. Чисельне інтегрування рівнянь (6) виконано методом Рунге - Кутта 4-го порядку.

    Координати підсупутникових точок в момент початку інтервалу вимірювання мали такі значення: ХХ = 0 °; ФГ = 60 ° N, Л2 = 120 °? ; ф2 = 60 °, Х3 = 120 ° W; Ф3 = 60<W. Висоти КА над поверхнею Землі: fy = 4 • 104 м; йз = 4,5 | 104 м; Лд = 4 • 104 м. Положення випромінювача на поверхні Землі належало випадковим в деякій зоні S, центр якої m0 = = 84 ° 58'W фя0 = 56 ° 30'ЛГ]. У розрахунках апріорне розподіл ймовірностей

    W "(x (0)) було вибрано рівномірним. Початкові значення оцінок х2 (0) = Хз0 і х3 (0) = ф50. Дисперсії цих оцінок залежать від кутових розмірів зони S, яка, в свою чергу, визначається шириною діаграм спрямованості ( ДН) приймальних антен, відстанню до Землі і кутами видимості на центр зони. Кутова ширина зони S по обох координатах прийнята рівною 3 °.

    Оцінка хх (0) формується за першими трьома величинам f * = ЛГ1 [рр (l); А.а0, фз0], обчисленим на основі спостережень z (l). Дослідження імовірнісних властивостей величин / 3 * показало, що вони мають близькі до гаусовим розподілу ймовірностей і коррелірованни між собою. Таким чином, в якості оцінки частоти fs (О) = (0) доцільно прийняти максимально правдоподібну оцінку f t. При виконанні чисельних розрахунків функція правдо-

    "До Г * it | ^

    подібності вибіркового вектора fs = | / * i f $ 2 fs3 I мала гаусів вид. Максимально правдоподібна оцінка fsl формувалася у вигляді лінійної комбінацією величин f * t з коефіцієнтами, залежними від елементів матриці ковариаций. Наближені значення коваріаційних моментів обчислювалися на основі Е / Т-перетворення (Unscented Transformation) [6, 7].

    Дослідження СКО Оф - оцінок широти і а ^ - СКО оцінок довготи виконано для o2fch = 10; 50; 100 Гц і afs = 50 Гц. Темп надходження даних в каналі вимірювання частоти 1 з.

    Результати отримані при статистичному усередненні по 200 незалежним реалізацій збурень в каналах спостережень і значень вектора х (0). На рис. 2, а, б для частоти

    fs = 10 ГГц показані залежності СКО оцінок від часу вимірювання.

    а б

    Мал. 2 - Точність оцінок координат ІІ: а - СКО оцінки широти; б - СКО оцінки довготи

    Потенційна точність методу

    Дослідження потенційної точності оцінок координат передбачає обчислення дисперсій оцінки (8). Відзначимо два способи вирішення цього завдання. Перший з них пов'язаний з рішенням рівняння Стратоновича [3, 7] для апостеріорної ПРВ з подальшим обчисленням матриці ковариаций оцінок. Другий орієнтований на застосування методу Монте-Карло для обчислення багатовимірного інтеграла (9). У зарубіжній технічній літературі останній метод отримання оцінок (метод фільтрації) називають «particle filter» (PF - фільтр частинок) [6, 8]. Цей метод фактично випливає з відомого у вітчизняній літературі методу «важливою» (суттєвою) вибірки. Суть методу в тому, що в якості оцінки апостеріорної ПРВ W (х (&) | Z * j (далі використовуємо позначення х (&) = Хк) допустимо прийняти емпіричну оцінку у вигляді

    = (17)

    де - нормовані вагові коефіцієнти; 8 ( ») - багатовимірна дельта-функція Діра-

    (0

    ка; - сукупність «точкових мас» частинок в п-вимірному просторі стані, які

    є незалежними вибірковими значеннями з істотною ПРВ | 'К). У іто-

    ге в якості оцінки (8) отримуємо середнє арифметичне у вигляді

    = = (18)

    г = 1

    Відомо [6, 8], що при виконанні досить природних вимог до суттєвої ПРВ (хк | г *) емпірична ПРВ (17) сходиться майже, напевно, до істинної ПРВ. При цьому (18) є асимптотично несмещенной оцінкою істинного статистичного середнього (8).

    При умові завдання суттєвої ПРВ у вигляді [6, 8]

    WI (xk \ Z *) = w (xk \ 4l1) ненормовані вагові коефіцієнти пов'язані рекурсивним співвідношенням

    де | х ^ 1 ') - ПРВ, що визначає правдоподібність спостережень для частинки .

    У таблиці наведено результати фільтрації координат для = 10 Гц і о / 5 = 50 Гц при

    використанні алгоритмів ЩУ? і РЕ. Показані СКО оцінок координат ІІ при усередненні по 100 незалежним реалізацій на тридцятому часовому кроці (Л = 30).

    Алгоритм про ^, кут. хв а, -, кут. хв Ао Примітка

    UKF 1,739 2,644 а = 0,1; р = 2; / = 1

    PF 2,631 3,383 Число «частинок» N = 24 000

    PF 1,827 2,664 Число «частинок» N =: 98 000

    PF 1,782 2,641 Число «частинок» N = 250 000

    Результати моделювання показують збіжність розрахункових значень СКО оцінок координат при збільшенні кількості частинок в алгоритмі PF. Відзначимо, що рішення задачі з двома КА дає в кінцевій точці інтервалу спостереження (А = 100). СКО оцінок для обох координат приблизно 24 'і одного КА - близько 50'.

    висновок

    Запропонований в роботі алгоритм обробки вимірювань поточної частоти радіосигналів, прийнятих на борту одного чи кількох КА, що знаходяться на високоеліптичних орбітах, забезпечує формування стійких оцінок координат наземного ІРІ з невідомої частотою випромінювання. Точність оцінок координат при СКО вимірювання частоти ofch = 10 Гц і Ofs = 50 Гц фактично відповідає потенційно досяжної точності Байесови оцінок при квадратичної функції втрат.

    Виконання всіх обчислень на поточному часовому кроці вимагає 1,5 мс при використанні ЕОМ з процесором Athlon 64 3000+.

    література

    1. Жданюк Б.Ф. Основи статистичної обробки траєкторних вимірювань / Б.Ф. Джу-нюк. - М.: Сов. радіо, 1978.

    2. Сейдж Е. Теорія оцінювання та її застосування в зв'язку і управлінні / Е. Сейдж, Дж. Меле; пер. з англ. - М.: Зв'язок, 1976.

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

    4. Julier, S.J. The Scaled Unscented Transformation / S.J. Julier // Proceedings of the American Control Conference. - V. 6, Anchorage, AK, USA, May 200. - P. 4555-4559.

    5. Julier, S.J. A new approach for filtering nonlinear systems / S.J. Julier, Uhlmann and Durrant-Whyte // Proceedings of the American Control Conference. - 1995. - P. 1628-1632.

    6. Geweke, J. Bayesian Inference in Econometrics Models using Monte Carlo Integration / J. Geweke. - Econometrica. - 1989. - V. 57.

    7. Стратонович P.JI. Умовні марковские процеси та їх застосування до теорії оптимального управління / P.JI. Стратонович. - М.: МГУ, 1966.

    8. Doucet, A. Monte Carlo Methods for Bayesian Estimation of Hidden Markov Models / A. Doucet // Application to Radiation Signals, PhD. Thesis, Univ. Paris-Sud, Orsay, 1997..

    9. Сильвестров С.Д. Точність вимірювання параметрів руху космічних апаратів / С.Д. Сильвестров [и др.]; під ред. С.Д. Сильвестрова. - М.: Сов. Радіо, 1970.

    Савін Олександр Олександрович Аспірант кафедри. радіотехнічних систем ТУСУРа Телефон: (3822) 41 38 92 Ел. пошта: Ця електронна адреса захищена від спам-ботів. Вам потрібно увімкнути JavaScript, щоб побачити її.

    Тисленко Володимир Ілліч

    Канд. техн. наук, доцент каф. радіотехнічних систем ТУСУРа

    Телефон: (3822) 41 38 92

    Ел. пошта: Ця електронна адреса захищена від спам-ботів. Вам потрібно увімкнути JavaScript, щоб побачити її.

    Savin А.А., Tislenko V.I.

    A quasi-optimal estimation of the coordinates of a ground radiation source in space-based system measuring radio signals frequency

    A relative openness of space links used for information transferring gives a real possibility

    for various malefactors to have an unauthorized access to these links. A radical way of struggling against

    such «pirates» assumes determination of a ground radiation source coordinates with the use of signals

    which can be received by one or more satellites. Calculation results of root-mean-square errors (RMS)

    of joint estimates for a ground radiation source's coordinates (latitude, longitude) and frequency

    are presented in the paper. The estimates are generated at the output of a quasi-optimal filter,

    which input signals are measurements of HF carriers of the ground radiation source received

    by onboard receivers of one or more satellites at high-elliptic orbits. RMS values ​​of the quasi-optimal

    estimates are compared with the values ​​of optimal Bayes estimates provided quadratic loss.


    Завантажити оригінал статті:

    Завантажити