王昱博，陳 清

（合肥工業(yè)大學(xué) 土木與水利工程學(xué)院，安徽 合肥 230009）

1 引言

巖體所受外部荷載情況復(fù)雜，荷載造成的內(nèi)部損傷的起裂部位、擴(kuò)展路徑隱藏于巖體內(nèi)部，難以直接觀測(cè)到損傷區(qū)域發(fā)展過(guò)程。因此，巖體損傷演化路徑和裂隙空間形態(tài)的表征對(duì)于研究巖體損傷演化的力學(xué)行為具有重要意義［1］?，F(xiàn)階段針對(duì)巖體（石）的三維重建大體上分為兩種：一類是利用數(shù)字圖像相關(guān)（DIC）、掃描電子顯微鏡（SEM）、計(jì)算機(jī)斷層掃描（CT）、聲發(fā)射（AE）等無(wú)損技術(shù)獲取巖石真實(shí)細(xì)觀結(jié)構(gòu)，經(jīng)處理在軟件中進(jìn)行三維重建還原，這種方法精細(xì)度高，但對(duì)設(shè)施、場(chǎng)地要求比較高，同時(shí)往往需要處理大量數(shù)據(jù)；另一類的方法從數(shù)值模擬的角度出發(fā)，通過(guò)程序算法模擬巖體（石）細(xì)觀損傷，利用二維位圖重構(gòu)三維形態(tài)，這種方法成本較低、具有較好的可重復(fù)性。

國(guó)內(nèi)外學(xué)者對(duì)巖體（石）細(xì)觀結(jié)構(gòu)的三維重建開(kāi)展了相關(guān)研究，取得了豐碩的成果。黃家國(guó)等采用掃描電鏡、核磁共振和計(jì)算機(jī)斷層掃描3種方法重建了孔隙頁(yè)巖的高密度三維模型［2］；熊自明等基于ASIFT及RANSAC算法，將隧道掌子面的數(shù)字?jǐn)z影進(jìn)行處理，實(shí)現(xiàn)非接觸式隧道掌子面裂隙紋理的重建［3］；Zhou X P等利用CT對(duì)含孔隙砂巖進(jìn)行研究，采用分層退火算法處理層位圖像重建砂巖的不規(guī)則特征［4］；張亞等利用Monte-carlo隨機(jī)模擬方法將水電站壩址區(qū)裂隙幾何參數(shù)導(dǎo)入模擬軟件中形成三維裂隙網(wǎng)格模型，對(duì)裂隙巖體滲流的研究帶來(lái)一定幫助［5］；何延龍等通過(guò)退火算法與分形算法結(jié)合的Hoshen-Kopelman算法將黏土礦物基團(tuán)的主要結(jié)構(gòu)特征分類，構(gòu)建了含黏土礦物分布的三維重建多孔介質(zhì)模型，較好地反映了主要黏土礦物的分布及產(chǎn)狀［6］；張挺等基于多孔介質(zhì)二維圖像和多點(diǎn)統(tǒng)計(jì)方法，通過(guò)訓(xùn)練圖像形成三維空間中的體素，重建出多孔介質(zhì)的重建模型［7］；趙延林等利用FLAC3D實(shí)現(xiàn)三維裂隙的隨機(jī)重建，并在重建模型中進(jìn)行細(xì)觀滲流計(jì)算，分析孔隙結(jié)構(gòu)對(duì)滲流的影響［8］。

學(xué)者們研究成果具有重要的應(yīng)用價(jià)值，但對(duì)于工程巖體而言，所研究對(duì)象基本為巖石試件，與工程尺度的巖體存在較大數(shù)量級(jí)差距，且模型重建過(guò)程相對(duì)復(fù)雜，方法較難遵循。本文在前人研究的基礎(chǔ)上，基于簡(jiǎn)便易行的圖像三維重建方法，結(jié)合數(shù)字圖像處理技術(shù)、數(shù)值模擬，于FLAC3D中建立含軟弱層的隧道數(shù)值模擬模型，對(duì)隧道模型掘進(jìn)過(guò)程中造成的損傷進(jìn)行三維重建。

2 圖像重建方法

