stk與matlab的接口及其在復(fù)雜航天任務(wù)仿真分析中的應(yīng)用

1、STK 與 Matlab 的接口及其在復(fù)雜航天任務(wù)仿真分析中的應(yīng)用,STK實(shí)驗(yàn)室二零一一年 八月,主要內(nèi)容,1 STK與Matlab接口簡(jiǎn)介,2 在復(fù)雜航天任務(wù)仿真分析中的應(yīng)用,使用STK與Matlab接口應(yīng)用場(chǎng)合與必要性,STK與Matlab互聯(lián)的設(shè)置與接口函數(shù)庫(kù),常用接口函數(shù)的使用,循環(huán)——星座對(duì)地觀測(cè)特性分析,循環(huán)+嵌套迭代——交會(huì)對(duì)接對(duì)遠(yuǎn)距離導(dǎo)引精度分析,迭代——交會(huì)對(duì)接目標(biāo)飛行器調(diào)相分析,1 STK與Matlab接口簡(jiǎn)介—

2、必要性,STK與Matlab互聯(lián)——兩個(gè)強(qiáng)大成熟商業(yè)軟件的強(qiáng)強(qiáng)聯(lián)合,,1 STK與Matlab接口簡(jiǎn)介—必要性,,有效的工程應(yīng)用例子——交會(huì)對(duì)接遠(yuǎn)距離導(dǎo)引精度論證,1 STK與Matlab接口簡(jiǎn)介—必要性,遠(yuǎn)距離導(dǎo)引精度仿真分析流程圖,,,,,,,,,,,,,,1 STK與Matlab接口簡(jiǎn)介—必要性,應(yīng)用效果,,,高效（快）,建模與修改快捷,可靠（好）,輔助實(shí)現(xiàn)關(guān)鍵技術(shù)環(huán)節(jié)突破,高效+可靠->技術(shù)優(yōu)勢(shì)->話語(yǔ)權(quán)->地

3、位,1 STK與Matlab接口簡(jiǎn)介—互聯(lián)與接口函數(shù),,,STK與Matlab互聯(lián)的初始設(shè)置,參考《STK在計(jì)算機(jī)仿真中的應(yīng)用 》,- 復(fù)制文件“agiCleanPath.m”、“agiInit.m”和“stkInit.m”到指定目錄,- 在Matlab中運(yùn)行指令“stkinit”,- 選擇M文件的位置,- 選擇Mex文件的位置,- 選擇默認(rèn)的其他設(shè)置,1 STK與Matlab接口簡(jiǎn)介—互聯(lián)與接口函數(shù),,,接口函數(shù),1 STK與Matl

4、ab接口簡(jiǎn)介—互聯(lián)與接口函數(shù),,,1 STK與Matlab接口簡(jiǎn)介—互聯(lián)與接口函數(shù),,,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,,1）初始化需要的接口函數(shù),stkinit—— 完成STK和Matlab的互聯(lián)和初始化,conid=stkOpen(stkDefaultHost); —— 返回互聯(lián)成功的主機(jī)端口的連接句柄,2）初始窗口管理,if stkValidScen == 1 stkUnload('/*')

5、end—— 如果已經(jīng)有打開的場(chǎng)景，則關(guān)閉場(chǎng)景。,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,,也可以通過交互提示進(jìn)行窗口管理，如：scen_open = stkValidScen;if scen_open == 1 rtn = questdlg('Close the current scenario?'); if ~strcmp(rtn,'Yes') stkClos

6、e(conid) return else stkUnload('/*') endend—— 如果已經(jīng)有打開的場(chǎng)景，則彈出提示對(duì)話框，詢問是否要關(guān)閉現(xiàn)有場(chǎng)景。,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,,3）建立場(chǎng)景,stkNewObj('/','Scenario','場(chǎng)景名稱');—— 建立給定名稱的場(chǎng)景。stkSet

