(原创文章,转载请注明来源,谢谢!)

前言

(前方多图高能预警,万字长文,希望你能耐心看完!)

首先告诉大家一个好消息。Eric Pottier教授的PolSARpro已更新到V6.0版本,各位可以到网址( https://www.ietr.fr/polsarpro-bio/ )自行下载。但是该版本仍在开发中,经过本人的测试,存在不少bug(例如没有帮助文档),各位根据实际情况使用。该版本较V5版本的重大改进是在生物量估算模块方面,该版本的简单介绍文档见百度云盘链接:链接:https://pan.baidu.com/s/1FJkhxkhkwtiEoICDRD19-A
提取码:h9b1 。感谢Eric Pottier教授,也希望该版本在后面能逐步完善。

距离上次PolSARpro处理Sentinel 1A的文章已近三月,期间许多人问我相关的问题,其中问得较多的是感兴趣区域的提取(或者说是掩膜提取),本文将对该操作进行介绍,并介绍PolSARpro中的Wishart监督和SVM监督分类过程(上次因为时间原因,没有详细介绍监督分类步骤,这次补上),希望可以给PolSAR相关研究者带来帮助。根据我的经验,SVM监督分类是PolSARpro中监督分类效果比较好的一种,也是目前遥感图像分类中常用成熟、稳定的算法。这个分类PolSARpro V5.1.3版本也能做,但是两个版本都有点bug,都不完美。

两个版本的PolSARpro可以共存,像博主两个版本都在使用,推荐保留两个版本,因为两个版本都有一些bug,协同使用可以对比一下。

后面,博主打算介绍欧空局SNAP处理软件以及SNAP1
的Python模块------snappy处理,以促进该开源处理遥感处理平台的推广使用及开发。

博主创建了一个关于欧空局SNAP处理交流的QQ群:欧空局SNAP处理交流群:665903216(这个群已满人),欧空SNAP处理交流群(二):1102493346。欢迎感兴趣者加入,共同交流,共同进步!

数据源

本次使用的数据与我本人的第一篇博客数据是相同的,链接网址:https://blog.csdn.net/lidahuilidahui/article/details/86485372
你可以选用你自己的数据。

数据预处理

利用桌面的PolSARproV6.0快捷方式启动PolSARproV6.0:
图1
图2
点击Enter,再点击OK,进入模块选择界面:
图3
点击Sentinel-1A/B mission(哨兵1A/1B卫星)模块(PolSARpro Bio模块也是可以的)点击Enter,即可进入以往的PolSARpro V5.1.3类似的界面(只是界面变成了灰色色调了而已,看到这个界面,竟想起许嵩的歌曲《灰色头像》):

图4
图5

数据读入、多视与地理编码

请参照我的第一篇博客(https://blog.csdn.net/lidahuilidahui/article/details/86485372 )完成数据读入、多视以及地理编码(含地理编码)操作,经过我的测试,如果感兴趣区(Region of Interest, ROI)提取操作早于地理编码前,后续将无法做地理编码,因为元数据文件metadata.xml文件(该文件含有控制点参数)会丢失。因此,感兴趣区的提取操作应在地理编码之后进行。这里选择的感兴趣区为上海市崇明区长兴岛和横沙岛,Swath带为IW2。在这个过程请保持耐心。地理编码成功后的界面如下(注意观察博主红色框标注部分):

图6

叠加到Google Earth

上次这个操作没有演示,这次展示一下,以后你可以用于PPT展示(It’s your show time! Wonderful!)。
关闭GIMP弹出的两个图片,查看生成的文件夹目录:
图7
其中,config_acquisition.txt和GEARTH_POLY.kml是叠加至Google Earth所必需的文件,丢失了这两个文件将无法叠加至Google Earth。提取感兴趣区操作这两个文件会丢失,因此,博主在感兴趣区的提取步骤前展示这个叠加效果。

进入C2文件夹,看到的文件应该是这样的:
图8
(如果你的界面中Display等菜单图标没有点亮,请退出后重新设置目录环境后进入,这时Display等菜单图标会点亮。估计是有点bug的缘故。)
选择Display—>Create KML File:
图9
input BMP File框中,选择:
图10
input GEARTH_POLY File框中,选择:
图11
(Reduction Factor:衰减因子,Transparency:透明度,保持默认设置即可)

点击Run,等待一下,见证奇迹的时刻到了:
图12
可以看到PauliRGB.bmp叠加至Google Earth上面了,并且吻合得很好(这表明两者的坐标系是完全一致的)。惊不惊喜,意不意外!

感兴趣区的提取

关闭Google Earth,退出Dispaly工具栏,进入下一操作。
选择Tools—>Data Set Management—>ROI Extraction(注意该操作仅能选一个感兴趣区):
图13
进入ROI Extraction 界面:
图14
点击Graphic Editor,点击Yes,进入感兴趣区编辑界面:
图15
选择PauliRGB.bmp作为勾选感兴趣区的图像,点击打开,进入如下界面:
图16
在弹出的图像窗口中点击鼠标右键,显示编辑栏菜单,其中:
Select area,选择感兴趣区(可以是多边形);
Delete area,删除感兴趣区;
Save configuration,将所勾选的感兴趣区保存为.txt文件;

该编辑工具快捷键(注意该工具与以往的图像处理软件快捷键不一样,例如ENVI,ArcGIS等)如下:
平移:按住鼠标左键拖动
放大:鼠标滚轮往后滚动
缩小:鼠标滚轮往前滚动

点击Select area,选择感兴趣区:
移动鼠标在合适的位置,点击鼠标左键放下一个多边形端点(锚点),移动鼠标,到另一个端点处,点击放下一个多边形端点,反复操作几次,直至勾选出整个感兴趣区(如果你多边形端点勾选错了,并不能仅仅修改或删除该个端点,仅能通过删除整个多边形,重新勾选多边形来实现修改),如图:
图17
那么如何结束感兴趣的绘制呢?
点击鼠标右键,选择Save configuration,结束绘制感兴趣区并保存下来。选择如果选择Select area和Delete area,则可以重新绘制感兴趣区。再强调一次,该操作仅能选择一个感兴趣区

点击Run,点击Yes,创建感兴趣区文件:
图19
成功时,会弹出一个掩膜文件和生成一个IW2_SNAP_ROI文件夹:
图20
进入IW2_SNAP_ROI文件夹,可以看到仅有一个C2文件夹,上一个步骤的仍能看到的config_acquisition.txt和GEARTH_POLY.kml两个文件均已丢失,已不能再叠加至Google Earth 上。进入C2文件夹中可以看到如下文件:
图21

极化滤波

将目录环境设置为IW2_SNAP_ROI:
图22点击Save&Exit,保存目录环境设置。

极化滤波操作,请参考我的第一篇博客(https://blog.csdn.net/lidahuilidahui/article/details/86485372 ),这里不再复述。
滤波界面:
图23

特征提取

实际上H-A-α分解就是提取极化散射特征的过程,分解得到的参数可以看作是一个分类特征。

双极化SAR数据的H-A-α分解操作,请参考我的第一篇博客(https://blog.csdn.net/lidahuilidahui/article/details/86485372 ),这里不再复述。

双极化SAR数据H-A-α分解参数界面:
图24
博主这里选择的是H(Entropy),A(Anisotropy),α(Alpha)三个特征(极化参数),你可以选择其它的特征,个数也可以多余3个或者少于3个,多少是你的自由,但是分类精度和运行速度会受到影响。
这两个步骤之后,进入IW2_SNAP_ROI_LEE的子目录C2文件夹,可以看到如下文件:
图25
可以看到C2文件夹中没有RGB文件或者PauliRGB文件,下面利用C2矩阵文件创建一个RGB假彩色图像用于选择训练样本。

图26
其中矩形框为生成的RGB文件名,我们将用这个RGB1.bmp(如下图所示)来选择训练样本。
图26
关闭生成的RGB1.bmp文件,点Exit退出Display工具栏,进入下一个步骤:监督分类。

监督分类

好了,快结束了。还能保持耐心吗?接下来会介绍Wishart、SVM监督分类,这里不会详细介绍其原理,相关的文献有许多,具体请自己查找相关文献(因为双极化SAR监督分类目前仅支持这两种方法的监督分类,不信的话,你可以测试下)

Wishart监督分类

这是早期到现在极化SAR数据常用的监督分类,该分类是基于复Wishart分布的极大似然分类,但是该分布仅在均匀区域中的建模中效果比较好,也就是说,Wishart监督分类一般只能在均匀区域地物分类中取得较好的效果,对于其它不均匀区域地物分类,往往分类效果不好。对于不均匀区域的多极化通道SAR图像的建模,至今仍是一个热点和难点所在,近几年,提出的K、Go分布、W分布都有其局限性。对其进行精确的建模,对于SAR滤波、目标识别、分类都具有非常重要的意义,尤其是在军事上的目标识别(例如舰艇、坦克等)。

Wishart分类器

需要指出的是:该分类器针对的是C2矩阵元素(即C11,C12,C21,C22,也就是这四个元素是可以看作服从复Wishart分布的)进行的。
开始操作吧:
选择Process—>Polarmetric Segmentation—>Wishart Supervised Classification。
图27
参数设置如下:
图28
其中,Window Size大小为boxcar filter器窗口的大小,实际Wishart等监督分类时往往会加一个滤波操作,如果你不想再滤波,只需设置窗口大小为1×1;

Color Map16只是定义样本标注的颜色,并不是说要选16个类别。

你也可以勾选Reject Class(拒绝类),Reject Radio:5(默认值),实际上是一个距离度量。我也不太懂,具体你可以看下帮助文档。一般都不勾选,你也可以勾选和不勾选了对比一下结果。

创建训练样本

选择创建训练样本的图像,这里使用的是上一步创建的RGB1.bmp图像。图29
打开后,如下所示:
图30
实际,在创建训练样本时,应该参考光学图像或者到地面实际调查测量。如果你想更精确地选择训练样本,你可以打开Google Earth 中查看相关信息。但是Google Earth 的影像不一定是最新的,最可靠的还是地面实际调查测量。

训练样本的创建:
在窗口中点击鼠标右键,会弹出一个菜单栏,其操作如下:
Add a new class:添加一个新类
Select area:选择一个样本区域
Delete area:删除一个样本区域
Save configuration:保存建立的样本区域,存储为.txt文件(勾选好全部样本点击这个选项,注意这个按键不是结束一个训练样本的绘制,而是直接结束了整个训练样本训练的程序)
图31
该窗口工具与ROI Extract工具快捷键有点类似,但是有所不同,刚开始操作往往不习惯,需要慢慢适应(如下:
平移:按住鼠标左键拖动
放大:鼠标滚轮往后滚动
缩小:鼠标滚轮往前滚动

在窗口中鼠标右键,创建类别点击Add a new class(每次点击,Class Num都会加1,但是类别名(数字序号)无法编辑(例如改为字符串的形式“water”),因此,操作时你需要清楚你的选择的类别是什么):
图32
这个时候,你会发现,点击鼠标坐标,仍然选择不了样本。继续,鼠标右键,点击Select Area,这时,你会发现终于可以选择样本了。
图33

移动鼠标到某个位置,放大(缩小)图像,鼠标左键放下一个多边形端点(锚点),再移动鼠标到某个位置,鼠标左键放下一个多边形端点,反复操作若干次,直至勾选出该个样本区域。(如果你多边形端点勾选错了,并不能仅仅修改或删除该个端点,仅能通过删除整个多边形,重新勾选多边形来实现修改)绘制多边形边界示意图如下:
图34
那么如何结束绘制一个样本的操作本?答案是:鼠标右键,选择Select Area
图35
早期的版本中value和样本数一致的,很遗憾,这几个版本都没。

然后你会发现,多边形绘制的确结束了,然而却多了一个烦恼,有一个端点在旁边,导致你无法移动到别处选择样本区域。怎么办?鼠标右键,选择Delete Area。OK,你会发现该端点消失
图36
博主这里选了4个类别:建筑、水体、裸土(含耕地)、植被(森林、草地、农作物等),每个类别勾选4-5个。(仅用于实验目的,就不要计较精度了)。
移动鼠标至合适位置,继续勾选样本区域(鼠标右键,点击Select Area)。(重复几次,取决于你这个类别的训练样本数量)直至勾选完第一个类别的所有训练样本。

鼠标右键,新建一个类别(Add a new class),创建新建类别的训练样本区域(Select area),勾选训练样本区域边界端点。反复勾选几次。(上述操作重复几次)

直到勾选完4个类别的全部训练样本。(红色:建筑;绿色:水体;紫色:裸土;黄色:植被)
图37
窗口中右键,选择Save configuration,保存样本区域数据。会发现Run Training Process 图标被点亮。
图38

训练样本

点击Run Training Process,训练样本,完成后如下图所示:
图39会生成类别中心集(矩形框所示)

分类

点击Run,执行分类。
图40
会弹出3个文件,其中,wishart_training_cluster_set.bmp为训练样本分布图;wishart_supervised_classs_3x3.bmp为分类结果图(其中3x3为滤波窗口大小),wishart_supervised_classs_3x3_RGB1,是对原来的选择的RBG1.bmp进行分类的结果,用处不大。我们需要的wishart_supervised_classs_3x3.bmp。

精度评价

我们把这三个图像关掉,点击CM Editor,查看混淆矩阵。
图41
可以看到混淆矩阵如下:
图41
可以看到分类的精度还是比较高(主要是类别和训练样本太少,不过建筑和水体的分类精度和光学影像进行对比,精度还是可以的)。
点击文本框的Exit退出混淆矩阵窗口。
点击Exit退出Wishart监督分类窗口。

查看分类结果

回到文件夹IW2_SNAP_ROI_LEE/C2中查看分类结果文件:
图42
注意博主标记的三个bmp文件,文件1是样本区域的分类图(相当从从分类结果中对样本区域进行ROI(感兴趣区)提取);文件2为对C2矩阵分类结果。文件3:创建的样本区域的分布图;

对这3幅图像进行联动操作,选择Display—>Display N BMP Files—>Open:

图43
依次打开文件1,2,3(Open点击3次),结果如下图所示:
图44
按住鼠标左键,移动至合适位置查看细节变化,鼠标滚轮控制放大缩小,可以看到三个图像有三个“十字丝”联动:
图45
分类结果查看操作结束,关闭图像,点Exit退出,点Exit退出。(分类结果制图等还是需要ArcGIS、QGIS等来制作)。

SVM监督分类

很遗憾,PolSARpro V6.0的SVM监督分类有缺陷,具体你可以试下。博主这个SVM监督分类,使用PolSARpro V5.1.3版本来做的。

设置目录环境:
图46
在分类之前,我们再做一个操作,将C2矩阵由线性单位转为分贝值,同时提取总功率Span分贝值用作分类特征(操作见下图):
图47
点击Run运行即可。

其中,Modulus表示取模(即复数的振幅,线性单位);10log(Modulus)表示将模转化为分贝值;
Span表示总功率,其中 S p a n = ∣ S 1 ∣ 2 + ∣ S 2 ∣ 2 Span=|S_1|^2+|S_2|^2 Span=S12+S22,对于C2矩阵 S p a n = ∣ C 11 ∣ + ∣ C 22 ∣ Span=|C_{11}|+|C_{22}| Span=C11+C22

成功后,会弹出4个bmp图像窗口,分别是C11_db,C12_db, C22_db, span_db(db表示分贝值得意思),如图所示:
图48
完成后关闭所有图像窗口,点击Exit,退出矩阵元素提取窗口。

SVM分类器

SVM(Support Vector Machine,支持向量机)算法具有小样本训练、支持高维特征空间的特点,可以很好地避免“维数灾难”和过学习问题。该算法成熟稳定,效果通常较好,目前在语音识别、自然语言处理、计算机视觉、遥感图像分类等多个领域仍在大量使用。

开始操作吧:
选择Process—>Polarmetric Segmentation—>SVM Supervised Classification。
图49
进入页面:
图50
可以看到图中GUI界面中有六个步骤,Step 1-Step 6,接下来完成这6步就可以。

创建训练样本

Step 1-Training Area(操作见下图):
图49
博客这里选择是的上一步Wishart监督分类创建的感兴趣区域文件,你也可以选择Graphic Editor按钮打开RGB图像重新选择训练样本区域。为了和Wishart监督分类对比,博主使用了同样的样本区域。

训练样本

Step 2-Training Area:
保持默认即可。

Step 3-Color Maps(操作见下图):
图51

SVM分类特征

Step 4-SVM Parameter Setting(操作见下图):
图52
如果选择C2,仅有将C2的元素作为分类特征。楼主选择的是other,分类特征选择的是:

特征变量alphaanisotropyC11_dbC12_dbC22_dbentropyspan_db
物理量平均散射角α反熵A C 11 C_{11} C11的模(分贝值) C 12 C_{12} C12的模(分贝值) C 22 C_{22} C22的模(分贝值)极化熵H总功率Span(分贝值)

操作5,也可以不选,仅是一个类别概率图像。其它的参数可以自行选择。

SVM参数优化

Step 5-Kernel Parameter(操作如下图):
图53
该步骤设置的是SVM核函数的类型(见矩形框所示):
RBF(Radial Basis Function, 径向基核函数,最常用的核函数),Polynomial(多项式核函数,其下的Degree设置的是多项式函数的阶数),Linear(线性核函数,即一次函数);(如果你对SVM算法很熟悉的话,应该一眼就能看明白)

操作1:RECOMMANDED(推荐参数)
(除非你SVM调参很厉害,否则建议你勾选这个选项)

操作2:点击Setup and Run,会弹出右边的小窗口;

操作3:点击Run RBF Kernel Paramete Optimisation:执行径向基核函数参数优化。这里使用的基于网格法对训练样本进行交叉验证(Cross Validation,CV)从而确定最优参数,主要确定参数C和G,其实C是前面的Cost(惩罚系数c),G是Gamma(核函数宽度g)。(其具体原理请查看相关文献)
执行后会弹出一个命令行窗口。
图54矩形框可以查看进度(百分比,%)。
完成后,右边小窗口可以看到找到的参数C和G的一个组合值。
图55
点击Exit and Save CV Parameters,将会看到Cost和Gamma值发生了变化,正是寻优的结果值256,0.50000。
图56

SVM分类

Step 6-Run Classification(执行分类):
执行后会弹出一个命令行窗口(接下来又是一个漫长的等待过程,请保持耐心),如下图所示:
图57
矩形框显示的是SVM分类的进度(百分比,%)

事实上,PolSARpro的SVM分类算法很好,比ENVI等软件的SVM分类快得多,并且可以自动寻找参数优化值。(可以说是甩ENVI等几条街,也不过分!)

完成后会弹出一个bmp图像(但这不是分类结果,这是之前勾选的概率图像)。
图58
可以把弹出的bmp图像关闭。

精度评价

很遗憾的是,这个版本的SVM分类仍然有点bug,看不到混淆矩阵。(据说早期版本是可以的,另外,Linux系统版本的bug会少很多,有条件的话可以装Linux系统版本的PolSARpro)。如果你点击 CM Editor,想查看混淆矩阵,很遗憾,找不到SVM的confusion_matrix.txt文件。所以,直接点Exit,退出SVM分类窗口就可以。

SVM分类结果

我们回到目录IW2_SNAP_ROI_LEE/C2,如图:
分类结果文件名为svm_classification_file:
图59
请参考wishart监督分类的分类结果查看操作:
图60

SVM与Wishart监督分类结果比较

我们依次打开下图所示的三个文件(文件1为SVM监督分类结果,文件2是SVM训练样本集,文件3是Wishart监督分类结果):
图61
对比结果如下:
图62
可以看到SVM分类更细碎些(斑点较多),总体上看两者的分类精度应该差距不大。(对比完可以退出所有的窗口了!)

后语

好了,终于结束了!(你看完此文和操作,花了很长时间,同样地,我写完此文也花了很长时间,但愿你能有所收获!)

请允许我牢骚多几句!关于PolSARpro软件的功能还有许多,很难详尽地介绍清楚,要想更好地使用和开发PolSARpro,需要对极化SAR的理论知识有较深的了解(学习还是很重要的!)。很遗憾这里使用的是欧空局的Sentinel-1双极化SAR数据(使用这个数据的主要原因是这个数据可以免费下载),不是全极化SAR数据,大量的实验表明,全极化SAR数据的分类精度高于双极化SAR数据、单极化数据(当然,博客的方法对于全极化SAR数据也是使用)。

事实上,全极化SAR数据可以做更多的操作,无论是极化分解、分类、极化干涉测量等都有更多应用。其实,博主更希望的使用国产数据来做些实验的,然而,由于种种原因(数据不开放,国产数据质量、兼容性等),国产卫星数据的利用率远远低于国外卫星数据的利用率,不禁使人感到担忧,没有使用,很难有反馈,也就很难有改进。

目前而言,国产卫星数据相比国外卫星数据(欧空局、美国、日本等),仍有很大的差距!不管怎样,还是希望各位有条件的话,还是多多使用国产卫星数据。无论好坏,用自己国家的数据,多少更放心些。博主后面也打算也打算分享一下国产卫星数据的应用,各位如果有时间的话,也可以写点国产卫星数据的应用的文章,促进国产卫星数据的推广使用。

如果博客内容有错的话,请在评论区批评指正!谢谢!

再打一次广告,博主创建的欧空局SNAP处理交流群:665903216(这个群已满人),欧空SNAP处理交流群(二):1102493346。欢迎感兴趣者加入,共同交流与进步!

参考文献

[1] (美)李仲森(Lee,J.),(法)鲍狄埃(Pottier,E.),洪文等. 极化雷达成像基础与应用[M]. 北京:电子工业出版社, 2013.06
[2] 尤淑撑,行敏锋,何彬彬等. 极化SAR土地资源调查监测技术与应用[M]. 北京:科学出版社,2018.07
[3] 刘涛,崔浩贵,谢凯等. 极化合成孔径雷达图像解译技术[M]. 北京:国防工业出版社,2017.12
[4] 马建文等.遥感数据智能处理方法与程序设计 (第二版)[M]. 北京:科学出版社,2010.01
[5]微信公众号:遥感事. 推送文章:Dragon系列-PolSARpro之SVM监督分类(https://mp.weixin.qq.com/s?__biz=MzI0MDQ2NTU2Ng==&mid=2247483736&idx=1&sn=2abef8d7b789b80ce851bd24f42071cf&chksm=e91b22a9de6cabbf8772db08e0d6aa82e8fc4ce7f8ff705134c709e9ac9fe5247030d090d44c&mpshare=1&scene=23&srcid=#rd)


  1. SNAP(Sentinel Application Platform)是欧空局为处理Sentinel卫星数据而开发的一个开源处理平台,功能强大,对欧空局相关卫星、美国、日本卫星等数据有很好的处理效果。 ↩︎

Logo

开放原子开发者工作坊旨在鼓励更多人参与开源活动,与志同道合的开发者们相互交流开发经验、分享开发心得、获取前沿技术趋势。工作坊有多种形式的开发者活动,如meetup、训练营等,主打技术交流,干货满满,真诚地邀请各位开发者共同参与!

更多推荐