Добавил:
Upload Опубликованный материал нарушает ваши авторские права? Сообщите нам.
Вуз: Предмет: Файл:

Вища геодезія для Mathcad

.pdf
Скачиваний:
20
Добавлен:
19.02.2016
Размер:
1.64 Mб
Скачать

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

2

2

 

 

 

 

 

 

 

 

 

 

 

 

 

d2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

d2R

 

 

R simplify

o 48 t2 24 t 2

26

 

 

 

 

 

 

 

 

2

 

 

 

 

 

 

 

 

 

 

 

 

 

dt

 

 

2

11

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

12 t 8

 

 

 

 

 

 

 

 

ȼɟɤɬɨɪ ɝɨɥɨɜɧɨʀ ɧɨɪɦɚɥɿ ɩɪɢ t = 0.5

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

532

 

 

 

 

 

 

 

 

N

(dR u d2R) u dR

1515.5

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

882

 

 

 

 

 

 

 

 

Ʉɪɢɜɢɡɧɚ ɬɚ ɜɟɤɬɨɪ ɤɪɢɜɢɡɧɢ ɩɪɢ t = 0.5

 

 

 

 

 

 

 

 

 

 

 

 

dR u d2R

 

 

 

 

 

 

 

0.038697769149986

k

 

 

 

0.13328905309374

k

 

N

 

0.110237723960157

 

 

 

dR

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

3

 

 

 

 

 

 

N

0.064156827801292

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ɇɨɪɦɚɥɶɧɚ ɬɚ ɝɟɨɞɟɡɢɱɧɚ ɤɪɢɜɢɡɧɚ

 

 

 

 

 

 

 

 

 

 

 

 

k

 

 

 

k

 

 

N n

 

 

0.049049462270012

k

 

 

 

 

k2 k

2

0.123935958969335

 

 

 

 

 

 

 

 

n

 

 

N

 

n

g

 

n

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

3.9. Геодезична лiнiя

Геодезичною лiнiєю називається така крива на поверхнi, у якої геодезична кривизна дорiвнює нулю в усiх її точках. Серед геометричних властивостей геодезичної лiнiї найважливiшою є те, що вона є найкоротшим шляхом мiж двома точками на криволiнiйнiй поверхнi.

Для рiвностi нулю геодезичної кривизни лiнiї на поверхнi необхiдно i

#„

достатньо, щоб вектор кривизни k кривої був в усiх точках направлений по

#„

нормалi до поверхнi, а тому достатньо, щоб вектор N був перпендикулярний до поверхнi всюди вздовж кривої.

Пiд натуральним параметром s для даної точки кривої розумiють довжину дуги кривої мiж даною точкою i деякою точкою кривої, прийнятою

з початкову.

 

 

 

 

 

Позначимо через θ

кут

мiж кривою на поверхнi та лiнiєю v = const

(з рис. 3.9 видно, що

θ —

це кут мiж векторами

˙

та

ru). Очевидно,

r

 

 

 

#„

 

#„

при перемiщеннi точки вздовж кривої змiнюються не тiльки її криволiнiйнi координати u, v, а може змiнюватись також i кут θ, тобто для кривої на поверхнi величини u, v, θ є функцiями натурального параметру s.

На поверхнi, координатна сiтка якої є ортогональною, геодезична лiнiя

51

описується наступною системою диференцiальних рiвнянь:

 

du

 

cos θ

,

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

= #„

 

 

 

 

 

 

 

 

 

 

 

=

 

#„

 

 

,

 

 

 

 

 

 

 

ds

 

| ru|

 

 

 

 

 

 

 

 

 

 

 

 

|

 

 

v

|

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

dv

 

sin

θ

 

 

 

 

 

 

 

 

 

 

 

 

 

 

r

 

 

 

 

 

 

 

 

 

(3.34)

ds

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

#„

 

 

 

 

#„

 

 

 

 

 

 

 

 

 

 

 

 

 

 

dθ

 

∂| ru|

cos θ

 

∂| rv|

sin θ

 

 

 

 

 

 

 

 

 

 

|

 

|·|

 

|

 

 

 

 

 

=

 

 

 

v

 

 

 

#„

 

#„ u

.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ds

 

 

 

 

 

 

 

 

ru

 

 

rv

 

 

 

