空中飛行器無源定位

1、評 分 專 用 頁編號：學員評閱記錄：評閱隊號評分備注教員評閱記錄：空中飛行器無源定位摘 要目標定位技術是導航與制導技術的重要基礎。在現有的導航與制導技術中，衛(wèi)星定位技術是精度最高的，也是較為理想的導航與制導技術。本文研究了利用同步衛(wèi)星確定空中飛行器的位置參數及如何提高定位精度及定位效率的優(yōu)化選擇問題。針對問題一，根據測向陣列方向和地球同步衛(wèi)星夾角與飛行器位置的關系，建立了基于最小二乘法的測量、計算誤差平方和目標函數模型，運用最小二乘法去逼近真實參數，當目標函數最小時，所得數據為最接近真實參數的數據。在此基礎上，構造了以余弦夾角法分析的聚類模型，剔除誤差較大的數據，以剩余數據為有效解，得到了更

2、加接近真實參數的數據，是模型進一步得到優(yōu)化。針對問題二，運用問題一構建的數學模型，求得五個不同時刻飛行器的位置坐標，接著采用最小二乘法曲線擬合的方法擬合出值、值、值與時間的函數關系并計算其擬合曲線的皮爾遜相關系數，驗證擬合曲線可靠性，進而解得時飛行器的空間位置坐標。然后將不同時刻的擬合值與計算所得值之間的空間距離作為原始非負時間序列，建立模型，對距離進行灰色預測，得到其距離誤差的預測值，最后得出時飛行器的距離誤差值作為其預測的可靠度。針對問題三，本文以幾何精度因子作為衛(wèi)星分布對定位精度影響的指標。通過分析衛(wèi)星數目的影響,得到了隨衛(wèi)星數目增加單調遞減,但遞減幅度逐漸變小的變化規(guī)律。綜合考慮衛(wèi)星定

3、位精度和定位效率認為6星組合方案最為適宜。計算所有的6星組合方案，找出其中最小的衛(wèi)星組合作為最終的優(yōu)選方案。在測量角度存在誤差限的情況下，我們對附表中的數據采用加入方差為正態(tài)分布的方法，然后進行500次仿真實驗，觀察加入方向誤差后產生的定位誤差。統(tǒng)計定位誤差，發(fā)現加入噪聲后有的定位結果誤差在以內，即定位精度可以認為。在本文的最后，針對每個問題對其結果進行了分析、對每個問題解決方法的優(yōu)缺點進行了分析，并提出了相應的改進方案。本文的特色在于運用灰色預測模型對飛行器位置的定位距離誤差進行科學的分析，將誤差的變化考慮進模型，進而提高了預測的精度。關鍵詞：無源定位 最小二乘法 灰色模型 幾何精度因子1

4、問題的重述目標定位技術是導航與制導技術的重要基礎。在現有的導航與制導技術中，衛(wèi)星定位技術是精度最高的，也是較為理想的導航與制導技術。目前，較為成熟的衛(wèi)星導航系統(tǒng)有GPS系統(tǒng)、Galileo系統(tǒng)等。衛(wèi)星定位的基本原理是目標接收機通過接收多顆衛(wèi)星的信號測量出目標距各衛(wèi)星的距離（偽距），再通過一定的計算確定出目標的位置。 對于空中飛行器，在其飛行過程中很容易接收到太空衛(wèi)星的信號?，F在考慮通過測量飛行器與地球同步衛(wèi)星的方向角來實現空中飛行器的自定位。在球心坐標系下，空中飛行器P 的空間坐標記為(x,y,z)，不妨設它同時能接收到N 顆同步衛(wèi)星的信號，其N 顆同步衛(wèi)星Xi的空間坐標分別記為(xi,yi,

5、zi)(i=1,2,N)。為了方便檢測與同步衛(wèi)星的方向角，在空中飛行器上固定安裝了兩個相互垂直的測向陣列，它們的指向分別為d1(d1x,d1y,d1z)和d2(d2x,d2y,d2z)。地球同步衛(wèi)星Xi 與空中飛行器P的位置關系示意圖如圖所示，i，i分別表示空中飛行器P的測向陣列方向d1，d2與地球同步衛(wèi)星Xi(i=1,2,N)的夾角?，F在請你們建立數學模型研究解決下面的問題：（1）通過測量空中飛行器測向陣列方向d1和d2與多顆地球同步衛(wèi)星的夾角i和i，建立空中飛行器定位的數學模型；對于附表1所給出的9顆同步衛(wèi)星的數據，試確定空中飛行器P的位置參數。（2）在某些特殊情況下，空中飛行器能直接檢測

