復(fù)合測量誤差分布的蒙特卡羅仿真

王輝

摘 要：針對用傳統(tǒng)方法計算復(fù)合測量誤差分布十分困難的情況，提出了基于蒙特卡羅法的復(fù)合測量量的誤差分布的統(tǒng)計處理方法。仿真分析表明，采用蒙特卡羅和傳統(tǒng)的方法統(tǒng)計的復(fù)合測量函數(shù)的均方差基本相同，相對于傳統(tǒng)的統(tǒng)計方法，蒙特卡羅法不用考慮復(fù)合測量函數(shù)的逆變換，可以利用計算機直接仿真具有特定分布的隨機誤差，能夠很快統(tǒng)計出復(fù)合測量隨機誤差的分布特性，具有一定的推廣價值。

關(guān)鍵詞：復(fù)合測量； 蒙特卡羅仿真； 均方差； Matlab

中圖分類號：TN911-34

文獻標識碼：A

文章編號：1004-373X(2011)09-0025-03

Monte-Carlo Simulation for Composite Measurement Error Distribution

WANG Hui

(The 93.Detachment, Unit 92941 of PLA, Huludao 125000, China)

Abstract： Monte Carlo simulation method is put forward for solving the difficulty of traditional composite measurement error distribution calculation methods. Simulation results prove that their root-mean square deviation are same, but the Monte-Carlo simulation does not need to calculate reverse transform. It can calculate random error distribution by computer and has generalization.

Keywords： composite measurement; Monte-Carlo simulation; root-mean square deviation; Matlab

0 引 言

在外測中，絕大多數(shù)問題是復(fù)合測量。通常，將需多個測量量才能解算得到所需要的量，稱為復(fù)合測量［1］。例如，在對外測數(shù)據(jù)進行事后處理時，利用雷達測量數(shù)據(jù)斜距R、方位角α和高低角β，通過坐標轉(zhuǎn)換得到目標在測量坐標系中的位置X=(x,y,z)琓［1］。復(fù)合測量誤差通常由直接測量與被測量之間的函數(shù)關(guān)系來計算，因此又稱為函數(shù)誤差。研究隨機變量的函數(shù)誤差的概率分布問題，傳統(tǒng)的做法為［2］：根據(jù)概率論的知識，設(shè)測量量x1,x2,…,xn的聯(lián)合密度函數(shù)為fx(x1,x2,…,xn)，而:

Yi=gi(x1,x2,…,xn), i=1,2,…,n

構(gòu)成x1,x2,…,xn到y(tǒng)1,y2,…,yn的一一對應(yīng)變換，則其存在逆變換，設(shè)其逆變換的雅可比行列式為:

J(y1,y2,…,yn)=

礸－11/祔1…礸－11/祔n

螵鰳

礸－1n/祔1…礸－1n/祔n

則得到y(tǒng)1,y2,…,yn的聯(lián)合概率密度函數(shù)為：

fy(y1,y2,…,yn)=J(y1,y2,…,yn)*

fx[g－11(y1,y2,…,yn),g－12(y1,y2,…,yn),…,

g－1n(y1,y2,…,yn)]

(1)

但是,當復(fù)合測量函數(shù)為非線性函數(shù)時,要求得 x1,x2,…,xn到y(tǒng)1,y2,…,yn的一一對應(yīng)變換的逆變換是十分困難的，所以采用上述解析方法是不可能的。

為解決上述問題,本文結(jié)合單脈沖雷達坐標轉(zhuǎn)換模型的隨機誤差分布的統(tǒng)計處理問題，介紹了一種使用簡單、計算精度高的方法——蒙特卡羅法，并利用Matlab軟件編程對復(fù)合測量量的統(tǒng)計特性進行分析，不僅可以模擬復(fù)合測量隨機誤差的分布，還可以計算標準差等數(shù)字特征，從而全面、完整地評定復(fù)合測量的結(jié)果。

1 蒙特卡羅仿真方法的一般原理

蒙特卡羅仿真的實質(zhì)是利用服從某種分布隨機數(shù)來模擬現(xiàn)實系統(tǒng)中可能出現(xiàn)的隨機現(xiàn)象，由于每次仿真試驗僅能描述所考察系統(tǒng)出現(xiàn)的一種可能狀態(tài)，故若能進行大量次數(shù)的仿真試驗，就能得到與現(xiàn)實所期望的情況相一致的統(tǒng)計結(jié)果。

