MATLAB曲線擬合在井筒中心測(cè)量中的應(yīng)用

袁會(huì)等

【摘 要】在煤礦生產(chǎn)中，常常需要測(cè)量立井井筒中心坐標(biāo)。介紹了用MATLAB最小二乘曲線擬合的方法測(cè)量井筒中心坐標(biāo)，編寫了擬合程序，并在李糧店煤礦副井井筒井口中心偏移測(cè)量中進(jìn)行了應(yīng)用，取得了較好的效果。

【關(guān)鍵詞】MATLAB；最小二乘法；超定方程組；曲線擬合；井筒中心測(cè)量

0 前言

MATLAB是由美國(guó)MathWorks公司開發(fā)的一套以矩陣計(jì)算為基礎(chǔ)的工程計(jì)算軟件。它集成數(shù)值計(jì)算、可視化和編程功能且便于使用，擁有強(qiáng)大的繪圖功能和為解決各種特殊的科學(xué)和工程計(jì)算問題的工具箱。它具有計(jì)算功能強(qiáng)、編程效率高、使用簡(jiǎn)便、易于擴(kuò)充等特點(diǎn)。MATLAB已成為計(jì)算機(jī)輔助設(shè)計(jì)、數(shù)據(jù)處理分析、算法研究與應(yīng)用開發(fā)的首選工具[1]。

在煤礦建井初期井筒中心位置是根據(jù)近井點(diǎn)用極坐標(biāo)法標(biāo)定的[2]。井筒位置用大木樁固定，并在木樁上釘小釘作為井筒中心標(biāo)志，標(biāo)定后按地面一級(jí)導(dǎo)線的精度對(duì)井筒中心坐標(biāo)進(jìn)行實(shí)測(cè)[3]，作為最后的測(cè)量成果。在煤礦生產(chǎn)中，由于提升系統(tǒng)改造、立井延深等原因，需要使用井筒中心坐標(biāo)，常常由于建井期間井筒中心坐標(biāo)丟失或現(xiàn)有井筒發(fā)生變形，需要重新測(cè)量井筒中心坐標(biāo)。較常用的方法是在井筒上測(cè)量任意三點(diǎn)，用解析法解算井筒中心坐標(biāo)[4]，這種方法只適用于井筒形狀較規(guī)則的豎井；采用最小二乘曲線擬合確定井筒中心坐標(biāo)[5]適用范圍廣、理論依據(jù)充分、精度高。

1 最小二乘曲線擬合的基本原理及MATLAB實(shí)現(xiàn)[6-9]

所謂曲線擬合是指設(shè)法找出某條光滑的曲線能最佳地?cái)M合數(shù)據(jù)。在曲線擬合時(shí)，并不要求擬合曲線一定經(jīng)過(guò)每一個(gè)數(shù)據(jù)點(diǎn)。其主要思想是使擬合后的曲線能反映這些離散數(shù)據(jù)的變化趨勢(shì)，使數(shù)據(jù)點(diǎn)的誤差平方和最小。也就是已知一組測(cè)定的數(shù)據(jù)（例如N個(gè)點(diǎn)（xi，yi））去求得自變量x和因變量y的一個(gè)近似解析表達(dá)式y(tǒng)=f（x）。若記誤差δi=f（xi）-yi，i=1，2，…，N，則要使誤差的平方和最小，即要求（1）式中d最小，這就是常用的最小二乘法原理。

用最小二乘法求擬合曲線時(shí)，首先要確定f（x）的形式。當(dāng)人們對(duì)研究對(duì)象的內(nèi)在特性和各因素間的關(guān)系有比較充分的認(rèn)識(shí)時(shí)，一般用機(jī)理分析的方法建立描述y=f（x）的數(shù)學(xué)模型，再用曲線擬合的方法確定模型中的參數(shù)。但如果由于客觀事物內(nèi)部規(guī)律的復(fù)雜性及人們認(rèn)識(shí)程度的限制，無(wú)法建立合乎機(jī)理規(guī)律的數(shù)學(xué)模型時(shí)，只有對(duì)觀測(cè)數(shù)據(jù)（xi，yi）進(jìn)行分析，繪制散點(diǎn)圖，先猜測(cè)y=f（x）的形式，并通過(guò)上機(jī)實(shí)驗(yàn)和誤差分析，不斷對(duì)比和計(jì)算，然后選出擬合數(shù)據(jù)較好的函數(shù)類型。最小二乘曲線擬合分為線性和非線性最小二乘擬合。

2.1 線性最小二乘擬合

利用最小二乘法進(jìn)行曲線擬合時(shí)，要用求偏導(dǎo)數(shù)的方法確定擬合系數(shù)，人工計(jì)算時(shí)計(jì)算量大且精度不高，利用MATLAB進(jìn)行線性擬合可用以下方法：

（1）多項(xiàng)式擬合用polyfit（x，y，n）命令，其使用格式為A=polyfit（x，y，n），其中x，y為已知的擬合數(shù)據(jù)，n為多項(xiàng)式次數(shù)。

（2）在最小二乘意義下，解超定方程組，利用A=R＼y命令，直接求得待求參數(shù)。

