林海星，程三友，王曦，陳靜，辜平陽，莊玉軍，趙欣怡，馬剛

1) 長安大學(xué)地球科學(xué)與資源學(xué)院，西安，710054；2) 中國地質(zhì)調(diào)查局西安地質(zhì)調(diào)查中心，西安，710054

內(nèi)容提要: 賽什騰地區(qū)隸屬柴北緣構(gòu)造帶，巖性復(fù)雜多變，前人對研究區(qū)巖性界線的劃分較為籠統(tǒng)，筆者等選擇Landsat-8 OLI、ASTER和Sentinel-2A為數(shù)據(jù)源，采用最佳波段指數(shù)確定各影像的波段組合，凸顯不同巖石的邊界；將ASTER短波紅外波段與Landsat-8 OLI、Sentinel-2A可見光—近紅外波段協(xié)同處理，構(gòu)成Landsat-8+ASTER(LA)數(shù)據(jù)和Sentinel-2A+ASTER(SA)數(shù)據(jù)，分析重采樣巖石標(biāo)準(zhǔn)光譜信息，擬定不同巖石波段運(yùn)算公式，基于多重分形理論選定不同巖石類型的閾值范圍，獲取主要巖性的分布；根據(jù)重采樣的黑云母標(biāo)準(zhǔn)光譜曲線，選取SA數(shù)據(jù)2262 nm波段和2336 nm波段進(jìn)行定向主成分分析，采用Crosta法閾值分割第二主成分，劃分黑云母異常等級，將其與巖性分布相關(guān)聯(lián)，識別出研究區(qū)主要巖性分布。通過巖石實測光譜分析、薄片鏡下鑒定與野外地質(zhì)調(diào)查相結(jié)合的方法完善解譯結(jié)果。巖性提取結(jié)果顯示，小賽什騰山東側(cè)新發(fā)現(xiàn)“U”型條帶，為輝長巖—英云閃長巖—二長花崗巖，重新圈定達(dá)肯大坂群第三巖組和第四巖組，灘間山群一組，花崗巖，二長花崗巖，黑云母花崗巖，流紋巖，似斑狀石英閃長巖，石英閃長巖，英云閃長巖，輝長閃長巖和輝長巖等巖石的邊界。此次基于多源遙感數(shù)據(jù)巖性識別方法研究對青海賽什騰地區(qū)野外地質(zhì)調(diào)查工作具有指導(dǎo)意義，可為高山峽谷區(qū)地質(zhì)填圖提供技術(shù)參考。

遙感技術(shù)具有覆蓋面積大、獲取信息速度快、不受地形的限制等優(yōu)點，隨著空間信息技術(shù)的發(fā)展及影像獲取快捷便利，已經(jīng)成為區(qū)域地質(zhì)調(diào)查工作的重要技術(shù)手段之一(薛重生，1997；陳星等，2014；何鵬等，2016；史俊波等，2016；張志軍等，2016；潘明等，2019；劉小雨等，2020；李娜等，2021)。賽什騰地區(qū)隸屬柴北緣構(gòu)造帶，地勢陡峻、氣候干旱、植被覆蓋稀少，基巖出露集中(程三友等，2021；申燕玲等，2021)。近年來，主要集中在魚卡—大柴旦—德令哈一帶進(jìn)行礦化蝕變信息提取研究(康高峰等，2007；鞠崎等，2009；王亞紅等，2009；徐廣東等，2013；魯立輝，2019；張焜等，2019；賴華亮等，2020)；關(guān)于賽什騰地區(qū)遙感巖性信息提取的研究程度較低，李根軍等(2021)基于ZY1-02D遙感數(shù)據(jù)提取賽什騰地區(qū)大理巖與二長花崗巖的巖性分布，但并未對其他巖性進(jìn)行分析與識別。因此，筆者等以ASTER，Landsat-8 OLI和Sentinel-2A為數(shù)據(jù)源，采用約翰霍普金斯大學(xué)標(biāo)準(zhǔn)光譜庫(JHU)巖石光譜信息，利用最佳波段指數(shù)、波段運(yùn)算和定向主成分分析等信息增強(qiáng)方法(張玉君等，2003；王建梅等，2005；劉磊等，2009；吳志春等，2018)，對賽什騰地區(qū)進(jìn)行巖性識別。

