卫星导航实验报告.doc

______________________________________________________________________________________________________________

电子科技大学通信与信息工程学院

本科教学
卫星与导航系列实验
标准实验报告

课程名称：
l 定位与导航原理与应用
l 定位与导航工程

电子科技大学教务处制表

电子科技大学通信与信息工程学院

标准实验报告

实验名称：导航信号传输模型仿真

电子科技大学教务处制表

电子科技大学
实 验 报 告
学生姓名：侯玉皓 学号：2012019030016 提交日期：2015.6.24
得分

项目名称：

实时卫星位置、速度和时间解算（PVT解算）及结果分析
【实验目的】
1) 理解实时卫星位置解算在卫星导航解算过程中所起的作用，了解为完成卫星位置解算所需的条件；
2) 了解 GPS 时间、卫星的额定轨道周期的含义，了解星历的构成、周期及应用条件；
3) 了解 Doppler 频移的成因、作用以及根据已知条件预测 Doppler 频移的方法；
4) 了解 Doppler 频移的变化范围及其与卫星仰角之间的关系；
5) 能够根据实验数据编写求解 Doppler 频移的相关程序。
【实验原理】
实时卫星位置解算在整个导航解算过程中具有举足轻重的作用，通常我们为了获得接收机的地理位置，需要对卫星发射导航电文时的时间及运行速度有所了解，所以可以说，卫星的实时速度和时间是解算卫星实时位置的基础，而卫星的实时位置又是解算接收机三维位置坐标的基础。可见卫星实时位置、时间及速度在整个定位过程中的重要地位。一般来说要确定接收机的三维位置，需要同时解算出至少四颗卫星的实时位置。
卫星某一时刻发出的信号可以分为三部分：载波（L1）、测距码（C\A）、导航电文。
对GPS某颗卫星进行实时位置的解算，需要已知这颗卫星的星历和周内时，这些信息都包含在速率为 50bps 的导航电文中（图3.1中的数据码）。导航电文通过测距码（C/A 码）进行扩频，然后用扩频的信号去调制频率为 L1的正弦波载波，然后卫星将调制后的载波信号播发出去。其模型可以用如下公式表示：
数据码
测距码
载波

(3.1)

其中和是调制幅度，是精码和粗码，它们都是对数据码的扩频码，数据码经过扩频后分为两路进行调制。在本地接收机收到卫星信号后，通过剥离载波L1，还原其扩频之后的信号，然后按照导航电文的格式最终将数据码编译成导航电文（数据码）。
它可分为5个部分：
l 遥测字 导航电文前导+授权用户信息
l 交接字 由于P码数据长，不易捕获，需要C\A码的辅助来捕获
l 数据块1 钟修正参数、期龄、星期编号等
l 数据块2 星历（在子帧2、3中），用于计算卫星位置
l 数据块3 历书（子帧4、5中）提供卫星布局，健康状况等信息
本实验的一个重点在于通过导航电文来获取其卫星发射时间和星历从而得到卫星的实时位置。
GPS卫星在空间中的位置是时间t的函数，要计算卫星的位置首先要收集齐时钟和星历参数，然后需要确定卫星的发射时刻。导航信号的发射时刻可以通过导航电文在每一子帧的的交接字中的周内时计数器（Z计数）得知，通过该计数器可以得到估算的发射时刻。同时在子帧1中包含钟改正参数来对估算的时钟进行修正。
导航电文中的数据块2是卫星星历信息，数据块3是卫星的历书信息。星历主要向用户提供有关计算卫星运行位置的信息，而历书主要向用户提供GPS卫星的概略星历及卫星的工作状态等。
在数据块2（子帧2和3）中包含有许多重要星历参数。星历数据参数如下表所示（一颗卫星）:

表3.1 星历数据参数
参数（子帧2中）
意义
参数（子帧3中）
意义

轨道半径的正弦调和修正幅度

倾斜角的余弦调和修正项幅度
Δn
计算值的平均移动误差
Ωe
在每星期历元轨道平面上升点的经度

