Аналитический расчет координат места судна при избыточном числе измерений

ФЕДЕРАЛЬНОЕ АГЕНСТВО МОРСКОГО И РЕЧНОГО ТРАНСПОРТА

ФЕДЕРАЛЬНОЕ ГОСУДАРСТВЕННОЕ ОБРАЗОВАТЕЛЬНОЕ  УЧРЕЖДЕНИЕ ВЫСШЕГО ПРОФЕССИОНАЛЬНОГО  ОБРАЗОВАНИЯ

МОРСКОЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ

 имени  адмирала Г.И.Невельского 

КАФЕДРА СУДОВОЖДЕНИЯ 

КУРСОВОЙ  ПРОЕКТ

ТЕМА: АНАЛИТИЧЕСКИЙ РАСЧЕТ КООРДИНАТ МЕСТА СУДНА ПРИ ИЗБЫТОЧНОМ ЧИСЛЕ ИЗМЕРЕНИЙ. 
 
 
 
 

                                                                              Выполнил: студент

                                                                                                 Лицай А. В.

                                                                                    Проверил: преподаватель

                                                                                                Домбинский А. П.   
 
 

Владивосток

2011

1. Составление уравнений линий  положения. Уравнение  поправок.

 

Счислимое место судна определяется координатами φ_с и λ_с. Если выполнены измерения навигационного параметра, то задача будет сводиться к вычислению поправок ∆φ и ∆ω (а затем и ∆λ) для перехода к обсервованному месту О с координатами φ₀ и λ₀. При любом способе обсервации можно воспользоваться азимутальным приемом перехода от счислимого места к обсервованному, сущность которого сводится к определению величины смещения ∆n и направления τ переноса из счислимого места к той линии положения, на которой ищется место судна (рис.1).

 

                       φ

                                             I

                   φ_o                                              O(φ_o,λ_o)

                                τ                                               g

                                           ∆n                                                       Рис. 1

                             τ

                                                                        I

               C(φ_c,λ_c)                     w_o                                     w

Пусть для определения  места судна измерен навигационный  параметр и₀, являющийся функцией координат.

и₀₌ f (φ₀, λ₀)  (1)

Наибольшее изменение  этой функции в данной точке характеризуется  градиентом, величина которого g определяется формулой

g=  (2)

или приближенно

g,  (3) 

где и₀ – измеренный и исправленный поправками навигационный параметр;

      и_с – рассчитанный по счислимым координатам «счислимый» навигационный параметр;

      v – ошибка навигационного параметра и₀, не известная при измерениях, как правило, малая величина.

Как видно из рис.1, линия положения ( I – I), на которой  находится обсервованное место О, определяется уравнением

            (4)

Точку К на линии положения принято называть определяющей точкой. С учетом уравнения (3) из уравнения (4) получаем

            (5)

откуда 

            (6) 

Это уравнение  иногда называют уравнением поправок (ошибок).

Уравнению поправок можно придать и другой вид. Для  этого разложим приращение функции  ∆u=u_o+v-u_c в ряд и, пологая, что счислимое место достаточно близко к обсервованному, ограничимся первыми членами разложения

      (7)

или 

Сравнивая уравнения (7) и (6), получим 

                                 (8)

отсюда

    .        (9)

Выражения (9) дают общий путь нахождения градиентов для  различных функций. Однако иногда проще  находить эти выражения из геометрических соображений.