圖1 柴達(dá)木盆地北緣小賽什騰山1∶50000地質(zhì)簡圖(據(jù)潘志明等，2016?修改)Fig.1 1∶50000 geological map of Xiaosaishiteng mountain in the northern margin of Qaidam Basin (modified according to Pan Zhiming et al., 2016# ?) Qhal—全新世沖積物；Qhdpl—全新世坡洪積物；Qp3pl—晚更新世洪積物；Pt1D2—達(dá)肯大坂群第二巖組；Pt1D1—達(dá)肯大坂群第一巖組；ηγ—二長花崗巖；πηγ—似斑狀二長花崗巖；πηο—似斑狀石英二長巖；δο—石英閃長巖；γδο—英云閃長巖；νδ—輝長閃長巖；ν—輝長巖Qhal—Holocene alluvium; Qhdpl—Holocene diluvium; Qp3pl—Late Pleistocene diluvium; Pt1D2—The second rock Formation of Daken Daban Group; Pt1D1—The first rock Formation of Daken Daban Group; ηγ—monzogranite; πηγ— porphyritic monzogranite; πηο— porphyritic quartz monzonite; δο—quartz diorite; γδο—tonalite; νδ—gabbro diorite; ν—gabbro

最佳波段指數(shù)確定各影像的波段組合，以此凸顯不同巖石的邊界，在小賽什騰山東側(cè)發(fā)現(xiàn)“U”型條帶。將ASTER短波紅外波段與Landsat-8 OLI、Sentinel-2A可見光—近紅外波段協(xié)同處理(戈文艷，2018)，構(gòu)成Landsat-8+ASTER(LA)數(shù)據(jù)和Sentinel-2A+ASTER(SA)數(shù)據(jù)，分析重采樣巖石標(biāo)準(zhǔn)光譜信息，擬定不同巖石波段運(yùn)算公式，基于多重分形方法選定不同巖石類型的閾值范圍(Mandelbrot，1974；成秋明，2006；李路，2020)，獲取主要巖性的分布。根據(jù)重采樣的黑云母標(biāo)準(zhǔn)光譜曲線，選取SA數(shù)據(jù)2262 nm波段和2336 nm波段進(jìn)行定向主成分分析(吳志春等，2018)，采用Crosta法閾值分割第二主成分(張玉君等，2003；劉磊等，2009)，劃分黑云母異常等級，將其與巖性分布相關(guān)聯(lián)。結(jié)合前人研究資料，對研究區(qū)主要巖性(二長巖、閃長巖、花崗巖、角閃石片麻巖和鈉長石片麻巖)進(jìn)行識別，綜合野外實地調(diào)查驗證、巖石樣品實測光譜分析和薄片鏡下鑒定，完善解譯結(jié)果，新發(fā)現(xiàn)輝長巖—英云閃長巖—二長花崗巖呈“U”型出露于小賽什騰山東側(cè)，研究區(qū)主要巖性有達(dá)肯大坂群第三巖組和第四巖組，灘間山群一組，花崗巖，二長花崗巖，黑云母花崗巖，流紋巖，似斑狀石英閃長巖，石英閃長巖，英云閃長巖，輝長閃長巖，輝長巖，并重新圈定各巖性的邊界，很好地指導(dǎo)了賽什騰山地區(qū)野外地質(zhì)調(diào)查工作。

1 研究區(qū)概況

研究區(qū)位于青藏高原柴達(dá)木盆地北緣(莊玉軍等，2020)，小賽什騰山呈北西走向，地勢北西高，南東低，切割強(qiáng)烈，地勢險峻，平均海拔3300 m，屬強(qiáng)烈剝蝕的構(gòu)造高山區(qū)(李根軍等，2021)。根據(jù)項目組成果(莊玉軍等，2020；王立軒等，2022)，研究區(qū)屬秦祁昆地層區(qū)的柴達(dá)木北緣小區(qū)，出露的地層為古元古代達(dá)肯大坂巖群(Pt1D)(陸松年等，2002)，巖體主要為石英閃長巖(δο)和黑云母花崗巖(γβ)(楊文軍等，2018；高萬里等，2019)(圖1)。區(qū)內(nèi)出露的古元古代達(dá)肯大坂巖群(Pt1D)為一套中深變質(zhì)的副變質(zhì)巖，劃分為混合片麻巖組(Pt1D1)和片麻巖組(Pt1D2)兩個巖組?；旌掀閹r巖組(Pt1D1)，主要為黑云斜長混合片麻巖(bpmg)夾混合巖化片麻巖；片麻巖巖組(Pt1D2)，主要為黑云斜長片麻巖(bpg)、黑云斜長角閃片麻巖(bphg)、二云母斜長片麻巖(dplg)等(郝國杰等，2004；陳世悅等，2016；王洪強(qiáng)等，2016)。研究區(qū)侵入巖主要呈大型巖基或巖株，東西向展布，侵入地層為達(dá)肯大坂巖群片麻巖組(Pt1D2)，與圍巖侵入界線清晰，主要為黑云母花崗巖(γβ)和石英閃長巖(δο)(楊文軍等，2018；高萬里等，2019)(圖1)。

2 數(shù)據(jù)與方法

2.1 數(shù)據(jù)源

本次研究采用的遙感數(shù)據(jù)ASTER、Landsat-8 OLI 和Sentinel-2A均來自于美國地質(zhì)調(diào)查局(https://earthexplorer.usgs.gov/)。在可見光—短波紅外范圍內(nèi)：ASTER數(shù)據(jù)有3個可見光—近紅外波段和6個短波紅外波段，空間分辨率分別是15 m，30 m(李進(jìn)波，2019)，時相為2000年10月16日，L1T級。Landsat-8 OLI數(shù)據(jù)有5個可見光—近紅外波段、2個短波紅外波段，空間分辨率30 m；1個全色波段，空間分辨率15 m，時相為2018年4月26日，L1T級(張玉君，2013)。Sentinel-2A是2015年6月23日發(fā)射，共計13個波段，波段2～4、8，空間分辨率10 m；波段5～7、8A、波段11～12，空間分辨率20 m；波段1、9～10，空間分辨率60 m(龔燃，2015)，本次研究采用的Sentinel-2A為 L1C級，時相為2019年9月29日。Sentinel-2A數(shù)據(jù)可見光—近紅外波段具有更高的光譜分辨率和空間分辨率，其可見光—近紅外波段的遙感巖性識別更具優(yōu)勢(王磊，2018)(圖2)。

圖2 ASTER、Landsat-8 OLI 和Sentinel-2A各波段空間分辨率、光譜分辨率比較(改自Masek et al., 2020)Fig.2 Comparison of spatial resolution and spectral resolution of ASTER、Landsat-8 OLI and Sentinel-2A (modified from Masek et al., 2020)

2.2 數(shù)據(jù)預(yù)處理

為了消除ASTER數(shù)據(jù)因傳感器之間的干擾而引起的誤差，需采用Cross-Talk軟件，對其進(jìn)行串?dāng)_校正(楊日紅等，2012)，輻射定標(biāo)串?dāng)_校正后的ASTER1～9波段，通過FLAASH模塊大氣校正輻射定標(biāo)后的ASTER數(shù)據(jù)(趙英時，2003)，選擇3次卷積法將大氣校正后的ASTER數(shù)據(jù)重采樣至30 m。Landsat-8 OLI和Sentinel-2A數(shù)據(jù)都需進(jìn)行輻射定標(biāo)和FLAASH大氣校正(龔燃，2015；張志杰等，2015；王磊，2018)。本次3種數(shù)據(jù)輻射定標(biāo)和大氣校正均在ENVI5.5.2軟件中實現(xiàn)。

