STK与matlab交互 Astrogator模块 （12）

一、背景介绍

本文结合控制的背景，采用小推力模型的条件下，实现每个时刻给卫星施加一个有限推力的脉冲，即发动机的采用恒推力发动机，同时考虑卫星的质量损耗，控制的变量为推力的方向以及发动机的推力。根据（11）中的数据作为文章的开头，首先获取初始化卫星，分析整个区间的时间序列。

clear;clc
uiApplication = actxGetRunningServer('STK12.application');
% Get our IAgStkObjectRoot interface
global root
root = uiApplication.Personality2;
checkempty = root.Children.Count;
if checkempty ~= 0
    root.CurrentScenario.Unload
    root.CloseScenario;
end
%% 根据你的需要设定场景的名称
root.NewScenario('chapter2');
StartTime = '29 Sep 2023 04:00:00.000';    % 场景开始时间
StopTime = '31 Oct 2023 04:00:00.000';     % 场景结束时间
root.ExecuteCommand(['SetAnalysisTimePeriod * "',StartTime,'" "',StopTime,'"']);
root.ExecuteCommand(' Animate * Reset');
%% ----------------------------------------------
%% A卫星
Sat_Name='A';
satellite=root.CurrentScenario.Children.New('eSatellite', Sat_Name);
satellite.Propagator;
sma=7071.393;
Ecc=0.001;
Inc=20;
w=0;
RAAN=0;
TA=359.999473339175;
satellite.SetPropagatorType('ePropagatorAstrogator'); 
satellite.Propagator;
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList Initial_State Propagate']);
%% 初始化卫星参数
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.CoordinateType Modified Keplerian']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Epoch ',StartTime,' UTCG']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.ElementType "Kozai-Izsak Mean"']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.sma ',num2str(sma),' km']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.ecc ',num2str(Ecc)]);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.inc ',num2str(Inc),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.w ',num2str(w),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.TA ',num2str(TA),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.RAAN ',num2str(RAAN),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' RunMCS']);
%% B卫星
Sat_Name2='B';
satellite2=root.CurrentScenario.Children.New('eSatellite', Sat_Name2);
satellite2.Propagator;
sma=7071.393;
Ecc=0.001;
Inc=20;
w=0;
RAAN=0;
TA2=359.998946678351;
satellite2.SetPropagatorType('ePropagatorAstrogator'); 
satellite2.Propagator;
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList Initial_State Propagate Maneuver']);
%% 初始化卫星参数
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.CoordinateType Modified Keplerian']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Epoch ',StartTime,' UTCG']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.ElementType "Kozai-Izsak Mean"']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.sma ',num2str(sma),' km']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.ecc ',num2str(Ecc)]);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.inc ',num2str(Inc),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.w ',num2str(w),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.TA ',num2str(TA2),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.RAAN ',num2str(RAAN),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' RunMCS']);
%% 生成LVLH坐标系 作为一个循环
MissionTime = '29 Sep 2023 04:05:00.000';    % 任务结束时间
data=root.ExecuteCommand(['Report_RM */Satellite/',Sat_Name,' Style "LVLH" TimePeriod "',StartTime,'" "',MissionTime,'" TimeStep 0.1']);%步长十分钟
LineNum=data.Count;
for i=1:LineNum-2
    struct=regexp(data.Item(i),',','split');
    Time{i}=struct{1};
end

注意最后的LVLH报告STK是没有自带的，需要自己定义，整个上面代码的目的是获得目标在LVLH坐标系下的时间序列，为下一步时间段的分析打好基础

二、任务仿真

在完成初始化后，本文采用控制的方式，每间隔0.1s改变一次发动机的推力和方向（该数值应该通过卫星的LVLH坐标系的位置和速度六个状态量计算出来），在这里直接给出猜想值。

%% 生成控制数据 （通过自带控制器实现）
Control=zeros(3,LineNum-3);% 第一项为控制推力的大小，第二项 控制推力的方向
for i=1:LineNum-3
    Control(1,i)=200;%恒推力模型
    Control(2,i)=30;
    Control(3,i)=20;
end

任务仿真开始，首先获得任务时间开始时的B卫星在A卫星LVLH坐标系下的位置和速度

 data2=root.ExecuteCommand(['Report_RM */Satellite/',Sat_Name,' Style "LVLH" TimePeriod "',Time{i},'" "',Time{i+1},'" TimeStep 0.1']);%步长十分钟ss=regexp(data2.Item(1),',','split');
x=str2double(ss{2});
y=str2double(ss{3});
z=str2double(ss{4});
vx=str2double(ss{5});
vy=str2double(ss{6});
vz=str2double(ss{7});
state{i}=[x;y;z;vx;vy;vz];%% 该变量记录下了每个时刻卫星的状态 

根据控制器计算，返回三个控制量，这里不考虑控制器的计算直接给出，然后将模型的参数进行更改，将发动机推力，有限机动的方向按照模型参数进行修改。

