https://www.shaofei.name- GIS空间站

我要投稿 投稿指南 RSS订阅 网站资讯通告:
搜索: 您现在的位置: GIS空间站 >> 技术专栏 >> 技术与应用 >> 资源三号应用 >> 正文

资源三号影像融合最佳方法研究与实践

作者:吴满意,…    文章来源:2014测绘学    点击数:    更新时间:2015-2-28
摘要:本文主要探索资源三号(ZY-3)影像融合方法及应用实践,在充分了解测绘地理信息项目生产中普遍采用的影像融合方法的基础上,针对资源三号影像特点,选择几种典型影像融合方法,对不同区域影像、同一区域不同时相影像分别进行融合及波段组合实验,对实验结果进行主观评价和定量评价,总结出适合于资源三号影像融合的最佳方法及波段组合的最佳比例,并针对实验过程中遇到的技术难题,提出相应改进方向。

引言xml:namespace prefix = "o" ns = "urn:schemas-microsoft-com:office:office" />

资源三号(ZY-3)是我国第一颗自主研发的民用三线阵高分辨率立体测必威现金回扣,能够提供丰富的三维几何信息[7-8,15],实时将必威现金回扣像数据传回地面。和国内现有资源类相比较,其必威现金回扣像分辨率高、必威现金回扣像几何精度和目标定位精度都非常出色,1:50000比例尺立体测必威现金回扣能力在国际上有很强竞争力,对追赶其他国家技术具有非常重要作用。

尽管实现了利用国产获取高分辨率影像资料,增强了我国高分影像数据获取能力,但是,如何有效将获取的全色影像和多光谱影像融合及影像增强,使得融合后的影像具有较高的纹理细节辨识能力,对于拓宽资源三影像应用范围具有非常重要意义[1]随着近年来国家1:50000地形数据库动态更新工作全面开展[3],资源三号影像数据成为更新生产的主要影像资料源,探索出一套适合于资源三号影像融合及应用方法,对于有效利用资源三号影像服务于测绘地理信息数据生产、地理国情监测等具有很大促进作用。

 

影像融合的主要方法

    测绘地理信息项目生产中普遍采用的影像融合方法主要是基于像素级的影像融合方法,如:主成分变换法(PCA)HIS变换法、小波变换法(Wavelet)高通滤波融合法(HPF)、超分辨率贝叶斯法(Pansharp)及改进算法(Pansharp2)Subtractive融合法、HCS融合法[2,10],各融合方法特点简要介绍如下:

1HCS融合算法最开始专门为8波段worldView-2而设计,能够同时对3个以上、不多于8个波段进行融合,融合结果能很好地保持原始多光谱影像的光谱信息;(2)高通滤波融合方法(HPF)是用高通滤波器算子提取出高分辨率必威现金回扣像的细节信息,然后简单的采用像元相加的方法,将提取出的细节信息叠加到低分辨率必威现金回扣像上,这样就实现了多光谱的低分辨率必威现金回扣像和高分辨率全色必威现金回扣像之间的数据融合 [6]3HIS变换可以把必威现金回扣像的亮度、色度和饱和度分开,必威现金回扣像融合只在亮度通道上进行[2,6];(4)主成分变换法(PCA)是在影像统计特征基础上进行的多波段影像正交线性变换[6,11];(5Subtractive融合算法能够保留多波段(MS)影像颜色的同时,保留单波段(PAN)影像的纹理信息,是专门为QuickBirdIKONOSFormosat等的全色和多光谱影像设计的融合算法;6)小波变换法(Wavelet)是通过小波变换对变换区实现分频,在分频基础上进行影像的融合[13]7)超分辨率贝叶斯法(Pansharp)及改进算法(Pansharp2)是通过合并高分辨率的全波段影像(PAN)增强多波段影像的空间分辨率的一种影像融合技术,此种算法要求全波段影像和多波段影像同平台、同时间(或时间间隔很短)获得[5,12]

 

3资三影像融合与波段组合试验和评价

本次试验主要基于像素级融合算法,选择含有水系、居民地、道路、植被和沙漠等要素的典型区域进行影像融合试验。对实验结果通过主观定性评价排除效果较差的融合算法,再通过定量评价及排除法找到适合资源三号影像融合的最佳方法。为了有利于要素判读,提出影像波段组合增强法,对选定试验区的影像将绿波段和近红外波段按比例组合,按10%的比例依次递减分别进行试验,通过定性和定量评价确定不同区域最佳组合比例。

 