7、TimePeriod('10 Apr 2003 00:00:00.0','12 Apr 2003 00:00:00.0','GREGUTC');—— 設(shè)置場(chǎng)景的起止時(shí)間和采用的時(shí)間系統(tǒng)。,stkSetEpoch('10 Apr 2003 00:00:00.0','GREGUTC'); —— 設(shè)置場(chǎng)景的歷元。注：在STK輸出的各類報(bào)告中，時(shí)間起點(diǎn)將以該歷

8、元作為時(shí)間零點(diǎn)。,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,,stkSyncEpoch;—— 同步aeroToolbox和STK場(chǎng)景歷元。注：aeroToolbox函數(shù)庫(kù)包含“安裝盤：\Program Files\AGI\STK\Matlab”目錄下各種atb開頭的內(nèi)部函數(shù)，可以與STK場(chǎng)景歷元保持獨(dú)立。STK場(chǎng)景時(shí)間發(fā)生改變時(shí)，通過調(diào)用該語(yǔ)句可保持二者時(shí)間同步。rtn = stkConnect(conid,'

9、Animate','Scenario/場(chǎng)景名稱','SetValues "10 Apr 2003 00:00:00.0" 60 0.1');rtn = stkConnect(conid,'Animate','Scenario/場(chǎng)景名稱','Reset');—— 設(shè)置STK場(chǎng)景動(dòng)畫歷元。,1 STK與Matlab接口簡(jiǎn)介—常用接

10、口函數(shù),,,4）建立航天器,stkNewObj('*/','Satellite','航天器名稱');—— 建立衛(wèi)星。,最基礎(chǔ)最常用的接口函數(shù)：stkSetPropClassical,語(yǔ)法：stkSetPropClassical('objPath', 'propagator', 'coordSystem', tStart, tStop, d

11、t, orbitEpoch, semimajorAxis, eccentricity,inclination, argOfPerigee, RAAN, meanAnomaly, coordEpoch),1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,,objPath – 有效路徑。Propagator – 選用的動(dòng)力學(xué)模型，可選擇為 ‘TwoBody’, ‘J2Perturbation’, ‘J4Perturbation’,

12、 ‘HPOP’ or ‘PODS’，分別表示“二體”、“J2攝動(dòng)”、“J4攝動(dòng)”、“HPOP模型”或“PODS模型”，區(qū)分大小寫。coordSystem – 選用的坐標(biāo)系，可選擇為 ‘Fixed’, ‘J2000’, ‘MeanOfDate’, ‘MeanOfEpoch’, ‘TrueOfDate’, ‘TrueOfEpoch’, ‘B1950’, ‘TEMEOfDate’, ‘TEMEOfEpoch’,

13、 ‘AlignmentAtEpoch’，區(qū)分大小寫。tStart, tStop – 起止時(shí)刻，以相對(duì)于歷元時(shí)刻的秒為單位。dt – 積分步長(zhǎng)，以秒為單位。orbitEpoch – 軌道歷元。,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,,semimajorAxis, eccentricity,inclination, argOfPerigee, RAAN, meanAnomaly– 6個(gè)軌道根數(shù)。其中，若坐

14、標(biāo)系選為‘AlignmentAtEpoch’，RAAN為升交點(diǎn)地理經(jīng)度，若坐標(biāo)系選為‘J2000’，RAAN為升交點(diǎn)赤經(jīng)。coordEpoch – 坐標(biāo)歷元, 所有 '...OfEpoch' 選項(xiàng)均需提供對(duì)應(yīng)的坐標(biāo)歷元。,同樣重要的接口函數(shù)：stkSetPropCart,語(yǔ)法：stkSetPropCart('objPath', 'propagator', 'coordSys

15、tem', tStart, tStop, dt, orbitEpoch, pos, vel, coordEpoch),1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,若動(dòng)力學(xué)模型選為“HPOP模型”，則還需要通過stkConnect函數(shù)設(shè)置各種力模型參數(shù)。,第一類：“HPOP force model” 設(shè)置重力、第三體引力、太陽(yáng)光壓、固體潮、海潮、輻射壓參數(shù)。,第二類：“HPOP drag model”

16、 設(shè)置大氣密度模型及參數(shù)。,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,第一類：“HPOP force model”,語(yǔ)法：HPOP Force {ForceModel} ,其中，中的參數(shù)表示需要用戶根據(jù)具體目錄安裝或參數(shù)使用情況設(shè)置的，{}中的參數(shù)表示STK固有的關(guān)鍵字，例如：stkConnect(conid,'HPOP','Scenario/場(chǎng)景名/Satellite/航天器

17、名','Force SolarRad Off')其中，為 'Scenario/場(chǎng)景名/Satellite/航天器名'，{ForceModel}為 SolarRad，為 Off。,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,{ForceModel} 和 的可選參數(shù)如表所示。,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,stkConnect(conid,'HPOP'

18、;,'Scenario/場(chǎng)景名稱/Satellite/航天器名稱','Force Gravity "C:\Program Files\AGI\STK\STKData\CentralBodies\Earth\JGM3.grv" 8 8');stkConnect(conid,'HPOP','Scenario/場(chǎng)景名稱/Satellite/航天器名稱',

19、'Force SolarRad Off');stkConnect(conid,'HPOP','Scenario/場(chǎng)景名稱/Satellite/航天器名稱','Force ThirdBodyGravity Sun Off');stkConnect(conid,'HPOP','Scenario/場(chǎng)景名稱/Satellite/航天器名稱'

20、;,'Force ThirdBodyGravity Moon Off');,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,第二類：“HPOP drag model ”,語(yǔ)法：HPOP Drag {DragState} {"AtmDensModel"},其中，與{}中的參數(shù)的含義不變。{DragState} 設(shè)置為On或Off。設(shè)置為Off則不需要輸入其他參數(shù)。設(shè)置為On時(shí)，為

21、-10.0~10.0之間的實(shí)數(shù)， 為0.0~9999.0之間的實(shí)數(shù)，{"AtmDensModel"}的有效輸入值和相應(yīng)的輸入值見表所示。,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,例如：stkConnect(conid,'HPOP','Scenario/場(chǎng)景名稱/Satellite/航天器名稱','Drag O

22、n 2.2 0.0022 "NRLMSISE 2000" Manual 150 150 2.667');stkConnect(conid,'HPOP','Scenario/場(chǎng)景名稱/Satellite/航天器名稱','Drag On 2.0 0.0022 "NRLMSISE 2000" File "C:\Program Files\A

23、GI\STK\DynamicEarthData\stkFluxGeoMag.fxm"');,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,,5）建立地面站,stkNewObj(' */ ','Facility','測(cè)站名稱')—— 建立地面站。stkSetFacPosLLA('Scenario/場(chǎng)景名稱/Facility/測(cè)站名稱', [phi;

24、 namda; h])—— 設(shè)置地面站的經(jīng)度、緯度和高度。stkConnect(conid,'SetConstraint','Scenario/場(chǎng)景名稱/Facility/測(cè)站名稱','ElevationAngle Min 最低仰角值')—— 設(shè)置地面站的最低仰角。stkConnect(conid,‘SetConstraint’,‘Scenario/場(chǎng)景名稱/Facility

25、/測(cè)站名稱’,‘Range Max 最大探測(cè)距離值');—— 設(shè)置地面站的最大探測(cè)距離。,1 STK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,,6）輸出計(jì)算結(jié)果,stkPropagate('*/Satellite/航天器名稱', t_start, t_stop)[T,r,v, cb] = stkEphemerisCBI('Scenario/場(chǎng)景名稱/Satellite/航天器名稱',dt,

26、t_start, t_stop)—— 運(yùn)行航天器并輸出給定起止時(shí)間內(nèi)的慣性坐標(biāo)系位置和速度。注：輸出航天器位置、速度之前必須運(yùn)行stkPropagate函數(shù)，否則力模型參數(shù)設(shè)置不生效。interval=stkAccess('Scenario/場(chǎng)景名稱/Satellite/航天器名稱','Scenario/場(chǎng)景名稱/Facility/測(cè)站名稱')—— 輸出航天器相對(duì)于測(cè)站的可見時(shí)間段。,1 S

27、TK與Matlab接口簡(jiǎn)介—常用接口函數(shù),,,7）關(guān)閉STK和Matlab連接,采用接口函數(shù)輸出的數(shù)據(jù)上雖然不如手工設(shè)置時(shí)形式多樣，但因?yàn)榭梢栽贛atlab環(huán)境下通過編程進(jìn)行數(shù)據(jù)處理，依然可以獲得多種需要的計(jì)算結(jié)果，并且可實(shí)現(xiàn)手工設(shè)置無(wú)法完成的循環(huán)計(jì)算和多重迭代功能。 上述1）~6）過程為利用STK/Matlab接口進(jìn)行計(jì)算分析的常用或基本函數(shù)。,stkClose(conid)注：STK和Matlab最多同時(shí)可創(chuàng)建2