6、到的同步衛(wèi)星數量較少，可以利用空中飛行器在勻速飛行過程中多次檢測的結果來實現定位。針對這種情況，試建立空中飛行器定位的數學模型；對附表2中給出的3顆同步衛(wèi)星的檢測數據，確定空中飛行器P在第70秒時的位置參數，并分析其可靠性。（3）當可用同步衛(wèi)星數量較多時，為了提高定位精度和定位效率，需要對可用的同步衛(wèi)星進行一定的優(yōu)選。試研究具體的優(yōu)選策略，并通過仿真，分析在檢測方向角誤差限為0.1°時空中飛行器的定位方法和精度。2 模型的假設（1）假設附表中所有的數據都是真實可靠的；（2）假設飛行器是一個質點，不考慮飛機的飛行姿態(tài)；（3）衛(wèi)星信號強度都足夠的強；（4）不考慮地球的實際曲率變化，認為地

7、球是一個均勻球體。（5）衛(wèi)星環(huán)繞運動以地心為中心的圓；（6）固定在飛行器上的兩個垂直測向陣列為單位向量；3 符號說明符號含義Xi(xi,yi,zi)同步衛(wèi)星的位置P(x,y,z)飛行器的位置PXi飛行器指向同步衛(wèi)星方向的矢量PXi與測向陣列d1(d1x,d1y,d1z)的夾角PXi與測向陣列d2(d2x,d2y,d2z)的夾角v(vx,vy,vz)飛行器的速度矢量j每個點到平均位置的距離與實際測量值之間的誤差與實際測量值之間的誤差d1d2測向陣列相互垂直關系確定的方程的解的誤差d1測向陣列模為1確定的方程的解的誤差d2測向陣列模為1確定的方程的解的誤差4 問題一模型的建立與求解4.1問題一的分

8、析同步衛(wèi)星處于地球赤道平面上，圍繞地球自轉軸旋轉，并且相對地球上一點靜止，因此可建立以地心為原點的三維坐標系，衛(wèi)星位置參數可確定。分析測向陣列方向和地球同步衛(wèi)星夾角與飛行器位置的關系，2個測向陣列以及飛行器位置參數屬于未知參數，即至少需要9個包含上述參數的方程才能確定一個飛行器位置參數以及2個測向陣列。由于方向角以及測向陣列存在誤差，故可求得多組未知參數解，因此，精度難以滿足實際需求。因此，利用方向角測量、計算誤差，以及飛行器高度為約束條件，得到使得計算誤差平方和最小的目標函數，將9顆衛(wèi)星同時代入關系式，利用最小二乘法，當目標函數最小時，可得到最逼近真實參數的參數解。在上述方法考慮了所有飛行器

9、參數解，即所有解都為有效解，并對最終結果產生影響。因此可進一步優(yōu)化，可考慮按夾角余弦法聚類，將相近的參數解聚成一類，將個體或者數量明顯較少的類別剔除，剩下類別的參數解作為有效解，并以坐標均值求解最終參數解。4.2基于最小二乘法的模型建立建立以地球中心為原點，以赤道平面為基本平面，X軸在基本平面內由地心向外指向經線與赤道的交點，Z軸過地心且垂直基本平面，Y軸與X，Z軸組成右手系的三維直角坐標系(見圖 1 )：圖 1 球心坐標系在此坐標系下，同步衛(wèi)星的位置是，飛行器的位置是，如圖 2 所示：圖 2 飛行器與衛(wèi)星位置關系圖建立同步衛(wèi)星位置參數模型： （地球同步衛(wèi)星所處經度角）將附表一的數據代入求得同