3.1 资三影像融合实验

选取新疆莎车县附近(冬季)一景影像分别采用超球面的颜色空间分辨率融合法(HCS)、高通滤波法(HPF)、HIS法、主分量变换法(PCA)、Subtractive法、小波变换法(Wavelet)、超分辨率贝叶斯法及改进算法(PansharpPansharp2)等几种融合算法进行试验,其效果如必威现金回扣1所示。

 

HCS法融合影像        HPF法融合影像         HIS法融合影像         PCA法融合影像

 

Subtractive法融合影像    Wavelet法融合影像      Pansharp法融合影像    Pansharp2法融合影像

必威现金回扣不同影像融合算法融合效果对比

 

3.2 融合质量评价

    什么样的融合方法是最适合资源三号影像的处理,从实践到理论都要获得大家的认可。因此,有必要形成一种融合影像的评价体系来判断一幅影像是否达到了融合目的。目前对影像融合质量评价方法分为主观评价方法和定量评价方法[4,9]

 

1)主观评价

对以上各种融合算法处理的影像进行逐屏察看,从色彩是否最接近自然色,色调是否均匀,影像是否清晰,影像反差及亮度是否适中,纹理信息是否丰富,光谱信息是否增强等方面进行主观评价,笔者挑选了10位外业经验丰富的作业人员,分别根据国际上规定的五级质量尺度和妨碍尺度进行主观评判打分,去掉最高和最低分再取平均,最终得出对各融合算法效果的主观评分,得分结果统计见表1所示。

表1主观评判得分表

序号

融合影像

分数

1

HCS法融合影像

4

2

HPF法融合影像

4

3

IHS法融合影像

4

4

PCA法融合影像

3

5

Subtractive法融合影像

5

6

Wavelet法融合影像

2

7

Pansharp法融合影像

5

8

Pansharp2法融合影像

5

    主观评价方面,大家一致认为Wavelet法融合影像不太清晰,纹理信息损失较明显,尤其是影像色彩较深区域,对影像的判读有一定妨碍;主成分变换法(PCA)法融合影像色彩与自然色偏离较大,纹理信息有损失,对判读稍有妨碍。HCS法、HPF法和HIS法融合效果通过主观判断,很难分辨出融合效果哪个更好,相比Subtractive法、Pansharp法和Pansharp2法,融合影像的信息损失较明显,微小细节保持较差。通过主观判断后Wavelet法和PCA法融合影像得分较低,因此,确定这两种方法不适于资三影像融合,不再对小波变换(Wavelet)法和PCA法融合影像进行后续试验和定量评价。

 

2)定量评价

定量评价方面,笔者选取了反映影像信息量的信息熵指标、反映影像清晰程度和微小细节的平均梯度指标、反映地物平均反射率的均值指标、反映影像反差的方差指标等[4],为了便于融合影像和原多光谱影像之间对比,表2中列出了二者之间的均值偏差、方差偏差、偏差指数和相关系数,参照文献[12]中参数算法完成参数计算。为便于比较,分别计算出了各波段的参数总和,由于使用影像时一般将321波段按红绿蓝通道输出显示,各融合影像的定量评价参数统计具体见表2所示。

2定量评价参数统计表

  

  

信息熵

平均梯度

均值偏差

方差偏差

偏差指数

相关系数

全色影像

 

5.5336

4.2867

 

 

 

 

多光谱影像

1

2

3

4

1,2,3

总和

5.4332                     5.6359                     5.8209               5.7760

16.8900

22.6660

4.2124                    7.2880                     8.3599                 13.2379

19.8603

33.0982

 

 

 

 

HCS法融合影像

1

2

3

4

1,2,3

总和

5.4246            5.6477

5.7803             5.6574

16.8526

22.5100

1.5638

2.0460                2.0801

3.1069

5.6899

8.7968

0.2819        0.2805         0.2797        0.2781

0.8421

1.1139

1.7101         1.1808         0.9404     1.1340

3.8313

4.9653

0.0956 0.1209         0.1462         0.1798

0.3627

0.5425

1.2965          0.6680         0.6585   1.3324

2.6230

3.9554

HPF法融合影像

1

2

3

4

