星敏感器姿态确定仿真综合实验.docx
- 文档编号:25918785
- 上传时间:2023-06-16
- 格式:DOCX
- 页数:21
- 大小:895.19KB
星敏感器姿态确定仿真综合实验.docx
《星敏感器姿态确定仿真综合实验.docx》由会员分享,可在线阅读,更多相关《星敏感器姿态确定仿真综合实验.docx(21页珍藏版)》请在冰豆网上搜索。
星敏感器姿态确定仿真综合实验
北京航空航天大学
专业综合实验报告
学院宇航学院班级111514
学号11151146姓名高荣荣
指导老师王海涌
2015年1月3日
星敏感器姿态确定仿真综合实验
摘要:
通过对电子星图模拟器和星敏感器PC仿真平台的操作,实现星敏感器姿态确定,以及借助星象天文馆,来实现电子星图模拟器对星的标定。
完成星敏感器系统仿真。
关键词:
星敏感器定姿星图模拟星图姿态矩阵
一实验目的
通过电子星图模拟器(ESS)和星敏感器PC仿真平台的操作,熟悉星图模拟和星敏感器姿态基本流程及各模块功能,完成星敏感器系统仿真。
二实验原理
1.星图模拟原理
(1)星图模拟系统是一种近似模拟星空的仿真系统。
为星敏感器算法调试、星敏感器产品测试及天文导航半物理仿真系统运行提供标准的星图输入,并提供已知参考星光矢量及星像中心的理想映射坐标。
分光学物理星模、电子星模和计算机软件星模三种。
(2)星图模拟系统的实现
旋转关系:
OZ是光轴指向:
星图模拟是根据光轴指向及旋角(
γ)所确定的既定视场,将其范围内的星空目标映射到星敏感器CCD面阵上并模拟出目标图像的过程。
其中,需要确定第二赤道坐标系、航天器本体坐标系、星敏感器坐标系的转换矩阵。
星敏感器固联(安装矩阵为常数阵),那么只考虑第二赤道坐标系和星敏感器坐标系之间的转换关系。
令O-UVW为第二赤道坐标系,令O’-XYZ为星敏感器坐标系。
那么星光矢量在两个坐标系下的分量列阵的关系可以表示为:
[X,Y,Z]T=Tsi[U,V,W]T,其中Tsi为转换矩阵。
2.星敏感器定姿基本原理
OsXsYsZs—星敏感器坐标系
Ouv—CCD成像面坐标系
OsO之间距离f为光学透镜的焦距
由图中的几何关系可得:
第n颗星的单位矢量在星敏感器坐标系中的分量列阵:
或,根据星像点质心坐标直接计算得单位矢量:
式中
为星敏感器测量误差矢量。
星敏视场内n颗恒星在星敏感器坐标系Sb中的单位矢量坐标分别为(X1,Y1,Z1),(X2,Y2,Z2),…,(Xn,Yn,Zn)。
经过星图匹配,获知这n颗恒星在惯性空间i系中单位矢量坐标:
[U1,V1,W1],[U2,V2,W2],…,[Un,Vn,Wn],是由星表中赤经赤纬经过直角坐标转换求得。
则存在以下矩阵转换关系式:
上式各矩阵赋予命名,对应关系式:
S=CTis
当观测星数n等于3颗时,星历矩阵C3×3为非奇异矩阵,那么转移矩阵Tis=C-1S。
n>3时,可以采用最小二乘法得到计算式:
Tis=(CTC)-1CTS
姿态矩阵Tsi=TsbTblTli,则载体相对于发射点惯性系l系的姿态矩阵Tbl=(Tsb)-1(Tis)-1(Tli)-1=(TliTisTsb)-1,这里Tsb为星敏感器固联常量安装矩阵,
发射点惯性坐标系l系相对于赤道惯性系i系的转移矩阵:
其中:
A为轴xl的方位角;
S为发射时刻的格林尼治恒星时,或发射时刻的春分点的格林时角GHAΥ;
(λ,φ)为发射点经纬度。
如果安装矩阵与火箭本体坐标系方向重合,则安装矩阵Tbs=单位阵I
那么:
Tbl=(TliTis)-1
欧拉角法表示火箭姿态
发惯系l,原点在发射点L,轴yl在当地铅垂线向上,轴xl在当地水平,且在名义射击平面内。
发射点惯性坐标系Lxlylzl转动3次转到本体坐标系oxbybzb。
转动顺序:
那么,从Sl到Sb的坐标变换矩阵:
矩阵元素对照得主值:
姿态角θ和ψ的取值范围都在[-90°,90°];
φ的取值范围都在[-180°,180°]
然后基于主值再获得真值。
三实验步骤
1.设定观星地址和观星时间,采用静基座天顶观测方式,解算出相机光轴的赤经、赤纬和旋角。
(1)设定观星地点为北京市延庆县马匹营村,北纬40°30′14.571″,东经116°04′20.815″
(2)设定北京时间为2011年8月20日22h22min222s
(3)2011年8月20日当天世界时0时对应的恒星时为12h16m40.6407s;
=26.3156
由此得RA=300.2372
Dec=40.5041
指定γ=0
2开启ESS,进行参数设置,选择“星图模拟”“单帧重复”菜单命令,将(
,
,
)计算值填入ESS的人机界面,观察“发送的星图”,此时ESS持续不断地将模拟星图数字视频传输给星敏感器PC仿真平台。
设置:
指定仿真星等阈值:
6
峰值饱和灰度对应星等:
8
得到模拟星图:
其中,星模软件使用方法见附录;
3.在PC仿真平台,存储接收的模拟星图,并选定六颗颗亮星。
结果如图所示。
得到亮星的坐标,星等,赤经和赤纬:
310.3580418
45.28035837
1.3
243.0664099
873.2034322
305.5571248
40.25669994
2.3
514.3430786
718.3319065
311.5528582
33.97025572
2.6
810.0597353
993.0520476
296.2438958
45.13073903
3.0
272.523317
369.3284183
303.4079165
46.74135835
4.0
190.4634936
621.7853499
299.0766043
35.08341932
4.0
783.1009519
466.7343707
4.运行虚拟天文馆软件stellarium.exe,选定
(1)中地址和观星时间,目
视对比模拟星图和stellarium界面,找出这六颗选定亮星的匹配星,记录其星号、赤经、赤纬和星等属性信息,计算在春分点赤道惯性系i系下的恒星方位矢量。
运行虚拟天文馆结果如图所示:
根据经纬度关系,并借助主观肉眼判断,匹配情况如下:
星名/星等
赤经/(J2000)
赤纬
天津四/1.25(脉动变星)
310.25/20h41m25s
45.16
天津一/2.20(恒星)
305.52/20h22m13.7s
40.17
天津九/2.45(恒星)
311.50/20h16m41s
34.01
δCyg-18Cyg-HIP97165A/2.9(双星)
296.25/19h45m20s
45.07
O1Cyg-31Cyg-V0695/3.8(食双星系统)
303.404/20h13m37.9s
46.46
ηCyg-21Cyg-HIP98110/3.85(恒星)
299.011/19h56m45s
35.06
注意事项:
(1)在运用虚拟天文馆软件时,要注意调节观察视角,一般为20度左右事宜。
本次实验,调节观察视角读数为22.8°。
(2)根据纬度关系显示,将虚拟天文馆的视场调节到指定的范围内,以便容易找到星星。
本次实验,大约根据纬度关系调节到45°的范围内。
(3)在记录赤经、赤纬的信息时,有两种表示的方法:
一种是J2000,另一种是J2011.2。
两者之间有非常细微的差距。
因为地轴的进动,使得春分点不断西退,而赤经赤纬又是依靠春分点的位置来确定,所以随着春分点位置的不断变化,赤经赤纬也会不断变化。
J2000版指的是以2000年春分点作为基准天体的坐标,日期版就是以现在的春分点作为基准的坐标。
差距是非常细微的,对一般的爱好者来说,忽略这个微小的差距也无妨,但是对于专业的观测来说,这点误差是要命的。
在本实验中,记录的是J2000下的赤经赤纬信息。
记录星等信息时,有目视星等和绝对星等两种星等信息。
绝对星等是假定把恒星放在距地球10秒差距(32.6光年)的地方测得的恒星的亮度,用以区别于视星等。
它反映天体的真实发光本领。
由于本实验中记录的是目视星等。
5.给定该相机焦距f的标定值,查找星敏感器PC平台下的日志文件获取6颗选定亮星的星象质心,是在星图预处理环节由灰度重心法计算获得的,计算星敏感器本体系s系下的星光矢量。
(1).焦距f=2889.2643,单位:
pixel,主点O的值(508.3,513.3)pixel。
(2).阵列平面坐标系和星敏感器本体坐标系的关系如图所示:
其中星敏感器本体系s系下的星光矢量由下式计算得到
6.基于i系下的恒星方位矢量和s系下的星光矢量,计算出观星时刻星敏感器相对于i系的姿态矩阵,并换算出星敏感器的光轴赤经赤纬和旋角,对比星敏感器PC仿真平台的计算结果。
Matlab程序见附录;
得到赤经赤纬及旋角为:
RA=300.2372
Dec=40.5041
四实验结论及误差分析:
1.以观测者所在位置,运用查表法得到的赤经赤纬为Ra=300.2372,Dec=40.5041,Gamma=0;再运用模拟星图ESS和虚拟天文馆得到的赤经赤纬为RaDeter=300.2385,DecDeter=40.5180,GammaDeter=-0.3402;所以,经运算得到绝对误差和相对误差如下:
绝对误差为:
|Ra-RaDeter|=0.0013
|Dec-DecDeter|=0.0139
|Gamma-GammaDeter|=0.3402
相对误差为:
η(Ra)=0.004%
η(Dec)=0.031%
η(Gamma)=0.340%
2.在星图模拟和虚拟天文馆的对比找星过程中,也存在一定的误差,主要是,星等误差,赤经和赤纬对比误差等等。
误差定量计算如下:
(1)赤经绝对误差:
绝对误差
310.25/20h41m25s
310.3580418
0.100
305.52/20h22m13.7s
305.5571248
0.030
311.50/20h16m41s
311.5528582
0.050
296.25/19h45m20s
296.2438958
0.010
303.404/20h13m37.9s
303.4079165
0.003
299.011/19h56m45s
299.0766043
0.065
分析:
最大误差为0.1°,因此,在合理误差范围之内,说明星星确认正确;
(2)赤纬绝对误差:
绝对误差
45.16
45.28035837
0.12
40.17
40.25669994
0.08
34.01
33.97025572
0.04
45.07
45.13073903
0.06
46.46
46.74135835
0.28
35.06
35.08341932
0.02
分析:
最大误差为0.28,计算相对误差为
0.6%,误差总体较小,因此也在合理范围之内。
说明确认星星正确;
(3)星等绝对误差
绝对误差
1.25
1.3
0.05
2.20
2.3
0.10
2.45
2.6
0.05
2.9
3.0
0.10
3.8
4.0
0.20
3.85
4.0
0.15
分析:
最大误差为0.20,计算相对误差一定比较大,但是不影响星星判定,星等的误差,大多还是由于虚拟天文馆的制作造成的,因此,从整体绝对误差上看,星星的误差还是在合理范围之内的。
因此,认为确定星星正确。
综上所述:
经过赤经,赤纬,星等的误差对比,确认模拟星图所找的星星在虚拟天文馆中得到一一对应。
误差可视为合理误差。
五实验总结
为期一周的综合实验到这里结束了,首先,要感谢王海涌老师的悉心教导,在王老师的帮助指导下,将天文导航这门课程的精髓在实验中得到了应用,收获颇丰。
本实验,主要就图模拟原理和星敏感器定姿原理展开了实验上的运用。
我们先利用观星地址和观星时间解算出光轴的赤经、赤纬和旋角。
再计算出由春分点确定的子午线和格林尼治子午线的时角。
第一步是根据观星时间查紫金山天文台出版的《天文年历》确定当日格林尼治零时的春分点时角,第二步由观星时间确定观星时刻春分点所在子午线和格林尼治子午线的时角。
进而计算由春分点第二赤道坐标系i系和地球坐标系e系间的基元旋转矩阵,由假定测着所在地理位置计算地理坐标系t系向地球坐标系e系的转移矩阵,从而确定地理坐标系t系向第二赤道惯性坐标系i系的转移矩阵,即可算出光轴的赤经、赤纬和旋角。
通过计算观星时刻的世界时计算儒略日从而计算地球转角θ、TT儒略世纪数t,最后得到格林尼治平恒星时GMST,再用上测者的位置信息,即可求得t系向i系的转移矩阵。
这种算法的好处是避免了查《天文年历》,比较方便。
在运行虚拟天文馆时,一定要注意光轴是和地面垂直这一条件。
按住键盘的向上的按键一直到头,则界面显示的就是正上方的星空。
同时也要注意调整视角的大小,因为星图模拟器的视角是20×20,要调整虚拟天文馆的FOV至合适值,大约是22.8°。
在使用虚拟天文馆时,要注意用已知信息和实际信息进行匹配,如,运用在模拟星图中的赤经和赤纬的信息在虚拟天文馆中进行信息对比,如果位置正确,注意调节视角等再进行匹配。
直至匹配完成。
附录:
Matlab程序:
1电算法
%%假定在北京古观象台(北纬40°30′14.571″,东经116°04′20.815″)进行观星
Lat=(39+54/60)*pi/180;
Lon=(116+18/60)*pi/180;
%%电算法解算格林尼治平恒星时
%时间为北京时2011年8月20日22点22分22.2秒
formatlong;
year=2011;%按UTC和UT1近似相等处理
month=8;
day=22;
hour=22;
minute=22;
second=22.2;
%估计预测DUT1
a=floor((14-month)/12);
y=year+4800-a;
m=month+12*a-3;
JDN=day+floor((153*m+2)/5)+floor((365+97/400)*y)-32045.5;
%求地球自转角θ
UTC_TAI=-34;%需查2011年的bulletina-xxiv-017.txt,文档中已给
JdUT1=JDN+(hour+minute/60+second/3600)/24;
Tu=JdUT1-2451545;
Theta=2*pi*(0.7790572732640+1.00273781191135448*Tu)*180/pi*3600;
%计算力学时TT,求GMST
secTT=second-UTC_TAI+32.184;%力学时0hTT
JdTT=JDN+(hour+minute/60+secTT/3600)/24;
t=(JdTT-2451545);
t=t/36525;
GMST=0.014506+Theta+4612.15739966*t+1.39667721*t.*t-0.00009344*t.^3+0.00001882*t.^4;%注意:
单位是角度秒,arcsecond
GMSTarc=deg2rad(mod(GMST/3600,360))%得到格林尼治平恒星时
temp=rem(GMST/15,24*3600);%GMST/15则为时间秒
H=floor(temp/3600);%恒星时“h”
M=floor(mod(temp,3600)/60);%恒星时“m”
S=mod(temp,3600)-M*60;%恒星时“s”
%%计算天顶方向在惯性系i系下的经纬度
%旋转矩阵
Fai=pi/180;
Pfai=Fai;
Thet=Fai;
L1=[100;
0cos(Fai)sin(Fai);
0-sin(Fai)cos(Fai)];
L2=[cos(Pfai)0-sin(Pfai);
010;
sin(Pfai)0cos(Pfai)];
L3=[cos(Thet)sin(Thet)0;
-sin(Thet)cos(Thet)0;
001];
%转移矩阵
Cts=L1*L2*L3;
Cet=[-sin(Lon)-sin(Lat)*cos(Lon)cos(Lat)*cos(Lon);
cos(Lon)-sin(Lat)*sin(Lon)cos(Lat)*sin(Lon);
0cos(Lat)sin(Lat)];
Cie=[cos(GMSTarc)-sin(GMSTarc)0;
sin(GMSTarc)cos(GMSTarc)0;
001];
R=Cie*Cet*Cts*[0;0;1];
%解算赤经赤纬
Alpha=degrees2dms(rad2deg(atan(R(2,:
)/R(1,:
)))),
Delta=degrees2dms(rad2deg(asin(R(3,:
)))),
2.matlab下打开raw文件;
clc;clearall;
fid=fopen('0.raw','r');
mk=1024;
nk=1024;
z=uint16(zeros(mk,nk));
fori=1:
mk
forj=1:
nk
z(i,j)=fread(fid,1,'uint16');
end
end
g=double(max(max(z)));
temp=double(255/g);
fori=1:
mk
forj=1:
nk
A(i,j)=z(i,j)*temp;
end
end
A=uint8(A);
figure;
imshow(A)
(3)计算赤经赤纬:
u0=508.3;
v0=513;
f=2889.2643;
star1=[113.649975131.88848347];
star2=[134.801883248.04160012];
star3=[135.906371147.15658602];
star4=[112.277795631.78433636];
star5=[135.159962441.78273636];
star6=[125.708962243.18811131];
Vstar1=[cosd(star1
(1))*cosd(star1
(2)),sind(star1
(1))*cosd(star1
(2)),sind((star1
(2)))];
Vstar2=[cosd(star2
(1))*cosd(star2
(2)),sind(star2
(1))*cosd(star2
(2)),sind((star2
(2)))];
Vstar3=[cosd(star3
(1))*cosd(star3
(2)),sind(star3
(1))*cosd(star3
(2)),sind((star3
(2)))];
Vstar4=[cosd(star4
(1))*cosd(star4
(2)),sind(star4
(1))*cosd(star4
(2)),sind((star4
(2)))];
Vstar5=[cosd(star5
(1))*cosd(star5
(2)),sind(star5
(1))*cosd(star5
(2)),sind((star5
(2)))];
Vstar6=[cosd(star6
(1))*cosd(star6
(2)),sind(star6
(1))*cosd(star6
(2)),sind((star6
(2)))];
uv1=[895.207618674.50844256];
uv2=[66.40130886885.4437202];
uv3=[105.8599846929.8961347];
uv4=[894.681476713.83111223];
uv5=[382.7642809940.7513796];
uv6=[341.0580967581.6978552];
VsensorBody1=[uv1
(2)-v0,-(uv1
(1)-u0),f]/sqrt((uv1
(1)-u0)^2+(uv1
(2)-v0)^2+f*f);
VsensorBody2=[uv2
(2)-v0,-(uv2
(1)-u0),f]/sqrt((uv2
(1)-u0)^2+(uv2
(2)-v0)^2+f*f);
VsensorBody3=[uv3
(2)-v0,-(uv3
(1)-u0),f]/sqrt((uv3
(1)-u0)^2+(uv3
(2)-v0)^2+f*f);
VsensorBody4=[uv4
(2)-v0,-(uv4
(1)-u0),f]/sqrt((uv4
(1)-u0)^2+(uv4
(2)-v0)^2+f*f);
VsensorBody5=[uv5
(2)-v0,-(uv5
(1)-u0),f]/sqrt((uv5
(1)-u0)^2+(uv5
(2)-v0)^2+f*f);
VsensorBody6=[uv6
(2)-v0,-(uv6
(1)-u0),f]/sqrt((uv6
(1)-u0)^2+(uv6
(2)-v0)^2+f*f);
S=[VsensorBody1;VsensorBody2;VsensorBody3;VsensorBody4;VsensorBody5;VsensorBody6];
C=[Vstar1;Vstar2;Vstar3;Vstar4;Vstar5;Vstar6];
Tsi=(S'*S)\S'*C;%采用最小二乘法
%*****浮点计算条件下,解算星敏光轴的赤经赤纬和旋角(α0,δ0,γ)******
DecDeter=asin(Tsi(3,3));%毫无疑问!
真值就是主值,单位:
弧度
ref=sind(89.99);
ifabs(Tsi(3,3))>ref%赤纬接近±90°
RaDeter=0;%赤纬位于天南极或天北极,赤经可以是任意值,直接令其指向春分点吧!
else
ifTsi(3,2)/cos(DecDeter)>ref%判断sin(α0),α0接近90°
RaDeter=90;
elseifTsi(3,2)/cos(DecDeter)<-ref%α0接近270°
RaDeter=270;
else
RaDeter=atand(Tsi(3,2)/Tsi(3,1))%主值,单位:
“°”
ifTsi(3,1)<0
RaDeter=RaDeter+180%真值,单位:
“°”
else
ifTsi(3,2)<0
RaDeter=RaDeter+360%真值,单位:
“°”
end
end
end
end
%计算旋角γ的主值和真值
ifTsi(1,3)/cos(DecDeter)>ref
GammaDeter
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 敏感 姿态 确定 仿真 综合 实验