2.3 融合增強(qiáng)處理

為有效利用多源遙感數(shù)據(jù)在光譜分辨率的優(yōu)勢，使不同遙感數(shù)據(jù)在巖性信息提取過程中能夠協(xié)同運(yùn)算，因此將不同空間分辨率的多源遙感影像，通過圖像融合增強(qiáng)處理，使其空間分辨率提高并達(dá)到一致。

通過對比Gram—Schmidt(GS)變換、Nearest Neighbor Diffuse(NND)變換和Principle Component Analysis(PCA)變換3種融合算法的結(jié)果，發(fā)現(xiàn)GS變換能夠保持影像的紋理和光譜信息；NND變換雖然保留了空間紋理信息，但是光譜信息損失較大，與原始影像相比，色彩出現(xiàn)失真；PCA變換色彩嚴(yán)重失真，色調(diào)變化較大(唐超等，2018；陳琪等，2021；李樹濤等，2021)。因此本研究選擇GS變換方法進(jìn)行融合，將預(yù)處理后的ASTER數(shù)據(jù)分別與Landsat-8 OLI全色波段和Sentinel-2A 10 m波段進(jìn)行融合；預(yù)處理后的Landsat-8 OLI多光譜波段與全色波段進(jìn)行融合；對于Sentinel-2A數(shù)據(jù)，由于其自身沒有全色波段，為了降低圖像融合對光譜的影響，通過計算10 m分辨率波段與其他波段的相關(guān)系數(shù)，分別選取與各波段相關(guān)系數(shù)最大的高分辨率波段進(jìn)行融合(杜挺，2015；梁麗娟等，2019)，最終得到4種融合數(shù)據(jù)為ASTER 10 m、ASTER 15 m、Landsat-8 15 m和Sentinel-2A 10 m，為后續(xù)參與巖性信息增強(qiáng)以及巖石光譜分析提供基礎(chǔ)數(shù)據(jù)。

3 巖性信息提取

3.1 圖像彩色增強(qiáng)

圖像色彩是最直觀區(qū)分地物邊界的依據(jù)，在巖性信息提取中選擇合適的波段組合，突顯研究區(qū)內(nèi)部不同巖性的色彩，為后續(xù)巖性邊界的劃分提供參考依據(jù)。根據(jù)最佳指數(shù)因子(OIF)(Chavez et al.，1982)，對數(shù)據(jù)進(jìn)行計算，OIF值越大說明組合波段所涵蓋的信息越多，冗余越小(王建梅等，2005；馬世斌等，2015)。計算公式如下：

(1)

式中Si為波段的標(biāo)準(zhǔn)差，Rij為波段之間的相關(guān)系數(shù)，計算得到數(shù)據(jù)的OIF指數(shù)。

由計算結(jié)果可知Landsat-8 15 m(756波段組合)，Sentinel-2A 10 m(938波段組合)，ASTER-15 m(724波段組合)和ASTER-10 m(924波段組合)的OIF值最高，信息量最大(表1)。因此分析目標(biāo)地質(zhì)體的診斷性光譜及各巖性之間的差異，結(jié)合假彩色合成圖像的地物色彩信息，確定了Landsat-8 15 m RGB756、Sentinel-2A 10 m RGB983、ASTER-15 m RGB724、ASTER-10 m RGB924 4種圖像，凸顯了不同巖性之間邊界(圖3)。

圖3 不同數(shù)據(jù)假彩色合成圖Fig.3 False color composite of different data(a)Landsat-8 數(shù)據(jù) 15 m，RGB756；(b)Sentinel-2A 數(shù)據(jù) 10 m，RGB938；(c)ASTER 數(shù)據(jù) 15 m，RGB724；(d)ASTER 數(shù)據(jù) 10 m，RGB924；(e)Landsat-8數(shù)據(jù) 15 m，“U”型構(gòu)造，RGB765；(f)Sentinel-2A數(shù)據(jù) 10 m，“U”型構(gòu)造，RGB983；(g)ASTER 數(shù)據(jù) 15 m，“U”型構(gòu)造，RGB742；(h)ASTER 數(shù)據(jù) 10 m，“U”型構(gòu)造，RGB942(a)Landsat-8 Data 15 m, RGB756；(b)Sentinel-2A Data 10 m, RGB938；(c)ASTER Data 15 m, RGB724；(d)ASTER Data 10 m, RGB924；(e)Landsat-8 Data 15 m,“U”structure, RGB765；(f)Sentinel-2A Data 10 m,“U”structure，RGB983；(g)ASTER Data 15 m,“U”structure，RGB742；(h)ASTER Data 10 m,“U”structure，RGB942

表 1 各數(shù)據(jù)RGB組合Table 1 RGB combination of data

通過觀察，可以發(fā)現(xiàn)在可見光波段具有較高光譜分辨率的Landsat-8和Sentinel-2A數(shù)據(jù)，各自的假彩色合成圖(圖3a和3b)與ASTER假彩色合成圖像(圖3c和3d)相比，顏色對比更加鮮明，可以清楚地看出巖石顏色主要為紫色，淺灰綠色，藍(lán)色和褐色4種色調(diào)，其中紫色巖石、淺灰綠色巖石和褐色巖石與其他巖石的邊界非常清晰。而在ASTER數(shù)據(jù)假彩色合成圖(圖3c和3d)中紫色整體顏色較淺，不易與周圍的淺灰綠色巖石區(qū)分，圖3c和3d右上角的淺灰綠色地物和褐色地物邊界模糊，無法進(jìn)行有效區(qū)分。結(jié)果表明，Landsat-8數(shù)據(jù)和Sentinel-2A(圖3a和3b)，要比相同分辨率的ASTER數(shù)據(jù)在巖石識別中更具有優(yōu)勢，巖石的邊界更加清晰，后期巖性識別過程中，有助于劃分巖性界線；且在該3種數(shù)據(jù)的假彩色合成圖像上均在小賽什騰山東側(cè)新發(fā)現(xiàn)半環(huán)形構(gòu)造(圖3e—3h)。