10、步衛(wèi)星在坐標系中的位置參數如下：表 1 同步衛(wèi)星的位置參數衛(wèi)星編號xyz11020140914027354216003-144213962404-241853454105-271043230106-303322929207-332282596008-403241232809-4175658680測向陣列方向矢量為和，矢量PXi和d1、d2的夾角分別為i、i，則利用矢量點乘關系利用一顆衛(wèi)星可以得到如下非線性方程組：由于兩個測向陣列本身相互垂直且模為1，但考慮到測量、計算過程中存在誤差，因此可得式中，、分別為計算誤差。由上面的分析可知，方向角測量、計算等過程中誤差，故可得、為計算誤差,聯立9顆衛(wèi)星

11、和測向陣列關系，建立21個方程組，使之滿足如下函數優(yōu)化問題：式中R為地球半徑，表示飛行器到地球質心的距離大于地球半徑，利用最小二乘法就可以求解出飛行器位置，測向陣列和這9個未知參數，此時，所得解處于最逼近滿足21個方程組的情況。4.3基于聚類分析的模型優(yōu)化 由于存在9個未知參數，因此需要9個方程組才能解出，觀察關系式可知，3顆衛(wèi)星反饋的的信號，利用最小二乘法，即能確定飛行器的最佳位置參數。因為實際過程中存在測向陣列和地球同步衛(wèi)星夾角的測量誤差，故極大可能存在每次選取3顆衛(wèi)星所確定的飛行器位置參數是不同的，根據9顆衛(wèi)星可得個飛行器位置參數，對84個參數解按夾角余弦法聚類，將數據進行標準化處理，因

12、此，可以得到標準化的飛行器位置參數（，），故其與地心遠點構成的向量，設論域，的觀測值為，即數據矩陣利用夾角余弦法，可得相似系數：根據構造模糊相似矩陣，在矩陣為基礎，用平方法求出的傳遞包閉，模糊等價矩陣即等于。然后，由大到小取一組，確定相應的截矩陣，可以進行分類。根據分類結果，剔除個體解或數目較少的類，將其視作測量、計算誤差較大引起偏離真實位置的情況，利用所剩下的參數解，求得其坐標均值，該坐標值與真實位置參數基本吻合，其中為剩余參數解個數。然后代入利用測向陣列和地球同步衛(wèi)星夾角得出的矢量關系式，求解出測向陣列和。4.4模型的求解利用基于最小二乘法建立的模型，將9顆衛(wèi)星數據同時帶進方程組，利用MA

13、TLAB求解出飛行器的位置為 ,測向陣列和。4.5模型的結果評價基于最小二乘法的求解模型，求解出比較接近真實參數的數據，充分利用所有衛(wèi)星提供的參數，全面考慮測量、計算存在的誤差，方法簡明易懂，結果較為合理；再進一步優(yōu)化的以夾角余弦法為基礎的聚類分析模型，篩選出明顯誤差的數據并剔除，在此基礎上求解，使結果更加合理，符合實際。5 問題二模型的建立與求解5.1問題二的分析由問題一分析可知，三個衛(wèi)星就可以確定飛行器的位置，但精度存在問題。所以當同步衛(wèi)星數量較少時，可利用飛行器在勻速飛行過程中多次檢測的結果來提高定位的精度。利用問題二所給出的五組數據，可以分別確定五個時刻的位置坐標，然后分別對五個時刻位

14、置坐標分別進行曲線擬合，得到飛行器的運動軌跡曲線方程，則可得衛(wèi)星在時的位置參數。再利用此運動軌跡對前五組數據進行檢驗，即可得到一組殘差值，將此值帶入灰色預測模型，得到其殘差預測值與擬合值相加即為所要預測飛行器更為精確的位置參數。通過計算前五個殘差值的方差與飛行器在相應階段內飛行距離的比值，即可算的飛行器位置參數的可靠性。5.2.1利用最小二乘法進行擬合最小二乘法曲線擬合原理：給定平面上的點，求，使在處函數值與實驗數據的偏差的平方和最小，即吻合度最高。Step1：利用問題二不同時刻檢測到的地球同步衛(wèi)星的相關數據，運用問題一所建立的數學模型，分別解得時刻空間飛行器的位置坐標（見表 2 ）：表 2