2. Основы метода наименьших квадратов

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

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

  Так, для оценки приращения  координат Δφ  и Δλ на плоскости могут быть измерены три навигационных параметра, т.е. n=3. В этой ситуации говорят. Что избыточность равна единице. Система уравнений линий положений будет иметь вид:

      g ₁ cos τ₁^.Δφ+ g ₁ sin τ₁^.Δω = u₀₁-u_c1

      g ₂ cos τ₂^.Δφ+ g ₂ sin τ₂^.Δω = u₀₂-u_c2 (10)

      g ₃ cos τ₃^.Δφ+ g ₃ sin τ₃^.Δω = u₀₃-u_c3

  Из-за наличия погрешностей измерений линии положения образуют фигуру погрешностей – треугольник (рис 2)

  Решение  любых двух уравнений из трёх  даст нам положение вершин  этого треугольника относительно  начала координат (счислимого места). Это означает, что решение любой пары уравнений не обращает в тождество оставшееся уравнение. Такая система не совместна, т.е. решение любой пары несовместно с третьим уравнением. Анализ размеров фигуры погрешностей и её поведение в последовательности измерений дают полезные сведенья о качестве измерительной навигационной информации.

 φ II

 I 

 

 III  

      I

 II III

      C(φ с,λ _с) λ 

  Для того, чтобы получить согласованное решение необходимо ввести дополнительное условия, которые можно получить, если более детально представить систему (10). Необходимо в окрестности фигуры погрешностей выбрать точку(предполагаемое решение), относительно которой и можно было бы сформулировать дополнительное условие.

  Рис.2 дает  представление о фигуре погрешностей при n=3. Отрезки u ^/1 u ^/₂ u ^/₃  называются невязками. В специальной литературе также встречается термин «поправка» или «ошибка линий положения» в зависимости от знака. В такой ситуации именно невязки определяют решение относительно фигуры погрешностей. Множество возможных сочетаний невязок определяет множество решений и задача заключается в подборе наиболее простого и физически интерполируемого условия. если обозначить Δu₁= u_0i- u_ci, то с учетом невязок систему (10) можно записать так:  

 g ₁ cos τ₁^.Δφ+ g ₁ sin τ₁^.Δω-Δu₁ = u₀₁-u_c1

      g ₂ cos τ₂^.Δφ+ g ₂ sin τ₂^.Δω-Δu₂ = u₀₂-u_c2 (11)

      g ₃ cos τ₃^.Δφ+ g ₃ sin τ₃^.Δω-Δu₃= u₀₃-u_c3

  Здесь величины u₁ u₂ u₃ – невязки, выраженные в единицах измерений навигационного параметра, что более удобно для дальнейших выкладок.

  Для согласования системы  (11) с рис.2 необходимо выразить  невязки u_i в линейных единицах т.е. разделить все челены на модули соответствующих градиентов:

      g ₁ cos τ₁^.Δφ+ g ₁ sin τ₁^.Δω-Δn₁ = u₀₁-u_c1

      g ₂ cos τ₂^.Δφ+ g ₂ sin τ₂^.Δω-Δn₂ = u₀₂-u_c2 (12)

      g ₃ cos τ₃^.Δφ+ g ₃ sin τ₃^.Δω-Δn₃= u₀₃-u_c3 

  Избыточная  система уравнений(10) превратилась  в систему (12) с недостаточным  числом уравнений, так как невязки  также неизвестны. Формально наиболее  простое решение такой системы  можно получить, если принять  следующее условие:

S= [u²]= u²₁+ u²₂ +u²₃=min,                  ₍₁₄₎ 

  При увеличении  числа измерений (n>3) возникает более сложная фигура погрешностей с числом вершин , но условие (13) даёт однозначное решение системы уравнений линий положения при любом значении n. Его суть сводится к поиску «центра тяжести» фигуры погрешностей измерений, т.е. к определению некоторого среднего значения координат из множества координат вершины фигуры. Решение системы уравнений называется оптимальным в смысле выполнения условия S=min, т.е. при выполнении критерия минимума суммы квадратов невязок.

  Избыточность  позволяет нам получить информацию о некоторых средних значениях координат,  а поэтому важным является утверждение, что оптимальная точка будет всегда находиться в нутрии фигуры погрешностей, если систематические погрешности δ_i=0. Метод наименьших квадратов является наиболее универсальным средством обработки избыточной навигационной информации. Его основы были разработаны Лежандром и Гауссом с 1795 по 1805 год.

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

3. Cоставление и решение нормальных уравнений

Обозначим                       g_icosτ_i→a_i

        g_icosτ_i→b_i

         Δu_i→c_i

Тогда система  уравнений (11) запишется следующим образом:

 a ₁^.Δφ+ b ₁^.Δω+l₁ = u₁

      a ₂^.Δφ+ b ₂^.Δω+l₂ = u_2 (15)

      a ₃^.Δφ+ b ₃^.Δω+l₃ = u₃ 

 Для решения  системы (15) воспользуемся критерием S=min и запишем его выражения через левые части уравнений