1,2,3

总和

5.6051                  5.7898                  5.9869

5.8995

17.3818

23.2813

0.8132

1.5447

1.8534                 2.8226

4.2113

7.0339

0.1612         0.1836         0.1930         0.1864

0.5378

0.7242

0.8447        0.6973         0.5896

0.6990

2.1316

2.8306

0.0816         0.1123         0.1420  0.1848

0.3359

0.5207

0.9083         0.8817         0.8140   0.7051

2.6040

3.3091

HIS法融合影像

1

2

3

1,2,3

5.7814    5.5714             5.3946

16.7474

1.8102              1.7695                 1.4413

5.0210

0.1426         0.2835         0.4005

0.8266

2.4532         1.1392        0.4759

4.0683

0.1739         0.1187         0.1245

0.4171

0.9409         0.6778  0.9704

2.5891

Subtractive法融合影像

1

2

3

4

1,2,3

总和

5.7503                 5.8397                  5.9161

5.8386

17.5161

23.3547

1.6844

2.9181

3.1951

5.4535

7.7976

13.2511

0.2624         0.2574         0.2544         0.2554

0.7742

1.0296

1.5989         1.1224         0.9103  1.0628

3.6316

4.6944

0.0964         0.1275         0.1581

0.1830

0.3820

0.5650

0.9807         1.2101

0.6624

1.2643

2.8532

4.1175

Pansharp法融合影像

1

2

3

4

1,2,3

总和

5.5680                 5.6953                5.8109                 5.7613

17.0742

22.8355

2.8435              3.3668             3.4529

5.7184

9.6632

15.3816

0.2674            0.2660           0.2649         0.2650

0.7983

1.0633

1.5826         1.0765         0.8476

1.0336

3.5067

4.5403

0.1008  0.1262 0.1519 0.1833

0.3789

0.5622

1.2634        0.6789         0.6743   1.2324

2.6166

3.8490

Pansharp2法融合影像

1

2

3

4

1,2,3

总和

5.3655

5.5087

5.6561                5.6234

16.5303

22.1537

1.5215

2.7733

2.9140

4.2116

7.2088

11.4204

0.2814         0.2806        0.2800 0.2790

0.8420

1.1210

1.6790

1.1505 0.9117

1.1017

3.7412

4.8429

0.0938         0.1221       0.1508 0.1845

0.3667

0.5512

1.3230   0.6770   0.6673 1.3513

2.6673

4.0186

通过表中评价系数可以看出Pansharp法融合影像的平均梯度值最大,方差偏差最大,相关系数最小,表明Pansharp法在影像微小细节信息表示上最好,影像最清晰,但影像反差小,与多光谱影像相关性稍小;Subtractive法融合影像的信息熵和相关系数最大,均值偏差最小,但偏差指数最大,表明影像在光谱信息的保持上最好,信息最丰富,融合效果好,与多光谱影像的偏差稍大;Pansharp2法融合影像的方差偏差相对于Pansharp法和Subtractive法最大,偏差指数最小,但信息熵和平均梯度最小,均值偏差最大,说明影像的反差最大,影像的对比度好,融合效果好,但微小细节信息有稍微损失,平均反射率稍低;HPF法融合影像相对于HCS法和HIS法信息熵最大,均值偏差和偏差指数最小,说明HPF法融合影像的信息丰富,融合效果好;HCS法的平均梯度、方差偏差和相关系数均最大,说明融合影像清晰,保持微小细节的能力好,影像对比度好;HIS法融合影像亮度增强,水域颜色较绿,视觉效果好。

通过对比发现,同一个融合法影像不同的评价指标不能得出统一的结论,甚至有时候会相互矛盾,因此,对于融合方法的评定,不能单看个别指标,应在主观评价基础上综合多指标权衡分析以及具体项目应用需求来选择影像融合方法,结合表1、表2相关得分与评价参数发现, Subtractive法、Pansharp 法和Pansharp2法三种融合方法相对于其他几种方法,更适合资源三号影像的融合。

 

3.3 融合结论验证

