999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

基于峰值權重的巖心高光譜礦化蝕變信息提取

2015-12-25 07:13:08張杰林趙學勝
自然資源遙感 2015年2期
關鍵詞:特征方法

張 媛,張杰林,趙學勝,袁 博

(1.中國礦業大學(北京)地球科學與測繪工程學院,北京 100083;2.核工業北京地質研究院遙感重點實驗室,北京 100029)

0 引言

高光譜遙感將傳統圖像的空間維與光譜維信息融合為一體,其高光譜分辨率的特點使得本來在寬波段多光譜遙感圖像中不可探測的地物在高光譜遙感圖像中能夠被探測和區分[1]。鉆探是地質勘查的重要技術手段之一,巖心是鉆探獲取的重要成果。傳統的巖心礦物識別主要是通過人工方法判別,這種方法對人的經驗依賴性較高,容易出現錯判或漏判等現象,且需耗費大量人力和時間。礦化蝕變信息是找礦的重要標志,鉆孔巖心作為鈾礦勘探成果的主要載體記錄著地質體垂向變化信息[2],通過對其中與鈾成礦有關的蝕變信息的提取,可對縮小找礦范圍和提高找礦準確度發揮重要作用。利用高光譜遙感技術對巖心進行測量,獲取巖心的光譜數據,并利用其光譜連續性強、分辨率高和圖譜合一等特點[3]進行巖心光譜數據分析和巖心礦化蝕變信息提取,是目前礦物識別與填圖及巖心編錄的重要發展方向之一。

光譜匹配是一種高光譜礦化蝕變信息提取的重要方法。光譜角填圖(spectral angle mapping,SAM)法[4-5]是光譜匹配中最重要的方法之一,其簡單高效、標量乘不變的特點使其得到了廣泛的應用,是目前巖心礦化蝕變信息提取中常用的方法。但由于使用SAM方法時巖心的全部光譜都參與運算[6],導致大部分次要光譜的相似性會壓制少數特征光譜之間的差異,給后續處理中的閾值設定帶來困難[7]。另外,由于巖心光譜的同質多象和類質同象等現象[8](主要表現在吸收譜帶漂移、吸收深度以及波形寬度的不同),使得光譜特征匹配存在很多不確定性,導致SAM的實際應用效果有限。對于礦物種類不同而光譜特征相似的“異物同譜”現象,為了加大光譜相似度之間的差異、更好地區分不同礦物,一些學者通過采取在差異性較大的波譜區間設置權重[9]、分組計算和歸一化計算[10]等方法改進SAM,提高了礦物識別和區分精度。目前SAM仍是大規模巖心礦化蝕變信息提取的重要方法[11-12]。在采用該方法提取巖心礦化蝕變信息時,一般都希望將屬于同一蝕變種類而光譜特征不同的蝕變礦物分為一類,利用其診斷性光譜吸收特征來判斷不同類別的蝕變礦物。但實際操作中,有時會將本屬同一類的蝕變礦物判斷為不同類別的蝕變礦物而造成漏分、錯分現象;盡管可以同時取多條波譜曲線作為參考端元來改善這種情況,但端元數量越多,誤差累積越大,對閾值設置的要求越高,且在短時間內很難找全代表該類蝕變礦物的全部光譜曲線。

針對上述問題,本文采取在同類蝕變礦物光譜曲線中尋找特征區間、并對診斷性吸收峰位置設置權重的方法,對傳統的SAM法進行改進,旨在提高典型礦化蝕變礦物的高光譜異常信息識別精度及提取效率。

1 基于峰值權重的SAM改進

傳統的SAM法是通過計算光譜矢量之間的夾角來判斷參考波譜與像元波譜之間的相似性。設參考波譜為x,圖像中的像元波譜為y,則它們之間的光譜角為

夾角α越小,說明越相似。參考光譜可以是實驗室光譜或野外測定光譜,也可以是從圖像中提取的像元光譜[13]。光譜角的大小只與光譜矢量的方向有關,與輻射亮度無關,從而減弱了照度和地形的影響[6]。然而,SAM是基于圖像整體光譜特征的匹配方法,即對圖像所有波段的信息同等處理;對于吸收峰相同而整體曲線有較大差異的同質光譜,則會增加光譜角計算的誤差,閾值較小時會出現漏判,閾值較大時會將其他類別錯分到此類。盡管通過選取同一種蝕變的多條典型光譜曲線作為參考端元可以改善這一缺陷,但這要求對參考端元的提取要全面而準確,需要用足夠的時間和精力去尋找更好的方法來提取有“同類異譜”現象的參考端元。通常,從不同深度和不同空間位置取得的巖心中,同一種礦化蝕變的波譜形態都不盡相同,即這一段巖心中的參考端元未必適應于下一段巖心。因此,采用遍歷所有巖心圖像去全面、準確、快捷地尋找參考端元的方法缺乏可操作性而不可取。同時,由于端元波譜識別的過程是人工判別的過程,其數量的增多也會提升誤差出現的概率。