28、個(gè)連接，因此每次程序運(yùn)行結(jié)束通常需要關(guān)閉連接，以免下次運(yùn)行重新建立連接時(shí)個(gè)數(shù)超過限制。,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—循環(huán),,循環(huán)——星座對(duì)地觀測(cè)特性分析,設(shè)某一對(duì)地觀測(cè)衛(wèi)星星座含多顆星，各星之間通過保持一定的幾何構(gòu)型，實(shí)現(xiàn)對(duì)某一特定區(qū)域的重點(diǎn)觀測(cè)。在系統(tǒng)設(shè)計(jì)階段，需要分析對(duì)某種給定的構(gòu)型，星座對(duì)某一特定區(qū)域進(jìn)行觀測(cè)時(shí)，各星觀測(cè)該區(qū)域的 時(shí)間間隔特性，如該特定區(qū)域重訪時(shí)間間隔平均值、最大值等。 對(duì)某一給定點(diǎn)

29、，通過手工設(shè)置STK中各類參數(shù)可計(jì)算和輸出星座中各星訪問該地點(diǎn)的時(shí)間間隔，但對(duì)于某一給定區(qū)域，需要在該區(qū)域內(nèi)密集采樣，采用手工輸入輸出方式已不再可能，因此需要通過Matlab作為計(jì)算引擎調(diào)用STK循環(huán)計(jì)算。 以下給出計(jì)算過程和結(jié)果。,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—循環(huán),,1）初始化并關(guān)閉所有打開的窗口,% 初始化stkinitremMachine = stkDefaultHost;% 關(guān)閉所有

