2.1 数据预处理
DMSP/OLS夜间灯光数据的预处理包含了重投影、重采样、研究区域裁剪和原始影像相互校正、不连续校正两大部分:(1)影像重投影、影像重采样与裁剪。将原始影像的投影坐标系转换为兰伯特等角圆锥投影(Lambert Conformal Conic),并将影像像元重采样为1 km,裁剪研究区年际影像。(2)相互校正与不连续校正。DMSP卫星传感器在获取地表数据的过程中受到多种因素的影响(如大气层的吸收和散射、太阳高度角、地形起伏度、传感器校准等)[33],导致不同传感器获取的同一年份的影像存在差异。同时,不同的OLS传感器在获取影像时并未进行星上辐射校正[34],造成同一卫星传感器获取的连续不同年份影像间相同位置的亮值像元DN值的异常波动。
2.1.1 相互校正
根据曹子阳的研究[32],选取中国黑龙江省社会经济发展较为稳定的鹤岗市作为不变区域标定区,采用不变参考区方法[34]对夜间灯光影像相互校正。以F162006中定标区域的像元DN值为标准,将其他年份的影像与参考数据F162006作回归曲线估计计算得到回归方程模型参数。
DNc=a×DN2+b×DN+c(1)
式中:a,b,c为二次回归模型参数; DNc与DN分别表示校正后与校正前影像的像元DN值。利用SPSS软件将待校正影像像元DN值与F16卫星传感器的2006年的辐射定标的夜间灯光影像像元DN值进行二次回归。SPSS软件计算的二次回归模型参数见表1。
2.1.2 不连续校正
经过相互校正后的夜间灯光影像仍然还存在像元DN值异常问题,因此,还需要进行影像的不连续校正。假设前一年的夜间灯光影像数据的亮值像元在后一年灯光影像中同一位置的像元也应保持为亮值像元,前一年夜间灯光影像中的亮值像元DN值,应小于或等于后一年影像中同一位置的亮值像元DN值[32],即满足公式2。
式中:DN(n-1,i),DN(n,i),DN(n+1,i)分别表示第n-1年,第n年和第n+1年通过相互校正和不连续校正后的夜间灯光影像的i像元的DN值。经过校正后的图像消除了像元DN值异常和“饱和”现象的问题(图1)。
2.2 基于夜间灯光数据的城市建成区提取
基于DMSP/OLS夜光数据提取城镇建成区的方法主要有4种:城市轮廓突变检测法[35]、经验阈值法[36]、基于较高分辨率影像数据空间比较法[37]、基于统计数据比较法[38]。通过比较,本研究选取“基于统计数据比较法”作为滇中城市群的城镇信息提取方法。借助ArcGIS 10.5软件,运用二分迭代,快速获取阈值,获取每个城镇的地理空间结构及空间信息。二分法迭代快速计算获取阈值,算法如下[39]。
DNt=int[(DNmax-DNmin)/2](3)
式中:DNmax和DNmin分别为某年某地区灯光影像像元最大灰度值和像元最小灰度值; DNt为该区域灯光影像的潜在阈值。
将DNt提取的城镇用地面积和统计资料的建成区面积进行比较,提取面积最接近统计年鉴值的则DNt设置为最优阈值。同时利用高分辨率数据(Google Earth影像)提取1993年、1998年、2003年、2008年、2013年的滇中城市群典型城镇边界,两者反复进行对比验证,进而互相精度校正,以获取最精确的提取结果。
利用最接近统计数据的阈值重分类(式4)夜间灯光影像,进而重建滇中城市群城市空间信息,得到滇中城市群基于DMSP/OLS夜间灯光遥感数据提取的建成区空间分布图(图2)。
式中:DN为像元灰度值; DNt为最优阈值; 0代表未建成区,1表示建成区。
2.3 城市扩张指数
为比较不同研究时段内城市扩张的快慢或强弱,采用城市扩张速率与城市扩张强度指数进行探究。针对所提取的建成区面积,可得面积增长率、城市扩张速率数学模型如下:
式中:At为面积增长率; Aa,Ab分别为前后年份提取的建成区面积; a为增长面积; Vt为扩张速率; t为时间间隔。
图2 1993年、1998年、2003年、2008年、2013年夜间灯光数据建成区空间分布
城市扩展强度指数是城市扩张空间变化的一个重要指标,通过分析城市扩展强度指数可定量地比较城市扩张的程度及速度。城市扩张强度数学模型如下:
式中:Aa,Ab分别为前后年份提取的建成区面积; t为时间间隔; Ca为扩张强度指数。