1.1 均勻隨機數(shù)發(fā)生器

大多數(shù)計算機的軟件庫中都包含著一個均勻隨機數(shù)發(fā)生器,它是以等概率產(chǎn)生在0和1之間的一個數(shù),這個隨機數(shù)發(fā)生器的輸出是一個隨機變量,其取值為0≤A≤1,從實際應(yīng)用角度看,可近似地認為計算機在(0,1)內(nèi)能輸出的數(shù)位足夠大,以至于它可以產(chǎn)生(0,1)間的任何連續(xù)值。

1.2 給定概率分布函數(shù)(PDF)的隨機變量的產(chǎn)生方法

均勻隨機數(shù)發(fā)生器產(chǎn)生隨機變量的范圍為0≤A≤1,其概率密度函數(shù)pA(a)=1。由于隨機變量A的概率分布函數(shù)(CDF)取值范圍為0≤F(a≤A)≤1,故可以用此隨機數(shù)發(fā)生器均勻地生成任意CDF的值。在已知該PDF的情況下,反求該隨機過程的隨機變量(即樣本值)。如某隨機過程,其隨機變量X的PDF為:

pX(x)=(1/b0)砮－x/b0, x≥0

(2)

那么,它的CDF為:

F(X)=p(x≤X)=∫琗－∞pxdx=1－e－X/b0,X≥0

(3)

由于F(X)∈［0,1］,其概率取值顯然決定于X的取值,而X取值概率又應(yīng)服從pX(x)的指數(shù)分布,這樣很自然會想到用均勻隨機發(fā)生器在［0,1］間隨機生成一個數(shù),以該數(shù)作為pX(x)的概率,從而得到隨機變量X。故有A=F(X)=1－e－X/b0,即:

X=b0*ln[1/(1－A)]

因此,(0,1)內(nèi)均勻分布的隨機變量可用來產(chǎn)生具有其他概率分布函數(shù)的隨機變量,一般方法是:假設(shè)想產(chǎn)生一個隨機變量C,其概率分布函數(shù)為F(C),由于0≤F(C)≤1,故可事先在(0,1)內(nèi)均勻產(chǎn)生一個隨機變量A,F(C)=A求反函數(shù)得C=F－1(A),這樣就解出具有F(C)=A的C值,重復(fù)這一步驟,求得具有CDF為F(C)的一系列新的隨機變量［3］C。

2 MC仿真測量系統(tǒng)的一般步驟

對于測量隨機誤差統(tǒng)計問題，仿真的對象為測量系統(tǒng)，所關(guān)心的問題是測量諸因素的變化對測量值的影響。

用蒙特卡羅仿真方法求解復(fù)合測量誤差分布的一般步驟如下［4］：

(1) 輸入各直接測量量R1i,R2i,…,Rni；α1i,α2i,…,αni和β1i,β2i,…,βni,剔除野值及系統(tǒng)誤差。

(2) 利用Matlab函數(shù)mean,std計算出算數(shù)平均值,,和標準偏差σR,σα,σβ。

(3) 根據(jù)測量過程中諸直接測量量的概率分布特征，仿真測量過程中諸直接測量量的大樣本偽隨機數(shù)，并繪制描述各輸入直接測量量誤差分布的統(tǒng)計直方圖。

(4) 將直接測量量的平均值,,加上諸直接測量量的仿真樣本，得出測量過程的“偽測量數(shù)據(jù)”，此時“偽測量數(shù)據(jù)”可以看作是實際過程數(shù)據(jù)。

(5) 按函數(shù)測量模型計算復(fù)合量y，并繪制函數(shù)隨機誤差分布的統(tǒng)計直方圖。

(6) 統(tǒng)計并輸出該復(fù)合量的最佳估值，標準差s。

3 實例仿真

本文的數(shù)據(jù)仿真程序選用Matlab作為軟件開發(fā)平臺。Matlab是基于向量或矩陣的數(shù)學軟件,具有非常強大的數(shù)值計算能力和卓越的數(shù)據(jù)可視化能力,Matlab的這種強大功能為誤差數(shù)據(jù)處理提供了便利,可以迅速編出科學高效的計算程序,大大提高了效率［5］。下面以實例說明之。

