C语言计算GPS卫星位置

“GPS广播星历计算卫星位置和速度”及“GPS伪距定位”计算试验1．试验内容及上交成果1.1 试验内容应用C语言按预定格式（函数、输入输出变量之名称、类型）编写“GPS广播星历计算卫星位置和速度”函数SatPos_Vel( )、“GPS伪距定位”函数Positioning( )。

将此两个函数组成文件F2.cpp，并包含于文件GPS_Positioning.cpp中。

编译、连接并运行文件GPS_Positioning.cpp，逐一时刻读取广播星历（Ephemeris.dat）、观测时间及伪距、卫星号（Observation.dat）信息，计算WGS84坐标系中观测时刻相应的卫星位置、速度以及载体位置，结果保存于文件Position.dat中。

1.2 上交成果磁盘文件F2.cpp、Position.dat，并存于“学号作者中文姓名”目录中。

2．函数说明2.1 星历文件读取函数void EFileReading(Efile)功能：读取星历文件，给星历数据结构体Ephemeris赋值。

输入变量：EFile 字符串，文件名。

2.2 观测数据读取函数int ObsReading(fp_Obs,Time,Rho,Mark)功能：从文件Observation.dat中读取某一时刻的伪距、卫星号。

读取成功函数值返回“1”，失败返回“－1”（读错，或至文件尾）。

输入变量：fp_Obs 文件指针；输出变量：Time double，时间（秒）；Rho double[12]，伪距（米）；Mark int[12]，卫星号，“－1”表示此通道无卫星、无伪距。

2.3 最小二乘估计函数int LeastSquareEstimation(Y,A,P,m,n, X)功能：最小二乘方法求解观测方程Y=AX+ε，其中观测值方差阵的逆阵为P（也称为权阵），得未知参数X。

成功返回“1”，失败返回“－1”（亏秩）。

输入变量：Y double[m]，观测方程自由项（米）；A double[m×n]，系数阵（无量纲），按第1行第1、2……n元素，第2行……顺序存放；P double[m]，权矩阵对角线元素，0或1（无量纲）；m int，观测值个数；n int，未知参数个数；输出变量：X double[n]，未知参数（米）。

坐标转换源代码--GPS定位程序(C++)GPS数据处理中为了满足不同的需要，处理的数据要进行坐标转换，得到在不同坐标系统下的结果，下面是笛卡尔坐标系，大地坐标系，站心地平坐标系（线型和极坐标形式）之间的转换源代码：头文件：#ifndef _COORDCOVERT_H#define _COORDCOVERT_H#include "stdlib.h"//WGS-84椭球体参数const double a=6378137.0;//长半轴const double flattening=1/298.257223563;//扁率const double delta=0.0000001;typedef struct tagCRDCARTESIAN{double x;double y;double z;}CRDCARTESIAN;typedef CRDCARTESIAN *PCRDCARTESIAN;//笛卡尔坐标系typedef struct tagCRDGEODETIC{double longitude;double latitude;double height;}CRDGEODETIC;typedef CRDGEODETIC *PCRDGEODETIC;//大地坐标系typedef struct tagCRDTOPOCENTRIC{double northing;double easting;double upping;}CRDTOPOCENTRIC;typedef CRDTOPOCENTRIC *PCRDTOPOCENTRIC; //站心地平坐标系（线坐标形式）typedef struct tagCRDTOPOCENTRICPOLAR{ double range;double azimuth;double elevation;}CRDTOPOCENTRICPOLAR;typedef CRDTOPOCENTRICPOLAR *PCRDTOPOCENTRICPOLAR;//站心地平坐标系（极坐标形式）//由笛卡尔坐标转换为大地坐标void CartesianToGeodetic (PCRDGEODETIC pcg, PCRDCARTESIAN pcc, double dSemiMajorAxis, double dFlattening);//pcg：指向所转换出的大地坐标的指针；//pcc：指向待转换的笛卡尔坐标的指针；//dSemiMajorAxis：参考椭球的长半轴；//dFlattening：参考椭球的扁率。

GPS卫星位置的计算利用C++编写了一段能计算单一瞬时卫星坐标的程序，在运行程序之前，需做部分准备工作：（1）在F盘下建立一名为“单一卫星广播星历”的txt文件。

（2）从“广播星历.txt”文件中拷贝从卫星PRN号开始的8行数据到“单一卫星广播星历.txt”中（3）在编辑选项中，将全部的“D”替换为“E”。

下面为我所选取的一个广播星历：18 06 8 25 6 0 0.0-2.472363412380E-04-1.023*********E-12 0.000000000000E+001.410000000000E+02-1.721875000000E+01 4.502687555010E-09 1.413760604187E+00-7.990747690201E-07 7.598234573379E-03 1.118145883083E-05 5.153709835052E+034.536000000000E+05-1.303851604462E-08-1.095067942661E-01 1.527369022369E-079.571235745530E-01 1.640000000000E+02-2.656176299285E+00-8.0374********E-09-5.193073455211E-10 1.000000000000E+00 1.389000000000E+03 0.000000000000E+002.000000000000E+00 0.000000000000E+00-1.024*********E-08 1.410000000000E+024.464490000000E+05 4.000000000000E+00程序设计部分：#include<stdio.h>#include<math.h>int main(){int i = 0;double n[50], n0, nn, t, tk, Mk, Ek, Vk, Yk, Gu, Gr, Gi, uk, rk, ik, xk, yk, zk, X, Y, Z, Lk, UT, yy, mm, JD, gpsz;FILE *fp;fp = fopen("F:\\单一卫星广播星历.txt", "r");if (fp == NULL){printf ("文件打开失败!\n");return 0;}while (! feof (fp)){fscanf(fp, "%lf", &n[i]);i++;}n0 = (sqrt(3986005E+8))/pow(n[17], 3);nn = n0 + n[12];/*计算卫星运行的平均角速度*/UT = n[4] + (n[5] / 60) + (n[7] / 3600);/*民用日的时分秒化为实数时*/if (n[1] >= 80)/*广播星历中年只有后两位，化为4位，参考1980年1月6日0点*/ {if (n[1] == 80 && n[2] == 1 && n[3] < 6){n[1] = n[1] + 2000;}n[1] = n[1] + 1900;}else{n[1] = n[1] + 2000;}if (n[2] <= 2){yy = n[1] - 1;mm = n[2] + 12;}if (n[2] > 2){yy = n[1];mm = n[2];}JD = (int)(365.25 * yy) + (int)(30.6001 * (mm + 1)) + n[3] + (UT / 24) + 1720981.5;/*化为儒略日*/gpsz = (int)((JD - 2444244.5) / 7);/*计算GPS周*/t = (JD - 2444244.5 - 7 * gpsz) * 24 * 3600;/*得出GPS秒*/tk = t - n[18];/*tk1为中间值，用以判断tk与正负302400的关系，然后返回到tk上*/ while (tk > 302400 || tk < -302400){if (tk > 302400){tk = tk - 604800;}else{tk = tk + 604800;}}/*计算归化观测时间*/Mk = n[13] + nn * tk;/*观测时刻的卫星平近点角*/Ek = Mk;Ek = Mk + n[15] * sin(Ek);Ek = Mk + n[15] * sin(Ek);/*迭代两次计算观测时刻的偏近点角*/Vk = atan(sqrt(1 - n[15] * n[15]) * sin(Ek)) / (cos(Ek) - n[15]);/*真近点角*/Yk = Vk + n[24];/*升交距角*/Gu = n[14] * cos(2 * Yk) + n[16] * sin(2 * Yk);Gr = n[23] * cos(2 * Yk) + n[11] * sin(2 * Yk);Gi = n[19] * cos(2 * Yk) + n[21] * sin(2 * Yk);/*摄动改正项*/uk = Yk + Gu;rk = n[17] * n[17] * (1 - n[15] * cos(Ek)) + Gr;ik = n[22] + Gi + n[26] * tk;/*经摄动改正后的升交距角、卫星矢径、轨道倾角*/xk = rk * cos(uk);yk = rk * sin(uk);zk = 0;/*卫星在轨道坐标系的坐标*/Lk = n[20] + (n[25] - 7.29211515E-5) * tk - 7.29211515E-5 * n[18];/*观测时刻t的升交点经度*/X = xk * cos(Lk) - yk * cos(ik) * sin(Lk);Y = xk * sin(Lk) + yk * cos(ik) * cos(Lk);Z = yk * sin(ik);/*卫星在WGS-84坐标系的坐标*/printf("该卫星在WGS-84坐标系中的坐标为：\nX = %lf m\nY = %lf m\nZ = %lf m\n", X, Y, Z);fclose(fp);return 0;}计算结果：该卫星在WGS-84坐标系中的坐标为：X = 9223153.692525 mY = 24133486.931401 mZ = 6032585.919385 m。