筆者通過對圖像像元進行分析發現:同類蝕變礦物的光譜曲線在某一波段區間的波形和吸收峰特征都十分相似;而在其他區間的差別則較大(圖1),這是造成SAM在同類異譜現象時產生較大誤差的主要原因,故將波形和吸收峰特征都十分相似的區間稱為特征區間。

圖1 同類異譜蝕變礦物典型波譜曲線Fig.1 Typical spectrum curves of altered minerals with same kind but different spectrum

為了減小非特征區間光譜波形差異對光譜角帶來的影響,可對每類蝕變礦物取特征區間[b1,b2](b1,b2分別代表特征區間的起、止波長),僅對像元在該特征區間內的光譜反射率進行處理(即局部光譜角分類)。由于每類蝕變的診斷性光譜吸收峰均包含在該類蝕變礦物的特征區間內,因而特征區間光譜也是區別不同類蝕變的重要局部信息。

因同一蝕變礦物的結構和組分不同,特征區間內吸收峰的形狀和寬度亦存在微小差異,且不同巖心圖像數據中同類蝕變礦物的診斷性光譜吸收峰的位置是確定的(有的可能會有1個像元內的微弱漂移),而不同物質的峰值位置則明顯不同。選取特征區間內的像元光譜進行的局部光譜角分類可以增強該區間波譜范圍內的光譜特征,凸顯光譜的吸收峰特性,拉大不同礦物在該區間的差異。考慮到峰值可能的漂移,本文在特征區間增加了光譜吸收峰位置附近同一波段范圍內光譜反射率的權重ω。對于取自遙感圖像的參考端元和圖像中的每一個像元,設反射率數組為 r=(x1,x2,…,xn),變換后反射率數組為r'=(y1,y2,…,yn)(其中n為波段總數);將有m個特征峰值的波段位置組成的數組表示為P=(a1,a2,…,am),其中aj表示第j個波峰的波段位置。i為第i個波段,在特征光譜區間[b1,b2]范圍內,對于每一個i∈[b1,b2],j∈[1,m]有

若i≠aj且i≠aj-1 且i≠aj+1,

則yi=xi;

否則

式中:ω為權重系數,應根據人工經驗進行選擇,ω=1時即為傳統的 SAM。由此再對特征區間[b1,b2]內的參考端元和圖像像元加權后的反射率數值計算光譜角,得到的結果可以拉大不同蝕變礦物在該光譜區間的差距,減小同類蝕變礦物光譜之間的差異,重點突出局部光譜信息;同時兼顧了峰值的位置,在閾值一定的情況下,能夠從遙感圖像中較好地提取蝕變礦物信息。此外,由于特征區間的縮小和吸收峰位置的強化,僅需1條代表該類型蝕變礦物的光譜曲線即可進行蝕變礦物信息提取,從而省去繁瑣的端元提取和人工判別、篩選等過程。

2 方法試驗與結果分析

本文采用的巖心高光譜圖像是利用高光譜測量系統對研究區內的巖心標本掃描得到的短波紅外(SWIR)數據,其技術指標及參數見表1。

表1 高光譜測量系統主要技術指標及參數Tab.1 Main technical performances of hyperspectral measurement system

2.1 數據預處理

為了更好地分析巖心高光譜數據的吸收特征和更準確地提取蝕變礦物信息,首先需要對截取的一小段巖心數據進行反射率轉換和噪聲去除。

本文采用統計學模型中的經驗線性法(empirical line,EL)對巖心數據計算反射率。EL定標方法假設圖像的DN值與反射率值R之間存在線性關系[14],即

式中:Gain為增益;Bias為偏移。在實際計算過程中,利用已測得的白板數據作為已知點的地面光譜,求得增益和偏移,進而計算巖心高光譜圖像反射率。