圖 4 重采樣巖石光譜曲線Fig.4 Resampled rock spectral curve(a)變質(zhì)巖重采樣光譜曲線： Ⅰ—鈉長石片麻巖—SA數(shù)據(jù)；Ⅱ—鈉長石片麻巖—LA數(shù)據(jù)；Ⅲ—角閃石片麻巖—SA數(shù)據(jù)；Ⅳ—角閃石片麻巖—LA數(shù)據(jù)；(b)巖漿巖重采樣光譜曲線：Ⅰ—閃長巖—LA數(shù)據(jù)；Ⅱ—閃長巖—SA數(shù)據(jù)；Ⅲ—二長巖—LA數(shù)據(jù)；Ⅳ—二長巖—SA數(shù)據(jù)；Ⅴ—花崗巖—LA數(shù)據(jù)；Ⅵ—花崗巖—SA數(shù)據(jù)；(c)Ⅰ—USGS標(biāo)準(zhǔn)光譜庫黑云母光譜曲線；Ⅱ—重采樣至SA數(shù)據(jù)重黑云母光譜曲線；Ⅲ—重采樣至LA 數(shù)據(jù)重黑云母光譜曲線(a) resampling spectral curves of metamorphic rocks: Ⅰ—albite gneiss SA data; Ⅱ—albite gneiss LA data; III—amphibole gneiss SA data; IV—amphibole gneiss LA data; (b) resampling spectral curves of magmatic rocks: I—diorite LA data; II—diorite SA data; III—monzonite LA data; IV—monzonite SA data; V—granite LA data; VI—granite SA data; (c) Ⅰ—spectral curve of biotite in USGS standard spectral library; Ⅱ—resample to SA data heavy biotite spectrum curve; Ⅲ—resample to LA data heavy biotite spectrum curve

3.2 波段運(yùn)算

巖石由于礦物成分的不同，其光譜吸收特征也存在一定的差異，因此采用波段運(yùn)算對研究區(qū)主要巖石進(jìn)行光譜分析，凸顯巖石信息。同時為充分利用不同數(shù)據(jù)的優(yōu)勢，筆者等將ASTER數(shù)據(jù)的短波紅外波段分別與Landsat-8 OLI和Sentinel-2A數(shù)據(jù)的可見光—近紅外波段進(jìn)行協(xié)同處理，最終構(gòu)成Landsat-8+ASTER(LA)數(shù)據(jù)和Sentinel-2A+ASTER(SA)數(shù)據(jù)(戈文艷，2018)。

3.2.1變質(zhì)巖信息提取

研究區(qū)出露的變質(zhì)巖主要為片麻巖(gn)，因此在約翰霍普金斯大學(xué)(JHU)標(biāo)準(zhǔn)波譜數(shù)據(jù)庫選擇與之相近的兩種變質(zhì)巖：角閃片麻巖和鈉長石片麻巖波譜曲線，將其重采樣至LA數(shù)據(jù)和SA數(shù)據(jù)(圖4a)。Sentinel-2A數(shù)據(jù)在可見光—近紅外波段中的光譜分辨率要優(yōu)于Landsat-8數(shù)據(jù)，鈉長石片麻巖在SA數(shù)據(jù)783 nm處有反射峰，704 nm有吸收谷，在LA數(shù)據(jù)則無法顯示這一光譜特征；在短波紅外波段上，角閃石片麻巖在2209 nm處有反射峰，2336 nm處有吸收谷，鈉長石片麻巖在2209 nm處有吸收谷，2262 nm處有反射峰(圖4a)。

選擇每種巖石的特征吸收波段，波段運(yùn)算方法凸顯巖性(趙英時，2003)。對于SA數(shù)據(jù)，采用角閃石片麻巖(B10-B12)/(B10+B12)，鈉長石片麻巖(B6-B5)/(B6+B5)；LA數(shù)據(jù)，采用角閃石片麻巖(B8-B10)/(B8+B10)，鈉長石片麻巖(B9-B8)/(B9+B8)。由于存在“異物同譜”、“同物異譜”的現(xiàn)象，因此波段運(yùn)算所得到的灰度圖像，其所包含的信息不僅有提取的目標(biāo)巖石，還可能包含其他巖石，而波段運(yùn)算又使得灰度圖像的信息高度集中，通過Crosta法對灰度圖進(jìn)行閾值選取，會使提取的巖性信息冗雜，所提取的巖性邊界無法與利用OIF指數(shù)獲得的假彩色合成圖像邊界對應(yīng)，而多重分形方法將灰度圖像的閾值提取問題轉(zhuǎn)換為散點圖直線擬合問題，在信息高度集中的情況下，放大了信息之間的差異，為巖性信息提取提供依據(jù)。因此進(jìn)行基于多重分形方法對每種巖性的比值灰度圖進(jìn)行閾值選定(Mandelbrot，1974；成秋明，2006；李路，2020)，得到常見變質(zhì)巖的巖性分布圖(圖5)。

3.2.2巖漿巖信息提取

研究區(qū)的主要侵入巖體為黑云母花崗巖和石英閃長巖，在約翰霍普金斯大學(xué)(JHU)標(biāo)準(zhǔn)波譜數(shù)據(jù)庫中，選擇相似的侵入巖如閃長巖、花崗巖和二長巖的波譜曲線，并將其重采樣至LA數(shù)據(jù)和SA數(shù)據(jù)(圖4b)。

