GPS卫星位置的计算(C程序计算)

1、GPS卫星位置的计算利用c卄编写了一段能计算单一瞬时卫星坐标的程序，在运行程序之前，需做部分准备工作：在F盘下建立一名为“单一卫星广播星历”的txt文件。从“广播星历.txt”文件中拷贝从卫星PRN号开始的8行数据到“单一卫星广播星历.txt”中在编辑选项中，将全部的“D”替换为“E”。下面为我所选取的一个广播星历：1806825600.0-2.472363412380E-04-1.023181539495E-120.000000000000E+001.410000000000E+02-1.721875000000E+014.502687555010E-091.413760604187E+00

2、-7.990747690201E-077.598234573379E-031.U8145883083E-055.153709835052E+034.536000000000E+05-1.303851604462E-08-1.095067942661E-011.527369022369E-079.571235745530E-011.640000000000E+02-2.656176299285E+00-&037477650349E-09-5.1930734552UE-101.000000000000E+001.389000000000E+030.000000000000E+002.0000000

3、00000E+000.000000000000E+00-1024454832077E-081.410000000000E+024.464490000000E+054.000000000000E+00程序设计部分：includemcludenitmam()mt1=0;doublen50?nO,nn,t,tk,Mk,Ek,Yk,Gu,Gi;Gi,uk,rk,ik,xk,yk,zk,X,Y,Z,Lk,UT,yy,rnm,JD,gpsz;FILE*fp;fp=fopen(F:单一卫星广播星历.txt,r);if(fp=NULL)pinitf(文件打开失败!n”)；return0;while(!feof

4、(fp)fscanf,H%lf&ni);i+;110=(sqrt(3986005E+8)/pow(nl7,3);mi=n0+n12;/*计算卫星运行的平均角速度*/UT=n4+(n5/60)+(n7)/3600)尸民用日的时分秒化为实数时*/if(nl=80)/*广播星历中年只有后两位，化为4位，参考1980年1月6口0点*/if(nl=80&n2=1&n36)nl=nl+2000;nl=nl+1900;elsenl=nl+2000;if(n22)yy=ni；mni=n2;JD=(int)(365.25*yy)+(int)(30.6001*(mill+1)+n3+(UT/24)+1720981

5、.5;/*化为儒略日*/gpsz=(int)(JD-2444244.5)/7);/*计算GPS周勺t=(JD-2444244.5-7*gpsz)*24*3600;/*得出GPS秒*7tk=t-n18;/*tk1为中间值，用以判断tk与正负302400的关系，然后返回到tk*/while(tk302400|tk302400)tk=tk-604800;elsetk=tk+604800;/*计算归化观测时间*/Mk=n13+mi*tk;/*观测时刻的卫星平近点角*/Ek=Mk;Ek=Mk+n15*sin(Ek);Ek=NIk+n15*sin(Ek);/*迭代两次计算观测时刻的偏近点角*/Vk=ata

6、n(sqrt(l-n15*n15)*sin(Ek)/(cos(Ek)-n15)y*真近点角*/Yk=Vk+n24;/*升交距角匕Gu=n14*cos(2*Yk)+n16*siii(2*Yk);Gr=n23*cos(2*Yk)+nll*sin(2*Yk);Gi=n19*cos(2*Yk)+n21*sm(2*Yk);/*摄动改止项日uk=Yk+Gu;rk=n17n17*(1-n15*cos(Ek)+Gr;ik=n22+Gi+n26*tk;/*经摄动改止后的升交距角、卫星欠径、轨道倾角*/xk=rk*cos(uk);yk=rk*sin(uk);zk=0;/*卫星在轨道坐标系的坐标*/Lk=n20+(n25-7.29211515E-5)*tk-7.29211515E-5*n18;/*观测时刻t的升交点经度*/X=xk*cos(Lk)-yk*cos(ik)*sm(Lk);Y=xk*siii(Lk)+yk*cos(ik)*cos(Lk);Z=yk*sm(ik);/*卫星在WGS-84坐标系的坐标日pnntf(n该卫星在WGS-84坐标系中的坐标为:nX=%lfmaiY=%lfmiiZ=