єдиний розв’язок

 

 

 

 

 

Система (3.34) має

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

u = u(s),

 

 

 

 

v = v(s),

 

 

θ = θ(s)

при заданих початкових умовах

 

 

 

 

 

 

 

 

 

 

 

 

u(0) = u0,

 

 

 

v(0) = v0,

 

 

θ(0) = θ0.

Першi двi початковi умови визначають криволiнiйнi координати (u0, v0) точки поверхнi, через яку проходить геодезична в початковiй точцi (при s = 0), а остання умова задає напрямок геодезичної в початковiй точцi. Єдинiсть розв’язку означає, що через точку поверхнi з координатами (u0, v0) в заданому напрямку проходить лише одна геодезична лiнiя.

52

4. Кривi на поверхнi елiпсоїда обертання

4.1. Криволiнiйнi координати на поверхнi елiпсоїда

Якщо в формулах (1.9) на стор.8 покласти H = 0, то отримаємо рiвняння

поверхнi елiпсоїда, записане в параметричнiй формi:

 

 

(4.1)

 

 

 

 

 

 

 

Y (B, L) = N cos B sin L,

 

 

 

 

 

 

 

 

 

 

 

X(B, L) = N cos B cos L,

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Z(B, L) = N(1 −e ) sin B,

 

 

 

 

 

 

 

 

 

 

a

 

 

 

 

 

 

 

 

 

 

 

 

 

 

де позначено N =

 

 

 

 

.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

1 −e2 sin2 B

 

 

 

 

 

 

 

 

 

 

 

Звичайно,

(4.1) можна переписати i у виглядi векторної функцiї

 

 

p

r (B, L) = N cos B sin L

.

 

 

(4.2)

 

 

 

 

 

 

 

 

 

 

 

N cos B cos L

 

 

 

 

 

 

 

 

 

 

 

#„

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

N(1 −e2) sin B

 

 

 

 

Результати з попереднiх пунктiв

застосуємо до елiпсоїда. В першу чергу

знайдемо координатнi вектори:

 

 

 

 

 

 

 

 

 

 

 

 

 

#„

 

 

M sin B cos L

 

 

 

 

 

#„

 

N cos B sin L

 

 

#„

r

 

 

M cos B

 

 

#„

 

r

 

 

 

 

(4.3)

rB

B =

2

rL

 

L =

0

 

M sin B sin L ,

 

 

 

N cos B cos L

де введено позначення M =

a(1 −e )

.

 

 

 

 

 

 

 

 

 

 

#„

 

 

 

 

 

(1 −e2 sin2 B)3/2

 

 

 

 

 

 

 

 

#„

Вектор rB

 

направлений по дотичнiй до меридiану, а вектор

rL — по

дотичнiй до паралелi.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

З (4.3) можна знайти, що

2 cos2 B,

 

 

r B ·

r L = 0,

 

(4.4)

 

 

 

 

 

rB2

= M2,

 

rL2 = N

 

 

 

 

 

 

 

 

#„

 

 

 

#„

 

 

 

 

 

 

 

#„

#„

 

 

причому остання рiвнiсть в (4.4) показує, що лiнiї меридiанiв та паралелей утворюють ортогональну координатну сiтку.

У вiдповiдностi з (3.23) та (4.4) запишемо для елiпсоїда першу квадра-

тичну форму

+ rL2 dL2

= M2 dB2 + N2 cos2 B dL2.

(4.5)

ds2 = rB2 dB2

#„

#„

 

 

4.2.Радiус кривизни меридiану, паралелi, нормального перерiзу на елiпсоїдi

Для знаходження радiуса кривизни меридiану покладемо L = const в

формулах (4.1). Тодi цi вирази можна розглядати як рiвняння, що описують

53

меридiан як задану в параметричнiй формi просторову криву:

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

X(B) =

 

 

a cos B cos L

,

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

1 −e2 sin2 B

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

a cos B sin L

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

p

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Y (B) =

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

,

 

 

 

 

 

 

 

 

 

(4.6)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

1

 

e

sin

 

B

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

a(1

e2) sin B

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Z(B) =

p

 

 

 

 

 

 

 

 

 

.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

1

 

e

2

sin

2

B

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

знайдемо за формулою (3.8), пiдставивши

Радiус кривизни

меридiану

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

p

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