利用圖像進(jìn)行三維重建的方法主要有三種，直接重建法、間接重建法和混合重建法［9］。其中直接重建法和間接重建法較為常用。

2.1 直接重建法

直接重建法大致分為四種算法，光線投射法、足跡表法、錯(cuò)切-變形法、體元投射法［10］。這類算法是一種根據(jù)圖片體素進(jìn)行投影從而重構(gòu)三維模型的技術(shù)，不生成中間幾何曲面。現(xiàn)今使用的直接重建法基本由光線投射法和足跡表法發(fā)展而來(lái)。光線投射法以一條光線穿過(guò)位圖序列構(gòu)成的體數(shù)據(jù)，在這條線上對(duì)體素進(jìn)行等距采樣，設(shè)置光照模型獲得、重構(gòu)模型，示意圖見(jiàn)圖1。而足跡表法則以構(gòu)成三維體數(shù)據(jù)的位圖存儲(chǔ)順序采集體素信息，讀取某種體素在每一幅位圖中的影響范圍，往一個(gè)平面投影體素，得到投影重構(gòu)模型。直接重構(gòu)法計(jì)算生成的三維圖像質(zhì)量好、穩(wěn)定，對(duì)于需要觀察內(nèi)部結(jié)構(gòu)信息的三維模型具有優(yōu)勢(shì)，缺點(diǎn)是計(jì)算量較大、計(jì)算時(shí)間長(zhǎng)。

圖1 光線投影法原理

2.2 間接重建法

間接重建法是CT技術(shù)進(jìn)行快速模型重構(gòu)常用的一種方法。這種方法從圖像構(gòu)成的數(shù)據(jù)集中根據(jù)閾值提取所需表面，并將此表面多邊形化，用不同的方法拼接多邊形，形成大量的中間曲面，疊加曲面以構(gòu)成三維表面模型［9］。間接重建法主要有移動(dòng)立方體法、輪廓拼接法、移動(dòng)四面體法、剖分立方體等算法［10］。典型的間接重建法是由Lorensen［11］提出的移動(dòng)立方體算法（MC），如圖2所示，該法通過(guò)閾值或極值等體素在立方體中剖分三角網(wǎng)格來(lái)近似等值面，將三角形網(wǎng)格組合拼接后，計(jì)算梯度向量為每個(gè)三角形生成法線，以便渲染陰影形成光照效果，構(gòu)建出物體的三維形態(tài)。這類算法簡(jiǎn)單、重構(gòu)速度快，但包括MC算法在內(nèi)的間接重構(gòu)法存在一些通病：有圖像前處理時(shí)間長(zhǎng)、模型精細(xì)度較低和內(nèi)部信息較難觀察等問(wèn)題。

圖2 移動(dòng)立方體算法的基本組合［11］

對(duì)于巖體的三維損傷重構(gòu)模型來(lái)說(shuō)，具有良好的可視性、精細(xì)度和易于觀察內(nèi)部損傷狀態(tài)等特性是十分必要的。matlab具有強(qiáng)大的可視化能力，間接重構(gòu)法和直接重構(gòu)法都可以在其中實(shí)現(xiàn)，但直接重構(gòu)法計(jì)算生成的三維模型質(zhì)量更好、更加精細(xì)，模型內(nèi)部信息均可呈現(xiàn)，用來(lái)展現(xiàn)巖體損傷的三維形態(tài)是一種較好的方法。因此，本文選擇在matlab中運(yùn)用直接重構(gòu)算法進(jìn)行基于圖像的巖體損傷重建。

3 重建圖像的獲取與損傷辨識(shí)

3.1 圖像切制

在計(jì)算機(jī)斷層掃描技術(shù)中，位圖的切制是按照一定順序或一定方向（或者理解為與位置有關(guān)）通過(guò)掃描或剖切研究對(duì)象得到的一系列二維平面圖。本方法中，首先在數(shù)值模擬軟件中編制損傷特征參量函數(shù)、圖片自動(dòng)切制函數(shù)，實(shí)現(xiàn)位圖的自動(dòng)切制。切制的方法類似于計(jì)算機(jī)斷層掃描層位圖的獲取方式，在特征參量指定量值范圍的云圖上沿模型高程方向，以一定的間距dd為距離切制多個(gè)平面，形成位圖圖序。在切制位圖時(shí)，減小間距dd以獲得更多位圖圖像，利于詳細(xì)描述損傷的重構(gòu)特征。圖3為損傷部位的提取與重建示意。

