基于3DMine的礦化域模型構建與資源量估算:以南美洲某淺成低溫熱液型金礦為例

王智綱, 楊晉升, 張瑞忠

(1. 招金有色礦業有限公司,山東煙臺265400; 2. 招金礦業股份有限公司,山東煙臺265400)

0 引 言

三維模型構建在地質資源儲量估算、找礦預測方面的應用越來越廣泛和規范(余牛奔等,2015;覃鵬等,2016;中華人民共和國自然資源部,2020;鄒進超等,2021)。傳統地質塊段法采用邊界品位、最低工業品位等指標圈定礦體及估算資源儲量,割裂了礦體的空間相關性,降低了資源估算的效率及找礦預測的準確性,削弱了資源儲量動態估算的優勢。目前,對于斑巖型、熱液脈型礦床,國內外采用低于邊界品位的礦化域模型進行資源儲量估算,不僅能大幅提高工作效率,而且能客觀反映礦化的空間分布規律,有利于找礦預測(杜菊民等,2022)。

通過對南美洲某淺成低溫熱液型金礦樣品進行統計分析,提出基于礦化域模型的資源量估算法和三維可視化表達,為資源儲量管理提供參考。

1 區域及礦床地質

安第斯成礦帶是全球重要的Cu、Au、Fe產區,根據地理位置可劃分為北、中、南安第斯成礦帶(盧民杰等,2016)。研究區地處北安第斯成礦帶北部的尼加拉瓜—哥斯達黎加金成礦帶,構造單元劃分屬于活動大陸邊緣的內巖漿弧,是環太平洋火山島弧的一部分,由科科斯(Cocos)板塊俯沖到加勒比板塊之下而形成,該火山弧帶自墨西哥南部經薩瓦爾多、危地馬拉、洪都拉斯、尼加拉瓜、哥斯達黎加至巴拉馬中部,長約2 000 km,產出眾多的斑巖型銅礦床和淺成低溫熱液型金礦床。

礦區出露地層為:中新—上新統Aguacate群安山質與玄武質火山巖、火山碎屑巖,Guacimal群中酸性侵入雜巖,上新—更新統Monteverde組含角閃石斑狀英安巖和流紋巖。金礦體賦存在Aguacate群中酸性火山巖中,主要巖性為英安巖與長英質陸源火山碎屑巖、凝灰巖、角礫巖(圖1)。

圖1 研究區區域地質圖(a)與礦床地質圖(b)1-第四系沖積物;2-更新—全新統安山質角礫巖;3-更新統熔巖流、凝灰巖和火山角礫巖;4-上新—更新統Monteverde組含角閃石斑狀英安巖和流紋巖;5-上新統Guacimal群中酸性侵入雜巖;6-中新—上新統Aguacate群安山質與玄武質火山巖、火山碎屑巖;7-花崗巖;8-火山泥流;9-礦化體;10-斷層線;11-斷層帶;12-熔巖流動方向;13-金礦床(點);14-研究區Fig. 1 Regional geological map (a) and deposit geological map (b) of study area

區域上,一系列平行排列的NW-SE向弧形斷裂(Liz斷層)的主要控礦斷裂整體傾向82°,傾角85°,斷面整體平直,局部形成S形彎曲,斷層帶內充填火山巖角礫。根據火山熔巖流動方向以及出現的火山角礫可知,礦區所處位置可能為火山活動中心,礦區侵入巖為細?；◢弾r、花崗閃長巖,多以脈狀產出,是與區內火山巖同源但不同形式的產物。

礦區礦化以金礦化為主,發育構造裂隙型熱液蝕變礦化和充填交代型細脈-浸染狀礦化,圍巖蝕變類型有絹云母化、硅化、青磐巖化、黃鐵礦化,硅化通常伴隨著貴金屬礦化,其兩側通常為伊利石-絹云母和黏土蝕變帶,形態上呈面狀、透鏡狀和囊狀,礦床成因類型為低硫化淺成低溫熱液型金礦床。

2 數據庫建立與樣品統計

2.1 鉆孔數據庫建立

