基于真實幾何特征的碎石道砟顆粒形態重構

王 威，張田爽，汪佳豪，相里海鵬，楊成忠，黃大維

(1.華東交通大學 軌道交通基礎設施性能監測與保障國家重點實驗室，江西 南昌 330013；2.華東交通大學 交通運輸工程學院，江西 南昌 330013)

有砟軌道是鐵路的傳統結構形式,主要由道床、軌枕、鋼軌及各種聯結件構成。道床作為有砟軌道的重要組成部分,是由形狀各異、大小不一的碎石道砟通過顆粒之間的咬合作用而形成的散粒體結構,其力學性能不僅與碎石材料的強度有關,且很大程度上受道砟顆粒形態特征的影響[1-3]。因此,在研究道床受力狀態,尤其是在分析其細觀力學行為時,構建與真實道砟形態特征相符的模型成為必要的環節之一[4-5]。

目前,道砟模型構建方法主要分為兩種,一種是基于圖形分析技術重構道砟顆粒形狀的方法,另一種是基于幾何理論構建道砟模型的方法。基于圖形分析技術重構道砟形狀的方法是對拍攝或掃描得到的數字圖形進行分析,直接獲取道砟顆粒模型,如Le等[6]、潘飛等[7]采用圖像分析系統獲取道砟顆粒的二維輪廓,井國慶等[8]、徐旸[9]通過激光掃描法重構了道砟三維模型。這些方法可以較為準確地獲取道砟模型,但實際運用中需要逐個采集道砟顆粒的數字圖形,工作量較大,重構效率較低。基于幾何理論構建道砟模型的方法是通過計算軟件批量構建道砟模型,如井國慶等[10]用圓形顆粒代替道砟顆粒對循環荷載作用下道砟破碎劣化進行了仿真,Liu等[11]、Lu等[12]利用隨機生成的不規則多邊形建立數值道砟模型進行仿真分析,Liu等[13]通過最優橢圓建立了道砟模型庫,張徐等[14]采用球形單元構建道砟模型。雖然這類方法提升了模型構建效率,但生成的道砟模型與真實道砟的形態特征存在一定的差距。近年來,傅里葉分析法被廣泛應用于顆粒材料表面細部特征的研究及其輪廓的重構,如Touiti等[15]利用傅里葉描述符分析了碳酸鹽砂、硅砂和鈣質砂表面形態特征的差異,Mollon等[16]、方浩等[17]基于傅里葉變換隨機重構了砂土顆粒模型,石崇等[18]、郜家奇等[19]通過傅里葉變換隨機重構了碎石顆粒模型。基于傅里葉變換法重構的顆粒模型在很大程度上還原了顆粒的細節形態特征。但現有的構建方法必須通過控制傅里葉描述符實現輪廓重構,而實際運用中傅里葉描述符的確定需要通過復雜的計算與處理流程才能獲得,不利于顆粒模型的快速準確構建。

為解決上述顆粒模型構建方法的不足,本文選取高速鐵路特級道砟,提取顆粒的平面幾何特征,采用傅里葉分析法,量化分析道砟顆粒形態特征,在此基礎上利用響應面法建立道砟幾何特征與傅里葉描述符之間的映射關系,進而提出使用道砟幾何特征代替傅里葉描述符的道砟顆粒模型快速重構方法,并通過比較道砟重構輪廓與實際道砟輪廓的相似程度,評估所提方法的可行性。

1 道砟形態特征分析

顆粒的形態特征是構建顆粒模型的關鍵,其中最常用的是平面圖形的輪廓特征。近年來,圖形(像)分析方法被廣泛地應用于評估顆粒狀材料的大小、形狀、棱角性以及表面質地系數等屬性。這些方法主要從兩個角度對顆粒形態進行評估:一是采用測量得到的幾何參數及通過這些參數構造的指標分析顆粒的幾何形態;二是采用數字圖像處理技術將圖形轉換為可以進行數值運算的指標再進行分析,如傅里葉變換法、Hough變換法、灰度計算法和小波變換法等,其中傅里葉變換法可以較好地反映顆粒的細節特征,本文將應用傅里葉分析法來表征道砟顆粒的形態特征。

1.1 基于測量值的道砟幾何形態特征分析