в неї (4.6):

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

v

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

dX

 

 

2

 

 

 

 

dY

2

 

 

 

dZ

2

 

3

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

+

 

 

 

+

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

dB

 

 

 

 

dB

 

 

dB

 

 

 

 

 

 

 

 

 

 

 

 

Rм =u

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

. (4.7)

 

 

2

 

 

2

 

 

 

 

 

 

 

2

 

 

 

 

 

 

 

2

 

 

 

 

 

 

2

 

2

 

 

u

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

2

 

 

u

 

 

dX d Y dY d X

 

+

 

dY d Z dZ d Y

 

 

+

 

dZ d X dX d Z

 

 

 

u

 

 

dB dB2

dB dB2

 

 

dB dB2

dB dB2

 

 

 

dB dB2

dB dB2

 

 

 

u

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

t

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Обчислення похiдних та спрощення отриманого виразу можна автоматизувати, якщо скористатись можливостями символьної математики системи MATHCAD, як це показано в наступному прикладi.

!"#"$%&'" ()'*"#+,- .'-,,-

 

 

a ! 0

L

! 0

e2 ! 0

B

! 0

a ! a

L

! L

e2 ! e2

B

! B

!"#$! %&' $#&()*+,-(. N(B), X(B), Y(B), Z(B). /#01,-2 3!"#$4 e2 54%(0+ 3!6+"!,-+343#-! ,!03+&2*4 $01**4 e2

a

N(B) !"

1 e2sin(B)2

X!" N(B)#cos(B)#cos(L)

Y!" N(B)#cos(B)#sin(L)

Z!" N(B)#(1 e2)#sin(B)

!"#$%&'$ (&)"*& %+, -.)/ 0" %)12 -$#3%! 0" *"!$4&'$ 5# ( *'3!!3

dX !

d

X

 

dY !

d

Y

 

dZ !

d

Z

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

dB

 

 

 

 

dB

 

 

 

dB

 

 

 

d2X !

d2

 

X

d2Y !

d2

 

Y

d2Z !

d2

 

Z

 

 

dB

2

 

 

 

dB

2

 

 

 

dB

2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

54

1 −e2 sin2 B
a

!"#$%&$ "'()*+,&&-: ".!,#" "'()*+/0#" 12 *3!"450#" ()*,+6&). 12 7&2#,&&). 8 3$9."!,&,8"#5 8)!27$ :"!#5+) 9+- "'()*+,&&- !29$5*5 .!)8)7&) #,!)9$2&5. ;2*1"*582&&-

.+/("8)< *+$8 factor (!"7.+292&&- &2 #&"%&).)) 12 simplify (*3!"4,&&-) 83+)82/16 &2 !,75+6121 *)#8"+6&)< "'()*+,&6; &,"'<$9&$*16 12 3"!-9". =< 8)."!)*12&&- 8)7&2(2/16 ,.*3,!)#,&12+6&" 8 <"9$ !"'"1).

p1 &' !dX2

 

 

dZ2"

3

 

factor

6

 

 

$ 1)

6

 

 

 

 

 

 

 

 

 

dY2

 

 

 

%

 

512#a #(e2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

2)9

 

 

 

 

 

 

 

 

 

 

 

 

 

 

simplify

(e2#cos(2#B) $ e2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

4

 

 

4

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

simplify

 

 

 

p2 &' (dX#d2Y $ d2X#dY)2 (((

 

 

 

 

 

 

64#a #(e2 $ 1)

 

 

 

 

 

 

 

 

 

 

 

%

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

2)6

 

 

 

 

 

 

 

 

2

 

 

 

2

 

 

factor

(e2#cos(2#B) $ e2

 

(dY#d2Z $ d2Y#dZ)

(dZ#d2X $ d2Z#dX)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

!"#"$%&$ '($)*+$ ,*-#' )./ -#)01!1 2-*,*'&* +3-*)0#. 4-* !5-$63&&0 ,*2$-*!"$,17+$

+$)*802#"$- assume, /2*9 )$',$./7 ,2#'#"* $:.#!"; )$51!"*+*( '%3&; '+0&&*( (<3

&3$:(0)&$ )./ 2$-32"&$=$ !5-$63&&/ !*+,$.;&*( ,*-#'0,).

 

 

 

 

 

 

 

 

 

 

 

 