鉆孔數據庫是按照地質數據結構來存儲和管理的庫,是礦山開采設計、資源量估算和找礦預測評價的基礎(周鄧等,2017)?；谌S地質建模技術,采用3DMine軟件建立鉆孔空間數據庫。①在保持原始數據完整性和準確性的基礎上,對收集到的地質圖、勘探線剖面圖、原始編錄數據進行整理,形成基礎地質資料。②依次提取工程定位表、測斜數據表、樣品分析結果表、地質要素表等。其中,工程定位表包含鉆孔的空間位置信息,測斜數據表包含鉆孔的方位、傾角、終止孔深等;樣品分析結果表包含樣品編號、采樣起止孔深、樣長、品位及其他;地質要素表主要包含巖性、斷層、蝕變等信息,是礦化域圈連與解釋的重要依據。通過鉆孔編號對上述文件建立關聯。③進行鉆孔數據庫的錄入、檢查和顯示。檢查內容主要為定位表與其他數據是否匹配、測斜表是否完整、樣品長度超限與采樣位置疊置等,數據檢查無誤后在圖形界面顯示(圖2)。

圖2 鉆孔三維空間位置示意圖Fig. 2 Schematic Diagram of 3D spatial location of drill holes

參與數據庫建庫的鉆孔共1 451個,孔深在36～500 m之間,采集Au化學分析樣品36 387件。

2.2 樣品統計分析

分別對原始樣品和礦化域樣品進行品位頻數和對數品位頻數統計,形成品位分布直方圖和曲線圖(圖3),對數直方圖總體呈正態分布。原始樣品Au品位最小值為0.01 g/t,最大值為182 g/t。品位平均值(1.303 71 g/t)與西舍爾估值(1.302 17 g/t)接近,但偏度(16.32)相對較大,需進行特異值處理。

圖3 原始與礦化域樣品品位分布直方圖和曲線圖(a) 原始樣品品位頻數;(b) 原始樣品對數品位頻數;(c) 礦化域樣品品位頻數;(d) 礦化域樣品對數品位頻數Fig. 3 Histograms and graphs of the grade distribution of samples in the original and mineralized domains(a) Grade frequency for original sample; (b) Logarithmic grade frequency for original sample; (c) Grade frequency for mineralized domain sample; (d) Logarithmic grade frequency for mineralized domain sample

2.3 礦化域邊界品位確定

為保證礦化域在形態上的完整性和連續性,需確定合理的邊界品位。根據原始與礦化域樣品品位分位數對比(表1),樣品在對數值-1.2左右存在微小拐點,對應的分位值接近20%,原始樣品品位接近0.07 g/t,大于該品位的樣品數占總數的80%,表明拐點之后的樣品金礦化連續性好,因此選取礦化域的邊界品位為0.07 g/t。

表1 原始與礦化域樣品品位分位數對比

3 礦化域模型構建

3.1 實體模型

礦化域實體建模的方案可概括為“平面+剖面”。根據地質規律和實際情況,結合鉆孔數據庫品位品級劃分與可視化表達結果,通過封閉的多段線在平面或剖面上圈定礦化域邊界,再將各剖面或平面之間的封閉多段線進行連接,建立礦化域實體模型。在各項連接過程中,采用NURBS圓滑工具對礦化域邊界進行平滑處理,驗證礦化域實體模型,確保模型無開放邊、自相交邊和無效邊,最終形成實體模型(圖4a)。

圖4 礦化域實體模型與塊體模型(a) 實體模型;(b) 塊體模型Fig. 4 Schematic diagram of entity model and block model of mineralized domain(a) Entity model; (b) Block model

3.2 塊體模型

塊體模型是品位賦值和資源儲量估算的基礎(馬東林等,2022)。單元塊體尺寸取決于礦體類型、規模、產狀、最小可采厚度、開采方式及勘探線網度等因素(滕正雙等,2015;劉曉寧等,2019;張磊等,2021),本次確立塊體尺寸為5 m×3 m×2 m,次分塊尺寸為2.5 m×1.5 m×1.0 m。為塊體模型建立相應關聯屬性,如資源類別、礦石密度、Au品位等。由于礦化形態呈面狀、透鏡狀,故無需建立旋轉塊體。圖4b為實體模型約束下最終形成的礦化域塊體模型,可以直觀地展示礦體在三維空間的形態及其變化特征,是礦化變化結構分析的基礎。

4 礦化變化結構分析

礦化變化結構分析是基于已建立的實體模型與塊體模型進行三維空間結構變化、品位分布及變化的分析,包括特異值的識別與處理、變異函數擬合。特異值指特高與特低品位,由于特高品位對品位均值、方差等影響較大,會引起品位高估,因此在變異函數擬合前需對特高品位進行識別與處理。

4.1 特異值識別與處理

特異值的存在會影響變量的分布特征,進而導致變異函數的平穩性下降。傳統識別方法有2種:①采用均方差的倍數確定特異值,即特高品位≥m+3σ(m為均值,σ為3倍均方差);②按照樣品品位變化系數來識別特異值。以分布密度曲線函數上拐點對應值作為特異值的下限(林吉飛等,2011),根據基本統計圖3c、d可知,樣品在高品位區間連續性相對較好,在對數值2附近出現拐點,該處對應的品位為91.85 g/t,因此將91.85 g/t作為特異值的下限。特異值確定后,觀察其所在空間位置,如果特異值樣品成帶分布,則將特異值樣品所在礦化域作為單獨的“富礦塊”對待。