圖3 損傷部位的提取與重建

3.2 損傷部位辨識(shí)

位圖中損傷特征參量的辨識(shí)，其本質(zhì)就是通過(guò)技術(shù)手段對(duì)損傷特征參量位圖中關(guān)鍵像素點(diǎn)的位置和顏色信息進(jìn)行辨識(shí)。為此，我們需要尋找合適的量化評(píng)價(jià)巖體裂隙宏細(xì)觀力學(xué)效應(yīng)的特征參量，獲取損傷參量的位圖。

YAI（yield approache index）［12-13］是一種基于Mohr-Coulomb強(qiáng)度破壞準(zhǔn)則的評(píng)價(jià)巖體危險(xiǎn)性程度的定量指標(biāo)。YAI估計(jì)了應(yīng)力誘發(fā)巖體屈服，出現(xiàn)破裂的風(fēng)險(xiǎn)程度，YAI值在0～1內(nèi)時(shí)巖體處于發(fā)生損傷的塑性狀態(tài)，裂紋逐漸蔓延，當(dāng)大于1時(shí)即宣告巖體的破壞。本文以此指標(biāo)量化巖體損傷程度并于FLAC3D中建立數(shù)值模擬模型，在數(shù)值模擬模型的YAI云圖上顯示巖體損傷的范圍及程度，計(jì)算完成后切制多張巖體損傷狀態(tài)的位圖，供辨識(shí)和提取。

由于位圖是由一系列帶有顏色特征像素點(diǎn)進(jìn)行不同的排列和染色所組成的，位圖中不同的顏色表示著相應(yīng)像素點(diǎn)帶有不同的顏色閾值，故進(jìn)行圖像處理時(shí)，位圖只需要顯示發(fā)生損傷部分。包含損傷信息的YAI位圖和大多數(shù)圖片一樣是RGB顏色格式，但實(shí)際上RGB不能很好地反應(yīng)圖像的形態(tài)特征，計(jì)算時(shí)也存在數(shù)據(jù)量較大的問(wèn)題，所以我們所獲得的位圖需要進(jìn)行灰度處理［14-16］?；叶忍幚淼姆椒ǚ譃榫€性處理和非線性處理兩大類，常用的線性處理算法是經(jīng)典心理學(xué)算法，這種算法通過(guò)確定三種顏色的權(quán)重將位圖灰度化，圖像處理過(guò)程中采用式（1）的RGB分量加權(quán)平均算法將彩色圖處理成灰度圖，此種RGB分權(quán)方式得出的灰度圖效果最為理想，可以將巖體損傷區(qū)域與其他區(qū)域清晰地區(qū)分出來(lái)。流程見(jiàn)圖4。

其中，f為處理后的灰度值，R為紅色灰度值，G為綠色灰度值，B為藍(lán)色灰度值。

圖4 圖像處理流程

對(duì)于一幅處理完成的灰度位圖，圖中像素點(diǎn)的灰度值為0～255之間。移動(dòng)鼠標(biāo)至matlab中灰度圖的某一區(qū)域時(shí)，界面上會(huì)顯示此區(qū)域的灰度值，此時(shí)，以此區(qū)域的灰度值將當(dāng)前區(qū)域與其他區(qū)域分割開(kāi)，這種成為圖像分割界限的灰度值被稱為閾值。因此，分割出位圖中損傷區(qū)域與非損傷區(qū)域的關(guān)鍵是找到損傷區(qū)域的閾值。如圖4中損傷區(qū)域的灰度值為76，以76為閾值灰度圖中的損傷區(qū)域分割出來(lái)，使用式（2）進(jìn)行邏輯計(jì)算，將位圖分割成由黑、白色的像素點(diǎn)所組成的二值圖。后期在重構(gòu)損傷模型的步驟中辨識(shí)和提取二值圖中的損傷部分即可。處理好的二值圖的顏色閾值只有單個(gè)變量，在后期識(shí)別過(guò)程中控制單個(gè)變量即可實(shí)現(xiàn)辨識(shí)過(guò)程，簡(jiǎn)化了辨識(shí)和提取的難度。