Sentinel-2A在可見光—近紅外波段上具有更高的光譜分辨率，二長巖在SA數(shù)據(jù)704.1 nm處有吸收谷，665 nm和783 nm處有反射峰，在LA數(shù)據(jù)中無法顯示這一光譜吸收特征；閃長巖在560 nm處有反射峰，665 nm處有吸收谷；短波紅外波段中，花崗巖在2209 nm處有吸收谷，2167 nm處達(dá)到反射峰值；二長巖2167 nm處達(dá)到反射最大值，2336 nm處之后呈上升趨勢，顯出一個極小的吸收谷，整體反射較低，選取每類巖石的特征波段進(jìn)行波段運(yùn)算。對于SA采用二長巖(B6-B4)/(B6+B4)，閃長巖(B2-B3)/(B2+B3)，花崗巖(B9-B10)/(B9+B10)；LA，二長巖(B8-B10)/(B8+B10)，閃長巖(B7-B10)/(B7+B10)，花崗巖(B7-B8)/(B7+B8)。同樣采用多重分形方法對每種巖性的比值灰度圖進(jìn)行閾值選定(Mandelbrot，1974；成秋明，2006；李路，2020)，得到常見巖漿巖的巖性分布圖(圖5)。

圖 5 巖性綜合提取(Sentinel-2A數(shù)據(jù)RGB432)Fig.5 Comprehensive lithology extraction (Sentinel-2A data RGB432)Qhal—全新世沖積物；Qhdpl—全新世坡洪積物；Qp3pl—晚更新世洪積物；Pt1D4—達(dá)肯大坂群第四巖組；Pt1D3—達(dá)肯大坂群第三巖組；O—ST1—灘間山群一組；γ—花崗巖；ηγ—二長花崗巖；γβ—黑云母花崗巖；λ—流紋巖；πδο—似斑狀石英閃長巖；δο—石英閃長巖；γδο—英云閃長巖；νδ—輝長閃長巖；ν—輝長巖Qhal—Holocene alluvium；Qhdpl—Holocene diluvium；Qp3pl—Late Pleistocene diluvium；Pt1D4—The fourth rock Formation of Daken Daban Group；Pt1D3—The third rock Formation of Daken Daban Group；O—ST1—Tanjianshan Group I；γ—granite；ηγ—monzogranite；γβ—biotite granite；λ—rhyolite；πδο— porphyritic quartz diorite；δο—quartz diorite；γδο—tonalite；νδ—gabbro diorite；ν—gabbro

圖6 黑云母異常分布(Sentinel-2A數(shù)據(jù)RGB432)Fig.6 Biotite anomaly distribution (Sentinel-2A Data RGB432)Qhal—全新世沖積物；Qhdpl—全新世坡洪積物；Qp3pl—晚更新世洪積物；Pt1D4—達(dá)肯大坂群第四巖組；Pt1D3—達(dá)肯大坂群第三巖組；O—ST1—灘間山群一組；γ—花崗巖；ηγ—二長花崗巖；γβ—黑云母花崗巖；λ—流紋巖；πδο—似斑狀石英閃長巖；δο—石英閃長巖；γδο—英云閃長巖；νδ—輝長閃長巖；ν—輝長巖Qhal—Holocene alluvium；Qhdpl—Holocene diluvium；Qp3pl—Late Pleistocene diluvium；Pt1D4—The fourth rock Formation of Daken Daban Group；Pt1D3—The third rock Formation of Daken Daban Group；O—ST1—Tanjianshan Group I；γ—granite；ηγ—monzogranite；γβ—biotite granite；λ—rhyolite；πδο— porphyritic quartz diorite；δο—quartz diorite；γδο—tonalite；νδ—gabbro diorite；ν—gabbro

3.3 主成分分析法(PCA)

筆者等選擇黑云母作為特征礦物對研究區(qū)內(nèi)部含黑云母的巖石分布進(jìn)行分析。將USGS礦物波譜曲線庫中的黑云母波譜曲線重采樣到協(xié)同數(shù)據(jù)(圖4c)，觀察可以發(fā)現(xiàn)，重采樣至SA數(shù)據(jù)的黑云母光譜曲線更加貼合USGS礦物波譜曲線庫中的黑云母波譜曲線，因此選擇SA數(shù)據(jù)進(jìn)行定向主成分分析(吳志春等，2018)。

根據(jù)重采樣至SA數(shù)據(jù)的黑云母光譜曲線，筆者等選擇其處于反射峰2262 nm波段和具有吸收谷的2336 nm波段做定向主成分變換(圖4c)，故包含黑云母礦物信息的主成分的特征系數(shù)貢獻(xiàn)值符號相反(吳志春等，2018)。定向主成分變換得到主成分特征向量矩陣如表2，發(fā)現(xiàn)PC2的特征系數(shù)貢獻(xiàn)值符號相反，因此選擇PC2進(jìn)行黑云母異常信息提取。在主成分變換中，信息量主要集中在PC1內(nèi)，PC2所包含的信息較少，且PC2所包含的信息像元值分布符合正態(tài)分布，多重分形方法更適用于信息高度集中的情況，因此選擇Crosta法，以標(biāo)準(zhǔn)差和平圴值計算Mean+2σ，Mean+2.5σ，Mean+3σ 3個值為界線提取黑云母礦物信息，獲得黑云母異常分布圖(張玉君等，2003；劉磊等，2009)(圖6)。

表 2 SA數(shù)據(jù)的黑云母主成分變換特征向量矩陣，主成分變換特征值以及黑云母異常等級劃分