参考时间的平近点角

倾斜角的正弦调和修正项的幅度

纬度辅角的余弦调和修正项幅度

参考时间的倾斜角

离心率

轨道半径的正弦调和修正项幅度

纬度辅角的正弦调和修正项幅度

近地点的幅度

长半轴平方根
Ω
赤径的速度

参考时间星历
IODE
数据与星历发布号

IDOT
倾斜角的速度

后续导航解算单元根据导航电文中相应的参数进行卫星位置解算、各种实时误差的消除、本地接收机位置解算以及定位精度因子（DOP）的计算等工作。也就是说，根据收到的导航电文，接收端就可以通过相关公式计算出发送电文时刻卫星的大致位置，这对于解算出接收机的地理位置尤为重要。
卫星的角速度和切向速度可以通过卫星轨道模型来进行估计。GPS卫星的额定轨道周期是半个恒星日，或者说 11 小时 58 分钟 2.05秒；各轨道接近于圆形，轨道半径（即从地球质心到卫星的额定距离）大约为26560km。由此可以计算得到卫星的运行的角速度和切向速度（如图3.3）：
ω=2π/(11*3600+58*60+2.05)≈0.0001458 rad/s (3.2)
然后通过角速度ω和已知的轨道半径rs（26560km）计算切向速度
vs=rs*ω≈26560km*0.0001458≈3874m/s （3.3）
本实验的另一个重要内容是 Doppler 频移的预测，即通过对卫星进行相隔时间为 1s 的多点测量（本实验给出了三点），进而估计 Doppler 频移。
l Doppler 频移产生原理：
由于卫星与接收机存在相对径向运动，因此信号接收频率会随这种相对运动而发生偏移，称为 Doppler频移，其直接表现是接收机接收到的卫星信号不在 L1 频点（1575.42MHz）上，而是在其基础上叠加了一个大约为-5KHz到+5KHz 左右的频率偏移。Doppler频移的存在给前端相关器进行频域搜索和卫星信号捕获带来了诸多困难。如果事先能够估算出Doppler频偏，就会降低捕获难度，缩短捕获时间，进而缩短接收机的启动时间。接收机的启动时间是衡量接收机性能好坏的重要参数之一，而实现卫星信号的快速捕获，缩短接收机的启动时间也是目前GNSS领域的热点问题。
l Doppler 频移的计算：
已知卫星位置和本地接收机的初始位置，根据空间两点间的距离公式，可以得出卫星到接收机的距离d。设某卫星在短时间 t 内经过点 S1、S2，接收机到点 S1、S2 的距离分别为d1、d2,。在 t 相对较小（本实验取 t=1s）的情况下，我们可近似认是卫星与接收机在 t 时间内的平均相对径向运动速度，再将此速度转换为频率的形式就可以得到 Doppler 频移的估计值。
设本地接收机的初始位置为 R（xr，yr，zr），卫星所经历的空间两点的坐标分别为 S1 （x1，y1，z1）、S2（x2，y2，z2），间隔时间为 t，卫星与接收机平均相对径向运动速度为 vd，光速为 c，Doppler 频移为 fd，则 Doppler 频移预测的具体公式如下所示：
d1
d2
22
L1
vd
22

图 3.1
其中，为一颗卫星不同时刻到接收机的伪距，为两个时刻之间的径向速度。
Doppler 频移与卫星的仰角有关，卫星仰角越大，Doppler 频移越小。当卫星的仰角为 90 度（即卫星在接收机正上方的天顶上）时，理论上 Doppler 频移为零。本实验根据卫星位置和本地接收机的初始位置算出卫星的仰角，来验证 Doppler 频移同卫星仰角之间的关系.