30、打開的窗口delete(get(0,'chiSatren'));conid=stkOpen(remMachine);scen_open = stkValidScen;if scen_open == 1 stkUnload('/*')end,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—循環(huán),,2) 建立場(chǎng)景、航天器和測(cè)站,% 建立場(chǎng)景并設(shè)置場(chǎng)景屬性stkNewObj('/','

31、Scenario','RemoteSensor');stkSetTimePeriod('1 Jan 2010 00:00:00.0','5 Jan 2010 00:00:0.0','GREGUTC');stkSetEpoch('1 Jan 2010 00:00:00.0','GREGUTC'); stkSyncEpoch;rt

32、n = stkConnect(conid,'Animate','Scenario/RemoteSensor','SetValues "1 Jan 2010 00:00:00.0" 60 0.1');rtn = stkConnect(conid,'Animate','Scenario/RemoteSensor','Reset&#

33、39;);% 建立航天器和測(cè)站stkNewObj('*/','Satellite','Sat1');stkNewObj('*/','Satellite','Sat2');stkNewObj('*/','Satellite','Sat3');stkNewObj('*/'

34、,'Satellite','Sat4');stkNewObj('*/','Facility','Station');,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—循環(huán),,3) 設(shè)置星座中各衛(wèi)星的屬性,t_start=0; t_stop=4*(24*3600); dt=60; orbitEpoch=t_start;a_sat=7000*1000;e_sat= 0;