simplify

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

R+ )*

 

p1

 

assume!"a

 

0!"e2 0!"e2 # 1!"B $

π

( $

 

2% 2%a%(e2 $ 1)

 

 

 

p2

 

 

 

 

 

3

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

rewrite!"sin

 

 

 

 

 

 

 

 

&2 $ 2%e2%sin(B)2'2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Виданий MATHCAD результат вже неважко спростити вручну:

 

 

 

 

 

 

 

 

 

 

 

Rм =

 

a(1 −e2)

 

 

,

 

 

 

 

 

(4.8)

 

 

 

 

 

 

 

 

 

 

(1 −e2 sin2 B)3/2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

що спiвпадає з введеною ранiше величиною M.

Аналогiчно, покладаючи B = const в (4.1) i розглядаючи X, Y , Z як функцiї скалярного параметру L, знайдемо радiус кривизни паралелi:

a cos B

 

Rп = N cos B = p1 −e2 sin2 B .

(4.9)

На рис. 4.1 та 4.2 показано залежнiсть радiусу кривизни меридiана та паралелi вiд геодезичної широти для Землi.

Через деяку точку Q на елiпсоїдi проведемо нормаль до поверхнi елiп-

соїда (вектор n на рис. 4.3). Площина, яка проходить через вектор n (так

#„

#„

звана нормальна площина) перетинає поверхню елiпсоїда вздовж лiнiї, яка називається нормальним перерiзом. Проводячи рiзнi площини через нормаль до поверхнi елiпсоїда, можна отримати безлiч нормальних перерiзiв. На рис. 4.3 зображено лише один з них — той, який утворений за допомогою нормальної площини, дотичної до паралелi, тобто розташованої пiд кутом 90до меридiану. Цей перерiз називається першим вертикалом. Введена ранiше величина N = p — це радiус кривизни першого

вертикалу, який проходить через точку з широтою B.

55

В загальному випадку нормальна площина може перетинатися з меридiаном пiд довiльним кутом A. В цьому випадку радiус кривизни нормального перерiзу буде виражатись формулою

RA =

 

 

 

 

1

 

 

 

.

(4.10)

 

 

 

 

 

 

 

 

2

 

 

2

 

 

 

 

cos

A

+

sin

A

 

 

 

M

 

 

N

 

 

 

 

 

 

 

 

 

 

 

 

 

Очевидно, пiдставляючи в (4.10) A= 0 та A= π/2, одержимо вiдповiдно радiуси кривизни меридiану та першого вертикалу.

При будь-якому розташуваннi нормального перерiзу, проведеного через точку Q, вiн буде залишатись плоскою кривою, а тому вектор головної нормалi до нього, побудований в точцi Q буде лежати в нормальнiй площинi, а значить, буде перпендикулярний до поверхнi елiпсоїда. Звiдси випливає, що в точцi Q геодезична кривизна нормального перерiзу дорiвнює 0.

Середнiм радiусом кривизни Rсер поверхнi в точцi називають границю середнього арифметичного значення радiусiв кривизни RA всiх можливих нормальних перерiзiв, проведених через дану точку поверхнi в рiзних напрямках:

 

 

 

 

 

 

 

= NM =

 

 

 

2

 

 

 

 

R

=

1

R dA=

1

 

 

dA

1

a

 

1 −e

 

.

 

(4.11)

сер

 

Z

 

A

 

Z

cos2 A + sin2 A

 

 

e2 sin2

B

 

 

 

 

 

0

 

 

 

 

0

M

 

N

 

 

 

 

 

 

 

 

 

 

Повною або гаусовою кривизною поверхнi називається добуток кри-

визн головних нормальних перерiзiв поверхнi в данiй точцi

 

 

 

 

 

 

 

 

 

 

 

 

 

 

K =

1 .

 

 

 

 

 

 

 

 

 

 

(4.12)

 

 

 

 

 

 

 

 

 

NM

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

0

 

 

 

 

 

 

 

 

 

 

 

 

0O

10O

20O

30O

40O

50O

60O

70O

80O

90O

0O

10O

20O

30O

 

40O

50O

 

60O

70O

80O

90O

Рис. 4.1. Залежнiсть радiусу кривизни

Рис. 4.2. Залежнiсть радiусу кривизни

меридiана вiд геодезичної широти B

 

паралелi вiд геодезичної широти B

 

 

 

 

 

 

 

 

 

 

56

 

 

 

 

 

 

 

 

 

 

 

 

 

4.3. Довжина дуги меридiана та паралелi

Оскiльки паралель являє собою коло, то довжина дуги паралелi отримується як добуток радiуса даної паралелi Rп = N cos B на рiзницю довгот L2 L1 крайнiх точок дуги:

 

 

 

a cos B

 

sп = (L2 L1) N cos B = (L2 L1)

 

 

.

(4.13)

 

 

 

1 −e2 sin2 B

Меридiан же на елiпсоїдi має форму

елiпса, а тому обчислення довжи-

 

p

 

 

 

ни дуги меридiану дещо складнiше. Нехай точка A1 на меридiональному перерiзi елiпсоїда (рис. 4.4) має геодезичну широту B. Вiзьмемо на нескiнченно малiй вiдстанi вiд точки A1 точку A2, яка має широту B + dB. Тодi рiзниця широт точок A1 i A2, якiй вiдповiдає дуга меридiана ds, буде dB. Розглядаючи елементарну дугу ds як дугу кола з радiусом M, будемо мати:

ds = M dB

де M — радiус кривизни меридiана в точцi з широтою B.

Довжина дуги меридiана на елiпсоїдi мiж точками з широтами B1 та B2 визначається iнтегралом

B2

B2

 

 

 

 

 

 

sм = Z

M dB = a(1 −e2) Z

 

 

dB

 

,

(4.14)

 

 

 

 

(1

e2 sin2

B)3/2

B1

B1

 

 

 

 

 

який не обчислюється в елементарних функцiях. Розкладемо пiдiнтеграль-

Z

#„

n

Q

Y

B

X

Рис. 4.3. До поняття першого вертикалу

57

ний вираз в (4.14) в ряд за степенями e

 

 

 

 

 

sм = a(1 −e2) ZB2

1 +

3

e2 sin2

B +

15

e4 sin4 B +

35

e6 sin6 B + . . . dB

4

8

16

B1

 

h

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

i

та виконаємо почленне iнтегрування:

(B2 B1)−

sм = a(1 −e2)"

1 + 43 e2 + 6445 e4

+ 256175 e6 + . . .

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

3

e2 +

15

e4 +

525

e6 + . . .

(sin 2B

sin 2B )+

 

 

 

 

8

 

 

 

 

 

32

 

 

 

 

1024

 

 

 

2

1

 

15

 

 

 

 

105

 

 

 

 

 

(4.15)

+

 

e4 +

 

e6 + . . . (sin 4B2 −sin 4B1)−

256

1024

35

e6 + . . . (sin 6B2 −sin 6B1) + . . .#

3072

Файл GeodeticLib.xmcdz мiстить функцiю MeridianArcLength, призначену для обчислення довжини дуги меридiана на елiпсоїдi. Формат її виклику наступний:

sм := MeridianArcLength(B1, B2, a, e)

Тут B1, B2 — геодезичнi широти початкової та кiнцевої точок дуги (в радiанах); a, e — велика пiввiсь та ексцентриситет елiпсоїда. Довжина дуги повертається в метрах. Обчислення реалiзовано за формулою типу (4.15), в

Q

2 Q1

M dB

B

Рис. 4.4. До знаходження довжини дуги меридiана

58

якiй взято дещо бiльшу кiлькiсть доданкiв (врахованi степенi e8, e10, e12).

На практицi можна користуватися наступною наближеною формулою:

sм = Mm(B2 B1)

1 + 8 (B2 B1)2 cos 2Bm

,

(4.16)

 

 

e2

 

 

похибка якої складає не бiльше 0,2 мм при довжинi дуги S < 100 км. В (4.16) Mm означає радiус кривизни меридiана, обчислений в точцi з широ-

тою B

=

1

(B + B ) (B

i B вимiрюються в радiанах, B > B ).

 

m

2

1

2

1

2

2

1

При довжинi дуги S 6 40 км в (4.16) можна вiдкинути доданок, який

мiстить e2, i скористатись формулою

 

 

 

 

 

 

 

 

sм = Mm(B2 B1),

 

(4.17)

похибка якої не перевищує 1,5 мм.

Використовуючи формулу Сiмпсона чисельного iнтегрування, отримає-

мо iншу наближену формулу для обчислення Sm:

 

sм =

1

(B2 B1)(M1 + 4Mm + M2),

(4.18)

6

де M1, M2, Mm — радiуси кривизни, обчисленi в точках

з широтами

B1, B2, Bm вiдповiдно. Похибка формули (4.18) не перевищує 0,2 мм при довжинi дуги до 400 км.

Якщо довжина дуги меридiана перевищує 400–500 км та при використаннi формул (4.16), (4.18) необхiдно зберегти вказану точнiсть обчислення, то дугу слiд роздiлити на частини i до кожної частини окремо застосувати формулу (4.16) або (4.18).

Можна наближено розв’язати i обернену задачу: за вiдомою довжиною дуги меридiана знайти широту (передбачається, що дуга починається на

екваторi). По сутi, ця задача зводиться до розв’язання рiвняння

Z B

sм = M dB

0

вiдносно B. Будемо шукати розв’язок у виглядi ряду

 

B = C +C e +C e2

+C e3

+ . . . ,

(4.19)

0

1

2

3

 

 

де C0, C2, C2, C3, . . . — невiдомi коефiцiєнти. Пiдставимо B1 = 0 та B2 =(4.19)

в (4.15):

1

sм

2 4 0

4

0

0

3

 

2 1

2 1

0

 

 

a

0 +

 

 

sм

=C

C e+

C

 

1

C

3

sin C

cos C

e2+ C

+

1

C

3

C

cos2 C

e3 +. . .

 

 

 

 

 

 

 

Перенесемо

 

 

в

праву частину та прирiвняємо 0 вирази, що стоять при

a

 

59

степенях e. Отримаємо систему рiвнянь

 

 

C1

= 0,

 

 

 

 

 

 

 

 

C0

sм

= 0,

 

 

 

 

a

 

 

 

 

 

1

 

 

 

3

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

C

 

 

 

 

C

 

 

 

 

sin C cos C = 0,

 

 

4

 

4

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

2

 

 

 

 

0

 

 

 

 

0

0

 

C3 +

 

C1

 

C1 cos C0 = 0,

 

2

2

 

. . .

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

1

 

 

 

 

3

 

 

2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

з якої послiдовно

знайдемо C , C , C , C , . . . i остаточно отримаємо формулу,

 

 

0 2

 

2

 

3

 

 

 

 

 

 

 

яка наближено розв’язує поставлену задачу:

 

3sc + ξ

e2

+

21

sc3 +

3

c2ξ

3

 

 

5

ξ e42+

 

 

 

B = ξ +

 

 

 

 

 

+

 

 

sc

 

 

 

 

 

4

32

8

64

64

9

6

(4.20)

 

 

151 5

 

21 4

 

 

25

 

3

 

39 2

 

21 −72ξ

 

 

+

 

sc

+

 

c

 

ξ −

 

sc

 

 

c

ξ +

 

 

sc

 

ξ e

+ . . .

192

32

 

192

128

 

768

256

де позначено ξ = saм , s = sin ξ, c = cos ξ.

Файл GeodeticLib.xmcdz мiстить функцiю B_via_ArcLength, яка за вiдомою довжиною s дуги меридiана, що вiдкладається вiд екватора, на основi формули типу (4.20) повертає вiдповiдне значення геодезичної широти B. Формат її виклику наступний:

B:= B_via_ArcLength(s, a, e)

4.3.1.Приклад обчислення довжини дуги меридiана

Для елiпсоїда WGS–84 обчислити довжину дуги меридiана мiж точками з широтами

B = 453017,221′′,

B = 492958,938′′.

1

2

Розв’язання

Переводимо широти з градусної мiри в радiанну:

B1 = 0,7942082994, B2 = 0,8639328310, Bm = B1 + B2 = 0,8290705652. 2

Знаходимо радiуси кривизни меридiана в точках з широтами B1, B2, Bm:

a(1 −e2)

M1 = (1 −e2 sin2 B1)3/2 = 6 367 947,027 м

a(1 −e2)

M2 = (1 −e2 sin2 B2)3/2 = 6 372 402,675 м

a(1 −e2)

Mm = (1 −e2 sin2 Bm)3/2 = 6 370 181,006 м

60