{
if(x==4||x==7) {}
else dat[x]=ch[y++];
}
break;
case 5: //东西经
dat[0]=ch[0];
break;
case 6: //航速
for(x=0;x<6;x++)
while((RMC_Data[i]!=',')&&(RMC_Data[i]!='*')) //RMC_Data[i]==','或者=='*'时结束循环。
{i++;}
if(j+1!=i)
{
for(j=j+1,k=0;j<i;j++)
{
uchar x,y;
send_byte('d');
switch(ch1)
{
case 0: //时间
for(x=0,y=0;x<8;x++)
{
if(x==2||x==5){time[x]=0x3a;}
else time[x]=ch[y++];
break;
case '*':
i=strlen(RMC_Data); //处理完成，软件结束。
break;
default:break;
}
}
// send_byte('e');
}
void GPRMC_Handle(uchar *ch,uchar ch1) //GPRMC信息，与，之间的子信息处理，存储到相应的结构体中元素中

C语言计算GPS卫星位置GPS（全球定位系统）是一种利用卫星定位来确定全球地理位置的技术。

GPS系统由一组卫星和地面接收器组成，可以帮助我们在地球上的任何位置确定自己的准确位置。

在C语言中，我们可以利用一些数学和物理公式来计算GPS卫星位置。

首先，我们需要了解GPS系统的原理。

GPS卫星运行在地球轨道上，同时向地球发送信号。

地面接收器接收到卫星发出的信号，并通过计算信号的距离和时间来确定自身的位置。

计算GPS卫星位置的关键是测量信号的传播时间。

当接收器接收到来自至少4个不同卫星的信号时，我们可以使用以下公式来计算GPS卫星的位置：速度=距离/时间由于信号的传播速度恒定（大约为光速），我们可以知道信号的传播时间等于距离与速度的比。

在C语言中，我们可以使用以下代码计算GPS卫星的位置：```c#include <stdio.h>#include <math.h>typedef structdouble x;double y;double z;} Point3D;Point3D calculateSatellitePosition(double distance, double latitude, double longitude, double altitude)Point3D position;position.x = (altitude + distance) * cos(E) * cos(longitude);position.y = (altitude + distance) * cos(E) * sin(longitude);return position;int maidouble latitude = 37.7749; // 纬度（假设）double longitude = -122.4194; // 经度（假设）double altitude = 0.0; // 海拔（假设）Point3D satellitePosition =calculateSatellitePosition(distance, latitude, longitude, altitude);printf("GPS卫星位置：(x=%.2f, y=%.2f, z=%.2f)\n", satellitePosition.x, satellitePosition.y, satellitePosition.z);return 0;```在这段代码中，我们首先定义了光速的常量。

（ ２ ）命名 空间 ，主 要用来 作为 一种 区别其他
Ｘ Ｍ Ｌ文 件 标 签 的 标 志 ； （ ３ ）元 数 据 信 息 ，描 述 该 文 档 所 创 建 的时 间 以 及航点坐标范围：
ｄ ａ ｔ ａ ， ｉ ｍｐ ｌ ｅ ｍｅ ｎ ｔ ｓ ｔ ｈ ｅ ｗａ ｙ ｐ ｏ ｉ ｎ ｔ ｄ ａ ｔ ａ ｂ ａ ｔ ｃ ｈ ｉ ｎ ｐ ｕ ｔ ， ａ ｎ ｄ ｉ ｍｐ ｒ ｏ ｖ ｅ ｓ
ｈ
图 １包含航点 的最小 Ｇ Ｐ Ｘ文件
在这 五个部分 中，Ｘ Ｍ Ｌ 序言 、命 名空间以及元
数据 信 息 可 以 不 用 改 动 ， 直 接 参 照模 板 也 可 以 使 ． Ｅ ｆ ｊ 。
在航 点信息部分 ，这里 只关注航 点部 分，忽略掉航 迹 、航线等 。一 个航 点信 息包括 经纬度 （ Ｗ Ｇ Ｓ 一 ８ ４坐 标） 、航 点名称 、航 点符号、创 建时 间等，在文件末 尾 附上 结 束 符 （ ／ ｇ ｐ ｘ ） 表示文件结束 。
Ｓ一 ８ ４ Ｃ Ｏ Ｏｒ ｄｉ ｎａ ｔ ｅ
等信 息 ；
（ ５ ）结束符，文件的结束标志 。
引 言
目前，手 持 Ｇ Ｐ Ｓ接收机在测绘 、物探 等勘 察工 作 中得 到了广泛的应用 ，已经成为 了野外勘察工 作 的常用工 具。手持 Ｇ Ｐ Ｓ轻巧方便 ，其定位 精度达 到
２ ～５ ｍ 的水 平 ，在 定 位 、导 航 等 方 面 深 受 野 外 勘 察
为野 外勘察 人员的常用工具。本文通过 Ｃ语言编程来处理
Ｗ Ｇ Ｓ ８ ４坐标航 点数据 ，实现了航 点数据 的批 量输 入，提高
了手 持 Ｇ Ｐ Ｓ的使用效率。
式 ）的缩写。Ｇ Ｐ Ｘ交换格 式文件 是 Ｇ Ｐ Ｓ设备 中的地 图 信 息 格 式 ，主 要 用 来 描 写 和 保 存地 理 数 据 信 息 。

c语言 gps课程设计一、教学目标本课程的目标是让学生掌握C语言在GPS领域的应用。

通过本课程的学习，学生将能够：1.理解GPS系统的基本原理和工作流程。

2.掌握C语言的基本语法和编程技巧。

3.能够使用C语言进行GPS数据的采集、解析和处理。

4.能够利用C语言实现简单的GPS导航功能。

二、教学内容本课程的教学内容主要包括以下几个部分：1.GPS系统的基本原理和工作流程。

2.C语言的基本语法和编程技巧。

3.GPS数据的采集、解析和处理。

4.GPS导航功能的实现。

三、教学方法为了达到本课程的教学目标，我们将采用多种教学方法，包括：1.讲授法：用于讲解GPS系统的基本原理和工作流程，以及C语言的基本语法和编程技巧。

2.案例分析法：通过分析具体的GPS应用案例，让学生了解GPS技术的实际应用。

3.实验法：通过实验让学生亲手操作，掌握GPS数据的采集、解析和处理方法，以及实现简单的导航功能。

四、教学资源为了支持本课程的教学内容和教学方法的实施，我们将准备以下教学资源：1.教材：选用合适的C语言和GPS相关教材，为学生提供理论学习的参考。

2.参考书：提供相关的参考书籍，拓展学生的知识面。

3.多媒体资料：制作课件和视频资料，直观地展示GPS系统和C语言的相关内容。

4.实验设备：准备GPS接收器、计算机等实验设备，让学生进行实际操作。

五、教学评估为了全面、客观地评估学生在C语言GPS课程中的学习成果，我们将采取以下评估方式：1.平时表现：通过学生在课堂上的参与度、提问回答、小组讨论等表现来评估其学习态度和理解程度。

2.作业：布置与课程内容相关的编程作业，评估学生对C语言编程和GPS应用的掌握情况。

3.考试：定期进行理论知识考试和编程实践考试，全面评估学生的知识掌握和实际应用能力。

4.项目报告：要求学生完成一个GPS相关的编程项目，通过项目报告评估学生的综合应用能力。

六、教学安排本课程的教学安排如下：1.教学进度：按照教材和大纲进行，确保覆盖所有重要知识点。

2.3.3
这是一种使用没有在标准下预定义的特殊内容的方法。它通常的格式为：$PmmmA,df1,df2,....,[CR][LF]P说明是属性内容，mmm定义为厂商信息代码，A(A-Z）表明信息类型。
NMEA-0183协议定义的语句非常的多，主要有$GPGLL、 $GPGGA、$GPVTG、$GPZDA等。因为本次课题我们需要使用GPS-GPZDA进行数据解析所以我们只介绍GPZDA标准格式的样式。
2.3.2
查询的标准格式为：$ttllQ,sss,[CR][LF]
前面两个字符表示请求者的信息源的识别码，后面两个字符是被查询设备的识别码，最后一个字符说明这是一个查询的信息。紧跟着的字符（SSS）包含了三个字的被查询内容的记忆信息。例如： $CCGPQ,ZDA[CR][LF] “CC”说明计算机正从“GP”这个设备查询ZDA中的内容。GPS将会每隔1秒发送内容直到有新的请求。
char chhour[5]={'\0'};
char chmin[5]={'\0'};
char chsec[8]={'\0'};
char chsum[10]={'\0'};
char afterch[128]={'\0'};
int i,j,n,m; /*定义并初始化一些数组，指针，整型变量等用于后续程序使用*/
通过这次课程设计我们可以了解，如何用GPZDA解析GPS的时间信号，在计算机上输出普通用户所能理解的语句，让用户能够在时间信息这方面更加精确地得到信息，给用户提供更可靠的更高水准的服务。
第一章
GPS全球卫星定位系统由21颗工作卫星和3颗备用卫星组成。它们分布在6个等间距的轨道平面上，每个轨道面上有4颗工作卫星。GPS能覆盖全球，能够连续、实时、隐蔽地定位，一次定位时间仅几秒到十几秒，用户不发射任何电磁信号，只要接受卫星导航信号即可定位，所以可全天候昼夜作业。

GPS编码格式及C语⾔解码有关磁偏⾓和地图定位的问题：地图的⽅向：上北、下南、左西、右东是⼤多数地图的⽅向，但这可不是通⽤原则，如果地图上有⽅向标，可以通过⽅向标了解到这些。

地磁极是接近南极和北极的，但并不和南极、北极重合，⼀个约在北纬72°、 西经96°处；⼀个约在南纬70°、东经150°处。

磁北极距地理北极⼤约相差1500km。

现在磁北极的位置在加拿⼤北⽅，在⼀天中磁北极的位置也是不停的变动，它的轨迹⼤致为⼀椭圆形，磁北极平均每天向北以40m。

⽬前磁北极正在逐渐离开加拿⼤，⼤约于2005年进⼊俄罗斯境内。

东经25度地区，磁偏⾓在1-2度之间；北纬25度以上地区，磁偏⾓⼤于2度；若在西经低纬度地区，磁偏⾓是5-20度；西经45度以上，磁偏⾓为25-50度。

在我国，正常情况下，磁偏⾓最⼤可达6度，⼀般情况为2-3度。

在我国除部分磁⼒异常的地⽅外，⼀般磁偏⾓都是西偏。

磁偏⾓还是不断有规律变化的，地图上的磁偏⾓只是测图时的磁偏⾓（磁北⽐真北偏左，加上磁偏⾓；磁北⽐真北偏右，减去磁偏⾓；在我国⼀般是加上）。

使⽤地图本⾝所注的磁偏⾓要注意出版年限，地图太⽼误差较⼤。

GPS 0183协议GGA、GLL、GSA、GSV、RMC、VTG解释$GPGGA例：$GPGGA,092204.999,4250.5589,S,14718.5084,E,1,04,24.4,19.7,M,,,,0000*1F字段0：$GPGGA，语句ID，表明该语句为Global Positioning System Fix Data（GGA）GPS定位信息字段1：UTC 时间，hhmmss.sss，时分秒格式字段2：纬度ddmm.mmmm，度分格式（前导位数不⾜则补0）字段3：纬度N（北纬）或S（南纬）字段4：经度dddmm.mmmm，度分格式（前导位数不⾜则补0）字段5：经度E（东经）或W（西经）字段6：GPS状态，0=未定位，1=⾮差分定位，2=差分定位，3=⽆效PPS，6=正在估算字段7：正在使⽤的卫星数量（00 - 12）（前导位数不⾜则补0）字段8：HDOP⽔平精度因⼦（0.5 - 99.9）字段9：海拔⾼度（-9999.9 - 99999.9）字段10：地球椭球⾯相对⼤地⽔准⾯的⾼度字段11：差分时间（从最近⼀次接收到差分信号开始的秒数，如果不是差分定位将为空）字段12：差分站ID号0000 - 1023（前导位数不⾜则补0，如果不是差分定位将为空）字段13：校验值$GPGLL例：$GPGLL,4250.5589,S,14718.5084,E,092204.999,A*2D字段0：$GPGLL，语句ID，表明该语句为Geographic Position（GLL）地理定位信息字段1：纬度ddmm.mmmm，度分格式（前导位数不⾜则补0）字段2：纬度N（北纬）或S（南纬）字段3：经度dddmm.mmmm，度分格式（前导位数不⾜则补0）字段4：经度E（东经）或W（西经）字段5：UTC时间，hhmmss.sss格式字段6：状态，A=定位，V=未定位字段7：校验值$GPGSA例：$GPGSA,A,3,01,20,19,13,,,,,,,,,40.4,24.4,32.2*0A字段0：$GPGSA，语句ID，表明该语句为GPS DOP and Active Satellites（GSA）当前卫星信息字段1：定位模式，A=⾃动⼿动2D/3D，M=⼿动2D/3D字段2：定位类型，1=未定位，2=2D定位，3=3D定位字段3：PRN码（伪随机噪声码），第1信道正在使⽤的卫星PRN码编号（00）（前导位数不⾜则补0）字段4：PRN码（伪随机噪声码），第2信道正在使⽤的卫星PRN码编号（00）（前导位数不⾜则补0）字段5：PRN码（伪随机噪声码），第3信道正在使⽤的卫星PRN码编号（00）（前导位数不⾜则补0）字段6：PRN码（伪随机噪声码），第4信道正在使⽤的卫星PRN码编号（00）（前导位数不⾜则补0）字段7：PRN码（伪随机噪声码），第5信道正在使⽤的卫星PRN码编号（00）（前导位数不⾜则补0）字段8：PRN码（伪随机噪声码），第6信道正在使⽤的卫星PRN码编号（00）（前导位数不⾜则补0）字段9：PRN码（伪随机噪声码），第7信道正在使⽤的卫星PRN码编号（00）（前导位数不⾜则补0）字段10：PRN码（伪随机噪声码），第8信道正在使⽤的卫星PRN码编号（00）（前导位数不⾜则补0）字段11：PRN码（伪随机噪声码），第9信道正在使⽤的卫星PRN码编号（00）（前导位数不⾜则补0）字段12：PRN码（伪随机噪声码），第10信道正在使⽤的卫星PRN码编号（00）（前导位数不⾜则补0）字段13：PRN码（伪随机噪声码），第11信道正在使⽤的卫星PRN码编号（00）（前导位数不⾜则补0）字段14：PRN码（伪随机噪声码），第12信道正在使⽤的卫星PRN码编号（00）（前导位数不⾜则补0）字段15：PDOP综合位置精度因⼦（0.5 - 99.9） 字段16：HDOP⽔平精度因⼦（0.5 - 99.9）字段17：VDOP垂直精度因⼦（0.5 - 99.9）字段18：校验值$GPGSV例：$GPGSV,3,1,10,20,78,331,45,01,59,235,47,22,41,069,,13,32,252,45*70字段0：$GPGSV，语句ID，表明该语句为GPS Satellites in View（GSV）可见卫星信息字段1：本次GSV语句的总数⽬（1 - 3）字段2：本条GSV语句是本次GSV语句的第⼏条（1 - 3）字段3：当前可见卫星总数（00 - 12）（前导位数不⾜则补0）字段4：PRN 码（伪随机噪声码）（01 - 32）（前导位数不⾜则补0）字段5：卫星仰⾓（00 - 90）度（前导位数不⾜则补0）字段6：卫星⽅位⾓（00 - 359）度（前导位数不⾜则补0）字段7：信噪⽐（00－99）dbHz字段8：PRN 码（伪随机噪声码）（01 - 32）（前导位数不⾜则补0）字段9：卫星仰⾓（00 - 90）度（前导位数不⾜则补0）字段10：卫星⽅位⾓（00 - 359）度（前导位数不⾜则补0）字段11：信噪⽐（00－99）dbHz字段12：PRN 码（伪随机噪声码）（01 - 32）（前导位数不⾜则补0）字段13：卫星仰⾓（00 - 90）度（前导位数不⾜则补0）字段14：卫星⽅位⾓（00 - 359）度（前导位数不⾜则补0）字段15：信噪⽐（00－99）dbHz字段16：校验值$GPRMC例：$GPRMC,024813.640,A,3158.4608,N,11848.3737,E,10.05,324.27,150706,,,A*50字段0：$GPRMC，语句ID，表明该语句为Recommended Minimum Specific GPS/TRANSIT Data（RMC）推荐最⼩定位信息 字段1：UTC时间，hhmmss.sss格式字段2：状态，A=定位，V=未定位字段3：纬度ddmm.mmmm，度分格式（前导位数不⾜则补0）字段4：纬度N（北纬）或S（南纬）字段5：经度dddmm.mmmm，度分格式（前导位数不⾜则补0）字段6：经度E（东经）或W（西经）字段7：速度，节，Knots字段8：⽅位⾓，度字段9：UTC⽇期，DDMMYY格式字段10：磁偏⾓，（000 - 180）度（前导位数不⾜则补0）字段11：磁偏⾓⽅向，E=东W=西字段16：校验值$GPVTG例：$GPVTG,89.68,T,,M,0.00,N,0.0,K*5F字段0：$GPVTG，语句ID，表明该语句为Track Made Good and Ground Speed（VTG）地⾯速度信息字段1：运动⾓度，000 - 359，（前导位数不⾜则补0）字段2：T=真北参照系字段3：运动⾓度，000 - 359，（前导位数不⾜则补0）字段4：M=磁北参照系字段5：⽔平运动速度（0.00）（前导位数不⾜则补0）字段6：N=节，Knots字段7：⽔平运动速度（0.00）（前导位数不⾜则补0）字段8：K=公⾥/时，km/h字段9：校验值GPS编码格式及C语⾔解码GPS接收机只要处于⼯作状态就会源源不断地把接收并计算出的GPS导航定位信息通过串⼝传送到计算机中。

#include <stdio.h>#include <math.h>#include <stdlib.h>#define bGM84 3.986005e14#define bOMEGAE84 7.2921151467e-5void main(){long double roota=0.515365263176E+04; //轨道长半轴的平方根（根号a）long double toe=0.720000000000E+04; //观测时刻toelong double m0=-0.290282040486E+00; //参考时刻toe的平近点角long double e=0.678421219345E-02; //轨道偏心率elong double delta_n=0.451411660250E-08;//卫星的摄动改正数△nlong double smallomega=-0.258419417299E+01;//近地点角距ωlong double cus=0.912137329578E-05;//纬度幅角正弦调和项改正的振幅（弧度）long double cuc=0.189989805222E-06;//纬度幅角余弦调和项改正的振幅（弧度）long double crs=0.406250000000E+01;//轨道半径的余弦调和项改正的振幅（m）long double crc=0.201875000000E+03;//轨道半径的正弦调和项改正的振幅（m）long double cis=0.949949026108E-07;//轨道倾角的余弦调和项改正的振幅（弧度）long double cic=0.130385160446E-07;//轨道倾角的正弦调和项改正的振幅（弧度）long double idot=-0.253939149013E-09;//轨道倾角变化率Ilong double i0=0.958512160302E+00; //轨道倾角（弧度）long double bigomega0=-0.137835982556E+01;//升交点赤经long double earthrate=bOMEGAE84; //地球自转的速率welong double bigomegadot=-0.856928551657e-08;long double t=0.720000000000E+04; //加入卫星钟差改正的归化时间long double A;long double n0=0,n,tk;long double mk,ek,tak,ik,omegak,phik,uk,rk;long double corr_u,corr_r,corr_i;long double xpk,ypk,xk,yk,zk;int i;printf("输入的数据：\n");printf("√a=%e \n",roota);printf("toe=%e \n",toe);printf("e=%e \n",e);printf("i0=%e \n",i0);printf("ω=%e \n",smallomega);printf("△n=%e \n",delta_n);printf("Ω0=%e \n",bigomega0);printf("I=%e \n",idot);printf("Cuc=%e \n",cuc);printf("Cus=%e \n",cus);printf("Crc=%e \n",crc);printf("Crs=%e \n",crs);printf("Cic=%e \n",cic);printf("Cis=%e \n",cis);printf("\n\n输出的结果为：\n",e);A=roota*roota;n0=sqrt(bGM84/(A*A*A));//平均角速度n0printf("n0=%.10lf \n",n0);tk=t-toe;//相对于参考时刻toe的归化时间tkprintf("tk=%.10lf \n",tk);n=n0+delta_n;//加摄动改正后的卫星平均角速度printf(" n=%.10lf \n",n);mk=m0+n*tk;//卫星平近点角printf("mk=%.10lf \n",mk);ek=mk;for(i=0;i<10;i++) ek=mk+e*sin(ek);//利用迭代法求偏近点角ekprintf("ek=%.10lf \n",ek);tak=atan2(sqrt(1.0-e*e)*sin(ek)/cos(ek)-e);//真近点角Vk的计算printf("Vk=%.10lf \n",tak);phik=tak+smallomega;//升交距角φk的计算printf("φk=%.10lf \n",phik);corr_u=cus*sin(2.0*phik)+cuc*cos(2.0*phik);//升交距角u的摄动改正δu printf("δu=%.10lf \n",corr_u);corr_r=crs*sin(2.0*phik)+crc*cos(2.0*phik);//卫星矢量r的摄动改正δr printf("δr=%.10lf \n",corr_r);corr_i=cis*sin(2.0*phik)+cic*cos(2.0*phik);//轨道倾角i的摄动改正δi printf("δi=%.10lf \n",corr_i);uk=phik+corr_u;//升交距角uprintf("uk=%.10lf \n",uk);rk=A*(1.0-e*cos(ek))+corr_r;//卫星矢量rprintf("rk=%.10lf \n",rk);ik=i0+idot*tk+corr_i;//轨道倾角iprintf("ik=%.10lf \n",ik);xpk=rk*cos(uk);//卫星在轨道平面坐标系的坐标ypk=rk*sin(uk);printf("xpk=%.10lf \n",xpk);printf("ypk=%.10lf \n",ypk);omegak=bigomega0+(bigomegadot-earthrate)*tk-earthrate*toe;//升交点经度Ωk的计算printf("Ωk=%.10lf \n\n",omegak);xk=xpk*cos(omegak)-ypk*sin(omegak)*cos(ik);//地心固定坐标系的直角坐标yk=xpk*sin(omegak)+ypk*cos(omegak)*cos(ik);zk=ypk*sin(ik);printf("Xk=%.4lf \n",xk);printf("Yk=%.4lf \n",yk);printf("Zk=%.4lf \n",zk);}。

GPS卫星位置的计算#include <stdio.h>#include <math.h>//GPS卫星数据结构typedef structdouble PRN; // PRN号double t; // 接收时间double Pseudorange; // 伪距double Azimuth; // 方位角double Elevation; // 俯仰角} Satellite;//计算GPS卫星的位置void calculateSatellitePosition(Satellite* satellite, double* X, double* Y, double* Z)double t = satellite->t;double Pseudorange = satellite->Pseudorange;double Azimuth = satellite->Azimuth;double Elevation = satellite->Elevation;double tau = c * (t - 0); // 信号传播时间//计算地心坐标系下的卫星位置double Xs = Pseudorange * cos(Elevation) * sin(Azimuth + omega_e * tau);double Ys = Pseudorange * cos(Elevation) * cos(Azimuth + omega_e * tau);double Zs = Pseudorange * sin(Elevation);//从地心坐标系转换到站心坐标系double latitude = *X; // 纬度double longitude = *Y; // 经度double sinLat = sin(latitude * PI / 180);double cosLat = cos(latitude * PI / 180);double sinLon = sin(longitude * PI / 180);double cosLon = cos(longitude * PI / 180);double X = -Xs * sinLon - Ys * sinLat * cosLon + Zs * cosLat * cosLon;double Y = Xs * cosLon - Ys * sinLat * sinLon + Zs * cosLat * sinLon;double Z = Ys * cosLat + Zs * sinLat;*X=X;*Y=Y;*Z=Z;int maiSatellite satellite;satellite.PRN = 1;satellite.t = 0;satellite.Azimuth = 45;satellite.Elevation = 30;double X = 0; // 站心坐标系的X坐标double Y = 0; // 站心坐标系的Y坐标double Z = 0; // 站心坐标系的Z坐标calculateSatellitePosition(&satellite, &X, &Y, &Z);printf("Satellite position: (%f, %f, %f)\n", X, Y, Z);return 0;此示例中，我们定义了一个Satellite结构体来保存GPS卫星的相关数据，包括PRN号、接收时间、伪距、方位角和俯仰角。

GPS卫星位置的计算利用C++编写了一段能计算单一瞬时卫星坐标的程序，在运行程序之前，需做部分准备工作：（1）在F盘下建立一名为“单一卫星广播星历”的txt文件。

（2）从“广播星历.txt”文件中拷贝从卫星PRN号开始的8行数据到“单一卫星广播星历.txt”中（3）在编辑选项中，将全部的“D”替换为“E”。

下面为我所选取的一个广播星历：18 06 8 25 6 0 0.0-2.472363412380E-04-1.023*********E-12 0.000000000000E+001.410000000000E+02-1.721875000000E+01 4.502687555010E-09 1.413760604187E+00-7.990747690201E-07 7.598234573379E-03 1.118145883083E-05 5.153709835052E+034.536000000000E+05-1.303851604462E-08-1.095067942661E-01 1.527369022369E-079.571235745530E-01 1.640000000000E+02-2.656176299285E+00-8.0374********E-09-5.193073455211E-10 1.000000000000E+00 1.389000000000E+03 0.000000000000E+002.000000000000E+00 0.000000000000E+00-1.024*********E-08 1.410000000000E+024.464490000000E+05 4.000000000000E+00程序设计部分：#include<stdio.h>#include<math.h>int main(){int i = 0;double n[50], n0, nn, t, tk, Mk, Ek, Vk, Yk, Gu, Gr, Gi, uk, rk, ik, xk, yk, zk, X, Y, Z, Lk, UT, yy, mm, JD, gpsz;FILE *fp;fp = fopen("F:\\单一卫星广播星历.txt", "r");if (fp == NULL){printf ("文件打开失败!\n");return 0;}while (! feof (fp)){fscanf(fp, "%lf", &n[i]);i++;}n0 = (sqrt(3986005E+8))/pow(n[17], 3);nn = n0 + n[12];/*计算卫星运行的平均角速度*/UT = n[4] + (n[5] / 60) + (n[7] / 3600);/*民用日的时分秒化为实数时*/if (n[1] >= 80)/*广播星历中年只有后两位，化为4位，参考1980年1月6日0点*/ {if (n[1] == 80 && n[2] == 1 && n[3] < 6){n[1] = n[1] + 2000;}n[1] = n[1] + 1900;}else{n[1] = n[1] + 2000;}if (n[2] <= 2){yy = n[1] - 1;mm = n[2] + 12;}if (n[2] > 2){yy = n[1];mm = n[2];}JD = (int)(365.25 * yy) + (int)(30.6001 * (mm + 1)) + n[3] + (UT / 24) + 1720981.5;/*化为儒略日*/gpsz = (int)((JD - 2444244.5) / 7);/*计算GPS周*/t = (JD - 2444244.5 - 7 * gpsz) * 24 * 3600;/*得出GPS秒*/tk = t - n[18];/*tk1为中间值，用以判断tk与正负302400的关系，然后返回到tk上*/ while (tk > 302400 || tk < -302400){if (tk > 302400){tk = tk - 604800;}else{tk = tk + 604800;}}/*计算归化观测时间*/Mk = n[13] + nn * tk;/*观测时刻的卫星平近点角*/Ek = Mk;Ek = Mk + n[15] * sin(Ek);Ek = Mk + n[15] * sin(Ek);/*迭代两次计算观测时刻的偏近点角*/Vk = atan(sqrt(1 - n[15] * n[15]) * sin(Ek)) / (cos(Ek) - n[15]);/*真近点角*/Yk = Vk + n[24];/*升交距角*/Gu = n[14] * cos(2 * Yk) + n[16] * sin(2 * Yk);Gr = n[23] * cos(2 * Yk) + n[11] * sin(2 * Yk);Gi = n[19] * cos(2 * Yk) + n[21] * sin(2 * Yk);/*摄动改正项*/uk = Yk + Gu;rk = n[17] * n[17] * (1 - n[15] * cos(Ek)) + Gr;ik = n[22] + Gi + n[26] * tk;/*经摄动改正后的升交距角、卫星矢径、轨道倾角*/xk = rk * cos(uk);yk = rk * sin(uk);zk = 0;/*卫星在轨道坐标系的坐标*/Lk = n[20] + (n[25] - 7.29211515E-5) * tk - 7.29211515E-5 * n[18];/*观测时刻t的升交点经度*/X = xk * cos(Lk) - yk * cos(ik) * sin(Lk);Y = xk * sin(Lk) + yk * cos(ik) * cos(Lk);Z = yk * sin(ik);/*卫星在WGS-84坐标系的坐标*/printf("该卫星在WGS-84坐标系中的坐标为：\nX = %lf m\nY = %lf m\nZ = %lf m\n", X, Y, Z);fclose(fp);return 0;}计算结果：该卫星在WGS-84坐标系中的坐标为：X = 9223153.692525 mY = 24133486.931401 mZ = 6032585.919385 m。

C语言计算星历位置GPS广播星历计算卫星位置和速度GPS（全球定位系统）广播星历是一种用于计算卫星位置和速度的星历数据。

星历是一种描述卫星在天空中的位置和速度的数据，可以通过接收卫星广播的星历数据来计算卫星的当前位置和速度。

GPS广播星历数据通常由地面控制站维护和广播，卫星周期性地发送星历数据以更新接收器上的星历数据库。

接收器通过接收到的星历数据计算卫星的位置和速度，并使用这些信息来计算接收器的位置。

星历数据通常包括卫星的轨道参数和时间参数。

轨道参数包括卫星的半长轴、偏心率、轨道倾角、升交点经度、近地点幅角和平近点角速度。

时间参数包括卫星的时钟校正参数和广播时间。

接收器使用这些参数来计算卫星的位置和速度。

首先，接收器通过测量卫星信号的到达时间和广播时间来计算卫星信号的传播时间。

然后，接收器使用传播时间和卫星的时间参数来计算卫星的时间误差。

接下来，接收器使用卫星的轨道参数来计算卫星的真实位置和速度。

接收器使用卫星的时间误差来修正卫星的广播时间，并将其转换为GPS时间。

然后，接收器使用修正后的广播时间和卫星的轨道参数来计算卫星的位置和速度。

计算得到的卫星位置和速度可以用于定位接收器的位置。

接收器通过测量多个卫星的信号传播时间来计算卫星到接收器的距离。

然后，接收器使用卫星的位置和速度来计算接收器的位置。

GPS广播星历计算卫星位置和速度是GPS定位的核心技术之一、接收器通过接收到的星历数据来计算卫星的位置和速度，然后使用这些信息来计算接收器的位置。

这种计算过程是GPS定位的基础，可以用于估计接收器的位置和速度。

总结起来，GPS广播星历是一种用于计算卫星位置和速度的星历数据。

接收器通过接收到的星历数据和卫星信号的传播时间来计算卫星的位置和速度。

这些计算结果可以用于定位接收器的位置，并在导航和定位应用中发挥重要作用。

C语言计算星历位置GPS广播星历计算卫星位置和速度C语言是一种通用的高级编程语言，可以用于计算星历位置以及计算GPS卫星位置和速度。

下面将详细介绍如何使用C语言来实现这些计算。

首先，我们需要了解星历和GPS广播星历的概念。

星历是一种描述天体位置的方法，它包含了每个天体的位置坐标、速度以及其他相关的信息。

星历常用于天文学研究和导航系统中。

GPS广播星历是由GPS卫星广播的星历信息，它包含了GPS卫星所处的位置、速度等信息。

通过接收并解码广播星历，我们可以计算出卫星的位置和速度。

在C语言中，我们可以使用数学库和一些公式来计算星历位置和GPS 卫星位置以及速度。

首先，我们需要导入数学库，可以使用`#include <math.h>`导入。

数学库提供了一些常用的数学函数，如计算平方根、计算三角函数等。

然后，我们需要根据星历或广播星历的信息，计算出卫星的位置和速度。

对于星历位置的计算，可以使用开普勒方程来逼近天体的真实位置。

开普勒方程的计算公式如下：E - e * sin(E) = M其中，E为偏近点角，e为偏心率，M为平近点角。

通过迭代计算，可以得到E的近似值。

然后，利用半长轴、偏心率和E的值，可以计算出卫星在轨道平面上的坐标。

对于GPS卫星位置和速度的计算，可以使用广播星历中的卫星钟差、偏心率修正项等信息。

具体的计算公式较为复杂，需要使用专门的算法进行计算。

在计算过程中，我们还需要考虑坐标系的转换，以确保最终计算得到的是相对于地球的地心坐标系中的位置和速度。

最后，我们可以将计算得到的卫星位置和速度输出，以便进行后续的处理或导航操作。

总结来说，使用C语言计算星历位置和GPS卫星位置和速度需要导入数学库并使用开普勒方程以及其他相关的计算公式来进行计算。

同时，还需要考虑坐标系的转换和其他相关的因素。

这只是一个简单的介绍，具体的实现可能需要更多的代码和算法。

GPS卫星位置计算实验学校：工业大学学号：20104168专业班级：测绘工程10-1班学院：土木与水利工程学院指导教师：庭叶一、实验目的1、通过对GPS卫星位置的计算，增强我们对《GPS测量原理与应用》课程的理解，使我们牢固掌握GPS测量的基本原理和过程，熟悉GPS测量数据处理的基本技能和计算方法。

2、提高运用计算机语言编程开发能力，运用测量平差知识、数学知识和计算机知识，通过给定的程序算法，编制程序实现GPS卫星位置的计算过程。

二、实验容与要求1、通过课堂对GPS测量原理的学习，在课后自主完成GPS卫星位置的计算过程并按照课本上的步骤设计程序完成GPS卫星位置的计算过程；2、要求每位同学独立完成GPS卫星位置的计算过程，编写程序，调试程序，并编写程序设计文档。

要求过程和算确、程序运行正确、设计文档完备；三、课程设计工具运用自己熟悉的编程开发语言（C、C#、C++、VC、VS、VB、FORTRAN等）。

四、实验总结1、实验简单介绍运行后的主界面比以前做的程序要漂亮些，给界面增加了一副图片，让界面看起来还是比较的漂亮，但由于起始的参数较多，因此，界面整体看起来比较的拥挤，但这个不影响程序的计算过程。

整个程序的功能不是很多，但满足了基本的要求，能够进行卫星的位置计算，并且能够达到很好的精度，与课本上给出的结果相差很小，符合实验的要求。

为了避免繁琐的数据的输入，在本次程序中我增加了一个的功能按钮，点击后，程序自动给相关的起始数据赋值并显示在界面上，如图所示：这些数据是教材《GPS测量原理与应用》41页的卫星星历数据；用户也可以根据需要更改相应的数据，如上图所示，当点击按钮时，就可以在对话框相应的位置计算出卫星的坐标，如下图所示：其中，计算出的卫星的坐标为Xk = 4589210.3510074820，Yk = 25371005.6999580190，Zk = -5618292.2998269377，这是卫星在地固坐标系中的坐标，如果还知道极移参数就可以计算出卫星在协议地球坐标系中的坐标了。

GPS模拟C语言第一篇：GPS模拟C语言GPS模拟# include # include # includeusing namespace std;ifstream fin(“1.txt”);//ofstream fout(“data_out.txt”);// 1 double x[5] = {-1, 0, 300, 0, 300};double y[5] = {-1, 450, 450, 0, 0};double t[5];double A, B, D,c = 2982000, p, r;double pi = 3.141592653;double xx, yy;//三角形算法void cal(int sign, double t1, double t2, double x1, double y1, double x2, double y2, double px, double py){double tt1 = c * t1, tt2 = c * t2;double q, ac;A = x2 *(x1*x1 + y1*y1x1 *(x2*x2 + y2*y2tt1*tt1)tt2*tt2);D = tt1 *(x2*x2 + y2*y2tt2 *(x1*x1 + y1*y1ac;r =(x1*x1 + y1*y10)< 1e-6){switch(h){case 1:cal(-1, t[2], t[4], 300, 0, 300,-450, 0, 450);break;case 2: cal(-1, t[1], t[3],-300, 0,-300,-450,300, 450);break;case 3: cal(1, t[1], t[2], 0, 450, 300, 450, 0, 0);break;case 4: cal(1, t[1], t[2],-300, 450, 0, 450, 300, 0);break;}break;} } system(“pause”);return 0;声源定位#include #include using namespace std;#define Max 50 int N=10;double xi[10];double yi[10];double ti[10][3];double det(double a[Max][Max],int n){ int i,j,start;double k,temp,result=1;for(start=0;start<=n-2;start++){i=1;while(a[start][start]==0&&start+i<=n-1){for(j=start;j<=n-1;j++){temp=-1*a[start][j];a[start][j]=a[start+i][j];a[start+i][j]=temp;}i++;}if(start+i==n)continue;for(i=start+1;i<=n-1;i++){k=a[i][start]/a[start][start];for(j=start;j<=n-1;j++)a[i][j]=a[i][j]-k*a[start][j];} } for(i=0;i<=n-1;i++)result*=a[i][i];return result;} int equalation(double a[Max][Max+1],double r[Max],int n){ double det(double a[Max][Max],int n);double det0[Max][Max],det1[Max][Max];int i,j,k;for(i=0;i<=n-1;i++)for(j=0;j<=n-1;j++)det0[i][j]=a[i][j];if(det(det0,n)==0)return 0;for(k=0;k<=n-1;k++){for(i=0;i<=n-1;i++)for(j=0;j<=n-1;j++){det0[i][j]=a[i][j];det1[i][j]=a[i][j];}for(i=0;i<=n-1;i++)det1[i][k]=a[i][n];r[k]=det(det1,n)/det(det0,n);} return 1;} double sqr(double x){ return x*x;} double g(double x,double y,double c,double xi,double yi,double ti){ return sqr(xi-x)+sqr(yi-y)-sqr(c*ti);} double fx(double x,double y,double c,int k){ int i;double sum=0;for(i=0;i<=N-1;i++)sum+=g(x,y,c,xi[i],yi[i],ti[i][k-1])*(xi[i]-x);return sum;} double fy(double x,double y,double c,int k){ int i;double sum=0;for(i=0;i<=N-1;i++)sum+=g(x,y,c,xi[i],yi[i],ti[i][k-1])*(yi[i]-y);return sum;} double fc(double x,double y,double c,int k){ int i;double sum=0;for(i=0;i<=N-1;i++)sum+=g(x,y,c,xi[i],yi[i],ti[i][k-1])*sqr(ti[i][k-1]);return sum;} double fxx(double x,double y,double c,int k){ int i;double sum=0;for(i=0;i<=N-1;i++)sum=sum-2*sqr(xi[i]-x)-g(x,y,c,xi[i],yi[i],ti[i][k-1]);return sum;} double fxy(double x,double y,double c,int k){ int i;double sum=0;for(i=0;i<=N-1;i++)sum=sum-2*(yi[i]-y)*(xi[i]-x);return sum;} double fxc(doublex,double y,double c,int k){ int i;double sum=0;for(i=0;i<=N-1;i++) sum=sum-2*c*sqr(ti[i][k-1])*(xi[i]-x);return sum;} double fyx(double x,double y,double c,int k){ int i;double sum=0;for(i=0;i<=N-1;i++)sum=sum-2*(xi[i]-x)*(yi[i]-y);return sum;} double fyy(double x,double y,double c,int k){ int i;double sum=0;for(i=0;i<=N-1;i++)sum=sum-2*sqr(yi[i]-y)-g(x,y,c,xi[i],yi[i],ti[i][k-1]);return sum;} double fyc(double x,double y,double c,int k){ int i;double sum=0;for(i=0;i<=N-1;i++)sum=sum-2*c*sqr(ti[i][k-1])*(yi[i]-y);return sum;} double fcx(double x,double y,double c,int k){ int i;double sum=0;for(i=0;i<=N-1;i++)sum=sum-2*(xi[i]-x)*sqr(ti[i][k-1]);return sum;} double fcy(double x,double y,double c,int k){ int i;double sum=0;for(i=0;i<=N-1;i++)sum=sum-2*(yi[i]-y)*sqr(ti[i][k-1]);return sum;} double fcc(double x,double y,double c,int k){ int i;double sum=0;for(i=0;i<=N-1;i++)sum=sum-2*c*sqr(sqr(ti[i][k-1]));return sum;} int main(){ int i=0,k;double x,y,c,x0,y0,c0;double a[Max][Max+1],r[Max];freopen(“1.txt”, “r”, stdin);freopen(“1_out.txt”, “w”, stdout);for(i=0;i<=9;i++){} cin>>xi[i]>>yi[i];xi[i] /= 10;yi[i] /= 10;cin>>ti[i][0]>>ti[i][1]>>ti[i][2];for(k=1;k<=3;k++){ cin>>x> >y;c=0.2;do {x0=x;y0=y;c0=c;i++;a[0][0]=fxx(x,y,c,k);a[0][1]=fxy(x,y,c,k);a[0][2]=fxc(x,y,c,k);a[1][0]=fyx(x,y,c,k);a[1][1]=fyy(x,y,c,k);a[1][2]=fyc(x,y,c,k);a[2][0]=fcx(x, y,c,k);a[2][1]=fcy(x,y,c,k);a[2][2]=fcc(x,y,c,k);a[0][3]=x*a[0][0]+y*a[ 0][1]+c*a[0][2]-fx(x,y,c,k);a[1][3]=x*a[1][0]+y*a[1][1]+c*a[1][2]-fy(x,y,c,k);a[2][3]=x*a[2][0]+y*a[2][1]+c*a[2][2]-fc(x,y,c,k);if(equalation(a,r,3)==0){cout<return 0;}x=r[0];y=r[1];c=r[2];} while(abs(x-x0)+abs(y-y0)+abs(c-c0)>=1e-6);cout<<<<第二篇：C语言模拟ATM机一、实验目的通过设计一个ATM机模拟操作的程序，全面运用课程的主要知识点，巩固对模块化程序设计、文件操作的理解，提高软件编程能力。

#include <stdio.h>#include <stdlib.h>#include "rinex.h"#include "sv.h"void WriteSatPosFile(FILE *SatPosFile,int prn,XYZCoor *SVPos);void GetSVPos(SVText *Text,double t,double *X,double *Y,double *Z);double Get_atan(double z,double y);int GetGPSTime(int year,int month,int day,int hour,int minute,double second,double *gpstime);int read_RinexEPP( FILE *RinexEPP_file,SVText *snv,int *year, int *month,int *day,int *hour,int *minute,double *second,double *gpstime,int *wn);void main(){char RinexEPPName[20];FILE *RinexEPP_file, *SatPosFile;double SVPosX,SVPosY,SVPosZ;int prn,i;SVText snv;char temp[200];int wn,year,month,day,hour, minute;double second,gpstime;printf("Please Input Rinex Nav File Name: ");scanf("%s",RinexEPPName);if ((RinexEPP_file = fopen(RinexEPPName, "rt")) == NULL) {fprintf(stderr, "Cannot open input file.\n");exit(1); }if ((SatPosFile = fopen("satpos.out", "w"))== NULL) {fprintf(stderr, "Cannot open outut file.\n");exit(1); }rewind(RinexEPP_file);for(i=0;i<3;i++) {fgets(temp,200,RinexEPP_file); }do {if(read_RinexEPP(RinexEPP_file,&snv,&year,&month,&day,&hour,&minute,&second,&gpstime,&wn) ) break;/* compute satellite position */GetSVPos(&snv,gpstime,&SVPosX,&SVPosY,&SVPosZ);/* write satellite position file */fprintf(SatPosFile,"%3i%19.11e%19.11e%19.11e\n",snv.prn,SVPosX,SVPosY,SVPosZ);}while(1);fclose(RinexEPP_file);fclose(SatPosFile);}/*********************** read_EPP() *************************/int read_RinexEPP( FILE *RinexEPP_file,SVText *snv,int *year, int *month,int *day,int *hour,int *minute,double *second,double *gpstime,int *wn){int j;char t1[30],t2[30],t3[30],t4[30],t5[30];/* PRN /EPOCH /SV CLK */if(fscanf(RinexEPP_file,"%d%d%d%d%d%d%s%s%s%s\n",&snv->prn,year,month,day,hour,minute,t1,t2,t3,t4)==EOF) return 1;*second=atof(t1);snv->af0=atof(t2);snv->af1=atof(t3);snv->af2=atof(t4);*year+=1900;*wn=snv->wn = GetGPSTime(*year,*month,*day,*hour,*minute,*second, gpstime);snv->tow = (long) gpstime;/* BROADCAST ORBIT -1 */if(fscanf( RinexEPP_file,"%s%s%s%s\n",t1,t2,t3,t4)==EOF) return 1;snv->aode=atof(t1);snv->crs=atof(t2);snv->deltan=atof(t3);snv->m0=atof(t4);/* BROADCAST ORBIT -2 */if(fscanf( RinexEPP_file,"%s%s%s%s\n",t1,t2,t3,t4)==EOF) return 1;snv->cuc=atof(t1);snv->e=atof(t2);snv->cus=atof(t3);snv->roota=atof(t4);/* BROADCAST ORBIT -3 */if(fscanf( RinexEPP_file,"%s%s%s%s\n", t1,t2,t3,t4)==EOF) return 1;snv->toe=atof(t1);snv->cic=atof(t2);snv->omega0=atof(t3);snv->cis=atof(t4);/* BROADCAST ORBIT -4 */if(fscanf( RinexEPP_file,"%s%s%s%s\n",t1,t2,t3,t4)==EOF) return 1 ;snv->i0=atof(t1);snv->crc=atof(t2);snv->omega=atof(t3);snv->omegadot=atof(t4);/* BROADCAST ORBIT -5 */if(fscanf( RinexEPP_file,"%s%s%s%s\n",t1,t2,t3,t4)==EOF) return 1 ;snv->idot=atof(t1);snv->wn=atof(t3);/* BROADCAST ORBIT -6 */if(fscanf( RinexEPP_file,"%s%s%s%s\n",t1,t2,t3,t4)==EOF) return 1 ;snv->tgd=atof(t3);return(0);}/* Using the datas from SV navigation Text to compute the SV position in ECEFCoor */ void GetSVPos(SVText *Text,double t,double *X,double *Y,double *Z){double mu=3.986005e14; /* Earth constant , unit m3/s2 */double n0,n; /* mean rate of angle */double tk, mk, ek1, ek2, ek;double Error=1.0e-12;double vk, phik;double deltau,deltar,deltai;double uk, rk, ik;double xk, yk;double omegak;double omegae=7.292115147e-5; /* earth self round rate of angle *//* compute SV mean rate of angle */n0 = sqrt(mu) / (Text->roota * Text->roota * Text->roota);n = n0 + Text->deltan ;/* compute time tk */tk = t - Text->toe;if(tk > 302400) tk -= 604800;else if(tk < -302400) tk +=604800;/* compute mean anomoly at measure time */mk = Text->m0 + n * tk;/* computer ek */ek1 = mk;do{ek2 = mk + Text->e * sin(ek1);if(fabs(ek2-ek1)<=Error ) break;ek1=ek2;}while(1);ek=ek1;/* computer vk *//* vk = atan(sqrt(1.0-Text->e *Text->e)*sin(ek)) / (cos(ek)-Text->e);*/vk = Get_atan(cos(ek)-Text->e,sqrt(1.0-Text->e *Text->e)*sin(ek));/* compute phik */phik = vk + Text->omega ;/* compute deltau,deltar,deltai */deltau = Text->cuc * cos(2.0*phik) + Text->cus *sin(2.0*phik) ;deltar = Text->crc * cos(2.0*phik) + Text->crs *sin(2.0*phik) ;deltai = Text->cic * cos(2.0*phik) + Text->cis *sin(2.0*phik) ;/* computer uk, rk and ik */uk = phik +deltau;rk = Text->roota *Text->roota * (1.0-Text->e * cos(ek)) + deltar ;ik = Text->i0 + deltai + Text->idot * tk ;/* compute xk,yk */xk = rk * cos(uk);yk = rk * sin(uk);/* compute omegak */omegak = Text->omega0 + (Text->omegadot - omegae)*tk -omegae * Text->toe ;/* compute ECEF */*X = xk * cos(omegak) - yk * cos(ik) *sin(omegak) ;*Y = xk * sin(omegak) + yk * cos(ik) *cos(omegak) ;*Z = yk*sin(ik) ;} /* End of Function *//* Caculate x=atan(y/z) */double Get_atan(double z,double y){double x;if (z==0) x=M_PI/2;else{if (y==0) x=M_PI;else{x=atan(fabs(y/z));if ((y>0)&&(z<0)) x=M_PI-x;else if ((y<0)&&(z<0)) x=M_PI+x;else if((y<0)&&(z>0)) x=2*M_PI-x;}}return x;}/********************************************** Convert date and time to GPS time***********************************************/int GetGPSTime(int year,int month,int day,int hour,int minute,double second,double *gpstime) {int dayofw,dayofy, yr, ttlday, m, weekno;/* Check limits of day, month and year */if (year < 1981 || month < 1 || month > 12 || day < 1 || day > 31) weekno = 0;/* Convert day, month and year to day of year */if (month == 1)dayofy = day;else{dayofy = 0;for (m=1; m<=(month-1); m++){dayofy += dinmth[m];if ( m==2 ){if (year % 4 == 0 && year % 100 != 0 || year % 400 == 0) dayofy += 1;}}dayofy += day;}/* Convert day of year and year into week number and day of week */ ttlday = 360;for (yr=1981; yr<=(year-1); yr++){ttlday += 365;if (yr % 4 == 0 && yr % 100 != 0 || yr % 400 ==0)ttlday += 1;}ttlday += dayofy;weekno = ttlday/7;dayofw = ttlday - 7 * weekno;*gpstime = (hour * 3600 + minute * 60 + second + dayofw * 86400); return weekno;}。