15、不同時刻飛行器的坐標時刻x坐標y坐標z坐標-4185.024606.536070.90-4154.034621.586083.15-4110.794667.036102.62-4088.174687.826126.47-4036.974711.756145.26Step2：通過以上表格數據，以時間為橫坐標，分別以軸、軸、軸數據為縱坐標，建立平面坐標軸，運用MATLAB軟件分別擬合出值、值、值與時間的函數關系，如下圖所示：圖 3 的坐標和時間的擬合曲線圖其中的坐標值擬合方程：圖 4 的坐標和時間的擬合曲線圖其中的坐標值擬合方程：圖 5 的坐標和時間的擬合曲線圖其中的坐標值擬合方程：Step3：皮

16、爾遜相關系數是趨勢線擬合程度的指標，它的數值大小可以反映趨勢線的估計值與對應的實際數據之間的擬合程度，擬合程度越高，趨勢線的可靠性就越高。經計算：坐標的皮爾遜相關系數為0.9848，坐標的皮爾遜相關系數為0.9848，坐標的皮爾遜相關系數為0.9848，各坐標的皮爾遜相關系數都接近1，說明曲線擬合吻合度較高，可靠性強。Step4：根據以上函數關系，當時，空中飛行器的位置坐標為(-3934.03、4797.28、6201.73)5.2.2建立模型對距離進行灰色預測灰色預測模型是通過少量的、不完全的信息，建立數學模型并做出預測的一種預測方法。對于灰色量的處理不是尋求它的統(tǒng)計規(guī)律和概率分布，而是將無

17、規(guī)律的原始數據，通過一定的方法處理，變成比較有規(guī)律的時間序列數據。即以數找數的規(guī)律，再建立動態(tài)模型。因為客觀系統(tǒng)無論多么復雜，它總是聯系的，有序的，有整體功能的。所以作為系統(tǒng)行為特征的數據，總是蘊含著某種規(guī)律?；疑蓴盗校簩覕档奶幚碇饕抢脭祿幚矸椒ㄈで髷祿g的內在規(guī)律，通過對已知數據列中的數據進行處理而產生新的數據列，以此來研究尋找數據的規(guī)律性，這種方法稱為數據的生成。建立模型：設原始非負時間序列為，為累加生成序列，即：的白化微分方程為：為發(fā)展灰數；為內生控制灰數，亦稱為灰作用量，設待辨識向量，將上式寫成矩陣的形式：按最小二乘法得式中 于是可得到灰色預測的 離散時間相應函數為：為所

18、得累加的預測值，將預測值還原即為：殘差序列相對誤差序列相應的預測模型模擬根據灰色預測模型得到時距離殘差的預測值如表 3 所示：表 3 距離殘差預測值（單位：）時間（s）距離殘差03.711022.432025.953029.934034.705040.106046.457053.75由此可知，時，空中飛行器的可信位置范圍為以(-3934.03、4797.28、6201.73)為圓心，以距離殘差預測值為半徑的球體內，位置可靠度即為。表 4 灰色預測模型的檢驗：時間實際數據模擬數據殘差相對誤差22.6822.400.26 1.12%21.1420.640.482.27%37.8237.070.78

19、2.55%31.4431.750.321.07% 平均相對誤差：說明灰色預測模型相對誤差較小，可靠性強，對于空中飛行器有良好的定位性。5.3問題二的結果分析問題二，利用5個時刻的測量數據，運用問題一中數學模型，解得5個時刻空間飛行器的空間位置坐標，在利用最小二乘法進行曲線擬合，得到其函數關系，進而得到所求時刻空間飛行器的位置坐標，最后再利用灰色模型求得距離誤差值，并將其作為可靠性指標，說明灰色模型可靠性強，對飛行器有較好的定位性。6 問題三模型的建立與求解6.1問題三的分析 當可用同步衛(wèi)星數量較多時，為了提高定位精度和定位效率，需要對可用的同步衛(wèi)星進行一定的優(yōu)選，因此需要建立具體的優(yōu)選策略，也