噪聲去除應盡量不改變光譜曲線在小范圍內的形態,即不改變吸收峰的位置。對反射率圖像中的每一個像元,利用其每奇數個波段求平均得到中間波段的反射率值,即

式中:b'i為處理后圖像的反射率值;bi為處理前圖像的反射率值;i為波段號;n為圖像的波段總數;m可根據去噪的程度選擇合適的值,噪聲較小時m可選擇1或2。“=bi”表示前m個波段和后m個波段的反射率值沒有參與計算,仍用原反射率值(如m=2時,表示第1和2波段與倒數第1和2波段的值沒有變化)。

該方法最大的優點是簡單易行,既有效地保留了波段間的相關性,又消除了噪聲對圖像光譜在小范圍內引起的劇變,提高了后續圖像處理和蝕變礦物信息提取的精度。

本文取m=2的方法可能會造成峰值在一個波段范圍內位置的改變(此情況出現的概率很小),但由于后續的光譜特征分析和識別、光譜端元提取以及對圖像的匹配分類均是在應用該方法去噪的基礎上進行的,因而不會影響圖像匹配分類的精度。

2.2 蝕變礦物光譜特征分析

蝕變礦物中不同離子(或基團)(如Fe2+,Fe3+,C和OH-等)的電子躍升、振動和轉動[15]使其在遙感圖像中具有診斷性光譜吸收特征,這些特征是識別蝕變礦物類別的基礎。本文中的巖心數據所在礦區圍巖蝕變發育,種類較多,與鈾礦有關的圍巖蝕變主要有赤鐵礦化、伊利石化、綠泥石化、碳酸鹽化、螢石化和堿交代等。綜合分析研究區主要熱液蝕變與鈾礦化的物譜關系,以便為深部鈾礦找礦提供基礎信息[2]。為便于分析,截取部分巖心數據進行試驗。巖礦鑒定結果顯示,該巖心段主要的蝕變為伊利石化(Al-OH)、綠泥石化(Mg-OH)和少量的碳酸鹽化(C)等3種類型,蝕變礦物幾乎覆蓋了整個巖心表面。

由于同種巖石或礦物因發育過程和發育狀態不同,其成分、結構及光譜特征就會產生一定的差異;同一種礦化蝕變也會因區域地質背景不同而使其波譜曲線呈現一定差異(如診斷性光譜吸收特征峰的漂移)。因此,為提高光譜匹配的精度,筆者選擇在圖像中確定波譜端元。

本文的巖心高光譜數據中,伊利石化、綠泥石化和碳酸鹽化礦物的診斷性光譜吸收特征[16]主要表現在短波紅外(SWIR,1 300~2 500 nm)波段:①伊利石化礦物由于鋁陽離子的電子躍遷以及Al-OH基團伸縮振動的合頻與倍頻作用,在短波紅外光譜段產生光譜吸收特征[17],特別是在2 210 nm和2 354 nm處形成明顯的吸收峰“二元結構”。在該巖心段的B80(R)B200(G)B48(B)假彩色合成圖像中顯示為亮白色;②綠泥石化礦物診斷性光譜吸收特征主要體現在2 258 nm(Mg-OH基團)和2 360 nm處,整體反射率較低,在該巖心段假彩色合成圖像中呈暗綠色點狀分布;③碳酸鹽化礦物由于受C影響,在2 342 nm處形成強烈的吸收峰(且左寬右窄),據此特征可區別于其他礦物。除上述特征外,3種蝕變礦物均在1 417 nm和1 910 nm處存在水汽吸收峰,且同一種物質的吸收峰均可能存在一個波段范圍內漂移。3種蝕變礦物在巖心高光譜圖像中的診斷光譜標志如表2所示。

表2 3種蝕變礦物在圖像中的診斷光譜特征Tab.2 Diagnostic spectrum features of three kinds of typical alteration minerals

2.3 特征區間選擇與吸收峰位置確定

根據伊利石化、綠泥石化及碳酸鹽化礦物的光譜特征,從巖心高光譜圖像中分別對3種蝕變礦物提取出多條端元光譜曲線(圖1)。