S= [u²]= u²₁+ u²₂ +u²₃=( a ₁^.Δφ+ b ₁^.Δω+l₁)² +( a ₂^.Δφ+ b ₂^.Δω+l₂)²+( a ₃^.Δφ+ b ₃^.Δω+l₃)²    (16)

Для определения  приращенияΔφ и Δλ, соответствующих минимуму критерия S, используем традиционный метод поиска экстремальных значений функций, взяв частные производные от S по Δφ  и Δω:

   и  .

Тогда

       2a₁(a₁Δφ+b₁Δω+l₁)+ 2a₂(a₂Δφ+b₂Δω+l₂)+2a₃(a₃Δφ+b₃Δω+l₃)=0 ₍₁₇₎

2b₁(a₁Δφ+b₁Δω+l₁)+ 2b₂(a₂Δφ+b₂Δω+l₂)+2b₃(a₃Δφ+b₃Δω+l₃)=0 

Выполнив сложение, получаем систему 2-х нормальных уравнений  в обозначениях Гаусса:                       

             [aa]Δφ+[ab]Δω+[al]=0 ₍₁₈₎

[ab]Δφ+[bb]Δω+[bl]=0

Так как  u_i=a_iΔφ+b_iΔω+l_i, то

[au]=0

[bu]=0

Матрица коэффициентов  системы (18)

Имеем следующие  свойства:

      -коэффициенты  на главной диагонали продолжительны;

      -коэффициенты, расположенные симметрично относительно  главной диагонали, равны. Это  симметричная система.

Универсальная запись нормальных уравнений в матричном  виде выглядит следующим образом                                  (A^TA)ΔX=A^TL

A^T- транспонированная  матрица А т.е., если

A=        ,то                 A^T= 

Формальная запись решения системы (19) при обращаемой матрице (А^TA) такова:

                                                                Х=(A^TA)^-1А^ТL

Расшифровка записи (20) приводит к формулам Крамера для решения системы. При определении места по трём линиям положения можно написать:

   =                       ;

      =          

   =          

     =            =                          (21) 

   Таким  образом получено правило Крамера, где D-главный определитель системы, а D_Δφ   и D_Δω  - определители для Δφ и Δω   соответственно.

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

([aa]+[ab])Δφ+([ab]+[bb])Δω+([al]+[bl])=0                 (22)

Способ решения  нормальных уравнений по правилу  Крамера при n>2 становится трудоёмким и не всегда устойчивым при малых значениях D.

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

     -способ  последовательного исключения искомых  величин;

     -способ  последовательных приближений (итерации)

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

  Способ итерации  легко реализуется на ЭВМ  , к недостатку стоит отнести итерационную процедуру, которая не даёт конечного решения, но быстродействие современных ЭВМ снимает этот вопрос. 

      4. Метод наименьших  квадратов при  неравноточных измерениях.

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

      Поэтому естественно считать, что поправки v_i, применяемые для формирования критерия S, должны быть обратно пропорциональны соответствующим СКП измерений, т.е. величинам m_i. Тогда критерий МНК будет рассчитан так:

       (23)

Величина р_i называется весом измерения и характеризует относительную точность измерений, входящих в группу. Таким образом критерий S в соответствии с (23)

  (24)

Исходя из формулы (24) система уравнений поправок для  трёх линий положения будет иметь  вид:

         (25)

После приведения системы (25) к нормальным уравнениям для двух неизвестных ∆φ и ∆ω

     (26)

Более компактно 

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

   (27) 

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

5. Оценка точности определения места по методу наименьших квадратов.

Поправки ∆φ и ∆ω, найденные методом наименьших квадратов, характеризуются следующими показателями точности (для неравноточных измерений) 

            (28)

где - СКП измерения с весом, равным единице;

D - главный определитель системы.

  при n › 5 вычисляется апостериорно, т.е. после вычисления поправок ∆φ и ∆ω

           (29)

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