通過統(tǒng)計提取巖性與異常等級之間的重疊面積，由圖7a可以發(fā)現(xiàn)，鈉長石片麻巖和角閃石片麻巖同黑云母異常的重疊面積分別為0.0516 km2和0.0299 km2，而花崗巖、二長巖和閃長巖與黑云母異常的重疊面積為0.3747 km2，0.4698 km2，0.339 km2，鈉長石片麻巖與花崗巖、二長巖和閃長巖的面積比分別為1∶7.26、1∶9.10、1∶6.57；角閃石片麻巖與花崗巖、二長巖和閃長巖的面積比分別為1∶12.53、1∶15.71、1∶11.34，兩類巖石分別與黑云母異常的重疊面積有著明顯的差異；且花崗巖、二長巖和閃長巖在一級異常上的面積分別為0.1181 km2、0.053 km2和0.015 km2，二級異常為0.0795 km2、0.0994 km2和0.046 km2，三級異常為0.1771 km2、0.3174 km2和0.278 km2，可以看出一級異常中花崗巖重疊面積最多，高出二長巖0.0651 km2，而二級異常和三級異常則是二長巖，分別高出花崗巖0.0199 km2，閃長巖0.1403 km2。

圖 7 黑云母異常與巖性分布統(tǒng)計(a：不同巖性與黑云母異常的重疊面積；b：不同巖性與黑云母異常重疊面積占巖性的面積百分比)Fig.7 Biotite anomaly and lithology distribution statistics(a: overlapping area of different lithology and biotite anomaly；b: percentage of overlapping area of different lithology and biotite anomaly in lithology)

由圖7b可知，花崗巖、二長巖和閃長巖與黑云母異常重疊面積占巖性面積的百分比中，一級異常分別為1.43%、0.38%和0.029%，二級異常分別為0.96%、0.72%和0.088%，三級異常分別為2.15%、2.30%和0.53%，發(fā)現(xiàn)在一級和二級異常中花崗巖占比最高，但是花崗巖與二長巖在一級異常和二級異常中面積百分比的差距從1.05%下降至0.24%，而三級異常二長巖的面積占比升高并超過花崗巖0.15%，這一趨勢可以反映二長巖相比花崗巖更集中在黑云母的二三級異常中；綜合圖7a和圖7b可以發(fā)現(xiàn)重疊面積較少的鈉長石片麻巖卻分別有著0.64%、0.54%和1.88%的面積百分比，尤其是在三級異常中，其所占的面積百分比與花崗巖和二長巖的面積百分比相近，而與之相比角閃石片麻巖分別為0.029%、0.075%和0.33%，表明鈉長石片麻巖與角閃石片麻巖相比更加富集黑云母異常，這也足以顯示出黑云母異常在區(qū)分巖性中有著一定的成效。

4 結(jié)果分析與野外驗證

為了驗證四種波段假彩色合成圖像(圖3)的巖性邊界劃分的效果，半環(huán)形構(gòu)造是否真實存在，以及研究區(qū)常見巖性提取是否有效，黑云母礦物分布與巖性之間的關(guān)系。筆者等所在項目組于2020年和2021年夏對研究區(qū)進(jìn)行野外實地驗證，采集大量巖石樣品，采用Spectral Evolution公司的SR-3500地物光譜儀對16個典型巖石樣品進(jìn)行光譜測試，得到巖石樣品的連續(xù)光譜曲線。

圖 8 不同巖石光譜曲線Fig.8 Spectral curves of different rocks(a)JHU標(biāo)準(zhǔn)波譜庫光譜曲線： Ⅰ—鈉長石片麻巖，Ⅱ—角閃石片麻巖；(b)變質(zhì)巖實測光譜曲線：Ⅰ—黑云母片麻巖(新鮮面)，Ⅱ—黑云母片麻巖(風(fēng)化面)，Ⅲ—云母石英片巖(新鮮面)，Ⅳ—云母石英片巖(風(fēng)化面)；(c)JHU標(biāo)準(zhǔn)波譜庫光譜曲線：Ⅰ—花崗巖，Ⅱ—二長巖，Ⅲ—閃長巖；(d)巖漿巖實測光譜曲線：Ⅰ—英云閃長巖(風(fēng)化面)，Ⅱ—英云閃長巖(新鮮面)，Ⅲ—二長花崗巖(風(fēng)化面)，Ⅳ—二長花崗巖(新鮮面)，Ⅴ—輝長巖(風(fēng)化面)，Ⅵ—輝長巖(新鮮面)(a)JHU standard spectrum library spectral curve: I—albite gneiss; II—amphibole gneiss;(b)measured spectral curves of metamorphic rocks: I—biotite gneiss (fresh surface); II—biotite gneiss (weathered surface); III—mica quartz schist (fresh surface); IV—mica quartz schist (weathered surface);(c)JHU standard spectrum library spectral curve: I—granite; II—monzonite; III—diorite;(d)measured spectral curves of magmatic rocks: I—tonalite (weathered surface); II—tonalite (fresh surface); III—monzogranite (weathered surface); IV—monzogranite (fresh surface); V—gabbro (weathered surface); VI—gabbro (fresh surface)

4.1 巖石的實測光譜特征分析對比

筆者等挑選4個典型巖石樣品的實測光譜曲線(圖8)通過與標(biāo)準(zhǔn)光譜庫內(nèi)的光譜曲線對比分析，得出以下結(jié)果。黑云母片麻巖在可見光—近紅外波段上有一段反射率呈明顯上升，在1364 nm、1849 nm、2153 nm、2215 nm和2278 nm處有反射峰，1403 nm、1908 nm、2204 nm、2250 nm和2332 nm處有吸收谷；云母石英片巖在可見光—近紅外波段反射率平緩，在1032 nm處呈突然下降的趨勢，在1933 nm和2348 nm處與黑云母片麻巖的反射特征相反，均出現(xiàn)反射峰。通過與鈉長石片麻巖的光譜曲線對比，黑云母片麻巖的光譜反射特征與鈉長石片麻巖相近，無法區(qū)分；由于云母石英片巖在1933 nm和2348 nm處有反射峰，因此可與鈉長石片麻巖區(qū)分。角閃石片麻巖在2320 nm、2390 nm處有吸收谷，2370 nm處有小反射峰，可以利用云母石英片巖在2248 nm處的吸收谷與角閃石片麻巖區(qū)分(圖8a和圖8b)。

