1.一种城市不同土地覆盖类型遥感自动解译方法,其特征在于,其步骤如下:
S1、利用遥感云平台采集X年的N幅Landsat卫星影像,并对N幅Landsat卫星影像进行数据预处理和云层剔除处理,得到N幅地表反射率图像;其中,Landsat卫星影像包括蓝色波段、绿色波段、红色波段、近红外波段、第一短波红外波段和第二短波红外波段这6个波段影像;
S2、基于地面调查数据确定城市包含的土地覆盖类型,基于GPS定位工具获取土地覆盖类型的地理坐标,并使用光谱采集软件采集分析地表反射率图像中的土地覆盖类型的像元光谱样本,得到用于识别土地覆盖类型的有效光谱波段;其中,土地覆盖类型包括蓝色彩钢板房、红色不透水面、水泥不透水面、其他不透水面、植被、水体、未识别类型这7种类型;
S3、基于水体遥感指数和步骤S1中得到的N幅地表反射率图像计算得到N幅水体遥感指数图像,采用均值图像重构方法对N幅水体遥感指数图像进行图像重构,得到X年水体遥感指数均值图像,记为WI图像;基于步骤S2中地面调查数据中包含的水体样本统计分析水体在WI图像中的分布区间,取其最小值作为阈值α;
S4、基于步骤S1中得到的N幅地表反射率图像中的第二短波红外波段图像,采用四分位数图像重构方法对N幅第二短波红外波段图像进行图像重构,取其重构的第零个四分位数图像为X年第二短波红外波段最小值图像,记为SWIRmin图像,取其重构的第二个四分位数图像为X年第二短波红外波段中值图像,记为SWIRM图像;基于步骤S2中地面调查数据中包含的水体样本统计分析水体在SWIRmin图像中的分布区间,取其最大值作为阈值β;基于步骤S2中地面调查数据中包含的红色不透水面样本统计分析红色不透水面在SWIRM图像中的分布区间,取其最小值作为阈值κ;基于步骤S2中地面调查数据中包含的水泥不透水面样本统计分析水泥不透水面在SWIRM图像中的分布区间,取其最小值作为阈值ν;
S5、创建蓝色彩钢板房遥感指数,基于步骤S1中得到的N幅地表反射率图像计算得到N幅蓝色彩钢板房遥感指数图像,然后采用四分位数图像重构方法对N幅蓝色彩钢板房遥感指数图像进行图像重构,取其重构的第三个四分位数图像为X年蓝色彩钢板房遥感指数重构图像,记为BI图像;基于步骤S2中地面调查数据中包含的蓝色彩钢板房样本统计分析蓝色彩钢板房在BI图像中的分布区间,取其最小值作为阈值γ;
创建蓝色彩钢板房遥感指数的方法为:
BPIi=(ξblue,i‑ξred,i)/(ξblue,i+ξred,i);
其中,BPIi表示蓝色彩钢板房遥感指数图像中像元i上的蓝色彩钢板房遥感指数值,ξblue,i表示地表反射率图像中像元i上的蓝色波段的地表反射率,ξred,i表示地表反射率图像中像元i上的红色波段的地表反射率,i=1,2,…,n,n为地表反射率图像中像元的总数量;
S6、基于步骤S1中得到的N幅地表反射率图像中的蓝色波段图像,采用四分位数图像重构方法对N幅蓝色波段图像进行图像重构,取其重构的第二个四分位数图像为X年蓝色波段中值图像,记为BM图像;基于步骤S2中地面调查数据中包含的蓝色彩钢板房样本统计分析蓝色彩钢板房在BM图像中的分布区间,取其最小值作为阈值δ;基于步骤S2中地面调查数据中包含的水泥不透水面样本统计分析水泥不透水面在BM图像中的分布区间,取其最小值作为阈值μ;
S7、基于植被遥感指数和步骤S1中得到的N幅地表反射率图像计算得到N幅植被遥感指数图像,采用四分位数图像重构方法对N幅植被遥感指数图像进行图像重构,取其重构的第四个四分位数图像为X年植被遥感指数最大值图像,记为VI图像;基于步骤S2中地面调查数据中包含的植被样本统计分析纯净像元的植被在VI图像中的分布区间,取其最小值作为阈值ε;统计分析混合像元的植被在VI图像中的分布区间,取其最小值作为阈值ζ;基于步骤S2中地面调查数据中包含的其他不透水面样本统计分析纯净像元的其他不透水面在VI图像中的分布区间,取其最大值作为阈值π;
S8、基于步骤S1中得到的N幅地表反射率图像计算得到N幅蓝色波段与红色波段的和值图像,采用四分位数图像重构方法对这N幅蓝色波段与红色波段的和值图像进行图像重构,取其重构的第二个四分位数图像为X年蓝色波段与红色波段的和值重构图像,记为RBsum图像;基于步骤S2中地面调查数据中包含的植被样本统计分析植被在RBsum图像中的分布区间,取其最大值作为阈值η;基于步骤S2中地面调查数据中包含的其他不透水面样本统计分析其他不透水面在RBsum图像中的分布区间,取其最小值作为阈值σ;
S9、基于步骤S1中得到的N幅地表反射率图像中的近红外波段图像,采用四分位数图像重构方法对N幅近红外波段图像进行图像重构,取其重构的第四个四分位数图像为X年近红外波段最大值图像,记为NIRmax图像;基于步骤S2中地面调查数据中包含的植被样本统计分析植被在NIRmax图像中的分布区间,取其最小值作为阈值θ;基于步骤S2中地面调查数据中包含的其他不透水面样本统计分析其他不透水面在NIRmax图像中的分布区间,取其最大值作为阈值ο;
S10、创建红色不透水面遥感指数,基于步骤S1中得到的N幅地表反射率图像计算得到N幅红色不透水面遥感指数图像,然后采用四分位数图像重构方法对N幅红色不透水面遥感指数图像进行图像重构,取其重构的第三个四分位数图像为X年红色不透水面遥感指数重构图像,记为RI图像;基于步骤S2中地面调查数据中包含的红色不透水面样本统计分析红色不透水面在RI图像中的分布区间,取其最小值作为阈值ι;
创建红色不透水面遥感指数的方法为:
RRIi=(ξred,i‑ξgreen,i)/(ξred,i+ξgreen,i);
其中,RRIi表示红色不透水面遥感指数图像中像元i上的红色不透水面遥感指数值,ξred,i表示地表反射率图像中像元i上的红色波段的地表反射率,ξgreen,i表示地表反射率图像中像元i上的绿色波段的地表反射率,i=1,2,…,n,n为地表反射率图像中像元的总数量;
S11、基于步骤S1中得到的N幅地表反射率图像中的红色波段图像,采用四分位数图像重构方法对N幅红色波段图像进行图像重构,取其重构的第二个四分位数图像为X年红色波段中值图像,记为RM图像;基于步骤S2中地面调查数据中包含的水泥不透水面样本统计分析水泥不透水面在RM图像中的分布区间,取其最小值作为阈值ξ;
S12、判断像元i在步骤S3中得到的水体遥感指数均值图像中的像元值是否大于阈值α,若不是,执行步骤S13,若是,判断像元i在步骤S4中得到的第二短波红外波段最小值图像中的像元值是否小于阈值β,若是,像元i为水体,若不是,执行步骤S13;
S13、判断像元i在步骤S5中得到的蓝色彩钢板房遥感指数重构图像中的像元值是否大于阈值γ,若不是,执行步骤S14,若是,判断像元i在步骤S6中得到的蓝色波段中值图像中的像元值是否大于阈值δ,若是,像元i为蓝色彩钢板房,若不是,执行步骤S14;
S14、判断像元i在步骤S7中得到的植被遥感指数最大值图像中的像元值是否大于阈值ε,若不是,执行步骤S15,若是,像元i为植被;
S15、判断像元i在步骤S7中得到的植被遥感指数最大值图像中的像元值是否属于区间范围[ζ,ε],若不是,执行步骤S16,若是,判断像元i在步骤S8中得到的蓝色波段与红色波段的和值重构图像中的像元值是否小于阈值η,若不是,执行步骤S16,若是,判断像元i在步骤S9中得到的近红外波段最大值图像中的像元值是否大于阈值θ,若不是,执行步骤S16,若是,像元i为植被;
S16、判断像元i在步骤S10中得到的红色不透水面遥感指数重构图像中的像元值是否大于阈值ι,若不是,执行步骤S17,若是,判断像元i在步骤S4中得到的第二短波红外波段中值图像中的像元值是否大于阈值κ,若不是,执行步骤S17,若是,像元i为红色不透水面;
S17、判断像元i在步骤S6中得到的蓝色波段中值图像中的像元值是否大于阈值μ,若不是,执行步骤S18,若是,判断像元i在步骤S4中得到的第二短波红外波段中值图像中的像元值是否大于阈值ν,若不是,执行步骤S18,若是,判断像元i在步骤S11中得到的红色波段中值图像中的像元值是否大于阈值ξ,若不是,执行步骤S18,若是,像元i为水泥不透水面;
S18、判断像元i在步骤S9中得到的近红外波段最大值图像中的像元值是否小于阈值ο,若不是,像元i为未识别类型,若是,执行步骤S19;
S19、判断像元i在步骤S7中得到的植被遥感指数最大值图像中的像元值是否小于阈值π,若是,像元i为其他不透水面,若不是,执行步骤S20;
S20、判断像元i在步骤S8中得到的蓝色波段与红色波段的和值重构图像中的像元值是否大于阈值σ,若是,像元i为其他不透水面,若不是,像元i为未识别类型;
S21、循环执行步骤S12至步骤S20,直至遍历完图像中的所有像元位置,完成研究区域的城市不同土地覆盖类型的遥感识别分类;
S22、采用3×3运算窗口对步骤S21中得到的分类结果进行滤波处理,消除分类结果中的孤立像元,完成最终的城市不同土地覆盖类型的遥感识别分类。
2.根据权利要求1所述的城市不同土地覆盖类型遥感自动解译方法,其特征在于,步骤S1中对Landsat卫星影像进行数据预处理和云层剔除处理的方法为:对Landsat卫星影像中所有像元值乘以修正系数0.0000275,再减去0.2,得到Landsat卫星影像的地表反射率图像;在地表反射率图像中像元i的位置上,如果蓝色波段、绿色波段、红色波段上像元值的极差的绝对值小于0.08,并且绿色波段上的像元值大于0.25,那么像元i的位置上的地物类型为云层,则将像元i的位置上所有波段的像元值改写为空值,实现云层剔除处理。
3.根据权利要求1或2所述的城市不同土地覆盖类型遥感自动解译方法,其特征在于,所述采用均值图像重构方法对N幅水体遥感指数图像进行图像重构的方法为:首先创建1幅像元行列数与水体遥感指数图像的像元行列数一致的空值图像;然后将N幅水体遥感指数图像进行图层堆叠,像元i的位置上存在n个水体遥感指数值,计算得到这n个数值的平均值,并写入空值图像的相应位置,依次遍历所有像元位置,得到水体遥感指数均值图像。
4.根据权利要求1所述的城市不同土地覆盖类型遥感自动解译方法,其特征在于,所述采用四分位数图像重构方法对N幅第二短波红外波段图像进行图像重构的方法为:首先创建5幅像元行列数与第二短波红外波段图像的像元行列数一致的空值图像;然后将N幅第二短波红外波段图像进行图层堆叠,像元i的位置上存在n个第二短波红外波段像元值,对这n个值按照从小到大的顺序排列成数列,分别将数列中位于第零、一、二、三、四的四分位点的值写入空值图像中的相应像元位置,依次遍历所有像元位置,即分别得到第二短波红外波段图像的第零、一、二、三、四的四分位数图像;其中第零、二、四的四分位点分别表示数组的最小值、中值、最大值。
5.根据权利要求1所述的城市不同土地覆盖类型遥感自动解译方法,其特征在于,所述步骤S22中消除分类结果中孤立像元的方法为:对于像元(x,y),图像中最大行数为X,最大列数为Y,其中x=2,3,……,X‑1;y=2,3,……,Y‑1;统计以像元(x,y)为中心的3×3窗口范围内不同类别的像元个数,设置像元(x,y)的数量权重为5,与像元(x,y)类别一致的像元个数为5+m,m为3×3窗口中与像元(x,y)类别一致的像元个数,令J表示为3×3窗口中与像元(x,y)类别不一致的类别中的最大像元个数,如果5+m小于J,则像元(x,y)的类别属性改写为与J类别一致的类别属性,否则,像元(x,y)的类别属性不变;依次遍历图像中所有可用像元位置,迭代3次,完成分类结果中孤立像元的消除。