Controlmeter=Control(:,i);
%% 利用Component Browser生成发动机模型
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList Initial_State Maneuver']);
root.ExecuteCommand(['ComponentBrowser */ Duplicate "Engine Models" "Constant Thrust and Isp" MyEngine']);
root.ExecuteCommand(['ComponentBrowser */ SetValue "Engine Models" MyEngine Thrust ',num2str(Controlmeter(1)),' N' ]);
root.ExecuteCommand(['ComponentBrowser */ SetValue "Engine Models" MyEngine Isp 3000 s' ]);
%% 利用STK
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.MnvrType Finite']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.Propagator "Earth Point Mass"']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.AttitudeControl Thrust Vector']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.AttitudeUpdate Update during burn']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.ThrustAxes "Satellite/A LVLH.Axes"']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.Spherical.Azimuth ',num2str(Controlmeter(2)),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.Spherical.Elevation ',num2str(Controlmeter(3)),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.EngineModel MyEngine']);
FiniteMnvr=satellite2.Propagator.MainSequence.Item(1).Maneuver;
FiniteMnvr.Propagator.StoppingConditions.Item(0).Properties.Trip=0.1;
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' RunMCS']);

然后预报0.1s之后的状态作为第二次控制参数的初始输入

%% 预报0.1s后的时刻  将该时刻作为下一时刻的初始状态
data3=root.ExecuteCommand(['Report_RM */Satellite/',Sat_Name,' Style "LVLH" TimePeriod "',Time{i},'" "',Time{i+1},'" TimeStep 0.1']);
sss=regexp(data3.Item(2),',','split');
xe=str2double(sss{2});
ye=str2double(sss{3});
ze=str2double(sss{4});
vxe=str2double(sss{5});
vye=str2double(sss{6});
vze=str2double(sss{7});

由于STK目前还没有实时规划的功能，因此把之前设置好的卫星直接去除，重新生成一颗新的卫星，其历元开始时刻即为0.1s之后的时间，初始位置即按照LVLH的数据生成。

satellite2.Unload
Sat_Name2='B';
satellite2=root.CurrentScenario.Children.New('eSatellite', Sat_Name2);
satellite2.SetPropagatorType('ePropagatorAstrogator'); 
satellite2.Propagator;
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList Initial_State Maneuver']);
%% 初始化卫星参数
Initials=satellite2.Propagator.MainSequence.Item(0).InitialState;
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.CoordinateSystem "Satellite/A LVLH"']);
Initials.Epoch=Time{i+1};
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Cartesian.X ',num2str(xe),' m']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Cartesian.Y ',num2str(ye),' m']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Cartesian.Z ',num2str(ze),' m']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Cartesian.Vx ',num2str(vxe),' m/sec']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Cartesian.Vy ',num2str(vye),' m/sec']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Cartesian.Vz ',num2str(vze),' m/sec']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' RunMCS']);

之后便可以重复上述动作，这样得到每一个时刻的位置和速度。

代码的总程序如下图所示

clear;clc
uiApplication = actxGetRunningServer('STK12.application');
% Get our IAgStkObjectRoot interface
global root
root = uiApplication.Personality2;
checkempty = root.Children.Count;
if checkempty ~= 0
    root.CurrentScenario.Unload
    root.CloseScenario;
end
%% 根据你的需要设定场景的名称
root.NewScenario('chapter2');
StartTime = '29 Sep 2023 04:00:00.000';    % 场景开始时间
StopTime = '31 Oct 2023 04:00:00.000';     % 场景结束时间
root.ExecuteCommand(['SetAnalysisTimePeriod * "',StartTime,'" "',StopTime,'"']);
root.ExecuteCommand(' Animate * Reset');
%% ----------------------------------------------
%% A卫星
Sat_Name='A';
satellite=root.CurrentScenario.Children.New('eSatellite', Sat_Name);
satellite.Propagator;
sma=7071.393;
Ecc=0.001;
Inc=20;
w=0;
RAAN=0;
TA=359.999473339175;
satellite.SetPropagatorType('ePropagatorAstrogator'); 
satellite.Propagator;
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList Initial_State Propagate']);
%% 初始化卫星参数
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.CoordinateType Modified Keplerian']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Epoch ',StartTime,' UTCG']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.ElementType "Kozai-Izsak Mean"']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.sma ',num2str(sma),' km']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.ecc ',num2str(Ecc)]);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.inc ',num2str(Inc),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.w ',num2str(w),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.TA ',num2str(TA),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.RAAN ',num2str(RAAN),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name,' RunMCS']);
%% B卫星
Sat_Name2='B';
satellite2=root.CurrentScenario.Children.New('eSatellite', Sat_Name2);
satellite2.Propagator;
sma=7071.393;
Ecc=0.001;
Inc=20;
w=0;
RAAN=0;
TA2=359.998946678351;
satellite2.SetPropagatorType('ePropagatorAstrogator'); 
satellite2.Propagator;
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList Initial_State Propagate Maneuver']);
%% 初始化卫星参数
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.CoordinateType Modified Keplerian']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Epoch ',StartTime,' UTCG']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.ElementType "Kozai-Izsak Mean"']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.sma ',num2str(sma),' km']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.ecc ',num2str(Ecc)]);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.inc ',num2str(Inc),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.w ',num2str(w),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.TA ',num2str(TA2),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Keplerian.RAAN ',num2str(RAAN),' deg']);
root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' RunMCS']);
%% 生成LVLH坐标系 作为一个循环
MissionTime = '29 Sep 2023 04:05:00.000';    % 任务结束时间
data=root.ExecuteCommand(['Report_RM */Satellite/',Sat_Name,' Style "LVLH" TimePeriod "',StartTime,'" "',MissionTime,'" TimeStep 0.1']);%步长十分钟
LineNum=data.Count;
for i=1:LineNum-2
    struct=regexp(data.Item(i),',','split');
    Time{i}=struct{1};
