Теория математической обработки геодезических измерений (стр. 4 )

б) УУ полюса в точке пересечения диагоналей геодезического четырёхугольника с измеренными углами выглядит так:

;

в) УУ дирекционных углов в полигонометрическом ходе, имеющем угловую привязку, имеет вид:

φα = ∑βизм.– k*π–αкон.+αнач.= 0;

г) УУ фигур для измеренных углов в треугольнике сети триангуляции – это простейшая линейная функция:

φβ = β1+ β2 + β3 – π = 0;

д) УУ превышений в нивелирном ходе тоже линейно:

φh = h1+ h2 + … +hn – Hкон + Hнач= 0;

Система линейно-независимых УУС будет иметь вид;

Fr 1(YT1n; ZT1q) = 0r1 → (K.47)

Для гарантии линейной независимости очередного вновь вводимого УУ достаточно включить в него измерение, ранее не задействованное в предшествующих УУ. Тем не менее, может наступить момент, когда все измерения уже задействованы, но количество УУС ещё не достигло числа r. Здесь всё будет зависеть от опыта и внимательности составителя УУ. Он должен удовлетворить все выше перечисленные требования, предъявляемые к УУС.

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

Формирование ковариационной матрицы измерений Ky рекомендуется выполнять так, что бы её элементы были числами одного десятичного порядка. С этой целью можно прибегнуть к масштабированию СКО измерений, выражая их, допустим, не в радианной мере, а в угловой или не в метрах, а в сантиметрах и т. п. При этом необходимо помнить, что соответствующие «невязки» УУС должны быть выражены в единицах того же масштаба.

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

Ky = K = diag{m12, m22 … mn2} = diag{π1, π2 … πn}. (K.48)

Величины πi, входящие в ковариационную матрицу некоррелированных измерений (K.48), представляют собой обратные «веса»:

πi = 1/pi.

Покажем это. Веса определяются как величины, обратно пропорциональные дисперсиям. В качестве коэффициента пропорциональности выступает величина s02 º 1. Известно, что mi2 = s02 / pi. Отсюда окончательно получаем mi2 = 1 / pi = πi.

2. Линеаризация

Линеаризованные УУС получаются путём разложения в ряд Тейлора исходной системы (K.47):

Br n vn1 + Wr1 = 0r 1. (K.14)

Коэффициенты Br n, неизвестные vn1 и свободные члены Wr1 этих уравнений традиционно записываются так:

, (K.49)

v1nT = (v1, v2, … vn), (K.50)

W1rT = (w1, w2, … wr). (K.51)

В алгебраической форме линеаризованные УУС имеют вид:

. (K.52)

Важнейшим моментом второго этапа является определение допустимости свободных членов («невязок») линеаризованных УУС.

Допустимые значения невязок УУС.

Вектор невязок Wr1 вычисляется по данным измерений yn1 и координатам опорных точек zq1:

Изучив статистические свойства векторов-оценивателей коррелатной версии, мы установили, что: E(W) = 0; а KW = N. Из последнего факта следует, что , т. е. среднее квадратическое значение j-ой невязки – это диагональный элемент матрицы коэффициентов нормальных уравнений коррелат, находящийся на пересечении j-ой строки и j-ого столбца.

Установление допустимости невязки w j начинается с проверки нулевой гипотезы о равенстве нулю её МО:

H0 = {E(wj) = 0} (K.53)

против альтернативной гипотезы

HA = {E(wj) ≠ 0}. (K.54)

Нулевая гипотеза (K.53) проверяется с помощью теста

tЭ = | wj |/σw j, (K.55)

который сопоставляется на уровне значимости α с квантилью tT распределения Стьюдента, характеризующегося r = n – k степенями свободы:

tT = arg(FСт.(n-1)= (1–α/2)). (K.56)

Воспользовавшись неравенством альтернативной гипотезы (K.54), получаем такое выражение для допустимого значения невязки j-го УУС:

tЭ = | wj |/σwj ≤ tT → wjдоп = tT* σwj = tT*. (K.57)

Квантиль распределения Стьюдента, называемая «доверительным множителем», часто заменяется квантилью нормального распределения на том же уровне значимости. Нормальные доверительные множители, соответствующие уровням значимости 0.05, 0.01 и 0.003 равны:

t0.05 = 1.96 ≈ 2; t0.01 = 2.58 ≈ 2.6; t0.003 = 3. (K.58)

Если вычисленная невязка лежит в допустимых пределах, это не означает, что она найдена без ошибок. К сожалению, не существует прямых контрольных выражений для проверки вычисления невязки. Выход один – вычисления «во вторую руку».

3. Нормализация

Нормальные уравнения коррелат в матричной форме имеют вид:

Nr r*Λr1 – Wr1 = 0r1. (K.20)

Коэффициенты нормальных уравнений представляют собой произведение двух ранее введённых матриц:

Nr r = Br n Kn n BTn r. (K.21)

Выполняя умножение данных матриц с учётом (K.48) и (K.49), получаем традиционную алгебраическую форму записи коэффициентов нормальных уравнений коррелат для случая некоррелированных измерений:

. (K.59)

Свободные члены НУ коррелат (K.51) переносятся из линеаризованных УУС. Вектор неизвестных (коррелат), записанный в строку, имеет вид:

Λ1rT= (λ1 l2 … lr). (K.60)

Теперь, выполняя умножение матрицы (K.59) на вектор коррелат (K.60) и вычитая из произведения вектор «невязок» (K.51), получаем алгебраическую запись системы НУ коррелат:

. (K.61)

4. Решение НУ

Решение НУ коррелат, предлагаемое в теоретической части изложения материала, осуществляется методом обращения матрицы коэффициентов Nr r системы НУ (K.61):

Lr 1 = Nr r-1 Wr 1. (K.22)

Это сделано по следующим соображениям.

Во-первых, предполагая линейную независимость УУС, мы вправе считать, что матрица коэффициентов Brn имеет полный строчный ранг, т. е. rank(B) = r. В таком случае матрица коэффициентов НУ Nr r = Br n K BnrT будет иметь такой же ранг, и её определитель не будет равен нулю, т. е. det(Nr r) ≠ 0. Это означает, что существует N-1 и решение системы (K.61) в форме (K.22).

Во-вторых, использование обратной матрицы N-1 значительно упростило как уже выполненные теоретические выкладки, так и предстоящие преобразования.

В развёрнутом виде уравнения (K.22) можно записать, учтя тот факт, что обратная матрица коэффициентов нормальных уравнений является ковариационной матрицей коррелат, т. е. Nr r-1 = KΛ:

. (K.62)

Обратите внимание, что каждый диагональный элемент Kjj – это квадрат среднего квадратического значения j-ой коррелаты!

Отдельно j-ая строка системы (K.62) выглядит следующим образом:

l j = Kj 1* w1 + … +Kj j* wj + … + Kj r* wr. (K.63)

Выше уже отмечалось, что вычисление элементов Kij обратной матрицы N-1 через союзную матрицу на практике представляет собой трудоёмкую вычислительную задачу. Много проще реализовывать обращение с использованием различных алгоритмов решения СЛАУ.

В геодезии традиционно широко применяется алгоритм , называемый «методом последовательного исключения неизвестных», и его модификации.

Суть метода Гаусса состоит в эквивалентной замене матрицы коэффициентов нормальных уравнений Nrr верхней треугольной матрицей Urr (для сокращения записей обратные веса π опущены, что не влияет на сущность алгоритма). Эквивалентное преобразование реализуется путём последовательной замены строк исходной матрицы их линейными комбинациями друг с другом. Итак, матрица

(K.59)

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

. (K.64)

Получаемая эквивалентная СЛАУ

(K.65)

легко решается в обратном порядке:

. (K.66)

Коррелаты Λ1rT = (l1, l2, … l1r) вычисленные методом обращения или методом Гаусса, либо каким-то другим способом, должны быть проконтролированы. Такой контроль осуществляется путём подстановки найденных корней в исходную систему НУ коррелат (K.20):

Nr r*Lr1 = Wr1. (K.67)

Вернёмся к вопросу получения элементов обратной матрицы с помощью некоторого алгоритма решения СЛАУ. Очевидно, что такой алгоритм будет линейным преобразованием свободных членов по некоторому правилу Z:

Lr1 = Nr r-1 Wr 1 = Zr r Wr1. (K.68)

Воспользовавшись тождеством N-1 = N-1 *I, выразим его через алгоритм Z:

N-1 = Z*I. (K.69)

Развернём этот результат подробнее:

. (K.70)

Из соотношений (K.70) следует очевидное правило вычисления столбцов обратной матрицы путём r-кратного использования алгоритма Z:

. (K.71)

5. МНК-оценивание

Коррелаты, полученные и проконтролированные на 4ом этапе, позволяют вычислить МНК-поправки в измерения yn1:

. (K.19)

Более подробно, с учётом (K.48), данная система записывается так:

. (K.72)

Получена коррелатная система уравнений поправок. Алгебраический эквивалент i-ой строки данной системы имеет вид:

. (K.73)

Вычисление МНК-поправок в измерения контролируется путём подстановки их в линеаризованные УУС (K.13):

Br n = – Wr1. (K.74)

Допустимые значения МНК - поправок.

В статистических свойствах векторов-оценивателей коррелатной версии, установлено, что: E() = 0; а =KBTN-1BK. Из последнего факта следует, что , т. е. среднее квадратическое значение i-ой поправки – это диагональный элемент ковариационной матрицы МНК-поправок в измерения, находящийся на пересечении i-ой строки и i-ого столбца этой матрицы.

Установление допустимости поправки начинается с проверки нулевой гипотезы о равенстве нулю её МО:

H0 = {E() = 0} (K.75)

против альтернативной гипотезы

HA = {E() ≠ 0}. (K.76)

Нулевая гипотеза (K.75) проверяется с помощью теста

tЭ = | |/, (K.77)

который сопоставляется на уровне значимости α с квантилью tT распределения Стьюдента, характеризующегося r = n – k степенями свободы:

tT = arg(FСт.(r)= (1–α/2)). (K.78)

Воспользовавшись неравенством альтернативной гипотезы (K.76), получаем такое выражение для допустимого значения МНК-поправки в i-ое измерение:

tЭ = | |/ = tT → доп = tT* = tT*. (K.79)

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

t0.05 = 1.96 ≈ 2; t0.01 = 2.58 ≈ 2.6; t0.003 = 3. (K.58)

Если все невязки Wj были в допустимых пределах, то и МНК-поправки окажутся допустимыми. Тем не менее, наличие недопустимых невязок побуждает нас локализовать измерения, содержащие «промахи». В такой ситуации анализ МНК-поправок на допустимость может сыграть положительную роль.

6. МНК-оптимизация (уравнивание)

Данный этап является заключительным шагом процедуры уравнивания, выполняемой с целью нахождения НДЗ измерявшихся величин. С вычислительной точки зрения он элементарен:

. (K.23)

Каждое измерение получает соответствующую МНК-поправку:

. (K.80)

В качестве контроля необходимо проверить удовлетворение исходных УУС (K.3):

Fr 1(; ZT1q) = 0r1. (K.81)

Уравнения (K.81) должны быть удовлетворены с той точностью, с которой вычислялись невязки wj.

7. Показатель точности

Апостериорное значение МПТ

, (K.33)

содержит в числителе квадратичную форму, получение которой необходимо проконтролировать. Следующая цепочка матричных преобразований даёт необходимое контрольное выражение (использовано выражение (K.19)):

=. (K.82)

Для случая независимых наблюдений формулы (K.82) и (К.33) имеют алгебраические эквиваленты, представленные в символах Гаусса:

, . (К.83)

Вычисление величины m2 по формуле (K.33) или (К.83) не контролируется.

8-9. Ковариации и a priori

Вычисление ковариационной матрицы МНК-поправок в измерения

может быть проконтролировано тождеством

= r. (K.84)

Стоящая в левой части формулы (K.84) сумма отношений дисперсий является следом произведения двух матриц: и K-1, т. е.

tr(*K-1) = . (K.85)

Покажем это на примере трёх измерений.

Пусть

, а .

Произведение этих матриц равно:

,

где kij = Kij / σj2, а след матрицы – это и есть искомая сумма (K.84).

Определим интересующий нас след для алгоритма коррелатной версии:

= r. (K.86)