特異值處理時應以構建的礦化域模型內部樣品點為對象,一般采用平均品位代替特異值的樣品品位,此次采用品位分布直方圖上出現不連續拐點處的品位代替特異值品位。特異值處理后需再次檢查,若平均值與西舍爾估值較接近,且偏度大幅減小,則特異值處理合理(楊豐銘等,2022),此次特異值處理后的平均值為1.368 59,西舍爾估值為1.368 85,偏度降為10.51,相對較合理。

4.2 變異函數擬合

變異函數是表示礦化范圍內區域化變量的相關關系和空間結構的數學工具(王炯輝等,2013),是普通克里格法估值的基礎與前提,地質統計學研究首先要進行區域化變量的結構分析,需確定的參數主要有塊金值、基臺值和變程:塊金值代表樣品變化的隨機成分;基臺值是變異函數趨于平穩狀態時達到的值;變程是變異函數趨于穩定時對應向量的長度,即各向異性影響范圍(宣良瑞等,2022)。

變異函數擬合前需進行樣品組合。為避免局部樣品分割過大或過小導致產生額外地質信息或過度平滑的情況,需綜合考慮組合樣品長度,原始樣長、最小可采厚度、勘探類型與網度以及塊體模型尺寸等均會影響組合樣長的選擇(黃松等,2019)。礦化域內樣品長度基本統計結果顯示,樣品長度多為1、2、4 m,平均值1.82 m,選擇2 m為組合樣品長度。通過樣品組合促使參與變異函數擬合的樣品平均值與西舍爾估值更為接近。

綜合考慮礦化域模型形態和組合樣品點的統計分析特征,采用球狀模型,設置16個扇區,依次確定礦化域品位具有最好連續性的主軸、次軸、短軸方向。通過不斷調整參數進行擬合,最終確立塊金值、基臺值和最大變程等球狀模型參數,優選變異函數模型。變異函數曲線特征(圖5)顯示,Au元素在3個方向上的變異函數在原點附近表現為塊金效應,反映Au品位在各方向呈隨機性變化。在一定范圍內,變異函數隨滯后距的增加而呈正相關,說明Au品位在此區間內是連續、有規律的,具有結構特性;超過此范圍,變異函數值在極限方差附近表現出小范圍波動,并最終趨于穩定。從變異函數擬合曲線可知,Au品位沿走向、傾向和厚度方向的基臺值均為1.21,沿各方向變程有所差異,表現出基臺值相同而變程不同的幾何各向異性。

圖5 變異函數擬合曲線Fig. 5 Variogram fitting graph(a) Principal axis; (b) Secondary axis; (c) Short axis

4.3 交叉驗證

目前,校驗變異函數擬合較為廣泛的方法有離散方差檢驗、交叉驗證等(趙向東等,2021),在3DMine中提供利用地質統計學進行校驗的交叉驗證法,其實質是比較克里格估值與原始值的偏差(施寶生等,2021)。因此,為確保理論變異函數在礦化域模型的適應性和可靠程度,需對擬合的變異函數模型進行交叉驗證(馬恒等,2019;景永波等,2021)。

表2和圖6為此次驗證結果。表2數據顯示,誤差均值(-0.005 6)與方差(0.436 6)趨近0,符合變異函數交叉驗證的判別依據,而克里格方差(0.928 4)與均值(1.150 7)趨近1,說明Au元素在礦床內具有恒定的變異性,因此統計的變異函數參數對礦化域進行品位賦值是合理、無偏的。圖6顯示,Au品位殘差分布聚集基本趨于0,直方圖顯示呈典型的正態分布,置信區間高于所要求的90%,說明在變異函數計算中并未出現系統的估值誤差。

表2 各向異性交叉驗證結果

圖6 交叉驗證殘差圖(a)與克里格估值直方圖(b)Fig. 6 Cross validation residual graph (a) and Kriging estimation histogram (b)

5 資源量估算與可視化

5.1 普通克里格估值