通過直接測量R,α,β來求復(fù)合測量量X=Rcos βcos α

Rsin β

Rcos βsin α。剔除了野值及系統(tǒng)誤差的R,α,β測量數(shù)據(jù)見表1。

上述程序為Matlab的腳本程序，直接運行，得到的結(jié)果如圖1和圖2所示。

圖1 模擬隨機數(shù)的誤差分布圖

圖2 蒙特卡羅仿真法輸出函數(shù)y誤差分布圖

由圖2可直觀地得到函數(shù)誤差分布的密度函數(shù)，由仿真程序求得函數(shù)y的標準差為80.506 6 m。又因為直接測量量R，α與β彼此之間不相關(guān)，根據(jù)誤差傳播理論，復(fù)合測量函數(shù)的標準差表示為［7］：

σ=(礷/礡)2σ2R+(礷/鄲)2σ2α+(礷/鄲)2σ2β

(4)

式中:σR,σα和σβ表示R,α和β的標準差。按式(2)計算坐標轉(zhuǎn)換模型中y的標準差：σy=81.965 0 m。可見，由隨機模擬法求得的標準差與理論計算值相差1.359 0 m，說明隨機模擬法的計算精度是相當高的。

4 結(jié) 論

通過比較可以看出，采用傳統(tǒng)的方法統(tǒng)計復(fù)合測量函數(shù)誤差的標準差以及用蒙特卡羅法對復(fù)合測量函數(shù)的標準差進行模擬計算的結(jié)果基本相同。但蒙特卡羅法相對于傳統(tǒng)的方法有明顯的優(yōu)勢：

(1) 雖然利用概率密度函數(shù)可以全面準確地反映復(fù)合測量結(jié)果的分布,但當概率密度函數(shù)的表達式較為復(fù)雜時,采用解析的方法非常復(fù)雜甚至是不可能的。采用蒙特卡羅法，只要通過數(shù)學建模，選擇合適的算法和編程，即可由計算機輸出的統(tǒng)計直方圖直觀地求出該復(fù)合測量的隨機誤差分布，計算出誤差的標準差等信息，與采用定量解析公式計算的方法相比，該方法不局限于顯函數(shù)的數(shù)學模型，也不受方差傳播非線性嚴重的限制，是一種解決復(fù)雜測量模型誤差分析的有效途徑。

(2) 在重復(fù)性的評定過程中，蒙特卡羅法能夠根據(jù)大樣本的仿真數(shù)據(jù)更加真實地反映出測量結(jié)果的分布特征。

參考文獻

［1］劉利生.外側(cè)數(shù)據(jù)事后處理［M］.北京：國防工業(yè)出版社，2000.

［2］王國玉.電子系統(tǒng)小子樣試驗理論與方法［M］.北京：國防工業(yè)出版社，2003.

［3］莊銘杰.蒙特卡羅仿真在移動通信中的應(yīng)用研究［J］.計算機仿真,2004,21(11):116-120.

［4］陳光,任志良,孫海柱.基于Matlab的函數(shù)誤差分布模擬與計算［J］.計量技術(shù),2006(5):116-120.

［5］王成剛.測試性建模與分析中的故障概率獲取方法研究［J］.測試技術(shù)學報,2010,24(1):10-13.

［6］田社平.基于Matlab的間接測量數(shù)據(jù)處理［J］.計量技術(shù),2004(3):6-8.

［7］費業(yè)泰.誤差理論與數(shù)據(jù)處理［M］.北京：國防工業(yè)出版社，2008.

［8］魯培耿.海軍裝備試驗常用詞典［M］.北京：國防工業(yè)出版社，2007.

［9］曲寶忠.海軍戰(zhàn)術(shù)導(dǎo)彈試驗與鑒定［M］.北京：國防工業(yè)出版社，2005.

［10］陳彥發(fā),呂國強,蔡斐.Ka波段E面波導(dǎo)濾波器的誤差研究[J].現(xiàn)代電子技術(shù),2009,32(13):11-14.

注：本文中所涉及到的圖表、注解、公式等內(nèi)容請以PDF格式閱讀原文