二長花崗巖在1418 nm、1913 nm、2210 nm、2276 nm和2347 nm處有吸收谷，1367 nm、2078 nm、2231 nm、2277 nm和2400 nm處有反射峰，其中1367 nm和1418 nm處的光譜特征與英云閃長巖，輝長巖不同，可以用其作為區(qū)分巖石依據(jù)。英云閃長巖在2400 nm處有吸收谷，與二長花崗巖，輝長巖的反射特征相反，可以區(qū)分輝長巖與其他巖石(圖8c)。輝長巖在2101 nm處為反射峰，而二長花崗巖，英云閃長巖則存在小的吸收谷，可以用于區(qū)分輝長巖與其他巖石。閃長巖在2340 nm、2400 nm處有吸收谷，2308 nm處有小反射峰，與英云閃長巖的光譜吸收特征相近，無法區(qū)分(圖8d)。

黑云母礦物分布集中分布在研究區(qū)的中部和西部，經(jīng)野外實際調(diào)查發(fā)現(xiàn)，黑云母一級異常和三級異常所在地區(qū)的巖石類型為黑云母花崗巖(γβ)、石英閃長巖(δο)和二長花崗巖(ηγ)；而將二級異常與三級異常重疊時，所在地區(qū)的巖石類型為達(dá)肯大坂巖群(Pt1D)。

4. 2 巖礦鑒定

為確保巖性類型劃分的準(zhǔn)確性，通過顯微鏡下觀察樣品薄片的鏡下特征，進(jìn)一步確定研究區(qū)內(nèi)各巖石的類型。筆者等選取6個典型薄片樣本進(jìn)行分析闡述(圖9)，具體的鏡下特征如下：

中粗?；◢弾r(圖9a)：花崗結(jié)構(gòu)、似斑狀結(jié)構(gòu)，塊狀構(gòu)造，主要由斜長石、石英組成，少量正條紋長石，遭受后期剪切作用礦物邊緣發(fā)生碎裂具弱定向性，長石具波狀消光、發(fā)生了鈉黝簾石化。

中粒二云二長花崗巖(圖9b)：中粒半自形—他形粒狀結(jié)構(gòu)，主要為石英、斜長石、微斜長石、白云母、黑云母，局部見蠕蟲結(jié)構(gòu)。

圖 9 采樣巖石顯微鏡下照片F(xiàn)ig.9 Micrographs of sampled rocks(a)中粗?；◢弾r；(b)中粒二云二長花崗巖;(c)片麻狀黑云母花崗巖；(d)糜棱巖化絹云石英片巖(原巖為酸性凝灰?guī)r);(e)黑云斜長片麻巖；(f)二云石英片巖(原巖為中酸性凝灰?guī)r)(a)medium coarse grained granite; (b)medium grained eryun monzogranite; (c)gneissic biotite granite; (d)mylonitized sericite quartz schist (original rock is acid tuff); (e)biotite plagioclase gneiss; (f)biotite quartz schist (original rock is medium acid tuff)

圖 10 綜合信息圖(Sentinel-2A數(shù)據(jù)RGB432)Fig.10 Comprehensive information map (Sentinel-2A data RGB432)

表 3 野外驗證光譜測試點坐標(biāo)Table 3 Coordinates of spectral test points for field verification

圖 11 巖性邊界野外實地驗證圖Fig.11 Field verification diagram of lithologic boundary(a)D1172處黑云石英片巖夾大理巖宏觀照片；(b)D1172處大理巖與黑云石英片巖邊界；(c)D1173處黑云母花崗巖；(d)D1174處黑云母花崗巖宏觀照片；(e—f)D1174處黑云母花崗巖近照(a) macro photo of biotite quartz schist with marble at D1172; (b) boundary between marble and biotite quartz schist at D1172; (c) biotite granite at D1173; (d) macro photo of biotite granite at D1174; (e—f) recent photo of biotite granite at D1174

片麻狀黑云母花崗巖(圖9c)：片狀粒狀變晶結(jié)構(gòu)，片麻狀構(gòu)造，變余花崗結(jié)構(gòu)，具定向性。主要為石英、斜長石、黑云母，少量白云母和矽線石，部分黑云母被擠壓定向和塑性變形，部分斜長石具蠕蟲結(jié)構(gòu)，原巖可能為黑云母花崗巖。

表 4 小賽什騰山目視解譯標(biāo)志Table 4 Visual interpretation sign of Xiaosaishiteng mountain

糜棱巖化絹云石英片巖(原巖為酸性凝灰?guī)r)(圖9d)：弱變質(zhì)變形域具變余凝灰結(jié)構(gòu)，變余碎屑主要為酸性斜長石晶屑，少量石英晶屑和長英質(zhì)巖屑；變余基質(zhì)主要為細(xì)粒長英質(zhì)和絹云母，發(fā)生了剪切定向和塑性流變，長石發(fā)生了絹云母化。部分長石晶屑和長英質(zhì)巖屑呈旋轉(zhuǎn)斑出現(xiàn)(長石旋轉(zhuǎn)斑具右性剪切特征)，整個巖石發(fā)生了剪切片理化。

黑云斜長片麻巖(圖9e)：花崗變晶結(jié)構(gòu)(半自形粒狀變晶結(jié)構(gòu))，片麻狀構(gòu)造，主要由斜長石組成，次要礦物為石英和黑云母，具定向排列，斜長石發(fā)生絹云鈉長石化和鈉黝簾石化，具有淺色礦物和暗色礦物分異條帶。

二云石英片巖(原巖為中酸性凝灰?guī)r)(圖9f)：變余凝灰結(jié)構(gòu)，片狀構(gòu)造。變余碎屑主要為石英晶屑、長石晶屑和長英質(zhì)巖屑，變余基質(zhì)為微細(xì)粒定向排列的長英質(zhì)和絹云母，變余晶屑和巖屑呈旋轉(zhuǎn)斑出現(xiàn)。部分變余長石晶屑包含渾圓狀、港灣狀石英，整個巖石發(fā)生了剪切定向。