普通克里格法是資源估算最常用的方法之一(劉佶林等,2019),與距離冪次反比法相比,具有引入各向異性、充分考慮局部變異性的優點。經過詳細勘探,礦化域內控制程度較高,樣品基本統計結果具正態分布特征,可以采用普通克里格法(王珂等,2019)。搜索橢球體參數的設置,一部分采用搜索橢球體擬合的變異函數確定,其中主軸代表礦體方向,由礦床地質可知:礦體走向為NW-SE向,主軸方位159°,與礦體方向一致;次軸傾角90°,接近直立,與控礦斷裂產狀基本一致;主軸傾伏角為9.3°,礦體存在小角度傾伏。此外,主軸/次軸與主軸/短軸可參考變異函數取值(表3),結合實體模型沿礦體方向的測量長度,最終確定比例分別為2.56、2.14,反映礦化的分布呈帶狀特征。

表3 普通克里格估值參數

Au品位屬性賦值時,主軸搜索半徑需考慮工程間距和變程。最少與最多樣品數的設置可以減少樣品的聚集,提高估值精度,考慮到礦區工程控制程度相對較高、礦體厚度差異較大的特點,選擇最少樣品數為3,最多樣品數為12。品位賦值后進行塊體約束,將邊界品位(0.2 g/t)以上的塊體確定為礦體。

5.2 資源量可視化

資源量的三維可視化表達可為礦體的定量預測和評價提供指導(毛先成等,2016),資源量的分類通過塊體模型的屬性字段著色來實現。根據資源地質可靠程度分類,地質統計學資源儲量分類需綜合考慮單個礦種勘查規范中推薦的探礦工程間距,通過設置搜索橢球體的半徑來實現,可直觀地描述礦體內部資源類別、品位變化及不同品級的分布特征(牛聰聰等,2022)。

根據工程控制程度,將資源量分為控制和推斷兩類,不同級別的資源量對應的勘探間距不同。研究區低溫熱液型金礦表現為帶狀礦化,形態復雜,走向、傾向和厚度方向的變程值分別為185、86、102 m。通常工程間距在礦體變化最大的方向上布設較小,工程間距小于變程方可有效控制礦體的變化,因此將基本工程間距設置為50 m×50 m,低于該工程控制間距下的資源級別為控制資源量,而推斷資源的工程間距則相應放大1倍,即低于100 m×100 m的間距、高于50 m×50 m的控制間距。

Au品位賦值時,通常根據設置搜索半徑和參與估值的工程數對資源類別進行劃分。對控制資源量與推斷資源量按照資源類型級別由高到低分別進行,其對應的橢球體搜索半徑分別為50、100 m,對應的工程數分別為≥3、≥2,得出分類結果(圖7a)。結果顯示:靠近地表和淺部位置的工程控制程度較高,對應資源級別也較高;深部控制程度低,推斷級別資源占比較少。在走向上礦化連續性較高,這與構造裂隙型熱液蝕變金礦化和充填交代型細脈-浸染狀型金礦化的實際情況相吻合。

圖7 資源量分類結果示意圖(a) 按資源類別分類;(b) 按品位品級分類Fig. 7 Schematic diagram of resource classification results(a) Classification by resource category; (b) Classification by content grade

根據邊界品位、最低工業品位等工業指標對資源進行分類,劃分為低品位Au資源(0.2～0.5 g/t)、中等品位Au資源(0.5～1.5 g/t)和高品位Au資源(>1.5 g/t)。品位區間統計結果(圖7b)顯示,Au資源具有低品位連續性好而高品位連續性差的特征,存在明顯的聚集趨勢,可劃分為大小2個聚集區,可能為兩處不同熱液活動的中心。此外,根據資源分布可知三維空間存在小角度側伏,說明成礦作用呈南高北低的特征,礦體北側及傾向延伸方向為重要找礦靶區,上述認識有利于進行成礦預測和開采管理。

6 結 論

通過分析南美洲某金礦區礦床地質特征,結合三維可視化理論,建立淺成低溫熱液型金礦礦化域三維模型,在對樣品基本統計的基礎上,采用地質統計學方法對資源量進行了估算和可視化表達。

(1)基本統計分析是礦化域模型構建和資源量估算的前提和基礎,可為特異值處理、礦化域邊界品位的選取等關鍵步驟提供充分依據。

(2)礦化域模型的建立保證了礦化的完整性和連續性,展示了礦化的空間形態品位分布,有利于復雜礦化體的資源量估算。對礦化域模型生成的塊體模型按照資源級別或品位品級進行分類以及按不同要素進行分類,可根據可視化表達效果對資源量進行科學管理。

(3)通過變異函數的擬合和交叉驗證、特異值的優化等方法能夠確保普通克里格估算的合理性和準確性。相較于傳統方法,普通克里格估算可節省大量的估算時間,降低估算過程的誤差率,在數字化礦山建設中具有較好的示范作用和應用前景。