end
%% 生成控制数据 （通过自带控制器实现）
Control=zeros(3,LineNum-3);% 第一项为控制推力的大小，第二项 控制推力的方向
for i=1:LineNum-3
    Control(1,i)=200;%恒推力模型
    Control(2,i)=30;
    Control(3,i)=20;
end
%% 开始循环
for i=1:LineNum-3
    data2=root.ExecuteCommand(['Report_RM */Satellite/',Sat_Name,' Style "LVLH" TimePeriod "',Time{i},'" "',Time{i+1},'" TimeStep 0.1']);%步长十分钟
    ss=regexp(data2.Item(1),',','split');
    x=str2double(ss{2});
    y=str2double(ss{3});
    z=str2double(ss{4});
    vx=str2double(ss{5});
    vy=str2double(ss{6});
    vz=str2double(ss{7});
    state{i}=[x;y;z;vx;vy;vz];%% 该变量记录下了每个时刻卫星的状态 
    %% 通过程序计算得到了第一个时刻的力的大小和方向
    Controlmeter=Control(:,i);
    %% 利用Component Browser生成发动机模型
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList Initial_State Maneuver']);
    root.ExecuteCommand(['ComponentBrowser */ Duplicate "Engine Models" "Constant Thrust and Isp" MyEngine']);
    root.ExecuteCommand(['ComponentBrowser */ SetValue "Engine Models" MyEngine Thrust ',num2str(Controlmeter(1)),' N' ]);
    root.ExecuteCommand(['ComponentBrowser */ SetValue "Engine Models" MyEngine Isp 3000 s' ]);
    %% 利用STK
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.MnvrType Finite']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.Propagator "Earth Point Mass"']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.AttitudeControl Thrust Vector']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.AttitudeUpdate Update during burn']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.ThrustAxes "Satellite/A LVLH.Axes"']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.Spherical.Azimuth ',num2str(Controlmeter(2)),' deg']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.Spherical.Elevation ',num2str(Controlmeter(3)),' deg']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Maneuver.FiniteMnvr.EngineModel MyEngine']);
    FiniteMnvr=satellite2.Propagator.MainSequence.Item(1).Maneuver;
    FiniteMnvr.Propagator.StoppingConditions.Item(0).Properties.Trip=0.1;
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' RunMCS']);
    %% 预报0.1s后的时刻  将该时刻作为下一时刻的初始状态
    data3=root.ExecuteCommand(['Report_RM */Satellite/',Sat_Name,' Style "LVLH" TimePeriod "',Time{i},'" "',Time{i+1},'" TimeStep 0.1']);
    sss=regexp(data3.Item(2),',','split');
    xe=str2double(sss{2});
    ye=str2double(sss{3});
    ze=str2double(sss{4});
    vxe=str2double(sss{5});
    vye=str2double(sss{6});
    vze=str2double(sss{7});
    %% 可以在上面生成一个断点 观察
    satellite2.Unload
    Sat_Name2='B';
    satellite2=root.CurrentScenario.Children.New('eSatellite', Sat_Name2);
    satellite2.SetPropagatorType('ePropagatorAstrogator'); 
    satellite2.Propagator;
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList Initial_State Maneuver']);
    %% 初始化卫星参数
    Initials=satellite2.Propagator.MainSequence.Item(0).InitialState;
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.CoordinateSystem "Satellite/A LVLH"']);
    Initials.Epoch=Time{i+1};
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Cartesian.X ',num2str(xe),' m']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Cartesian.Y ',num2str(ye),' m']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Cartesian.Z ',num2str(ze),' m']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Cartesian.Vx ',num2str(vxe),' m/sec']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Cartesian.Vy ',num2str(vye),' m/sec']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' SetValue MainSequence.SegmentList.Initial_State.InitialState.Cartesian.Vz ',num2str(vze),' m/sec']);
    root.ExecuteCommand(['Astrogator */Satellite/',Sat_Name2,' RunMCS']);
end
data2=root.ExecuteCommand(['Report_RM */Satellite/',Sat_Name,' Style "LVLH" TimePeriod "',Time{i},'" "',Time{i+1},'" TimeStep 0.1']);%步长十分钟
ss=regexp(data2.Item(1),',','split');
x=str2double(ss{2});
y=str2double(ss{3});
z=str2double(ss{4});
vx=str2double(ss{5});
vy=str2double(ss{6});
vz=str2double(ss{7});
state{LineNum-2}=[x;y;z;vx;vy;vz];%% 得到每个时刻的状态