Итак, контрольное тождество (K.84) доказано.

В среднем, дисперсии МНК-поправок в измерения и дисперсии независимых измерений (они же дисперсии истинных поправок) относятся как число избыточных измерений «r» к числу всех измерений «n»:

/ n = r / n. (K.87)

Из последнего соотношения следует, что модули МНК-поправок в измерения в среднем короче модулей истинных поправок , т. е.

/ ≈ . (К.88)

Априорная ковариационная матрица уравненных измерений, как это показано в параграфе о статистических свойствах, представляет собой разность ковариационных матриц измерений yn1 и МНК-поправок к ним :

.

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

= k, (K.89)

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

Во-первых,

= tr(*K-1), (K.90)

а, во-вторых,

= n–r =k. (K.91)

10. Ковариации a posteriori

Апостериорная ковариационная матрица уравненных измерений получается путём умножения априорной матрицы на апостериорное значение МПТ m2. Диагональные элементы матрицы участвуют в контроле процесса построения предваряющих её ковариационных матриц.

Апостериорная ковариационная матрица уравненных измерений

, (K.92)

позволяет использовать контроль (K.90), модулированный делением следа произведения *K-1 на апостериорное значение МПТ μ2:

. (K.93)

Параметрическая версия МНК-оптимизации и оценки точности

данных (Параметрический способ уравнивания)

Процесс математической обработки результатов измерений в ГП на основе параметрической версии отличается от вышеизложенной коррелатной версии тем, что каждая измеряемая величина Yi представляется в виде явной функции Fi от вектора линейно независимых параметров Xk1, каковыми могут быть как некоторые измерявшиеся величины, так и другие параметры построения, например, координаты определяемых ГТ сети. При этом возникают те же три задачи, расширенные появлением дополнительных неизвестных – параметров. Первая – это нахождение НДЗ измеренных величин и параметров, вторая – оценка точности измерений и третья – оценка точности НДЗ измеренных величин и параметров, а так же функций от параметров.

Первая задача традиционно называется уравниванием измерений или, с учётом того, что «уравнивание» выполняется под условием минимума суммы квадратов уклонений результатов измерений от некоторых оптимальных значений, – МНК-оптимизацией.

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

Имеем «n» измеренных величин y1, y2, … , yn, вектор реальных значений которых выглядит следующим образом Y1nT = (Y1,…,Yn). Они образуют некоторую систему и связаны между собой посредством линейно независимых параметров X1kT = (X1,…,Xk). Такая система представляет собой математическую модель ГП и называется системой параметрических уравнений связи (ПУС):

Yn1 = Fn1(XT1k; ZT1q), (П.1)

где «k» – число линейно независимых параметров, равное числу необходимых измерений, определяемому целью создания ГП.

Результаты измерений yi (i = 1, 2, …, n), как это было и в коррелатной версии, характеризуются некоторой ковариационной матрицей

, (К.6)

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

= = (К.7)

Здесь – среднее арифметическое i – ого измерения, равное

= (К.8)

Корреляционные моменты (ковариации) Kij могут быть оценены по формуле

(K.9)

По-прежнему, предполагается, что:

1) результаты измерений yi необходимо отягощены только случайными ошибками и представляют собой элементы спектров соответствующих СВ Yi, т. е. yi Yi. Тот факт, что измерения свободны от постоянных систематических ошибок, моделируется условием совпадения математического ожидания СВ Yi и реального значения Yi измерявшейся величины:

E(Yi) = Yi ; (K.10)

2) координаты опорных точек zq1 = Zq1 рассматриваются как безошибочные константы;

3) приближённые значения параметров x1, x2, … , xk трактуются как неслучайные величины, на которые не действует оператор математического ожидания, т. е.

E(xj) = xj ; (П.2)

4) обратная ковариационная матрица таких параметров тождественно равна нулевой матрице:

. (П.3)

Резюмируя сказанное, сконцентрируем имеющуюся информацию.

Дано:

Проект ГП.

Схема, чертёж высотной или плановой геодезической сети.

Математическая модель ГП.

Yn1 = Fn1(XT1k; ZT1q) – ПУС, где Yi и Zt – реальные значения измеряемых величин и исходных данных, Xj реальные значения параметров;