通過分析圖1中3種蝕變礦物的典型波譜曲線可以看出,在波段[182,256]區間內,伊利石化礦物在2 054~2499 nm波形相似且包括了所有的診斷性光譜吸收特征,吸收峰位置分別為2 210 nm,2 354 nm和2 445 nm;綠泥石化礦物在1 916 nm處的吸收峰之后波形相似特征明顯,但由于后幾個波段噪聲大,故取1 952~2 451 nm,即波段[165,248]為特征區間,主要的吸收峰位置為2 258 nm和2 360 nm;碳酸鹽化礦物的主要吸收峰為2 342 nm,所以取2 246~2 499 nm,即特征區間為[214,256]。由于僅需要對添加了峰值權重的特征區間內的曲線部分進行光譜角匹配,便可將同一類蝕變礦物的多條“異譜”端元(圖1中僅為一部分)簡化成1條端元參與分類計算。

2.4 蝕變信息提取與驗證

分別取3種蝕變礦物的1條波譜端元作為參考光譜,根據其不同的特征區間和峰值位置,利用IDL編程語言實現上述蝕變信息提取算法。在閾值相同的情況下,權重系數分別取ω=1(普通SAM),ω=2,ω=3和ω=4進行試驗。圖2示出用傳統的SAM算法(圖2(a))和基于不同峰值權重的匹配算法(圖2(b)-(d))提取蝕變礦物的結果。

圖2 不同權重下蝕變礦物提取結果對比Fig.2 Comparison among results of altered mineral classification by setting different weights

從圖2的分類圖中可以直觀地看出,用基于峰值權重的匹配算法提取的蝕變礦物比用傳統的SAM方法提取的蝕變礦物更為全面而準確。在該巖心段圖像與成礦有關的蝕變中,伊利石化占主導,成面狀分布;伴有綠泥石化呈點狀分布;而碳酸鹽化較少。通過對蝕變礦物提取結果的分析和統計,可以較為精確地對蝕變信息進行定位,并促進定量化遙感在巖心礦化蝕變信息提取中的應用。

為了驗證本文礦化蝕變信息提取方法的有效性,分別對傳統SAM和在不同權重系數下基于峰值權重匹配方法提取的蝕變礦物分類結果進行總體精度和Kappa系數的計算和對比。總體精度反映了正確分類像元占總分類像元的比例;Kappa系數則從統計意義上表明被評價分類結果在多大程度上優于隨機分類結果。利用ERDAS軟件中的Accuracy Assessment模塊進行了分類精度評估,該方法是將專題分類圖像中的特定像元與已知類別的參考像元進行比較(即對分類數據與地面真值進行對比)。由于截取的區域較小,在圖像中隨機選取200個點(包括背景值),對每個隨機點查驗該點像元的光譜曲線,判斷分類結果是否與地面真值相符,得到分類精度(表3)。

表3 不同權重系數匹配方法的分類精度Tab.3 C lassification accuracy of differentm atching methods by setting different weights

從表3可以看出,由于強化了特征區間內吸收峰的特點,基于峰值權重的方法比傳統的SAM方法精度高。在ω=2時,分類總體精度可提高20%以上。但ω的取值并不是越大越好,而是隨著ω的增大,提取精度降低。在實際工作中,ω要根據圖像信噪比情況和同類物質吸收峰差異的大小來選擇。

3 結論

從巖心高光譜數據中提取礦化蝕變信息常用的SAM方法,由于自身基于整體光譜匹配的局限性,使其在信息提取時造成較大誤差。因此,本文通過對局部光譜特征區間中的峰值合理地添加權重,對傳統的SAM法進行改進,使其能夠克服上述缺陷。實驗結果證明:

1)用改進的SAM法提取蝕變礦物的總體精度從58%提高到70%~80%(取決于不同的權重系數)。

2)由于本文方法對端元的數量要求不高,每種蝕變礦物僅需1~2條端元光譜曲線即可參與光譜角計算,因而在提高分類精度的同時,也提高了端元選擇和提取的效率。特別是在大數據量、大規模的巖心蝕變礦化信息提取和復雜繁瑣的巖心編錄工作中,本文方法應用效果顯著。

3)通過對與鈾成礦關系密切的3種礦化蝕變信息的提取和分析,可為進一步找出鈾礦富集規律、縮小普查勘探的范圍、提高找礦效率提供依據。

下一步的工作是開展如何根據巖心高光譜圖像的信噪比及光譜吸收峰深度進一步確定峰值位置的權重系數,并從中尋找規律的研究。

[1] 張婷婷,唐菊興,郭 娜,等.高光譜遙感在斑巖礦床蝕變信息提取中的應用[J].礦物學報,2011(s1):922.Zhang T T,Tang JX,Guo N,etal.Application of hyperspectral remote sensing in alteration information extraction of porphyry deposits[J].Acta Minalogica Sinica,2011(s1):922.