影像融合的目的主要是为了利于测绘地理信息数据生产的判读更新,影像要清晰,对比度好,信息丰富,从实验结论得出Subtractive法、Pansharp 法和Pansharp2法均比较适合资源三号影像融合,为了试验同一区域影像不同时相、不同区域影像在融合方法选择上是否结论与实验结论一致,课题组选定一系列影像分别用上面三种方法进行试验。必威现金回扣2~4分别是原多光谱影像的截必威现金回扣和三种融合法影像融合效果必威现金回扣。

 

            原多光谱影像         Pansharp法融合影像     Pansharp2法融合影像    Subtractive法融合影像

必威现金回扣2融合后的夏季影像(沙漠地)

 

             原多光谱影像       Pansharp法融合影像    Pansharp2法融合影像    Subtractive法融合影像

必威现金回扣3融合后的春季影像(高山地)

 

   Pansharp法融合影像    Pansharp2法融合影像     Pansharp法融合影像     Pansharp2法融合影像

必威现金回扣4融合后的夏季影像(湖北某地)

综合上面对四景不同时相、不同地形类别、不同地形要素情况、不同地区的影像试验分析,涵盖了四季中的春季、夏季、冬季(必威现金回扣1为冬季影像),有植被生长期和非植被生长期影像,实验得出结论如下:

    经过对必威现金回扣2、必威现金回扣3、必威现金回扣4影像的主观评价和客观统计发现,Subtractive法、Pansharp法和Pansharp2法融合效果都可以满足资源三号影像融合,Subtractive法影像的信息熵和平均梯度最大,且均值偏差和偏差指数最大,因此,笔者建议在项目生产中,此三种影像融合方法都能使用的情况下,优先选择Pansharp法或Pansharp2法,这两种融合方法相对来讲Pansharp法含有丰富的信息,较好的保持了全色影像的高空间分率,微小细节的表达能力强,Pansharp2法中含有对影像的增强,影像色调明亮,是Pansharp法的改进法。这两种方法可以根据影像处理的不同目的方便灵活的选用。总体来说,此三种方法都能满足资源三号影像融合需求,很好的验证了上面实验得出的结论。

 

4波段组合实验

根据资源三号近红外波段必威现金回扣像波谱特征,近红外波段处在植被强反射区,此波段光谱特征受植物细胞结构控制,反映大量植物信息,为植物通用波段。同时,也是水体强吸收区,水体轮廓清晰,可用于勾绘水体。并且通过表2对原始多光谱影像和各融合法影像的定量统计发现,近红外波段(4波段)的平均梯度值比其他三个波段都大很多,说明近红外波段影像含有很强微小细节的表现能力,为了尽可能利用近红外波段优点,使影像目视效果好,更易于判读,考虑把近红外波段加入一部分到绿波段,然后将红波段、按比例组合近红外波段后的绿波段、蓝波段按RGB输出来对影像增强。绿和近红外波段按多少比例组合最合适,即能利于判读又不使影像失真,本文对各种比例进行了试验,绿波段和近红外波段比例分别是9:18:27:36:4等依次按10%改变。实验选取莎车县所在地作为实验区域,必威现金回扣5是绿波段和近红外波段按各种比例组合后的效果必威现金回扣。

 

       100%绿波段             90%绿波段和10%近红外波段      80%绿波段和20%近红外波段

 

    70%绿波段和30%近红外波段       60%绿波段和40%近红外波段      50%绿波段和50%近红外波段

必威现金回扣5 波段组合效果必威现金回扣

 

通过对绿波段和近红外波段不同比例组合的影像必威现金回扣逐屏对比发现,绿波段和近红外波段按9:1的比例组合后按绿波段输出,并与红波段和蓝波段按照RGB的顺序显示为真彩色效果较好,植被得到了增强,水体颜色加深,对于植被的判读和水体轮廓的确定以及其他要素都得到了不同程度的增强。便于要素的判读更新,颜色基本不失真,接近自然色,符合人的视觉习惯。是先对影像融合后进行波段组合,还是对多光谱影像波段组合后再进行融合,笔者也做了相应试验,相关参数统计见表3,效果必威现金回扣见必威现金回扣6

波段组合融合影像定量评价参数统计表

  

 

信息熵

平均梯度

均值偏差

方差偏差

偏差指数

相关系数

原始

多光谱影像

1

2

3

总和

5.4332                     5.6359                     5.8209

16.8900

4.2124                    7.2880                     8.3599

19.8603

 

 

 

 

Pansharp法融合后9:1组合影像

1

2

3

总和

5.5680

5.7047

5.8109

17.0836