35、i_sat= 60;w_sat= 0; Raan_sat=0;M_sat=0 ;stkSetPropClassical('*/Satellite/Sat1','J2Perturbation','J2000',t_start,t_stop,dt,orbitEpoch,a_sat,e_

36、sat,i_sat,w_sat,Raan_sat,M_sat);stkPropagate('*/Satellite/Sat1', t_start, t_stop)stkSetPropClassical('*/Satellite/Sat2','J2Perturbation','J2000',t_start,t_stop,dt,orbitEpoch,a_sat,e_sat,

37、i_sat,w_sat,Raan_sat+pi/4,M_sat+pi);stkPropagate('*/Satellite/Sat2', t_start, t_stop)stkSetPropClassical('*/Satellite/Sat3','J2Perturbation','J2000',t_start,t_stop,dt,orbitEpoch,a_sat,e_

38、sat,i_sat,w_sat,Raan_sat+pi/2,M_sat+2*pi);stkPropagate('*/Satellite/Sat3', t_start, t_stop)stkSetPropClassical('*/Satellite/Sat4','J2Perturbation','J2000',t_start,t_stop,dt,orbitEpoch,a_

39、sat,e_sat,i_sat,w_sat,Raan_sat+3*pi/4,M_sat+pi);stkPropagate('*/Satellite/Sat4', t_start, t_stop),2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—循環(huán),,4）循環(huán)計(jì)算區(qū)域內(nèi)各點(diǎn)的重訪時(shí)間間隔及最大值,namda_min=0;namda_max=2*pi;n_namda=40;deta_namda=(namda_max-namda_m

40、in)/n_namda;phi_min=-60*pi/180;phi_max=60*pi/180;n_phi=30;deta_phi=(phi_max-phi_min)/n_phi;for i=1:n_namda+1 for j=1:n_phi+1 remain=(n_namda+1)*(n_phi+1)-(i-1)*(n_phi+1)-j namda=namda_min+(i-1)

41、*deta_namda; X(i,j)=namda*180/pi; phi=phi_min+(j-1)*deta_phi; Y(i,j)=phi*180/pi; stkSetFacPosLLA('Scenario/RemoteSensor/Facility/Station', [phi; namda; 0]);,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—循環(huán),,stk