【实验步骤】
1) 查阅资料建立相应模型，在C/C++或者Matlab平台上根据星历数据及其定义实现对卫星实时位置的解算；
2) 运行主程序以取得可视卫星的实时导航数据（如 GPS 时间、各颗卫星的星历等）。将实验平台仪器的USB端口接入电脑，待驱动安装成功后，打开实验一程序；
3) 在“选择 GPS 时刻”列表框的下拉菜单中，任意选择一个GPS 时刻。（注：北斗和 GPS 系统由于存在系统时差而具有不同的周内时。这里的 GPS 时刻，对于GPS卫星指其系统周内时，对于北斗卫星则表示将北斗的周内时加上系统时差换算之后的 GPS 系统周内时）；
4) 在“所选时刻可视卫星星历”列表框中出现所选时刻天空中所有可视卫星的星历信息，如图 3.6 所示。选定一颗卫星，将“所选时刻可视卫星星历”中该卫星对应的参数输入到1）中的解算代码中，计算卫星位置。
5) 在“选择卫星号”列表框的下拉菜单中，出现所选时刻天空中所有可视 卫星的序号。北斗卫星的编号从 101 开始，即北斗 1 号星的编号为 101。选择与4）中对应的卫星序号，在“卫星位置信息”中会列出所选时刻该卫星的实时位置如图3.7。对比该位置与之前代码解算的结果。并将其记录在表格中（表格一）；
6) 在“卫星位置信息”列表框中同时会出现所选卫星在所选 GPS 时刻一秒和两秒后的所对应的 ECEF 坐标系下的三维坐标以及接收机在 ECEF 坐标系下的初始位置坐标，这些数据用于求解 Doppler 频移，根据附表记录其值（表格一）；
7) 在“卫星位置信息”列表框中还会出现该卫星在 11 小时 58 分后的ECEF 位置坐标，这是根据卫星在所选 GPS 时刻的星历数据推算出来的，用以验证卫星的额定轨道周期。根据附表记录其值（表格一）；
8) 根据步骤 6）记录的数据，在C/C++或Matlab环境下编写代码，实现对Doppler 频移估值的求解，将所得数据记录在附表中（表格一）；
9) 重复前面实验，记录并解算出所选时刻天空中所有可视卫星的相关数据，按附表格式将所得数据记录下来（表格二）；
10) 重复前面实验，比较并分析不同时刻同一卫星的仰角、ECEF 坐标系下的坐标以及Doppler频移的差异和同一时刻不同卫星仰角、坐标及Doppler频移差异；
11) 重复步骤 2 到步骤 11，选择不同时间段进行记录、求解、分析。
【核心程序代码】
%% Find satellite's position ----------------------------------------------

%Restore semi-major axis
a = eph(prn).sqrtA * eph(prn).sqrtA;

%Time correction
tk = check_t(time - eph(prn).t_oe);

%Initial mean motion
n0 = sqrt(GM / a^3);
%Mean motion
n = n0 + eph(prn).deltan;

%Mean anomaly
M = eph(prn).M_0 + n * tk;
%Reduce mean anomaly to between 0 and 360 deg
M = rem(M + 2*gpsPi, 2*gpsPi);

%Initial guess of eccentric anomaly
E = M;

%--- Iteratively compute eccentric anomaly ----------------------------
for ii = 1:10
E_old = E;
E = M + eph(prn).e * sin(E);
dE = rem(E - E_old, 2*gpsPi);
if abs(dE) < 1.e-12
% Necessary precision is reached, exit from the loop
break;
end
end
%Reduce eccentric anomaly to between 0 and 360 deg
E = rem(E + 2*gpsPi, 2*gpsPi);
%Compute relativistic correction term
dtr = F * eph(prn).e * eph(prn).sqrtA * sin(E); %Ïà¶ÔÂÛÐÞÕýÏî

%Calculate the true anomaly
nu = atan2(sqrt(1 - eph(prn).e^2) * sin(E), cos(E)-eph(prn).e);

%Compute angle phi
phi = nu + eph(prn).omega;
%Reduce phi to between 0 and 360 deg
phi = rem(phi, 2*gpsPi);