其中g(shù)（i,j）為邏輯運(yùn)算符，i、j為像素點(diǎn)坐標(biāo)，f（i,j）為像素點(diǎn)灰度值，T為閾值。

由此，可以利用Matlab精確地識(shí)別出二值圖中巖體的損傷區(qū)域。此外，在數(shù)值模擬軟件的自動(dòng)出圖階段應(yīng)對(duì)每一幅圖進(jìn)行有規(guī)律地命名，在Matlab中根據(jù)同樣命名規(guī)律運(yùn)用循環(huán)語(yǔ)句將選用的位圖序讀入到矩陣中，以實(shí)現(xiàn)圖像的逐幅提取，達(dá)到自動(dòng)處理圖像的功能，如圖5所示。

圖5 二值圖位圖圖序

4 基于位圖圖序的三維重構(gòu)

4.1 數(shù)值模型描述

內(nèi)含軟弱帶的非均質(zhì)隧道進(jìn)行損傷重建。馬蹄形隧道直徑10m，隧道左上方含一條軟弱帶。隧道軸線沿地應(yīng)力的大主應(yīng)力方向，模型的尺寸為長(zhǎng)30m，寬15m，高30m，模型單元數(shù)為338333，節(jié)點(diǎn)數(shù)為58282。數(shù)值模擬模型如圖6所示。

表1 巖石力學(xué)參數(shù)

圖6 開(kāi)挖隧道計(jì)算模型

該隧道軸線沿地應(yīng)力中間主應(yīng)力方向，大主應(yīng)力為0.6MPa；中間主應(yīng)力為側(cè)向壓力0.45MPa；最小主應(yīng)力為垂直方向0.26MPa。隧道的開(kāi)挖方式采用全斷面一次性開(kāi)挖，不考慮支護(hù)結(jié)構(gòu)，隧道進(jìn)尺為3m一個(gè)工況，共分5個(gè)工況開(kāi)挖。

4.2 損傷模型重建

依據(jù)開(kāi)挖工況及已介紹的方法截取重建損傷模型所需位圖，截取間距為0.05m，自下而上截取，獲取模型的YAI云圖特征參量位圖600幅，為便于辨識(shí)損傷區(qū)域，將位圖進(jìn)行二值化處理，損傷區(qū)域?yàn)槲粓D中灰度值為0（黑色）的部分，部分位圖見(jiàn)圖7及圖8。

圖7 部分灰度位圖

圖8 部分二值化位圖

利用二值化的位圖進(jìn)行巖體損傷區(qū)域的三維重建，直觀展現(xiàn)損傷區(qū)域的空間形態(tài)，具體步驟如下：

（1）封裝成三維體數(shù)據(jù)。

將由FLAC3D導(dǎo)出的n幅損傷位圖以命名序列讀入matlab中，利用cat函數(shù)對(duì)選取的N幅圖像數(shù)據(jù)進(jìn)行三維體數(shù)據(jù)集的構(gòu)造，得到X×Y×N的三維矩陣。

（2）優(yōu)化三維體數(shù)據(jù)字段量。

所形成的三維體數(shù)據(jù)由多幅圖像構(gòu)成，為了防止計(jì)算時(shí)內(nèi)存溢出，采用reducevolume和smooth函數(shù)減少數(shù)據(jù)字段和平滑處理，減少的數(shù)據(jù)量比例越大體數(shù)據(jù)集中數(shù)據(jù)越少，因此在保證運(yùn)算順利的基礎(chǔ)上應(yīng)盡量抽取更多的數(shù)據(jù)以保證三維重構(gòu)的完整性。

（3）轉(zhuǎn)換位圖坐標(biāo)。

三維體數(shù)據(jù)矩陣D的數(shù)據(jù)量減少時(shí)，模型像素點(diǎn)坐標(biāo)可能會(huì)發(fā)生改變，需要檢查模型高度和寬度方向所占的像素點(diǎn)坐標(biāo)與其在模型對(duì)應(yīng)點(diǎn)的實(shí)際坐標(biāo)的一致性，參考模型在實(shí)際坐標(biāo)中的位置對(duì)選取云圖的像素點(diǎn)矩陣x、y、z進(jìn)行坐標(biāo)轉(zhuǎn)換。

