Вычисление элементов орбиты ИСЗ по координатам и составляющим скорости

Вычисляем постоянные площадей
с1=y·z-y∙z
с1=5438458.535·402.86--754.84·11490078.419=10864105220.881м2с
с2=z·x-z∙x
с2=11490078.419·-5362.74-402.86·-11901.617=-61613556977.293м2с
с3=x·y-x∙y
с3=-11901.617·-754.84--5362.74·5438458.535=29174045904.430м2с
c=с12+с22+с32
c=10864105220.8812+-61613556977.2932+29174045904.4302=69031761820.761м2с
2. Вычисляем фокальный параметр
p=c2μ
p=69031761820.76123.9860044∙1014=11955290.717 м
3. Вычисляем угол наклона плоскости орбиты к плоскости экватора
i=arccosс3c
i=arccos29174045904.430 69031761820.761=65º.00΄00,1''
4. Вычисляем долготу восходящего узла Ω
sinΩ=c1c∙sini
sinΩ=11510864105220.881 69031761820.761·sin65º.00΄00,1''=0.173647755
cosΩ=c2c∙sini
cosΩ=-61613556977.293 69031761820.761·sin65º.00΄00,1''=0.984807828
Ω=9º.59΄59,9''
5 . Вычисляем квадрат скорости спутника и геоцентрическое расстояние до спутника
v2=x2+y2+z2
v2=-5362.742+-754.842+5402.862=29491104.728м2с2
r=x2+y2+z2
r=-11901.6172+3905438458.5352+11490078.41992=12712154.615 м
6.Вычисляем постоянную энергии
h=v2-2∙μr
h=29491104.728 -2∙3.9860044∙101412712154.615=-33220599.475м2с2
7. Вычисляем большую полуось орбиты спутника
a=-μh
a=-3.9860044∙1014-33220599.475=11998592.629 м
8. Вычисляем составляющие вектора Лапласа
D=x∙x+y∙y+z∙z
D=-11901.617·-5362.74+5438458.535·-754.84+11490078.419·402.86=
D=587541594.078 м
D=μr+h
D =3.9860044∙1014112712154.615+-33220599.475=-1864747.374
f1=D∙x-D∙x
f1=-1864747.374·-11901.617-587541594.078·-5362.74=
f1=3173028797443.880
f2=D∙y-D∙y
f2=-1864747.374·5438458.535-587541594.078·-754.84=
f2=-9697851286389.240
f3=D∙z-D∙z
f3=-1864747.374·11490078.419-587541594.078·402.86=
f3=-21662790056130.600
f=f12+f22+f32
f=3173028797443.8802+-9697851286389.2402+-21662790056130.6002=23945623907894.300 м
9