[2] 張杰林,黃艷菊,王俊虎,等.鈾礦勘查鉆孔巖心高光譜編錄及三維礦物填圖技術研究[J].鈾礦地質,2013,29(4):249-255.Zhang JL,Huang Y J,Wang JH,et al.Hyperspectral drilling core logging and 3Dmineralmapping technology for uranium exploration[J].Uranium Geology,2013,29(4):249-255.

[3] 甘甫平,王潤生.高光譜遙感技術在地質領域中的應用[J].國土資源遙感,2007,19(4):57-60.doi:10.6046/gtzyyg.2007.04.13.Gan F P,Wang R S.The application of the hyperspectral imaging technique to geological investigation[J].Remote Sensing for Land and Resources,2007,19(4):57-60.doi:10.6046/gtzyyg.2007.04.13.

[4] Kruse F A,Lefkoff A B,Boardman JW,et al.The spectral image processing system(SIPS)-interactive visualization and analysis of imaging spectrometer data[J].Remote Sensing of Environment,1993,44(2/3):145-163.

[5] Robila S.An investigation of spectralmetrics in hyperspectral image preprocessing for classification[C]//Geospatial Goes Global:from Your Neighborhood to theWhole Planet.ASPRSAnnual Conference,Baltimore,Maryland.2005:7-11.

[6] Hecker C,van der Meijde M,van derWerff H,et al.Assessing the influence of reference spectra on synthetic SAM classification results[J].IEEE Transactions on Geoscience and Remote Sensing,2008,46(12):4162-4172.

[7] 王 瀛,郭 雷,梁 楠.逐波段修正負相關的光譜角填圖算法[J].火力與指揮控制,2013,38(2):22-25.Wang Y,Guo L,Liang N.A spectral anglemapper algorithm modified negative correlation band by band[J].Fire Control& Command Control,2013,38(2):22-25.

[8] 甘甫平,王潤生.遙感巖礦信息提取基礎與技術方法研究[M].北京:地質出版社,2004.Gan F P,Wang R S.Studies on Basis and Technique Methods of Remote Sensing Minerals Information Extraction[M].Beijing:Geological Publishing House,2004.

[9] 何中海,何彬彬.基于權重光譜角制圖的高光譜礦物填圖方法[J].光譜學與光譜分析,2011,31(8):2200-2204.He Z H,He B B.Weight spectral angle mapper(WSAM)method for hyperspectral mineralmapping[J].Spectroscopy and Spectral Analysis,2011,31(8):2200-2204.

[10] 唐 宏,杜培軍,方 濤,等.光譜角制圖模型的誤差源分析與改進算法[J].光譜學與光譜分析,2005,25(8):1180-1183.Tang H,Du P J,Fang T,et al.The analysis of error sources for SAM and its improvement algorithms[J].Spectroscopy and Spectral Analysis,2005,25(8):1180-1183.

[11] 別小娟,張廷斌,孫傳敏,等.藏東羅布莎蛇綠巖遙感巖礦信息提取方法研究[J].國土資源遙感,2013,25(3):72-78.doi:10.6046/gtzyyg.2013.03.13.Bie X J,Zhang T B,Sun CM,etal.Study ofmethods for extraction of remote sensing information of rocks and altered minerals from Luobusha ophiolite in east Tibet[J].Remote Sensing for Land and Resources,2013,25(3):72- 78.doi:10.6046/gtzyyg.2013.03.13.

[12] 張志軍,甘甫平,李賢慶,等.基于ASTER數據的蝕變礦物信息提取——以哈密黃山銅鎳礦區為例[J].國土資源遙感,2012,24(2):85-91.doi:10.6046/gtzyyg.2012.02.16.Zhang Z J,Gan F P,Li X Q,etal.The extraction of altered mineral information based on ASTER data:A case study of the Huangshan copper- nickel ore district in Hami[J].Remote Sensing for Land and Resources,2012,24(2):85-91.doi:10.6046/gtzyyg.2012.02.16.

[13] 徐元進.面向找礦的高光譜遙感巖礦信息提取方法研究[D].武漢:中國地質大學,2009.Xu Y J.Research of Prospecting- Oriented Approaches to Information Extraction of Rocks and Minerals Using Hyperspectral Remote Sensing Data[D].Wuhan:China University of Geosciences,2009.

