untitled7.m

GPS 的地球卫星运转情况仿真 还可以生成GPS信号

%计算所有的卫星位置后
%也要同时知道哪一颗卫星对用户来说是可见的
%因为只有根据可见的卫星才能定位
clear;clc;close all;
a=26560;  %a is the long radii of a ellipse
e=0.02;
E=[0:0.1:2*pi];
x=a*(cos(E)-e);
y=a*sqrt(1-e^2)*sin(E);
z=0*E;
TimeNow=0;%
%
%
global SatellitePosition %定义全局变量
SatellitePosition=ones(24,4)
figure(1);
DrawEarth(0); %调用程序五
hold on;
DtoR=2*pi/360;
A1=[32.8 92.8 152.8 212.6 272.8 332.8];
DrewSatelliteOrbit;%调用程序四

%平均近地角 
Ctable=[10 50 160 260;80 180 220 320;10 130 250 340;50 150 170 300;100 210 310 340;140 150 240 350];
Wx=ones(1,1);Wy=ones(1,1);Wz=ones(1,1);
simple=1;
for k=1:6
   A=A1(k)*DtoR;  %升交点的经度 
    B=55*DtoR;   %轨道倾角  
    for m=1:4
        C=Ctable(k,m)*DtoR+TimeNow*2*pi/24;;  %近地点的幅角
x=a*(cos(C)-e);
y=a*sqrt(1-e^2)*sin(C);
z=0*C;                                       
R3=[cos(A) -sin(A) 0;sin(A) cos(A) 0;0 0 1];
R1=[1 0 0;0 cos(B) -sin(B);0 sin(B) cos(B)];
R2=[cos(C) -sin(C) 0;sin(C) cos(C) 0;0 0 1];                                          
L1=length(E);
R312=R3*R1*R2;
Ans=R312*[x;y;z;]; %将卫星轨道坐标系转换到地球赤道坐标  
Wx=[Wx Ans(1,:)];
Wy=[Wy Ans(2,:)];
Wz=[Wz Ans(3,:)];
x1=Ans(1,:);
y1=Ans(2,:);
z1=Ans(3,:);
 DrawSatellite(x1,y1,z1,k)  %调用程序三
 temp=k*m;
 SatellitePosition(temp,:)=[x1 y1 z1 1];

 hold on 
    end
end


%用余弦定理地心到用户和用户到卫星之间的夹角
%如果两者之间的夹角小于90度，就认为该卫星是不可见的
EarthCenterPos=[0 0 0];
UserPosition=[6400 0 0];
for k=1:24
    Temp=SatellitePosition(k,1:3)-UserPosition;
    Dist1=Temp*Temp';
   Temp=UserPosition-EarthCenterPos; 
   Dist2=Temp*Temp';
   Temp=SatellitePosition(k,1:3)-EarthCenterPos;
   Dist3=Temp*Temp';
   JiaJiao=acos((Dist1+Dist2-Dist3)/2/sqrt(Dist1)/sqrt(Dist2));
   if(JiaJiao<=pi/2)
      SatellitePosition(k,4) =0;
   end
end
figure(2)
DrawEarth(0); %调用程序五
hold on;
DrewSatelliteOrbit;%调用程序四
for k=1:24
    if ( SatellitePosition(k,4) ==1)
        TempX= SatellitePosition(k,1);
        TempY= SatellitePosition(k,2);
        TempZ= SatellitePosition(k,3);
       DrawSatellite2(TempX,TempY,TempZ,simple);
        %调用程序三绘制出用户可见的卫星
    end
end
TempX=UserPosition(1);
TempY=UserPosition(2);
TempZ=UserPosition(3);
Cube=50;
BoxPlot3(TempX,TempY,TempZ,Cube,Cube,3*Cube,1);
%标记出用户的位置
hold on
SatellitePosition=[SatellitePosition;UserPosition 0];
%用变量CalUserPos来记录迭代计算结果
CalUserPos=CalUserPosition(SatellitePosition);%调用程序九
[m,n]=size(CalUserPos);
for(k=1:m)
    TempX=CalUserPos(1);
    TempY=CalUserPos(2);
    TempZ=CalUserPos(3);
    Cube=50;
    BoxPlot3(TempX,TempY,TempZ,Cube,Cube,3*Cube,3);
    hold on;
end