42、Connect(conid,'SetConstraint','Scenario/RemoteSensor/Facility/Station','ElevationAngle Min 56.8164'); interval1=stkAccess('Scenario/RemoteSensor/Satellite/Sat1','Scenario/Remot

43、eSensor/Facility/Station'); interval2=stkAccess('Scenario/RemoteSensor/Satellite/Sat2','Scenario/RemoteSensor/Facility/Station'); interval3=stkAccess('Scenario/RemoteSensor/Satelli

44、te/Sat3','Scenario/RemoteSensor/Facility/Station'); interval4=stkAccess('Scenario/RemoteSensor/Satellite/Sat4','Scenario/RemoteSensor/Facility/Station');,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—循環(huán),,if size(in

45、terval1)==[0 0] & size(interval2)==[0 0] & size(interval3)==[0 0] & size(interval4)==[0 0] Z(i,j)=20; else temp1=[interval1.start,interval2.start,interval3.star

46、t,interval4.start,interval1.stop,interval2.stop,interval3.stop,interval4.stop]; temp2=sort(temp1); Z(i,j)=max(diff(temp2))/3600; end endend,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—循環(huán),,5）繪制計(jì)算結(jié)果并關(guān)閉連接,save('

47、;mat')set(gcf,'color',[1,1,1])figuresurf(X,Y,Z)xlabel('namda/deg')ylabel('phi/deg')stkClose(conid),2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—迭代,,迭代——交會(huì)對(duì)接目標(biāo)飛行器調(diào)相分析,交會(huì)對(duì)接任務(wù)中，目標(biāo)飛行器在飛船入軌前2~3個(gè)月先行入軌，目標(biāo)飛行器入軌后在一天時(shí)間內(nèi)通過兩次

48、變軌進(jìn)入某一高度的近圓軌道，之后通過大氣衰減和主動(dòng)變軌，使2~3個(gè)月后飛船入軌時(shí)刻目標(biāo)飛行器軌道的相位、高度和偏心率滿足要求。在任務(wù)分析階段，需要根據(jù)給定的空間環(huán)境設(shè)計(jì)目標(biāo)飛行器的圓化軌道高度，在后續(xù)2~3個(gè)月的時(shí)間內(nèi)盡可能通過大氣衰減作用滿足各種終端約束條件，減少主動(dòng)變軌的能量消耗。,第4圈變軌,,,第13圈變軌,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—迭代,,為了重點(diǎn)說(shuō)明STK/Matlab接口的使用，對(duì)上述設(shè)計(jì)進(jìn)行簡(jiǎn)化，假定不進(jìn)行調(diào)整相位

49、、高度和偏心率的變軌，設(shè)計(jì)過程簡(jiǎn)化為迭代搜索滿足飛船入軌時(shí)刻相位約束條件的目標(biāo)飛行器圓化軌道的高度，為單層迭代過程。 需要補(bǔ)充說(shuō)明的是，在STK中盡管Astrogator模式可實(shí)現(xiàn)一定功能的迭代，但在Astrogator模式下空間環(huán)境參數(shù)無(wú)法根據(jù)需要進(jìn)行設(shè)置，并且難以設(shè)定某些特殊的約束條件，這種情況下可以利用Matlab的編程能力通過STK/Matlab接口實(shí)現(xiàn)迭代求解。 以下給出計(jì)算過程。,2 在復(fù)

50、雜航天任務(wù)分析中的應(yīng)用—迭代,,1）初始化和并關(guān)閉所有打開的窗口,與前述過程相同，不再重復(fù)。2）建立場(chǎng)景和航天器目標(biāo)飛行器入軌后，在第4圈近地點(diǎn)變軌抬高遠(yuǎn)地點(diǎn)高度，在第13圈遠(yuǎn)地點(diǎn)變軌進(jìn)入圓軌道，之后進(jìn)入大氣衰減自由飛行階段。以下建立場(chǎng)景后，再建立3個(gè)航天器分別用于不同階段的軌道預(yù)報(bào)。 stkNewObj('/','Scenario','GM');stkSetTimePeriod

51、(‘6 Sep 2011 2:11:53', 6 Dec 2011 2:11:53 ','GREGUTC');stkSetEpoch(6 Sep 2011 2:11:53 ','GREGUTC'); stkSyncEpoch;rtn = stkConnect(conid,'Animate','Scenario/GM','SetValue

52、s " 6 Sep 2011 2:11:53 " 60 0.1');rtn = stkConnect(conid,'Animate','Scenario/GM','Reset');stkNewObj('*/','Satellite','MB_1');stkNewObj('*/','Sa

53、tellite','MB_2');stkNewObj('*/','Satellite','MB_3');,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—迭代,3）外推至第4圈近地點(diǎn),% 設(shè)置J2000坐標(biāo)系中初始軌道根數(shù)t_start_1=0; t_stop_1=0.187*(24*3600); dt=60;orbitEpoch=t_start_1;coordEpoch=

54、t_start_1;a_MB=6645999;e_MB= 0.01129;i_MB=42.747;w_MB=124.342; Raan_MB= 7.309;M_MB= 0 ;,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—迭代,,% 選用HPOP模型，設(shè)置大氣密度模型、重力場(chǎng)模型、光壓和第三體攝動(dòng)模型。stkSetPropClassi

55、cal('*/Satellite/MB_1','HPOP','AlignmentAtEpoch',t_start_1,t_stop_1,dt,orbitEpoch,a_MB,e_MB,i_MB,w_MB,Raan_MB,M_MB,coordEpoch);stkConnect(conid,'HPOP','Scenario/GM/Satellite/MB_1'

56、;,'Drag On 2.2 0.0039 "NRLMSISE 2000" Manual 150 150 2.667');stkConnect(conid,'HPOP','Scenario/GM/Satellite/MB_1','Force Gravity "C:\Program Files\AGI\STK\STKData\CentralBodie

57、s\Earth\JGM3.grv" 8 8');stkConnect(conid,'HPOP','Scenario/GM/Satellite/MB_1','Force SolarRad Off');stkConnect(conid,'HPOP','Scenario/GM/Satellite/MB_1','Force Third

58、BodyGravity Sun Off');stkConnect(conid,'HPOP','Scenario/GM/Satellite/MB_1','Force ThirdBodyGravity Moon Off');stkPropagate('*/Satellite/MB_1', t_start_1, t_stop_1)[T_MB_1, r_MB_1,

59、v_MB_1, cb] = stkEphemerisCBI('Scenario/GM/Satellite/MB_1',dt, t_start_1, t_stop_1);,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—迭代,,4）迭代求解第4圈近地點(diǎn)變軌速度增量（或求解圓化軌道高度）,% 采用牛頓迭代法求解，約束條件為第89天目標(biāo)飛行器軌道面與飛船入軌點(diǎn)共面時(shí)刻的相位在114°±10 °范圍內(nèi)。deta

60、_v_2=9;error=100;while abs(error)>1 Cal_remain; index_RG=find(T_interp>89); error=beta_interp(index_RG(1))-114 deta_v_2=deta_v_2+error/3600;End 其中，Cal_ remain為需要循環(huán)計(jì)算的用戶自定義子函數(shù)，用于計(jì)算給定一個(gè)第4圈近地

61、點(diǎn)變軌速度增量后，從第4圈近地點(diǎn)變軌外推至第13圈遠(yuǎn)地點(diǎn)并進(jìn)行軌道圓化變軌后的軌道，軌道外推一直到第89天目標(biāo)飛行器軌道面與飛船入軌點(diǎn)共面時(shí)刻。輸出結(jié)果外共面時(shí)刻的相位。,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—迭代,,Cal_ remain子函數(shù)為：% 第4圈變軌外推至第13圈變軌 pos_2=r_MB_1(:,length(r_MB_1));vel_2=v_MB_1(:,length(v_MB_1));ix=vel_2/norm(

62、vel_2);ir=pos_2/norm(pos_2);iz=CrossProduct(ir,ix)/norm(CrossProduct(ir,ix));iy=CrossProduct(ix,iz);thita_2=0;vel_2=vel_2+deta_v_2*(cos(thita_2)*ix+sin(thita_2)*iy);t_start_2=t_stop_1;t_stop_2=0.782*(24*3600);dt=

63、60;orbitEpoch=t_start_2; stkSetPropCart('*/Satellite/MB_2','HPOP','J2000',t_start_2, t_stop_2, dt, orbitEpoch, pos_2, vel_2)stkConnect(conid

64、,'HPOP','Scenario/GM/Satellite/MB_2','Drag On 2.2 0.0039 "NRLMSISE 2000" Manual 150 150 2.667');stkConnect(conid,'HPOP','Scenario/GM/Satellite/MB_2','Force Gravity &

65、quot;C:\Program Files\AGI\STK\STKData\CentralBodies\Earth\JGM3.grv" 8 8');,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—迭代,,stkConnect(conid,'HPOP','Scenario/GM/Satellite/MB_2','Force SolarRad Off');stkConnect(conid

66、,'HPOP','Scenario/GM/Satellite/MB_2','Force ThirdBodyGravity Sun Off');stkConnect(conid,'HPOP','Scenario/GM/Satellite/MB_2','Force ThirdBodyGravity Moon Off');stkPropagat

67、e('*/Satellite/MB_2', t_start_2, t_stop_2)[T_MB_2, r_MB_2, v_MB_2, cb] = stkEphemerisCBI('Scenario/GM/Satellite/MB_2',dt, t_start_2, t_stop_2);,%第13圈變軌外推至飛船入軌時(shí)刻pos_3=r_MB_2(:,length(r_MB_2));vel_3=v_M

68、B_2(:,length(v_MB_2));u_e=3.986005e+014;ix=vel_3/norm(vel_3);ir=pos_3/norm(pos_3);iz=CrossProduct(ir,ix)/norm(CrossProduct(ir,ix));iy=CrossProduct(ix,iz);it=CrossProduct(iz,ir);vel_3=sqrt(u_e/norm(pos_3))*it;,2 在復(fù)

69、雜航天任務(wù)分析中的應(yīng)用—迭代,t_start_3=T_MB_2(length(T_MB_2));t_stop_3=92*(24*3600);dt=60;orbitEpoch=t_start_3; stkSetPropCart('*/Satellite/MB_3','HPOP','J2000&#

70、39;,t_start_3, t_stop_3, dt, orbitEpoch, pos_3, vel_3)stkConnect(conid,'HPOP','Scenario/GM/Satellite/MB_3','Drag On 2.2 0.0039 "NRLMSISE 2000" Manual 150 150 2.667');stkConnect(conid,&

71、#39;HPOP','Scenario/GM/Satellite/MB_3','Force Gravity "C:\Program Files\AGI\STK\STKData\CentralBodies\Earth\JGM3.grv" 8 8');stkConnect(conid,'HPOP','Scenario/GM/Satellite/MB_3&#

72、39;,'Force SolarRad Off');stkConnect(conid,'HPOP','Scenario/GM/Satellite/MB_3','Force ThirdBodyGravity Sun Off');stkConnect(conid,'HPOP','Scenario/GM/Satellite/MB_3','

73、;Force ThirdBodyGravity Moon Off');stkPropagate('*/Satellite/MB_3', t_start_3, t_stop_3)[T_MB_3, r_MB_3, v_MB_3, cb] = stkEphemerisCBI('Scenario/GM/Satellite/MB_3',dt, t_start_3, t_stop_3);,2 在復(fù)雜航天任

74、務(wù)分析中的應(yīng)用—迭代,%計(jì)算共面時(shí)刻相位T_MB_all=[T_MB_1,T_MB_2(2:length(T_MB_2)),T_MB_3(2:length(T_MB_3))]'/(24*3600);R_MB_all=[r_MB_1,r_MB_2(:,2:length(r_MB_2)),r_MB_3(:,2:length(r_MB_3))]';V_MB_all=[v_MB_1,v_MB_2(:,2:length(v

75、_MB_2)),v_MB_3(:,2:length(v_MB_3))]';data_RGD=load('data_pos_RGD.m'); % 讀入入軌點(diǎn)在慣性系中的位置。 T_RGD_all=data_RGD(:,1);r_RGD_all=data_RGD(:,2:4);

76、step_window=1/60/24; t_start=0;t_end=92;T=[t_start:step_window:t_end]';index1=find(T_MB_all>t_start & T_MB_allt_start & T_RGD_all<t_end);R_RGD=[interp1(T_RGD_all(index2),r_RGD_all(index2,1),T),int

77、erp1(T_RGD_all(index2),r_RGD_all(index2,2),T),interp1(T_RGD_all(index2),r_RGD_all(index2,3),T)];,2 在復(fù)雜航天任務(wù)分析中的應(yīng)用—迭代,j=1;for i=1:length(T) AngMomentum_MB=CrossProduct(R_MB(i,:),V_MB(i,:))/norm(CrossProduct(R_MB(i,:),

78、V_MB(i,:))); alpha(i)=acos(Dotproduct(R_RGD(i,:)/norm(R_RGD(i,:)),AngMomentum_MB'))*180/pi; beta(i)=acos(Dotproduct(R_RGD(i,:)/norm(R_RGD(i,:)),R_MB(i,:)/norm(R_MB(i,:))))*180/pi; if DotProduct(CrossProdu

79、ct(R_RGD(i,:),R_MB(i,:)),AngMomentum_MB)1 & alpha(i)90 index(j)=i; T_interp(j)=interp1([alpha(i),alpha(i-1)],[T(i),T(i-1)],90); beta_interp(j)=interp1([T(i),T(i-1)],[beta(i),beta(i-1)],T_inter