i = 1, 2, …, n; t = 1, 2, …, q; j = 1, 2, …, k.

Числовые данные.

yT1n = (y1, y2, … yn) – вектор результатов измерений;

xT1k = (x1, x2, … xk) – вектор приближённых значений параметров;

zT1q = (z1, z2, … zt) – вектор координат опорных пунктов;

– ковариационная матрица измерений (стохастическое расширение математической модели).

Теоретические посылки.

Y1nT=(Y1, Y2, … Yn) – вектор СВ Yi, являющихся вероятностными моделями измерений;

E(Yn1) = Yn1 – условие отсутствия постоянных систематических ошибок в измерениях;

E(xj) = xj и – неслучайность приближённых значений параметров;

– априорное значение масштабного показателя точности измерений.

Найти:

1) НДЗ (уравненные значения) измеренных величин Ỹn1= и параметров ;

2) апостериорное значение показателя точности измерений – ;

3) апостериорные ковариационные матрицы НДЗ , и функций от НДЗ параметров : , и .

Решение.

1.  Нахождение наиболее достоверных значений измеренных

величин и параметров (Уравнивание).

Представим векторы реальных значений измеряемых величин и параметров в виде сумм исходных векторов и малых поправок к ним:

Yn1 = yn1 + vn1, (K.11)

Xk1 = xk1 + Xk1. (П.4)

Поправки к измерениям vn1 и приближённым значениям параметров Xk1 полагаем значительно мèньшими по модулю самих результатов измерений yn1 и приближённых параметров xk1 т. е. |vn1| << | yn1| и |Xk1| << | xk1|. Такое предположение основывается на том, что геодезические измерения yn1 выполняются с относительными СКО порядка не ниже 10-3÷10-4, а приближённые значения параметров xk1 вычисляют по тем же измерениям.

При подстановке (K.11) и (П.4) в ПУС (П.1) получаем такой результат:

yn1 + vn1 = Fn1(xT1k + XT1k; ZT1q), (П.5)

Разложив функцию (П.5) в ряд Тейлора в окрестностях точек yn1 и xk1 и ограничившись только линейными членами, получим линеаризованные ПУС (ЛПУС):

yn1 + vn1 = Fn1(xT1k;zT1q) + {∂Fi/∂Xj}*Xk1 . (П.6)

Частные производные функций Fi по параметрам Xj – это матрица коэффициентов линеаризованных ПУС Ank. Числовые значения производных находят по приближённым значениям параметров xk1. Разность векторов yn1 и Fn1(xT1k;zT1q) представляет собой вектор «свободных членов» ЛПУС:

Ln1 = yn1 – Fn1(xT1k zT1q). (П.7)

С учётом введённых обозначений, уравнения (П.6) запишутся таким образом:

An k Xk1 – Ln1 = vn 1. (П.8)

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

rank(An k) = k. (П.9)

В таком случае матрица An k будет матрицей полного столбцового ранга.

Система ЛПУС (П.8) содержит, кроме неизвестных поправок Xk1 к параметрам xk1, так же не известные поправки vn 1 в измерения yn1. Это означает, что она будет иметь бесчисленное множество решений.

Для нахождения единственных значений зависимых векторов и на систему (П.8) накладывается МНК-ограничение:

Ank * – Ln1 =

. (П.10)

Y =

Необходимым условием существования экстремума функции Y является равенство нулю её частных производных:

. (П.11)

Система (П.11) содержит «k» уравнений и «n» неизвестных , записанных в строку. Транспонируя эту систему, получаем «лемму Гаусса»:

. (П.12)

Эти неизвестные являются функциями (П.8) «k» параметров, что позволяет получить систему нормальных параметрических уравнений, число которых равно числу неизвестных:

. (П.13)

Здесь

Nkk = AT K-1 A – (П.14)

матрица коэффициентов, а

Gk1 = AT K-1 L – (П.15)

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

Предположив линейную независимость вектора параметров, мы установили, что матрица коэффициентов ЛПУС An k является матрицей полного столбцового ранга. В таком случае матрица коэффициентов нормальных уравнений Nkk = AT K-1 A будет квадратной матрицей полного ранга, т. е.