急求GPS卫星坐标计算C#程序!

急求GPS卫星坐标计算C#程序!
急求GPS卫星坐标计算C#程序!

急求GPS卫星坐标计算C#程序!
try
{
//参数赋值
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 = 0; i < 7; i++)
{
for (int j = 0; j < 4; j++)
{
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] = 31; dinmth[1] = 28; dinmth[2] = 31;
dinmth[3] = 30; dinmth[4] = 31; dinmth[5] = 30;
dinmth[6] = 31; dinmth[7] = 31; dinmth[8] = 30;
dinmth[9] = 31; dinmth[10] = 30; dinmth[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 = 0; i < 12; i++)
{
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 = 1981; yr < year; yr++)
{
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 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