20、就是如何選，選幾顆，選哪幾顆。而選擇幾何分布好的衛(wèi)星組合可以提高定位精度和定位效率。通常我們用幾何精度因子來表征用戶和可見衛(wèi)星在空間幾何分布的好壞。從顆衛(wèi)星中選擇顆衛(wèi)星，共有種選擇方案。從可見衛(wèi)星中選擇參加導航定位計算的衛(wèi)星數目不同, 的取值也不相同, 與衛(wèi)星數目之間有一定的變化規(guī)律，分析衛(wèi)星數目和之間的關系，從而找到一個合適的衛(wèi)星數目。同時，衛(wèi)星數目越多，意味著方程數量越多，計算量越大，定位效率也就越低。因此，為了兼顧定位的效率，應該在保證定位精度的前提下，盡可能減小定位衛(wèi)星的數量。根據對定位精度和定位效率的要求，確定一個合適的衛(wèi)星數目。確定衛(wèi)星數目為之后，還需對衛(wèi)星的幾何布局進行優(yōu)選，遍歷

21、從所有衛(wèi)星中選擇顆衛(wèi)星的所有方案，從中選擇最小的衛(wèi)星組合，即為優(yōu)選的結果。確定定位衛(wèi)星后，給測量角度加上最大為0.1°的高斯噪聲，進行多次的仿真實驗，統(tǒng)計定位精度。6.2基于幾何精度（）因子模型幾何精度GDOP：是衡量定位精度的很重要的一個系數，它代表GPS 測距誤差造成的接收機與空間衛(wèi)星間的距離矢量放大因子。實際表征參與定位解的從接收機至空間衛(wèi)星的單位矢量所勾勒的形體體積與成反比，故又稱為幾何精度因子。實際上，GDOP的數值越大，所代表的單位矢量形體體積越小，即接收機至空間衛(wèi)星的角度十分相似導致的結果，此時的GDOP會導致定位精度變差。好的GDOP, 是指其數值小，代表大的單位矢量

22、形體體積，導致高的定位精度。好的幾何因子實際上是指衛(wèi)星在空間分布不集中于一個區(qū)域，同時能在不同方位區(qū)域均勻分布。假設地球同步衛(wèi)星的數目為4個，定義下列向量和矩陣：，可將上式進一步簡化，得到：由于陣階數為4*3，故在上式兩側左乘的轉置矩陣就能獲得最小二乘解。若將用分量形式表示，即：定義幾何精度因子如下：圖 6 GDOP計算流程圖6.3優(yōu)選方法的確定由于在最優(yōu)定位衛(wèi)星的選擇中，問題可抽象為一個組合問題，即已知顆地球同步衛(wèi)星，從中選出顆組合用于定位計算，以作為目標函數，進而計算出每一種組合對應的值，即相同數量的衛(wèi)星條件下，有種組合?？紤]衛(wèi)星空間分布對精度的影響和方便計算的要求，在地球同步衛(wèi)星數量相同

23、的條件下，選取兩組分布差別較大的組合，運用問題（1）的方法聯立方程組即可求得飛行器的位置參數，并根據所構建的模型。計算出不同數量下的衛(wèi)星不同組合得到的值，(見表 5 )：表 5 不同衛(wèi)星組合位置參數對照表衛(wèi)星組合位置坐標(單位:km)（-5210.62，6608.68，3217.94）7.52(-5203.93， 6605.61， 3149.20)8.16(-5198.70， 6604.12，3122.95)4.69(-5202.86， 6605.96， 3166.78)4.13(-5198.26， 6605.17， 3124.96)0.76(-5201.33， 6605.27， 3157.3

24、6)0.72(-5200.26， 6605.29， 3139.00)0.36(-5202.51，6605.51，3147.92)0.31(-5202.40， 6605.72，3150.56)0.22(-5198.50， 6604.57，3121.25)0.19(-5188.74， 6605.39，3128.43)0.11根據上表數據計算出同一數量衛(wèi)星不同組合的平均值，見表6： 表 6 不同數量衛(wèi)星組合均值表衛(wèi)星數4567897.844.410.740.330.200.11圖 7 衛(wèi)星數量對GDOP的影響由上圖可知：值隨地球同步衛(wèi)星數目的增加而顯著減小，地球同步衛(wèi)星數目在4-6時，值顯著降低，下