得益于圖形和圖像分析方法的發展,道砟顆粒的幾何尺寸可以通過計算機準確獲取,避免了手工測量帶來的誤差。表1列出了顆粒材料常見的基本幾何尺寸及其含義,利用這些幾何尺寸可以進一步構造并量化道砟顆粒形態的指標,如長寬比、圓形度、球形度、棱角度、規則度以及粗糙度等[20-22]。

表1 顆粒基本幾何尺寸

1.2 基于傅里葉變換的道砟形態特征分析

在平面xOy坐標系下,道砟顆粒輪廓可以看作是一系列散點的集合,見圖1。假設任意一點Pi在道砟顆粒邊界上沿逆時針方向運動,則道砟顆粒的輪廓可以由一組坐標序列{(xn,yn|n= 0,1, 2, …,N-1}表示,其中,N為坐標序列的長度。

定義中心距離函數r(n)為道砟顆粒輪廓上的點到其幾何中心點的距離,即

(1)

式中:xc、yc為道砟顆粒幾何中心點的坐標,計算式為

(2)

(3)

根據離散傅里葉變換(DFT),中心距離函數可以表示為[23]

(4)

其逆變換表達式為

(5)

式中:k為數字化頻率,k= 0, 1, … ,N-1;F(k)為中心距離函數r(n)的傅里葉系數,這些系數在頻率上表達了輪廓的形狀特征,其中低頻系數表達了輪廓的宏觀特征,高頻系數表達了輪廓的細節信息;j為虛數單位。

由于中心距離函數只有實數序列,其傅里葉變換是對稱的,因此只有N/2個傅里葉系數表達了道砟輪廓的形狀特征。

在實際計算過程中,道砟的位置、縮放和旋轉信息對道砟形態的分析影響不大,為保證描述符具有平移、旋轉、尺度和起始點不變性,要對提取的描述符進行歸一化處理。已有研究表明[24-25],中心距離函數的傅里葉系數已具有平移不變性,其幅值|F(k)|具有平移、選擇和起始點不變性。由傅里葉變換公式(4)可知,F(0)代表的是中心距離的均值,因此,只需將原始得到的傅里葉系數幅值|F(k)|同除以|F(0)|便可實現尺度歸一化。經過上述處理,得到的具有平移、旋轉、尺度和起始點不變性的傅里葉描述符為

(6)

式中:|F(k)|為傅里葉系數幅值;Dn為傅里葉描述符,其典型頻譜圖見圖2。

圖2 傅里葉描述符頻譜

Mollon等[16]研究表明,D0反映顆粒的平均粒徑;D1反映中心函數中選擇的中心點與顆粒形心的偏移程度,當中心選擇在顆粒的形心時,D1為0;D2反映顆粒的伸長程度;D3～D7表征顆粒的宏觀特征,如三角形性、四邊形性等;D8～D64表征顆粒的微觀特征,如表面紋理、粗糙性等。此外,Mollon等[16]指出傅里葉描述符D2～D64呈對數線性遞減趨勢,并提出D3、D8線性組合方程近似表示3

Dn=2α·log2(n/3)+log2(D3)3

(7)

Dn=2β·log2(n/8)+log2(D8)8

(8)

式中:α、β均為線性系數。

2 基于真實幾何特征的碎石道砟顆粒形態重構方法

對道砟顆粒輪廓的構建是從細觀層面分析道床力學性能,揭示道床劣化機理的基礎工作之一。本節將通過傅里葉分析與響應面法建立道砟幾何尺寸與傅里葉描述符之間的映射關系,提出基于幾何形態特征快速重構道砟顆粒形態的方法。

2.1 道砟形狀樣本獲取

本文選取原礦為花崗巖的特級道砟作為試驗對象,其材質性能指標均符合相關規范。根據TB/T 2140—2008《鐵路碎石道砟》[26]對高速鐵路特級道砟粒徑的規定,本文選取粒徑范圍為22.4～31.5、31.5～40.0、40.0～50.0、50.0～63.0 mm四個粒徑組的道砟進行研究。為建立道砟幾何特征與傅里葉描述符之間的關系,試驗在每個粒徑組分別選取100個顆粒,共400個顆粒進行數字圖像采集與分析。為避免道砟上附著的細小顆粒對圖像采集產生不利影響,試驗前對選取的道砟進行了清洗。采用最大分辨率為1 920×1 080、有效像素為200萬的USB高清工業相機采集每個道砟顆粒圖像,進而運用數字圖形處理技術提取道砟輪廓信息。道砟輪廓信息獲取方法如圖3所示。首先對采集到的圖像進行直方圖均衡化處理,改善圖像灰度信息;其次采用雙邊濾波和形態學處理銳化道砟邊緣,提高分割的準確性;然后選取合適的閾值將道砟與背景分割開來;最后使用邊緣檢測函數尋找道砟輪廓,提取輪廓坐標序列。

圖3 道砟輪廓獲取方法

本研究中,單個道砟輪廓由長度超過3 000的坐標序列構成,直接采用離散傅里葉變換(DFT)會導致計算效率低下,因此本文采用快速傅里葉變換FFT代替DFT以提高計算效率。在使用FFT時,要求序列的長度N為2的指數倍數。為滿足該要求,利用邢小軍等[27]提出的重采樣方法,依次設置23、24、…、210、211個重采樣點,分別計算各重采樣圖形與原圖形面積的相對誤差,結果見圖4。由圖4可知,隨著采樣點數的增加,重采樣圖形的面積逐漸接近于原圖形的面積,當采樣點M的個數為26即64時,相對誤差為0.22%,當采樣點的個數增大到27即128時,相對誤差降低為0.000 11%,而當采樣點個數增大到28即256時,相對誤差為0.000 075 %。由此可見,當采樣點個數為128時,重采樣圖形便可較為準確地表示原圖形,與Das[28]研究結果一致。

圖4 采樣點個數對道砟形狀的影響

對道砟輪廓重采樣后,根據定義計算表1中道砟的幾何尺寸,同時采用FFT計算道砟輪廓形狀的傅里葉描述符,最終獲得400組道砟幾何尺寸和400組傅里葉描述符。方浩等[17]研究指出傅里葉描述符D2、D3、D8可以較全面地反映顆粒的形態特征,且其余描述符均可由這三個描述符通過線性組合方程近似表示[16]。因此,本文選取D2、D3、D8這三個傅里葉描述符作為表征道砟顆粒形狀的關鍵參數。圖5為本文道砟樣本的傅里葉描述符D2、D3、D8分布曲線。由圖5可知,不同粒徑下道砟顆粒的傅里葉描述符D2、D3、D8的分布曲線有著相似的變化規律,因此,可以推測道砟顆粒的幾何尺寸與其傅里葉描述符之間存在一定的聯系。

圖5 傅里葉描述符D2、D3、D8分布曲線

2.2 響應面模型建立

為建立道砟幾何尺寸與傅里葉描述符D2、D3、D8之間的映射關系,采用方差分析研究各幾何尺寸對傅里葉描述符的影響程度,并選取影響程度較顯著的因素建立響應面模型。在統計學中,自變量對因變量的影響程度可以通過顯著性水平p的大小反映,p越小,自變量對因變量的影響程度越大,當p<0.05時,表明自變量對因變量有顯著影響;當p<0.01時,表明自變量對因變量有極顯著影響。利用方差分析統計各幾何尺寸對傅里葉描述符D2、D3、D8的影響顯著性,結果見表2。

表2 幾何尺寸與D2、D3、D8的顯著性水平p

由表2可知,各幾何尺寸對D2均有顯著影響;除周長G和最小弗雷特直徑Fmin外,其他幾何尺寸對D3有顯著影響;僅長度L,最小外接圓半徑Rc,最大內切圓半徑Ri對D8有顯著影響。因此,為了在簡化模型的同時提高回歸模型的準確度,將幾何尺寸長度L,最小外接圓半徑Rc,最大內切圓半徑Ri作為自變量,D2、D3、D8分別作為因變量,采用響應面法建立幾何尺寸與傅里葉描述符之間的回歸模型,計算得D2、D3、D8的擬合回歸方程為

D2=0.111 835+0.010 545×L-0.008 854×Rc-

0.022 189×Ri+0.000 185×L×Rc-0.000 061×

L×Ri+0.000 109×Rc×Ri-0.000 080×L2-

(9)

D3=0.065 779-0.012 043×L+0.023 921×Rc-

0.004 476×Ri-0.000 130×L×Rc+0.000 098×

L×Ri-0.000 136×Rc×Ri+0.000 070×L2-

(10)

D8=0.006 736-0.000 482×L+0.001 748×Rc-

0.001 296×Ri-0.000 050×L×Rc-8.692 96×

10-8×L×Ri-8.524 71×10-8×Rc×Ri+

(11)

等高線圖可以直觀反映單一自變量對因變量的影響規律及兩個自變量之間的交互作用,當等高線呈橢圓形表明自變量交互作用顯著,當等高線呈圓形表明自變量交互作用不顯著。限于篇幅,本文僅說明長度L、最小外接圓半徑Rc及最大內切圓半徑Ri對D2的影響以及這三個幾何尺寸之間交互效應的強弱,結果見圖6,其中固定變量均取其變化范圍的中值。圖6(a)為長度L和最小外接圓半徑Rc對傅里葉描述符D2影響的等高線,其中最大內切圓半徑Ri為定值(Ri=20 mm)。由圖6(a)可知,在指定區間內,當Rc一定時,隨著L的增大,D2值呈先增大后減小的趨勢;當L較小時,隨著Rc的增大,D2逐漸減小,當L較大時,隨著Rc的增大,D2逐漸增大;圖中等高線呈橢圓形,表明在Ri不變的情況下,L和Rc的交互效應較強。圖6(b)為長度L和最大內切圓半徑Ri對D2影響的等高線,其中最小外接圓半徑Rc為定值(Rc=35 mm),由圖6(b)可知,當Ri一定時,隨著L的增大,D2呈先增大后減小的趨勢;當L一定時,隨著Ri的增大,D2逐漸減小;等高線呈橢圓形,表明在Rc不變的情況下,L和Ri交互效應較強。圖6(c)為最小外接圓半徑Rc和最大內切圓半徑Ri對D2影響的等高線,其中L為定值(L=90 mm)。由圖6(c)可知,當Ri一定時,隨著Rc的增大,D2呈先增大后減小的趨勢;當Rc一定時,隨著Ri的增大,D2逐漸減小;等高線近似圓形,表明在L不變的情況下,Rc和Ri交互效應較弱。

圖6 長度、最小外接圓半徑及最大內切圓半徑對D2的影響

2.3 相似性評估

道砟顆粒輪廓形狀復雜,僅僅依靠宏觀特征量不能對其相似程度進行全面評估,還需考慮道砟細節信息的相似程度。因此,參考周志宏等[29]相似度評估方法,結合道砟的宏觀特征量和細節信息提出一個綜合相似度因子用于評估道砟顆粒形態重構結果的可靠性。

在眾多道砟顆粒的宏觀特征量中,基本幾何尺寸不僅簡單易獲得,還可以有效反映道砟顆粒的形狀特征,因此本文選用表1所示的8個幾何尺寸作為描述道砟宏觀形狀的特征量,并引入百分比差異的概念比較兩個顆粒是否相似,其定義式為

(12)

式中:Δxi為基于表1列舉的幾何尺寸指標計算的重構道砟顆粒與真實道砟顆粒幾何特征量百分比差異;xi,F為重構道砟顆粒的幾何特征量;xi,R為真實道砟顆粒的幾何特征量。

(13)

(14)

幾何特征相似率可從宏觀特征上反映道砟的相似程度,但忽略了道砟顆粒之間細節信息的相似程度,因此引入道砟形狀差異率來評價道砟之間的形狀差異。道砟形狀差異率定義為兩個道砟非重疊部分的面積與平均面積的比值,見圖7。

圖7 道砟形狀差異率

在極坐標中,道砟顆粒和非重疊部分的面積可以采用積分方法求得。假設真實道砟與重構道砟的面積分別為AR、AF,非重疊部分的面積為Adif,則根據定義,道砟形狀差異率Rat可以表示為

Rat=Adif/Aave×100%

(15)

式中:Aave為真實道砟顆粒與重構道砟顆粒面積的平均值,即

Aave=(AR+AF)/2

(16)

要評估道砟重構輪廓與真實輪廓的相似度,需綜合考慮道砟幾何特征相似率和道砟形狀差異率,因此,采用加權的方法計算道砟顆粒的相似度,計算式為

Par=[r1·Gsr+r2·(1-Rat)]×100%

(17)

r1+r2=1

(18)

式中:Par為道砟顆粒相似度;r1、r2為加權因子。

綜上所述,通過式(5)～式(11)便可以根據道砟的長度L、最小外接圓半徑Rc以及最大內切圓半徑Ri重構道砟形態,并根據式(12)～式(18)計算重構道砟顆粒形態與真實道砟顆粒形態的相似程度,具體實現步驟見圖8。

圖8 道砟顆粒形態重構流程

3 基于真實幾何特征的碎石道砟顆粒形態重構方法的評價

為驗證本文提出的基于幾何真實形態的碎石道砟顆粒形態重構方法的可靠性,根據圖8所示流程重構了圖3中的道砟顆粒,結果見圖9。其中,圖9(a)為中心距離函數的對比,由圖可知重構道砟中心距離函數有著與真實道砟輪廓中心距離函數相似的變化規律;圖9(b)為輪廓形狀對比,由圖可知重構道砟的輪廓與真實道砟輪廓的形狀基本相同,并且很大程度上還原了真實道砟的細節特征;圖9(c)為重構道砟與真實道砟采樣點中心距離的對比,由圖可知重構值與真實值比較吻合,數據點落在直線y=x周圍,擬合優度R2達到了 0.982 1,表明重構道砟與真實道砟較為接近。

圖9 重構道砟輪廓形狀分析

表3所示為該道砟重構輪廓與真實輪廓的幾何特征量,據此計算重構道砟輪廓與真實道砟輪廓的幾何特征相似率Gsr為96.69%。根據面積的定義,在極坐標系下,計算道砟形狀差異率為4.81%。取r1=0.7,r2=0.3,根據式(17)計算圖9所示重構道砟輪廓與真實道砟輪廓的相似度為96.24%,表明重構道砟與真實道砟的幾何形態比較接近。

表3 幾何特征量計算結果

進一步地,將400個道砟樣本所有采樣點的重構值與真實值進行對比分析,結果見圖10。由圖10可見,所有樣本的重構輪廓與真實輪廓比較吻合,整體擬合優度為0.925 9。此外,采用所提方法對道砟重構輪廓與實際道砟輪廓的相似度進行分析,圖11為所有道砟重構輪廓與實際輪廓相似度的分布情況。由圖11可知,相似度的平均值為95.21%,其中80%以上道砟輪廓的重構值與真實值的相似度超過94%。結合圖10和圖11可知,本文提出的基于幾何形態的碎石道砟顆粒形態重構方法整體上可以較為準確的模擬真實道砟輪廓。然而,由于部分道砟的表面尖角較多,使用實數傅里葉變換會導致重構顆粒在尖角處發生稍許偏差,因此在今后的研究中,需注重道砟尖角的還原,以進一步提升道砟輪廓重構的準確性。

圖10 全部采樣點重構值與真實值對比分析

圖11 重構輪廓與真實輪廓相似度的統計分析

4 結論

本文利用數字圖像處理技術和快速傅里葉變換分析道砟顆粒表面形態特征,并采用響應面法建立碎石道砟幾何形態特征與傅里葉描述符的響應關系,提出基于真實幾何特征的碎石道砟顆粒形態重構方法,最后通過相似度分析法驗證了重構方法的可靠性,得出以下結論:

1)基于傅里葉變換法分析道砟形態特征時,當采樣點數取128時,可滿足道砟輪廓刻畫的精度要求,此時得到的傅里葉描述符D2、D3、D8可作為表征道砟輪廓的關鍵參數。

2)道砟顆粒長度L、最小外接圓Rc、最大內切圓半徑Ri對傅里葉描述符D2、D3、D8有顯著影響,可以采用響應面法建立兩者之間的映射關系,進而根據傅里葉逆變換重構道砟顆粒形態。

3)本文所提出的基于真實幾何特征的碎石道砟顆粒形態重構方法可以較為準確地刻畫道砟顆粒的細節特征,道砟重構輪廓與真實輪廓的平均相似度為95.21%,且80%以上道砟的相似度超過94%。