[14] 鄧書斌.ENVI遙感圖像處理方法[M].北京:科學出版社,2010.Deng S B.Processing Methods of ENVIRemote Sensing Imaging[M].Beijing:Science Press,2010.

[15] 地質部情報研究所.礦物巖石的可見-中紅外光譜及其應用遙感專輯(第一輯)[M].北京:地質出版社,1980.Information Research Institute of Geology Ministry.Remote Sensing Album of the Mineral Rock Visible-infrared Spectroscopy and Its Application(the First Series)[M].Beijing:Geological Publishing House,1980.

[16] 傅 錦,裴承凱,韓曉青,等.基于多元統計技術的鈾礦蝕變信息高光譜模型[J].世界核地質科學,2007,24(3):166-171.Fu J,Pei C K,Han X Q,etal.Hyperspectralmodel for information of uranium mineral alteration technique based onmathematical statistics analysis[J].World Nuclear Geoscience,2007,24(3):166-171.

[17] 雷天賜,祝明強,周萬蓬.高光譜數據挖掘在蝕變礦物識別與提取中的應用[J].資源環境與工程,2005,19(3):213-219.Lei TC,Zhu M Q,Zhou W P.Application on datamining of hyperspectrum to identification and extraction of alterationminerals[J].Resources Environment& Engineering,2005,19(3):213-219.

猜你喜歡
特征方法
抓住特征巧觀察
新型冠狀病毒及其流行病學特征認識
如何表達“特征”
不忠誠的四個特征
當代陜西(2019年10期)2019-06-03 10:12:04
學習方法
抓住特征巧觀察
用對方法才能瘦
Coco薇(2016年2期)2016-03-22 02:42:52
四大方法 教你不再“坐以待病”!
Coco薇(2015年1期)2015-08-13 02:47:34
賺錢方法
捕魚
主站蜘蛛池模板: 99re经典视频在线| 成年人国产网站| 日韩毛片免费| 精品久久综合1区2区3区激情| 中美日韩在线网免费毛片视频| 国产精品福利导航| 成人午夜视频在线| 青青操国产视频| 亚洲成人动漫在线观看 | 国产一区二区免费播放| 国产导航在线| 尤物亚洲最大AV无码网站| 欧美日韩亚洲国产| 麻豆精品在线| 亚洲人在线| 免费在线一区| 特级做a爰片毛片免费69| 婷婷色在线视频| 91精品人妻一区二区| 久草网视频在线| 亚洲欧美日韩另类在线一| 毛片视频网| 自拍欧美亚洲| 毛片免费在线视频| 欧美日韩第二页| 欧美第九页| 亚洲天堂网在线观看视频| 国产成熟女人性满足视频| 国产亚洲精品精品精品| 国产成人精品高清不卡在线| 91亚洲影院| 精品三级网站| 亚洲色欲色欲www在线观看| 婷婷综合在线观看丁香| 午夜精品一区二区蜜桃| 国产高清免费午夜在线视频| 中文成人无码国产亚洲| 人妻夜夜爽天天爽| 日韩在线欧美在线| lhav亚洲精品| 五月婷婷导航| 亚洲精品色AV无码看| 亚洲久悠悠色悠在线播放| av在线5g无码天天| 日韩在线欧美在线| 国产视频一二三区| 手机在线国产精品| 成年人福利视频| 国产97公开成人免费视频| 国产乱子伦视频在线播放| 青草午夜精品视频在线观看| 一级毛片免费不卡在线视频| 亚洲国内精品自在自线官| 亚洲大尺码专区影院| 亚洲色图综合在线| 制服丝袜一区二区三区在线| 亚洲欧美成aⅴ人在线观看| 欧美色综合网站| 国产一区二区三区日韩精品| 国产三区二区| 国产美女无遮挡免费视频| 亚洲中文字幕在线观看| 久久五月天国产自| 18禁不卡免费网站| 久青草国产高清在线视频| 国产精品一区二区无码免费看片| 97无码免费人妻超级碰碰碰| 超级碰免费视频91| 欧洲精品视频在线观看| 国产在线视频欧美亚综合| 玖玖精品视频在线观看| 中文字幕精品一区二区三区视频 | 国产成人在线小视频| 影音先锋丝袜制服| 一级毛片在线播放免费| 欧美视频在线播放观看免费福利资源 | 亚洲欧美自拍中文| 无码内射在线| 99久久精品视香蕉蕉| 欧美激情二区三区| 噜噜噜久久| 国产特级毛片|