25、降速率較大，而當定位衛(wèi)星的數量達到6顆以后，即衛(wèi)星數量為6-9顆時，其值變化不大，即此時增加衛(wèi)星數量，對的改善已經不明顯。因此綜合考慮定位精度和定位效率，6顆同步衛(wèi)星較為適宜。 從9顆衛(wèi)星選擇6顆衛(wèi)星，共有種選擇方案，根據問題一所構建的數學模型計算出每一種方案對應的值，我們選擇其中對應最小的一種選擇方案。經計算，當選擇第1,2,3,5,7,9顆衛(wèi)星時，其值最小。經度分別為，該組合對應的為。6.4角度誤差影響的仿真正態(tài)分布：是一個在數學、物理及工程等領域都非常重要的概率分布。若隨機變量服從一格數學期望、方差為的高斯分布，記為。其期望值決定了其位置，其標準差決定了分布的幅度。在多次測量中，設備的測

26、量誤差是正態(tài)分布的，而且均值為。問題中給出的精度可以認為是測量誤差的標準差。在這種理解下，用各自的誤差限對測量誤差進行無量綱化(也可以看成是一種加權法)處理是合理的，即求解的無約束優(yōu)化問題更合理。由6.3的分析可知，當六顆衛(wèi)星定位時，定位精度和效率最優(yōu)，故構建數學模型如下：運用計算機仿真技術進行500次仿真實驗得到實驗定位位置與未加入誤差時定位位置的距離誤差，統(tǒng)計這些距離誤差，得到其分布圖(見 圖 8 )：圖 8 距離誤差分布圖由上圖可以看出，絕大部分距離誤差都分布在250km以內,因此，認定定位精度為250km。6.4問題三的結果分析問題三，利用幾何精度因子來衡量衛(wèi)星幾何分布的好壞，遍歷所有

27、方案后，發(fā)現6星組合定位的平均精度和效率較好，其定位結果與9星定位的結果很接近，說明6星定位的組合是較為合理的。在給測量值加入最大為0.1°的高斯分布后，大量的仿真實驗顯示，加入的誤差給定位帶來了250km的偏差。說明模型建立的定位方法對于角度測量有著很高的要求，否則就有可能造成較大的誤差，即說明本模型對測量誤差有很高的要求，測量誤差是影響定位準確與否的重要因素。7 模型的優(yōu)缺點分析與改進方向7.1優(yōu)點：（1） 問題一基于最小二乘法的求解模型，求解出比較接近真實參數的數據，充分利用所有衛(wèi)星提供的參數，全面考慮測量、計算存在的誤差，方法簡明易懂，結果較為合理；再進一步優(yōu)化的以夾角余弦法

28、為基礎的聚類分析模型，篩選出明顯誤差的數據并剔除，在此基礎上求解，使結果更加合理，符合實際。（2） 問題二在問題一的基礎上利用同一時刻給出的三顆衛(wèi)星的參數算出各時刻對應的飛行器位置參數，然后利用最小二乘法擬合成曲線，推算出運動模型。并用模型預測擬合位置與實際位置距離，則其實際位置等于擬合值跟預測值的和，將誤差的變化考慮進模型，提高了定位的精度。（3） 在解決問題三時，考慮衛(wèi)星優(yōu)選策略時使用了精度幾何因子GDOP的方法，兼顧定位的精度和效率，選擇了合適的衛(wèi)星個數和衛(wèi)星分布，這種方法精度比較高。當測量角度發(fā)生偏差時，假設偏差符合正態(tài)分布，通過500次的仿真，得出距離誤差分布圖，更加形象直觀。 7.

29、2缺點：（1） 在問題一求解過程中沒有考慮到那些個別會引起測量偏差較大的衛(wèi)星對結果的影響，即沒有對衛(wèi)星進行優(yōu)選。（2） 在問題二中僅僅利用三顆衛(wèi)星提供的參數得出飛行器在不同時刻的位置參數本身就存在較大誤差，求解過程中飛行器的位置參數在軸上的數據殘差過大。（3） 在問題三中，本文僅僅只考慮了衛(wèi)星個數帶來的運算量對定位效率的影響 ，沒有沒有考慮其他因素的影響，這是不符合實際情況的。7.3改進方向：問題二求解過程中，將飛行器各軸上的位置參數分別進行數據擬合，最后得到時飛行器的位置在一個球體范圍內，難以定位，可綜合各軸位置參數進行擬合得出飛行器的空間軌跡隨時間得變化規(guī)律，可得預測到更為精確直觀的飛行器