%Correct argument of latitude
u = phi + ...
eph(prn).C_uc * cos(2*phi) + ...
eph(prn).C_us * sin(2*phi);
%Correct radius
r = a * (1 - eph(prn).e*cos(E)) + ...
eph(prn).C_rc * cos(2*phi) + ...
eph(prn).C_rs * sin(2*phi);
%Correct inclination
i = eph(prn).i_0 + eph(prn).iDot * tk + ...
eph(prn).C_ic * cos(2*phi) + ...
eph(prn).C_is * sin(2*phi);

%Compute the angle between the ascending node and the Greenwich meridian
Omega = eph(prn).omega_0 + (eph(prn).omegaDot - Omegae_dot)*tk - ...
Omegae_dot * eph(prn).t_oe;
%Reduce to between 0 and 360 deg
Omega = rem(Omega + 2*gpsPi, 2*gpsPi);

%--- Compute satellite coordinates ------------------------------------
satPositions(1, satNr) = cos(u)*r * cos(Omega) - sin(u)*r * cos(i)*sin(Omega);
satPositions(2, satNr) = cos(u)*r * sin(Omega) + sin(u)*r * cos(i)*cos(Omega);
satPositions(3, satNr) = sin(u)*r * sin(i);
【结果记录与分析】
1、 按附表格式整理实验数据，并整理所编程序。（以一颗卫星为例）
附表一
卫星时间
卫星序号:10

仰角
Doppler 频移

解算坐标
ECEF坐标
308640

-6898114.365083
-6898114.840231
69.366891°

4.836185542725942 KHz

10231973.056976
10231973.400231
11426817.711968
11426816.812402
308641
-6899460.241363
-6899460.726453
69. 366891°

4.836197123663138 KHz

10232948.980243
10232949.326035
11424268.230451
11424267.987032
308642

-6890805.435930
-6890805.874032
69. 366891°

4.836208705610145 KHz

10233924.760354
10233924.785343
11421718.243299
11421718.207653
11小时58 分后的坐标 ECEF 位置坐标
X:-6898413.749678
Y:10229384.563445
Z:11423157.630254
2、对同一时刻不同卫星的 Doppler 频移进行比较，根据实际数据得出卫星仰角Doppler频移之间的关系。
附表二
GPS 时间
可视卫星序号
仰角
Doppler 频移
308641
10
69. 366891°

4.836197123663138 KHz
308641
20
16.607849°
4.893789169568091 KHz
308641
32
44.365742°
4.869897124663638 KHz

由上表可得: 卫星的仰角越大，其多普勒频移越小。
3、比较并分析不同时刻同一卫星的仰角、ECEF 坐标以及 Doppler 频移的差异。
附表三
GPS世间
卫星序号
ECEF坐标
仰角
Doppler 频移

308641
10

-6898413.749678
69. 366891°

4.836197123663138 KHz

10229384.563445
11423157.630254
308642
10

-6898413.209384
69. 366891°

4.836208705610145 KHz
10229384.048932
11423157.103970
由上表分析得出什么结论？
答：在时间相差很近的两个时刻里，卫星的坐标，仰角以及多普勒频移变化都很小，几乎不变化。仰角越大，多普勒频移越小。

4、比较当前时刻卫星的 ECEF 坐标及由当前星历推算出的该卫星在 11 小时 58 分后的 ECEF 坐标，思考二者为何只是大致相同而不是绝对一致。
当前时刻与 11 小时 58 分后 2 号卫星的位置信息由下表给出：
GPS世间
卫星序号
ECEF坐标
仰角
308641
10

X:-6898413.749678
13.058204°

Y:10229384.563445
Z:11423157.630254
351721

10

X:-6898413.034178
13.058204°

Y:10229384.984645
Z:11423157.204354
分析其原因可能是：
可能是卫星的严格周期为11 小时 58 分2秒，而不是简单的11 小时 58 分；
可能是因为接收机的计时误差；可能因为电离层、对流层误差导致。
【意见与建议】
对本实验过程及方法、手段的改进建议：
希望可以把那款软件拷下来，这样方便学生在寝室调整代码。

Welcome To
Download !!!

欢迎您的下载，资料仅供参考！

精品资料