2.2 非線性最小二乘擬合

擬合函數(shù)y=f（x）的待定系數(shù)不能全部以線性形式出現(xiàn)即為非線性最小二乘曲線擬合。在大多數(shù)情況下，一般將非線性擬合函數(shù)（如指數(shù)、雙曲線函數(shù)等）轉(zhuǎn)化為線性擬合函數(shù)求解。不能轉(zhuǎn)化的，MATLAB常用以下命令實(shí)現(xiàn)。

（1）lsqcurvefit（）命令，其使用格式為x=lsqcurvefit（fun，x0，xdata，ydata），其中fun是要擬合的非線性函數(shù)，x0是初始參數(shù)，xdata，ydata是擬合點(diǎn)的數(shù)據(jù)，該函數(shù)最終返回系數(shù)矩陣。

（2）nlinfit（）命令，其應(yīng)用格式為beta=nlinfit（x，y，fun，beta0），其中x和y是擬合點(diǎn)數(shù)據(jù)，fun是回歸（擬合）的函數(shù)，beta0是初始函數(shù)。

（3）lsqnonlin（）命令，其應(yīng)用格式為x=lsqnonlin（fun，x0），其中fun為擬合函數(shù)，x0為初始參數(shù)，將輸出的系數(shù)結(jié)果放在變量x中。

2 基于MATLAB井筒中心測(cè)量的思路及程序[10-12]

2.1 井筒中心測(cè)量的思路

顯然，井筒從理論上應(yīng)該是一個(gè)圓形，它的數(shù)學(xué)模型可以用圓方程來(lái)表示，圓的一般方程為：

3 應(yīng)用實(shí)例

2013年11月李糧店煤礦副井井筒潰沙，造成地表下沉及井筒偏斜，需要測(cè)量井口中心偏移量。我們對(duì)上井口井壁測(cè)量了23組點(diǎn)坐標(biāo)（見表1）。

然后利用擬合程序進(jìn)行圓擬合，求得井筒上口中心坐標(biāo)x0=298.869、y0=655.061及半徑r=3.286，繪制了擬合圖形（圖1）。并與建井標(biāo)定坐標(biāo)相比較，求出了偏移值。

4 結(jié)束語(yǔ)

（1）曲線擬合求取井筒中心坐標(biāo)，是比較科學(xué)合理的方法，利用MATLAB強(qiáng)大矩陣運(yùn)算功能，通過(guò)求解超定方程組進(jìn)行圓曲線擬合，簡(jiǎn)化了程序編制，提高了運(yùn)行效率。

（2）最小二乘法曲線擬合求取井筒中心坐標(biāo)，是建立在所有測(cè)量值基礎(chǔ)上的，測(cè)量的數(shù)據(jù)量越大，取得的井筒中心坐標(biāo)越接近真值。

（3）各數(shù)據(jù)點(diǎn)的采集盡量在井筒水平斷面上均勻選取，并進(jìn)行高精度測(cè)量，以便提高井筒中心坐標(biāo)的求取精度。

【參考文獻(xiàn)】

[1]王正林，劉明，陳連貴，等編.精通MATLAB[M].3版.北京：電子工業(yè)出版社，2013.

[2]張國(guó)良，朱家鈺，顧和和，等編.礦山測(cè)量學(xué)[M].徐州：中國(guó)礦業(yè)大學(xué)出版社，2001.

[3]中華人民共和國(guó)能源部.煤礦測(cè)量規(guī)程[M].北京：煤炭工業(yè)出版社，1989.

[4]柏運(yùn)輝.確定井筒中心的一種方法.礦山測(cè)量[J].1987（03）：31-32.

[5]張全武.用最小二乘擬合方法確定井筒中心坐標(biāo)[J].礦山測(cè)量，1997（03）：37-38.

[6]張韻華，陳效群，奚梅成，編.數(shù)值計(jì)算方法與算法[M].2版.科學(xué)出版社，2006.

[7]陳光，任志良，孫海柱.最小二乘曲線擬合及MATLAB實(shí)現(xiàn)[J].兵工自動(dòng)化，2005，24（3）：107-108.

[8]何正風(fēng)，張德豐，丁偉雄，等編.MATLAB在數(shù)學(xué)方面的應(yīng)用[M].北京：清華大學(xué)出版社，2012.

[9]宣明，王新成，阮婧，等編.數(shù)學(xué)建模與數(shù)學(xué)實(shí)驗(yàn)[M].杭州：浙江大學(xué)出版社，2010.

[10]陳懷琛，高淑萍，楊威，編.工程線性代數(shù)（MATLAB版）[M].北京：電子工業(yè)出版社，2007.

[11]同濟(jì)大學(xué)數(shù)學(xué)系，編.工程數(shù)學(xué)線性代數(shù)[M].5版.北京：高等教育出版社，2007.

[12]周本虎，瞿勇，李薇，等編.MATLAB與數(shù)學(xué)實(shí)驗(yàn)[M].北京：中國(guó)林業(yè)出版社，2005.

[責(zé)任編輯：劉展]