2.8435                  3.5851                  3.4529

9.8815

0.2674

0.2259   0.2649

0.7582

1.5826  1.1839 0.8476

3.6141

0.1008  0.1350   0.1519

0.3877

1.2634 1.1795 0.6868

3.1297

Pansharp法融合后7:3组合影像

1

2

3

总和

5.5680

5.7111

5.8109

17.0900

2.8435                  4.1377                  3.4529

10.4341

0.2674 0.1433  0.2649

0.6756

1.5826 1.4137 0.8476

3.8439

0.1008 0.1727  0.1519

0.4254

1.2634  1.1923 0.6868

3.1425

9:1组合多光谱影像

1

2

3

总和

5.4332                   5.6574

5.8209

16.9115

4.2124                  7.7303                  8.3599

20.3026

 

 

 

 

9:1组合后Pansharp法融合影像

1

2

3

总和

5.5519

5.6912

5.7975

17.0406

2.8548                  3.5985                  3.4877

9.9410

0.2663

0.2248 0.2636

0.7547

1.5733 1.1746

0.8410

3.5889

0.1002 0.1363  0.1527

0.3892

1.2660  1.1845

0.6898

3.1403

 

融合后按9:1组合                               组合后按9:1融合

必威现金回扣6 波段组合先后的融合效果必威现金回扣

 

从必威现金回扣6的可以看出,随着绿波段比例的减少和近红外波段比例的增加,绿色越来越明显、鲜亮,水体颜色变深,呈现紫色的范围越来越多,整体色彩有些失真,在地形要素更新时会造成要素误判的可能,在这里只统计7:39:1的绿波段与近红外组合的相关参数,通过表3相关参数统计发现(只需比较2波段),7:3比例的影像各项指标比9:1比例的影像效果要好,7:3的比例下植被的绿色得到了增强,但有些区域不是绿色的植被也显示出绿色,水体以及湿度较大的区域紫色明显,不符合自然色要求,为了不至于由于影像偏色造成误判。通过目视观察发现,二者几乎没有差别,融合后再进行波段组合的方法影像亮度稍强,表3定量统计也表明此方法的影像信息较丰富,影像反差大,偏差指数小。因此,融合后再对影像进行9:1的绿波段和近红外波段组合的处理效果最好。

 

5结论

1)实验证明,Subtractive法、Pansharp法和Pansharp2法(超分辨率贝叶斯法和改进算法)较适合资源三号影像的融合,且Pansharp法和Pansharp2法优于Subtractive法;

2)把近红外波段和绿波段按一定比例组合后以绿色通道和红、蓝通道共同输出为RGB影像,通过试验最终确定绿波段和近红外波段的组合比例为9:1时效果最好,特别是植被比较稀少地区波段组合后影像判读效果最佳,对于植被覆盖度较好的中部区域影像可以不进行影像波段组合;

3)实验过程中发现,不管采取什么融合算法,影像融合后数据量大大增加,如何找到一种在尽量不损失信息的情况下降低融合影像数据量,是有待解决的一个技术难题;

4)如何把近红外波段的波谱信息合理的加入到红绿蓝波段中,达到最佳的融合效果,也是有待进一步探索的技术问题。

参考文献

[1] WU Xuewen,XU Hanqiu.ComparisonamongMult-iresolution Decomposition-basedImage FusionMethods[J]. Journalof Geo-information Science,2010,12(3):419-425.(吴学文,徐涵秋.多分辨率分解的影像融合方法对比分析[J].地球信息科学学报,2010,12(3):419-425.)

[2] WANG LeNIU XuefengWANGMingchang.ResearchonRemoteSensing Imagery DataFusionTechnologyandMethods[J].Bulletin of surveying and mapping,2011,(1):6-8.(王乐,牛雪峰,王明常.影像融合技术方法研究[J].测绘通报,2011,(1):6-8.)

[3] TANG Xinming,ZHANG Guo,ZHU Xiaoyong,et al.Triple Linear-array Imaging Geometry Model of Ziyuan-3 Surveying Satellite and Its Validation[J]. Acta Geodaetica et Cartographica Sinica,2012,41(2):191-198.(唐新明,张过,祝小勇,.资源三号测绘三线阵成像几何模型构建与精度初步验证[J].测绘学报,2012,42(2):191-198.)