При малом числе  измерений (n<5) оценка m₍₁₎ становится весьма приближенной, поэтому принимается m₍₁₎=1.

Элементы среднего квадратического эллипса погрешностей так же выражаются через коэффициенты нормальных условий и рассчитываются по следующим формулам: 

      (30) 

Если  величина положительная , то a, вычисленный по третьей формуле (30) определяет направление малой оси, если же отрицательная, то угол a определяет направление (относительно меридиана) большой оси эллипса.

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

Во всех формулах D - главный определитель системы нормальных уравнений. При n < 5 за величину принимается единица.

    Решить  задачу № 1/134, 

 
 

φ_с = 36° 20,0'S λ_с = 129° 30,0'Е

Аналитическое решение

Промежуточные результаты расчета сведем в таблицу: 

Контроль результатов  расчета дает

[раа] + [pab] + [pal] = [pas] =7,821

[pab] + [pbb] + [pbl] = [pbs] =5,026

Результаты  могут отличаться на 0,001 - 0,002 за счет погрешностей округления чисел.

Запишем систему нормальных уравнений с  рассчитанными коэффициентами

4,129∆φ+2,747∆W= -0,945

2,747∆φ+2,906∆W= 0,627

Рассчитаем определители системы уравнений

Д=4,453

Д_∆φ = -4,469

Д_∆ω=5,185

∆φ= Д_∆φ/Д= -1,004-1,0' к S

∆W=Д_∆W/Д= 1,164 1,2' к E

= = 1,4901,5' к E

φ_c  = 36º 20,0' S    λ_с= 129º 30,0' E 

∆φ=        1,0'  к S    ∆λ=          1,5'  к E

φ_o= 36º 21,0' S    λ_o= 129º 31,5' E

Рассчитаем вспомогательную  величину q и параметры эллипса погрешностей

q=5,628 
 
 

0,976

 

Следовательно угол 2Т_а лежит в 3-ей четверти

2Ta=257,6 º;   Та=128,8º 

Графоаналитическое  решение

Промежуточные результаты расчета занесем в  таблицу 

 

    По  параметрам τ и ∆n выполним прокладку линий положения в масштабе 1 см = 1 миля и проверим значения углов θ между линиями положения. 
 
 
 
 
 
 
 
 
 
 
 

Снимем с прокладки  приращения координат в принятом масштабе

∆φ₁₃= -0,8'     ∆W₁₃= -2,0'

∆φ₁₄= -2,2'     ∆W₁₄= 4,1' 

Рассчитаем  вероятнейшие координаты: 
 
 
 

==1,6141,6 ' к E 

φ_c  = 36º 20,0' S    λ_с= 129º 30,0' E 

∆φ=         -1,5'  к S    ∆λ=          1,6'  к E

φ_o=  36º 21,5' S    λ_o= 129º 31,6' E 

    Нанесем на рисунок обсервованную точку по полученным приращениям координат ∆φ и ∆W. Для оценки точности полученного вероятнейшего места все необходимые данные для расчета имеются в приведенной выше таблице.

    Построим  «полигон весов», последовательно складывая  векторы длиной Р_лп, направленные под углами 2τ. Масштаб построения «полигона весов» необходимо выбрать самостоятельно с тем, чтобы рисунок поместился примерно на половину тетрадного листа. В нашем примере можно выбрать масштаб 1 см = 1 единице веса. 
 
 
 
 
 
 
 

    Снятая  с рисунка величина замыкающего  вектора q – 5,7 см = 5,7 ед. веса.

    Построим  биссектрису угла между вектором q и направлением на N. В результате получим направление малой полуоси эллипса погрешностей Т_b. Перпендикулярно Т_b проведем направление большой полуоси Т_а и снимем транспортиром угол между полученным направлением Т_а и N. В нашем примере он равен: Т_а=129°.

Рассчитаем  веса эквивалентных линий положения

P_max=

P_min= 
 
 
 

    Нанесем эллипс погрешностей на рисунок, где выполнена прокладка линий положения. Построим его в принятом масштабе 1 см = 1 миля в месте расположения обсервованной точки, ориентируя направление Т_а относительно N.