（4）計(jì)算巖體損傷百分率。

通過(guò)識(shí)別每張位圖損傷區(qū)域的像素面積乘以高度值，計(jì)算損傷總體積，顯示于重建模型旁。

（5）提取等值面進(jìn)行三維構(gòu)造。

在完成模型像素點(diǎn)坐標(biāo)轉(zhuǎn)換后，調(diào)用函數(shù)isosurface從體數(shù)據(jù)提取等值數(shù)據(jù)繪制等值面。二值圖中損傷區(qū)域?yàn)楹谏谏幕叶戎禐?，提取出三維體數(shù)據(jù)中灰度值為0的網(wǎng)格點(diǎn)，連接仰角相同的點(diǎn)以形成等值面，在顯示平面累計(jì)投影，將等值面的頂點(diǎn)和面的信息傳遞給patch函數(shù)，定義和構(gòu)造所選區(qū)域圖像的顏色和光線。

（6）優(yōu)化顯示效果。

最后，基于三維體數(shù)據(jù)調(diào)用isonormal、lighting和view等函數(shù)及命令分別定義所選區(qū)域的光照效果、光線陰影以及觀察視角等屬性，通過(guò)繪制模型邊界輪廓圖顯示出破損三維體在模型中的具體位置，重構(gòu)模型如圖9所示。

圖9 YAI云圖及隧道損傷特征三維空間形態(tài)

圖9中的分圖分別為不同開(kāi)挖步對(duì)應(yīng)的巖體損傷特征三維空間形態(tài)，伴隨著開(kāi)挖的進(jìn)行，隧道各開(kāi)挖步造成的巖體損傷在三維空間中逐步累積，由于受到隧道左上方的軟弱層的影響，巖體損傷區(qū)域靠近軟弱層一側(cè)的發(fā)展和分布較多，損傷重建模型很好地反應(yīng)了這一點(diǎn)。結(jié)合圖9和圖10可以看出，初步開(kāi)挖3m時(shí)巖體的損傷體積(DR)損傷率(DV)分別為11.4183m3和0.090621%，至隧道開(kāi)挖完成，DR和DV分別為223.4181m3和1.7732%，巖體損傷程度得到了量化，損傷的累積變化近似看作線性關(guān)系。

圖10 巖體損傷體積和損傷率

因此，根據(jù)上述數(shù)值模擬和三維重建的結(jié)果，可以認(rèn)為該方法較好地量化了工程巖體的損傷程度、預(yù)測(cè)了工程作業(yè)過(guò)程可能造成的巖體損傷，形象地重現(xiàn)了模型中巖體損傷的形態(tài)和分布情況。

5 結(jié)論

本文基于圖像三維重建技術(shù)，結(jié)合數(shù)字圖像處理技術(shù)、數(shù)值模擬，進(jìn)行了巖體損傷的表征，得出以下結(jié)論：

（1）通過(guò)設(shè)置不同的特征參數(shù)，圖像重建方法能夠重建諸如位移、應(yīng)力、YAI、單元安全度等三維形態(tài)圖，有助于綜合分析巖體損傷機(jī)制，從而實(shí)現(xiàn)巖體內(nèi)部損傷三維形態(tài)的有效表征與直觀展示。

（2）隧道損傷重建案例表明：通過(guò)辨識(shí)和重構(gòu)的巖體內(nèi)部損傷空間結(jié)構(gòu)模型，能夠清晰地展現(xiàn)工程巖體損傷區(qū)的空間分布形態(tài)，有效分析、推測(cè)后續(xù)工況中工程巖體損傷的發(fā)生位置、擴(kuò)展路徑、演化特征等巖體損傷規(guī)律。

（3）圖像重建方法不受工程尺度效應(yīng)限制，可作為傳統(tǒng)巖體無(wú)損檢測(cè)方式的補(bǔ)充方案，輔助分析工程巖體穩(wěn)定性相關(guān)問(wèn)題，為工程災(zāi)害控制的優(yōu)化設(shè)計(jì)提供一定的參考，具有較好的適用性和工程價(jià)值。