30、位置參數。問題三中，考慮到選取衛(wèi)星組合的人為主觀因素，可盡量多的選取衛(wèi)星組合對飛行器進行定位，然后根據需要賦予不同權重，從而得到最優(yōu)解。參考文獻1 孫仲康,郭福成,馮道旺.單站無源定位跟蹤技術. 國防大學出版社. 2008.11第1版2 程銘東，劉利姣，黃光明.基于遺傳算法的多傳感器網絡中目標定位算法J.中國科技論文在線3 陳玲，無源定位與跟蹤算法研究D. 北京：北京航空航天大學，2005。4 陳小平，騰云龍，康榮雷等.幾何精度因子改進算法研究J.電子科技大學學報，2008，(6)5 卓金武.Matlab在數學建模中的應用. 北京航天航空大學出版社. 2011.4第1版6 鄧聚龍.灰色預測與決

31、策,華中工學院出版社，1986.附錄：附錄1.第一問的求解的主函數function F=f1(x)R=6367;h=35800;theta1=142/360*2*pi;theta2=163.0/360*2*pi;theta3=172.0/360*2*pi;theta4=136.0/360*2*pi;theta5=130.0/360*2*pi;theta6=125.0/360*2*pi;theta7=110.0/360*2*pi;theta8=89.0/360*2*pi;theta9=76.0/360*2*pi; a1=100.32/360*2*pi;b1=11.58/360*2*pi;a2=1