為了驗證運(yùn)用遙感手段提取的巖性位置和范圍與實際情況的是否一致，特進(jìn)行野外實地驗證，將驗證點坐標(biāo)(表3)運(yùn)用ArcGIS10. 6 疊置在綜合巖性提取分布圖(圖10)上，野外露頭照片見圖11。

通過野外實地踏勘露頭觀察，在研究區(qū)西部D1172處(圖10b)觀察有灰色黑云石英片巖，片理發(fā)育，風(fēng)化破碎嚴(yán)重，發(fā)育細(xì)小沖溝，內(nèi)部夾灰白色大理巖條帶(圖11a和圖11b)；D1173處(圖10b)深灰色片麻狀黑云母花崗巖，整體地形切割強(qiáng)烈，風(fēng)化破碎嚴(yán)重，內(nèi)部沖溝發(fā)育(圖11c)；D1174處(圖10b)深灰色片麻狀黑云母花崗巖，整體呈深灰色—土黃色，山頂尖棱，山脊線波狀起伏，地形起伏明顯，地形切割強(qiáng)烈，風(fēng)化破碎嚴(yán)重，內(nèi)部細(xì)小沖溝發(fā)育，無植被，中低山地貌(圖11d—11f)。

發(fā)現(xiàn)遙感圖上顯示的半環(huán)形構(gòu)造最外圍巖石為輝長巖，內(nèi)部分別是英云閃長巖和二長花崗巖。野外重點踏勘發(fā)現(xiàn)，閃長巖和鈉長石片麻巖與實際情況非常吻合(圖11)。綜合野外實地調(diào)查驗證、巖石樣品實測光譜分析、薄片鏡下鑒定以及項目組的相關(guān)成果(莊玉軍等，2020；王立軒等，2022)，建立小賽什騰山目視解譯標(biāo)志(表4)，獲得小賽什騰山1∶50000遙感解譯圖(圖12)，很好地指導(dǎo)了賽什騰山地區(qū)野外地質(zhì)調(diào)查工作，也可為鄰區(qū)地質(zhì)填圖工作提供參考。

圖 12 柴達(dá)木盆地北緣小賽什騰山1∶50000遙感解譯圖 Fig. 12 1∶50000 remote sensing interpretation map of Xiaosaishiteng mountain in the northern margin of Qaidam BasinQhal—全新世沖積物；Qhdpl—全新世坡洪積物；Qp3pl—晚更新世洪積物；Pt1D4—達(dá)肯大坂群第四巖組；Pt1D3—達(dá)肯大坂群第三巖組；O—ST1—灘間山群一組；γ—花崗巖；ηγ—二長花崗巖；γβ—黑云母花崗巖；λ—流紋巖；πδο—似斑狀石英閃長巖；δο—石英閃長巖；γδο—英云閃長巖；νδ—輝長閃長巖；ν—輝長巖Qhal—Holocene alluvium；Qhdpl—Holocene diluvium；Qp3pl—Late Pleistocene diluvium；Pt1D4—The fourth rock Formation of Daken Daban Group；Pt1D3—The third rock Formation of Daken Daban Group；O—ST1—Tanjianshan Group I；γ—granite；ηγ—monzogranite；γβ—biotite granite；λ—rhyolite；πδο— porphyritic quartz diorite；δο—quartz diorite；γδο—tonalite；νδ—gabbro diorite；ν—gabbro

5 結(jié)論

(1)基于Sentinel-2A和Landsat-8 OLI數(shù)據(jù)在可見光—近紅外波段的光譜分辨率優(yōu)勢，采用最佳波段組合，進(jìn)行圖像色彩增強(qiáng)，結(jié)果表明該方法獲得假彩色合成圖像地物顏色豐富，界線清晰，與相同分辨率的ASTER數(shù)據(jù)在地物邊界劃分中更具優(yōu)勢。

(2)將ASTER數(shù)據(jù)的短波紅外波段與Landsat-8 OLI數(shù)據(jù)和Sentinel-2A數(shù)據(jù)的可見光—近紅外波段進(jìn)行協(xié)同處理，得到的SA數(shù)據(jù)和LA數(shù)據(jù)，最大程度上保留巖石和礦物的光譜診斷特征。基于巖石和礦物光譜診斷特征的不同，采用波段運(yùn)算的方法突顯協(xié)同數(shù)據(jù)中的不同巖性信息，結(jié)合前期的野外工作進(jìn)行驗證，可以為下一步的野外地質(zhì)工作提供指導(dǎo)。

(3)依照黑云母礦物在短波紅外波段上的波譜特征，基于定向主成分分析法圈定特征礦物集中分布區(qū)域，結(jié)合野外調(diào)查可知：黑云母的異常等級與所在地區(qū)的巖石類型有關(guān)，當(dāng)黑云母異常等級集中在特定范圍的時候，可以指示達(dá)肯大坂巖群(Pt1D)所在的范圍，為鄰區(qū)的巖性填圖工作提供參考。

致謝：西安地質(zhì)調(diào)查中心的何世平研究員、王立軒碩士以及長安大學(xué)的王盼龍博士指導(dǎo)，幫助完成野外調(diào)查工作；巖石光譜分析由長安大學(xué)成礦作用及其動力學(xué)實驗室劉磊教授和梁永春博士指導(dǎo)；巖石薄片鏡下鑒定由西安地質(zhì)調(diào)查中心葉芳研究員完成，在此一并表示衷心感謝。

注 釋/Note

? 潘志明, 提振海, 徐旭明. 2016. 青海省冷湖行委野駱駝泉地區(qū)三幅1∶5萬區(qū)調(diào)報告. 河北區(qū)域地質(zhì)礦產(chǎn)調(diào)查研究所.