美国GPS卫星导航系统的GPS 卫星位置的计算步骤

1、计算卫星运行的平均速度n

2、时间tk计算

3．计算观测瞬间的卫星平近点角

4．计算偏近点角Ek

5、计算真近点角Vk

6．计算升交距角φk

7、计算摄动改正项δ u、δ r、δ i

8、计算经过摄动改正的升交距角uk、卫星矢径rk、和轨道倾角Ik

9．计算卫星在轨道平面上的位置

10、计算观测时刻的升交点经度f

11、计算卫星在地心坐标系中空间直角坐标

12.卫星在协议地球坐标系中的坐标计算

{

//参数赋值

double ss1 = 0.00000001, ss2 = 1000000,we=7.2921151467*ss1*1000

double pi = 3.1415926535

string[] str = textBox2.Text.Trim().Split(' ')

double[,] a = new double[7, 4]

int k = 0

for (int i = 0i <7i++)

{

for (int j = 0j <4j++)

{

a[i, j] = Convert.ToDouble(str[k])

k++

}

}

//规划时刻计算

int[] dinmth = new int[12]

double year, month, day, hour, minute, second, dayofy, ttlday, dayofw, gpstime, tk

int weekno

dinmth[0] = 31dinmth[1] = 28dinmth[2] = 31

dinmth[3] = 30dinmth[4] = 31dinmth[5] = 30

dinmth[6] = 31dinmth[7] = 31dinmth[8] = 30

dinmth[9] = 31dinmth[10] = 30dinmth[11] = 31

year = Convert.ToDouble(textBox21.Text)

month = Convert.ToDouble(textBox22.Text)

day = Convert.ToDouble(textBox23.Text)

hour = Convert.ToDouble(textBox24.Text)

minute = Convert.ToDouble(textBox25.Text)

second = Convert.ToDouble(textBox26.Text)

if (month == 1)

{

dayofy = day

}

else

{

dayofy = 0

for (int i = 0i <12i++)

{

dayofy += dinmth[i]

if (i == 1)

{

if (year % 4 == 0 &&year % 100 != 0 || year % 400 == 0)

dayofy += 1

}

}

dayofy += day

}

ttlday = 360

for (int yr = 1981yr <yearyr++)

{

ttlday += 365

if (yr % 4 == 0 &&yr % 100 != 0 || yr % 400 == 0)

{

ttlday += 1

}

}

ttlday += dayofy

weekno = (int)(ttlday / 7)

dayofw = ttlday - 7 * weekno

gpstime = hour * 3600 + minute * 60 + second + dayofw * 86400

tk = gpstime - a[3, 0] * ss2

if (tk >302400)

{

tk -= 604800

}

else if(tk<-302400)

{

tk += 602400

}

app.value1 = tk.ToString("e2")

//平均角速度计算

double n0, n

n0 = (Math.Sqrt(3986004.418) / Math.Pow(a[2, 3], 3)) * ss1

n = n0 + a[1, 2] * ss1

app.value2 = n.ToString("e2")

//平近点角计算

double M0, Mk

M0 = a[1, 3]

Mk = M0 * 10 + n * tk

if (Math.Abs(Mk) >2 * pi)

{

int w1 = (int)(Mk / (2 * pi))

Mk = Mk - 2 * pi * w1

}

app.value3 = Mk.ToString("e2")

//偏近点角计算

double Ek, E0

Ek = Mk

do

{

E0 = Ek

Ek = Mk + a[2, 1] * 0.1 * Math.Sin(E0)

} while (Math.Abs(Ek - E0) >0.00000001)

app.value4 = Ek.ToString("e2")

//真近点角

double angle1, angle2, Vk

angle1 = Math.Sqrt(1 - Math.Pow(a[2, 1] * 0.1, 2)) * Math.Sin(Ek)

angle2 = Math.Cos(Ek) - a[2, 1] * 0.1

Vk = Math.Atan(Math.Abs(angle1 / angle2))

if (angle1 >= 0 &&angle2 <= 0)

{

Vk = pi - Vk

}

else if (angle1 <= 0 &&angle2 <=0)

{

Vk = pi + Vk

}

else

{

Vk = 2 * pi - Vk

}

app.value5 = Vk.ToString("e")

//升交角距

double Qk

Qk = Vk + a[4, 2] * 10

if (Math.Abs(Qk) >2 * pi)

{

int w2 = (int)(Qk / (2 * pi))

Qk = Qk - 2 * pi * w2

}

app.value6 = Qk.ToString("e")

//升交角距改正,轨道向径改正，轨道倾角改正计算

double Fu, Fr, Fi

Fu = a[2, 0] * ss1 * 1000 * Math.Cos(2 * Qk) + a[2, 2] * ss1 * 1000 * Math.Sin(2 * Qk)

Fr = a[4, 1] * 1000 * Math.Cos(2 * Qk) + a[1, 1] * 100 * Math.Sin(2 * Qk)

Fi = a[3, 1] * ss1 * 100 * Math.Cos(2 * Qk) + a[3, 3] * ss1 * 10 * Math.Sin(2 * Qk)

app.value7 = Fu.ToString("e")

app.value8 = Fr.ToString("e")

app.value9 = Fi.ToString("e")

//改正后升交角距、轨道向径、轨道倾角计算

double Uk, Rk, Ik

Uk = Qk + Fu

Rk = Math.Pow(a[2, 3] * ss2 * 0.01, 2) * (1 - a[2, 1] * 0.1 * Math.Cos(Ek)) + Fr

Ik = a[4, 0] + Fi + a[5, 0] * ss1 * 0.1 * tk

app.value10 = Uk.ToString("e")

app.value11 = Rk.ToString("e")

app.value12 = Ik.ToString("e")

//升交点轨道直角坐标系坐标计算

double Xk, Yk

Xk = Rk * Math.Cos(Uk)

Yk = Rk * Math.Sin(Uk)

app.value13 = Xk.ToString("e")

app.value14 = Yk.ToString("e")

//升交点经度计算

double N

int w3

N = a[3, 2] * 10 + (a[4, 3] * ss1 - we) * tk - we * a[3, 0] * ss2

if (Math.Abs(N) >2 * pi)

{

w3 =(int)(N / (2 * pi))

N = N - 2 * pi * w3

}

app.value15 =N.ToString("e")

//地固坐标系下卫星空间直角坐标计算,Ik为ik

double X, Y, Z

X = Xk * Math.Cos(N) - Yk * Math.Cos(Ik) * Math.Sin(N)

Y = Xk * Math.Sin(N) + Yk * Math.Cos(Ik) * Math.Cos(N)

Z = Yk * Math.Sin(Ik)

app.value16 = X.ToString("e")

app.value17 = Y.ToString("e")

app.value18 = Z.ToString("e")

new Form2().ShowDialog()

}

catch

{

MessageBox.Show("请输入卫星时刻或又一次加载卫星参数导致出错")

return

}