32、25.58/360*2*pi;b2=36.04/360*2*pi;a3=135.93/360*2*pi;b3=46.34/360*2*pi;a4=92.91/360*2*pi;b4=5.99/360*2*pi;a5=85.46/360*2*pi;b5=6.95/360*2*pi;a6=79.24/360*2*pi;b6=11.99/360*2*pi;a7=60.73/360*2*pi;b7=29.81/360*2*pi;a8=35.82/360*2*pi;b8=54.63/360*2*pi;a9=21.33/360*2*pi;b9=69.25/360*2*pi; xm1=(R+h)*cos(t

33、heta1);ym1=(R+h)*sin(theta1);zm1=0;xm2=(R+h)*cos(theta2);ym2=(R+h)*sin(theta2);zm2=0;xm3=(R+h)*cos(theta3);ym3=(R+h)*sin(theta3);zm3=0;xm4=(R+h)*cos(theta4);ym4=(R+h)*sin(theta4);zm4=0;xm5=(R+h)*cos(theta5);ym5=(R+h)*sin(theta5);zm5=0;xm6=(R+h)*cos(theta6);ym6=(R+h)*sin(theta6);zm6=0;xm7=(R+h)*cos(t

34、heta7);ym7=(R+h)*sin(theta7);zm7=0;xm8=(R+h)*cos(theta8);ym8=(R+h)*sin(theta8);zm8=0;xm9=(R+h)*cos(theta9);ym9=(R+h)*sin(theta9);zm9=0; if(x(1)2+x(2)2+x(3)2<63672) F=inf;else F=(xm1-x(1)*x(4)+(ym1-x(2)*x(5)+(-x(3)*x(6)-(xm1-x(1)2+(ym1-x(2)2+(zm1-x(3)2)0.5*cos(a1)2+. (xm1-x(1)*x(7)+(ym1-x(2)*x(8)+

35、(-x(3)*x(9)-(xm1-x(1)2+(ym1-x(2)2+(zm1-x(3)2)0.5*cos(b1)2+. (xm2-x(1)*x(4)+(ym2-x(2)*x(5)+(-x(3)*x(6)-(xm2-x(1)2+(ym2-x(2)2+(zm2-x(3)2)0.5*cos(a2)2+. (xm2-x(1)*x(7)+(ym2-x(2)*x(8)+(-x(3)*x(9)-(xm2-x(1)2+(ym2-x(2)2+(zm2-x(3)2)0.5*cos(b2)2+. (xm3-x(1)*x(4)+(ym3-x(2)*x(5)+(-x(3)*x(6)-(xm3-x(1)2+(ym3-x(

36、2)2+(zm3-x(3)2)0.5*cos(a3)2+. (xm3-x(1)*x(7)+(ym3-x(2)*x(8)+(-x(3)*x(9)-(xm3-x(1)2+(ym3-x(2)2+(zm3-x(3)2)0.5*cos(b3)2+. (xm4-x(1)*x(4)+(ym4-x(2)*x(5)+(-x(3)*x(6)-(xm4-x(1)2+(ym4-x(2)2+(zm4-x(3)2)0.5*cos(a4)2+. (xm4-x(1)*x(7)+(ym4-x(2)*x(8)+(-x(3)*x(9)-(xm4-x(1)2+(ym4-x(2)2+(zm4-x(3)2)0.5*cos(b4)2+.

37、(xm5-x(1)*x(4)+(ym5-x(2)*x(5)+(-x(3)*x(6)-(xm5-x(1)2+(ym5-x(2)2+(zm5-x(3)2)0.5*cos(a5)2+. (xm5-x(1)*x(7)+(ym5-x(2)*x(8)+(-x(3)*x(9)-(xm5-x(1)2+(ym5-x(2)2+(zm5-x(3)2)0.5*cos(b5)2+. (xm6-x(1)*x(4)+(ym6-x(2)*x(5)+(-x(3)*x(6)-(xm6-x(1)2+(ym6-x(2)2+(zm6-x(3)2)0.5*cos(a6)2+. (xm6-x(1)*x(7)+(ym6-x(2)*x(8)+

38、(-x(3)*x(9)-(xm6-x(1)2+(ym6-x(2)2+(zm6-x(3)2)0.5*cos(b6)2+. (xm7-x(1)*x(4)+(ym7-x(2)*x(5)+(-x(3)*x(6)-(xm7-x(1)2+(ym7-x(2)2+(zm7-x(3)2)0.5*cos(a7)2+. (xm7-x(1)*x(7)+(ym7-x(2)*x(8)+(-x(3)*x(9)-(xm7-x(1)2+(ym7-x(2)2+(zm7-x(3)2)0.5*cos(b7)2+. (xm8-x(1)*x(4)+(ym8-x(2)*x(5)+(-x(3)*x(6)-(xm8-x(1)2+(ym8-x(

39、2)2+(zm8-x(3)2)0.5*cos(a8)2+. (xm8-x(1)*x(7)+(ym8-x(2)*x(8)+(-x(3)*x(9)-(xm8-x(1)2+(ym8-x(2)2+(zm8-x(3)2)0.5*cos(b8)2+. (xm9-x(1)*x(4)+(ym9-x(2)*x(5)+(-x(3)*x(6)-(xm9-x(1)2+(ym9-x(2)2+(zm9-x(3)2)0.5*cos(a9)2+. (xm9-x(1)*x(7)+(ym9-x(2)*x(8)+(-x(3)*x(9)-(xm9-x(1)2+(ym9-x(2)2+(zm9-x(3)2)0.5*cos(b9)2+.

40、(x(4)*x(7)+x(5)*x(8)+x(6)*x(9)2+. (x(4)2+x(5)2+x(6)2-1)2+. (x(7)2+x(8)2+x(9)2-1)2;endend附錄2.最小二乘法求解過程clc;clear; xx=-6000;-3000;-1000;yy=2000;4000;7000;zz=1000;3000;5000;vx=-30;0;30;vy=-30;0;30;vz=-30;0;30;xxx=;x0=-6000;2000;1000;1;0;1;2;1;6;30;30;30;mm=10;nn=x0;for i=1:3 for j=1:3 for k=1:3 for ii=1:3 for jj=1:3 for kk=1:3 aa1=xx(i);bb1=yy(j);cc1=zz(k);vx1=vx(ii);vy1=vy(jj);vz1=vz(kk); x0=aa1;bb1;cc1;1;0;1;2;1;6;vx1;vy1;vz1; options = optimset(optimset,'MaxFunEvals',100000); x,resnorm=lsqnonlin(f1,x0,options); xxx=xxx;x; if resnorm<mm nn=x; mm=resnorm; end end end en