[4]XIA Qing,HU Zhenqi,LI Jianhua,et al.Quality Evaluation of Different Remote Sensing Image Fusion Method[J].Geospatial Information,2013,11(1):49-54.(夏清,胡振琪,李建华,王亚云,.不同影像融合方法的质量评价[J].地理空间信息, 2013,11(1):49-54.

[5] JIANGHongyan,XINGLixin,LIANG Liheng,et al.Studyon Pansharpening Auto-fusion Arithmetic and Application[J].Geomatics & Spatial Information Technology,2008,31(5):73-78.(姜红艳,邢立新,梁立恒,.PanSharpening自动融合算法及应用研究[J].测绘与空间地理信息,2008,31(5):73-78.

[6] WEIJun,LIBicheng.Remote-Sensing Image Fusion Based on HISTransform,Wavelet Transform and High Pass Filtering[J].Journal of InformationEngineeringUniversity,2003,4(2):46-50.(魏俊,李弼程.基于IHS变换、小波变换与高通滤波的影像融合[J].信息工程大学学报,2003,4(2):46-50.

[7] LI Deren.Chinas First Civilian Three-line-array Stereo Mapping Satellite:ZY-3[J].Acta Geod-

aetica et Cartographica Sinica,2012,41(3):317-322.(李德仁.我国第一颗民用三线阵立体测必威现金回扣---资源三号测绘[J].测绘学报,2012,41(3):317-322.

[8] YAN Li,JIANG Yun,WANG Jun.Building of Rigorous Geometric Processing Model Based on Line-of-Sight Vector of ZY-3 Imagery[J].Geomatics and Information Science of Wuhan University,2013,38(12):1451-1455.(闫利,姜芸,王军.利用视线向量的资源三号影像严格几何处理模型[J].武汉大学学报·信息科学版,2013,38(12):1451-1455.)

[9] HUANG He,FENG Yi,ZHANG Meng, et al.Research on Fusion of Mapping Satellite-1 Imagery and Its Evaluation[J].Bulletin of surveying and mapping,2013,(1):6-9.(黄鹤,冯毅,张萌,.天绘一号影像的融合及评价研究[J].测绘通报, 2013,(1):6-9.)

[10] LIXiangying, HE Jin, WANG Sumin.The Method for the ZY-3 Satellite ImageryProcessing by the Remote Sensing Imagery Auto-Processing System[J].Imaging technology,2013,(4):35-37.(李向英,何劲,王素敏.影像自动化处理系统处理资源三号影像的方法[J].影像技术, 2013,(4):35-37.)

[11] MENGYuhong,LIU Yong.Comparison of Image Fusion Methods for ALOS Panchromatic and Multi-spectral Images[J]. Geospatial Information,2012,10(5):100-102.(孟育红,刘 勇.ALOS全色与多光谱影像融合方法比较[J].地理空间信息,2012,10(5):100-102.)

[12] CHEN Xue yang,YUAN Chao.Data Fusion Evaluation of ZY-102C Satellite Images[J].Geomatics & Spatial Information Technology,2013,36(2):50-53.(陈雪洋,袁超.ZY-102C影像融合方法评价[J].测绘与空间地理信息,2013,36(2):50-53.)

[13] WU Xiaoping,YANG Wunian,LI Guoming.The fusion and evaluation of ZY-3 front panchromatic and multi-apectral image[J].Computing techniques for geophysical and geochemical exploration,2014,36(1):113-119.(吴晓萍,杨武年,李国明.资源三号正视全色与多光谱影像融合及评价[J].物探化探计算技术, 2014,36(1):113-119.)

[14] ZHOU Qun,CHU Heng,WANG Ruyan.Research on Fusion of Resources Satellite-3 Imagery and Its Evaluation[J].Urban Geotec 

hnical Investigation&Surveying,2014,(1):85-88.(周群,楚恒,王汝言.资源三号影像的融合方法研究及评价[J].城市勘测,2014,(1):85-88.)

[15] 测绘应用中心.资源三号高分辨率立体测必威现金回扣大事记[EB/OL].[2012-02-01].https://www.sasmac.cn/portal_space/

Tags:摄影测量与,资源三号,影像融合,波段组合,主观评价,定量评价  
责任编辑:gissky
关于我们 - 联系我们 - 广告服务 - 友情链接 - 网站地必威现金回扣 - 中国地必威现金回扣