| 高分二号异轨立体数据的森林高度提取 |  [PDF全文] |
收稿日期: 2017-04-01
2. University of Chinese Academy of Sciences, Beijing 100049, China
3. Research Institute of Forest Resources Information Techniques of Chinese Academy of Forestry, Beijing 100091, China
4. Forestry Science and Technology Information Research Institute of Chinese Academy of Forestry, Beijing 100091, China
森林生态系统是陆地生态系统的主体，虽然其分布面积仅占陆地生态系统总面积的30%，但森林植被生物量占陆地植被生物量总量的80%(Kindermann 等，2008；Pan 等，2011)，其初级生产力(GPP)占陆地生态系统总初级生产力的75% (Beer 等，2010；Pan 等，2011)。因此，准确刻画森林植被碳储量的空间分布格局及其动态变化趋势是理解碳循环机制和陆地生态系统碳收支核算的基础(Lefsky 等，2005)。作为森林植被碳储量的主体，森林地上生物量主要取决于4个变量，即树高、胸径、林分密度和树种。树高和胸径的关系符合一定的生态学规律，可用生长方程进行描述。因此树高是森林地上生物量的关键指示因子，对森林高度的精确估算是提高森林地上生物量估算的关键。
近年来在森林生物量遥感估算方面取得的主要进展都源于对森林高度的直接探测，激光雷达被认为是估算森林地上生物量的最佳数据源，其原因在于它通过激光测距实现了对包括森林高度在内的森林垂直结构的直接测量(Drake 等，2002)；干涉雷达数据(InSAR)可与已知的林下地形作差分进行有效林分高估测(庞勇 等，2003) 。欧洲正在积极发展的极化干涉雷达技术，其关键在于通过极化组合变换实现散射相位中心的分离，从而形成对森林空间结构和地形的直接探测(Papathanassiou和Cloude，2001)。
除激光雷达、干涉雷达和极化干涉雷达外，光学摄影测量数据(或称立体观测数据)也具备对森林垂直结构的直接探测能力。摄影测量通过对不同角度观测引起的视差的测量来提取地物的高程(即垂直分布)。由于光学图像不具备对森林冠层的穿透能力，光学摄影测量数据获取的是森林冠层顶部的高程，在林下地形数据的配合下，具备对森林高度的探测能力。最近十几年，研究人员利用航空数据或超高分辨率星载摄影测量数据在提取森林结构参数方面做了有益的探索。Sheng等人(2001)在利用几何体对针叶树的冠形进行近似的基础上，利用高分辨率航空摄影测量数据对单木冠形进行了重建。Gong等人(2002)提出了一种利用高分辨率航空摄影测量数据对橡树冠层边缘进行精确提取的方法。St-Onge等人(2004)将航空摄影测量数据与激光雷达获取的地形数据相结合提取植被高度，结果表明所提取的树高精度可达0.59 m。Toutin (2004)利用IKONOS立体像对对地形进行了提取，并与激光雷达数据进行对比，结果发现，在裸露地表IKONOS立体像对对地形的估算精度可达1.5 m。St-Onge等人(2008)使用IKONOS立体像对联合激光雷达获取的林下地形数据获取了森林高度，结果表明其与激光雷达获取的森林高度的相关系数R2达到0.81。Hobi和Ginzler(2012)对利用星载的WorldView-2获取的立体像对提取DSM的能力进行了评价。
除了这些超高分辨率(<0.5 m)的航空/航天摄影测量外，专业的星载摄影测量系统不断发射升空，如日本2006年发射的ALOS卫星上搭载了摄影测量传感器PRISM，中国也发射了资源三号(ZY-3)星载摄影测量系统。基于星载摄影测量数据的森林结构参数提取逐渐得到重视；如Ni等人(2015a)对ASTER GDEM数据的分析表明，星载摄影测量数据中的确含有森林高度信息；对ALOS/PRISM数据的分析结果表明，在林下地形数据的支持下，高分辨率摄影测量数据对森林垂直结构的刻画能力与机载激光雷达(LVIS)相当，各种扰动造成的林班清晰可见，与机载大光斑激光雷达(LVIS)提取的重心高(RH50)的相关性(R2)可达0.74，RMSE为2.56 m (Ni 等，2014a)；基于ZY-3数据的研究结果表明，利用冬季的高分辨率摄影测量数据，可以有效提取落叶林区的林下地形，以机载小光斑激光雷达获取的林下地形为参考，对林下地形的提取精度可达2.67 m(Ni 等，2015b)。
虽然已有研究结果证实了星载摄影测量数据在森林高度提取方面的潜力，但光学遥感数据最大的问题受是云雨等天气因素的影响严重。区域森林地上生物量产品的生产需要充分挖掘潜在数据源，通过多源数据的融合降低天气因素的制约。近几年，以宽画幅高分辨率为特色的国产高分系列对地观测系统飞速发展，特别是GF-2的全色图像的空间分辨率达到了0.8 m。虽然GF-2不是专业的星载摄影测量系统，但具备±35°的的侧摆能力，理论上在同一地区，通过侧摆形成的重复观测可构成交轨立体观测。摄影测量对高程的测量精度主要取决于基高比(立体交会角)和图像空间分辨率。同等的观测精度条件下，提高图像空间分辨率可降低对基高比的要求。因此GF-2号数据具备立体观测的潜力。本研究在机载激光雷达(LiDAR)数据提供的林下地形(DEM)数据的支持下，以机载激光雷达数据获取的森林高度(CHM)为参考，重点分析GF-2通过侧摆构成的异轨立体观测数据用于对森林高度测量的可行性和精度。
2、试验区与数据 (2.1) 试验区的选择本文的研究区位于内蒙古自治区根河市大兴安岭林区，51°01′°N，121°21′°E。该研究区地处大兴安岭西北坡，与鄂伦春旗，额尔古纳市，牙克石市以及黑龙江省漠河县和塔河县毗邻。海拔高度范围在650—1269 m，坡度15°以内的缓坡在80%以上，地势起伏较为平缓。气候类型属于寒温带温润性森林气候，年平均气温为–5.4 ℃，历史最低气温为–52.6 ℃。雨季主要集中在6—8月份，年均降雨量在450—550 mm。该研究区植被覆盖率达95%以上，森林覆盖率为73.6%，主要树种为兴安落叶松和白桦，其次是樟子松、云杉等。林下植被主要有大叶杜鹃、东北赤杨和绣线菊等。
(2.2) 试验数据本文采用的GF-2数据分别获取于2015年6月20日和2016年7月19日，它们的轨道侧摆角(滚动角)分别为0.00118°和20.4984°，它们的空间范围分别如图1中的黄色和蓝色多边形所示；重复覆盖区域长19.1 km，宽14.5 km。
|
| 图 1 高分二号数据和机载激光雷达数据空间覆盖范围(背景图为Googleearth影像) Figure 1 Spatial coverage of GF-2 and the LiDAR data (The background is the image from Googleearth) |
图1中白色多边形所示的区域是机载激光雷达数据的覆盖范围。机载激光雷达数据在2012年8月30日—9月14日，由运-5飞机搭载莱卡ALS60小光斑激光雷达系统获取。系统工作波段为1064 nm，激光光束发散角为0.22 mrad，飞行相对高度为1800 m，数据获取频率为2—4个点/m2。激光测距辅以平台的定位测量数据(GPS)和惯性导航测量(INS)数据实现对激光雷达点云坐标的解算。与RTK采集的地面控制点的比较结果显示，点云的水平定位精度优于0.5 m，垂直定位精度优于15 cm。采用TerraScan软件和必要的人工编辑，将点云分为地面点，植被点和其他点。对地面点云栅格化和插值处理得到0.5 m分辨率的林下地形数据(DEM)；对地面点和植被点进行栅格化得到0.5 m分辨率的数字表面高程模型(DSM)，两者相减得到森林高度(CHM)。
3、试验方法 (3.1) GF-2的DSM提取利用摄影测量数据进行高程自动提取包含两步，即构建成像几何模型和密集点云匹配。成像几何模型用于描述图像坐标与地理坐标的几何关系。成像几何模型分为严格成像几何模型和通用成像几何模型。严格成像几何模型是通过严格追踪从卫星轨道、卫星姿态和相机成像等物理过程来建立的几何模型，其定位精度高，但由于其需要用到卫星平台和有效载荷的一些属于技术机密的关键参数，对于未公开的卫星，难以建立严格的成像几何模型。通用成像几何模型则是利用有理函数来建立图像坐标与地理坐标的几何关系。由于通用成像几何模型不涉及机密的关键参数，而且模型形式统一，不同的图像只需要提供有理函数的系数(RPC)，便于相关数据处理软件的开发和维护，因此受到广泛的欢迎。
目前的通用商业软件(如ENVI，ERDAS、PCI)都已经集成了较好的密集点云匹配算法，本文使用GF-2数据及同时提供的RPC参数，采用ENVI提供的立体观测数据机处理功能(Topographic菜单下的DEM Extraction wizard工具)，以2015年6月20日的全色图像为左片，2016年7月19日的全色影像为右片，进行立体影像处理得到摄影测量点云数据(Ni 等，2012)。处理过程包括连接点的选择和密集匹配参数设置两个关键步骤。软件提供了连接点自动选择工具，该工具需要设定所需的连接点个数、搜索窗口大小、移动窗口大小和最小的相关系数。本数据处理过程采用软件默认连接点数25，搜索窗口大小采用251，移动窗口采用65，最小相关系数采用0.85。自动选择工具会在云边缘或纹理信息欠丰富的地方产生错误，因此需要对自动选择的连接点进行目视检查。
密集匹配参数设置主要包括，最小相关系数、移动窗口大小，地面起伏状况(terrain relief)和地形精细程度(terrain detail)。本文采用的最小相关系数为0.85；移动窗口大小采用默认的5×5；地面起伏状况有low，Moderate和High共3个级别，本文选用的是High；地形精细程度根据图像分辨率和大小会有不同的级别，针对GF-2数据软件给出了7个级别的选择，本文采用最高的级别。在完成参数和相应文件名的设置后，软件进入自动化处理状态，直到最后得到摄影测量点云数据。
(3.2) GF-2的DSM数据纠正与分析综合考虑试验区内单木树冠的尺寸和GF-2数据的空间分辨率，以10 m的网格对GF-2摄影测量点云数据进行栅格化得到GF-2的DSM，以落入每个像元的所有点的平均值作为该像元的高程值。激光雷达的DEM和DSM通过平均值重采样的方法将其分辨率变为10 m。
通用成像几何模型易于生产和使用，但其精度比严格成像几何模型要低。此外，本研究采用激光雷达提供的DEM做为林下地形，不同源高程数据之间的高程基准会存在一定的差异。因此需要对GF-2的DSM与激光雷达DEM进行水平位置的配准和高程差异的纠正。水平位置的配准采用基于仿射变换模型的高精度DEM数据配准算法来实现(Ni 等，2014b)，所采用的变换模型如式(1)所示，其基本过程是对两个高程数据分别进行渲染，利用渲染后的图像进行同名点的选择，进一步利用同名点进行图像配准参数的计算，利用图像配准参数对GF-2 DSM进行重采样，实现两套数据的水平配准。由于数据覆盖范围较小，两个高程基准的差异可认为是常数，可以利用裸露地表处两个数据的高程差值对它们之间的高程差异进行纠正。
| $\left[ {\begin{array}{*{20}{c}}\!\!\!\! {{x_m}}\!\!\!\! \!\!\!\! {{y_m}}\!\!\!\! \end{array}} \right] = {{A}}\left[ {\begin{array}{*{20}{c}}\!\!\!\! {{x_s}}\!\!\!\! \!\!\!\! {{y_s}}\!\!\!\! \end{array}} \right] + \left[ {\begin{array}{*{20}{c}}\!\!\!\! {{x_0}}\!\!\!\! \!\!\!\! {{y_0}}\!\! \!\!\end{array}} \right],\;\;{{A}} = \left[ {\begin{array}{*{20}{c}}\!\!\!\! {\cos\varphi } & {\sin\varphi } \!\!\!\! \!\!\!\! { - \sin\varphi } & {\cos\varphi } \!\!\!\! \end{array}} \right]\left[ {\begin{array}{*{20}{c}}\!\!\!\! {{k_x}} & \omega \!\!\!\! \!\!\!\! 0 & {{k_y}}\!\!\!\! \end{array}} \right]$ | (1) |
式中，
经过纠正后的高分二号DSM与激光雷达DEM相减得到摄影测量CHM，通过与激光雷达CHM做像元对像元的散点分析，对摄影测量CHM的精度进行评价。
4、试验结果与分析图2(a)展示了利用GF-2核线影像合成的红蓝立体图像(左片与右片分别采用红色和蓝色通道合成，可用红蓝立体镜进行观察)，图2(b)为图2(a)中白色方框所示位置的局部放大图。由于GF-2数据的基高比(立体交会角)较小，立体观测能力主要依靠高图像分辨率来实现，因此图2(a)的立体感觉较弱，相比之下分辨率较高的图2(b)的立体感觉较强。
|
| 图 2 红蓝彩色合成的核线图像(红色通道为左片，蓝色和绿色通道为右片) Figure 2 Red-blue stereoscopic image (Red: Left epipolar image, Green and blue：Right epipolar image) |
|
| 图 3 GF-2立体观测数据处理结果 Figure 3 Stereoscopic processing results of GF-2 images |
图3展示了利用GF-2数据进行立体处理的结果。图3(a)为两台相机重复观测区域的左片全色影像，图3(b)为摄影测量点云数据栅格化成10 m分辨率DSM时，每个像元内的点数图像。理想情况下，每个像元内最多100个点。可以看出，大部分像元的点云密度较高，这是由于GF-2立体交会角较小，立体像对的图像纹理相似性较高。点云密度较低的像元主要分布在低矮灌草丛覆盖的平坦区域，这是因为立体像对的获取时间间隔长达一年，这些地方的图像纹理发生了较大变化。图3(c)为立体处理得到的DSM，图3(d)左图为图3(c)白色矩形框覆盖区域的局部放大图，图3(d)右图为激光雷达得到的对应的DSM，可以看出GF-2与激光雷达获取的DSM的空间格局非常相似，可以清楚地看到树林斑块和林间道路等细节信息。
图4展示了GF-2的DSM与激光雷达DEM重叠区域的渲染图像，共选择了同名点18个，X和Y方向的准精度分别为0.138个像元和0.145个像元，相应的配准参数如表1所示。
图5(a)以剖面线的形式展示了GF-2的DSM与激光雷达DEM高程差异，两者高程差异为156 m。图5(b)展示了消除高程差异后的剖面线。可以看出，消除高程差异后，GF-2的DSM与激光雷达DSM吻合较好。
|
| 图 4 GF-2的DSM与激光雷达DEM渲染图 Figure 4 Shaded image of GF-2 DSM and LiDAR DEM |
|
表 1 激光雷达与GF-2的DEM/DSM配准参数 Table 1 Parameters for the co-registration between GF-2 DSM and LiDAR DEM |
|
| 图 5 GF-2的DSM与激光雷达DEM高程差异纠正前后对比 Figure 5 Correction of the elevation divergence between GF-2 DSM and LiDAR DEM |
图6展示了GF-2的CHM与激光雷达CHM的对比。图6(a)为GF-2 的DSM与激光雷达DEM相减得到的GF-2的CHM图像；图6(b)为激光雷达CHM图像。由图6(a)和图6(b)明显可以看出两者的空间分布非常相似，表明GF-2的观测数据对森林冠层高度的刻画能力较好。图6(c)为图6(b)与图6(a)像元对像元的散点图。由于像元较多，图6(c)是点密度分割图，以0.05 m为栅格大小，图像灰度表示的是落入每个栅格的点数，右侧图例给出了落入每个栅格内的像元数与灰度的对应关系。可以看出GF-2的CHM与激光雷达CHM呈现明显的线性相关关系。线性拟合的显示，两者的相关性R2=0.51，均方根误差为3.6 m。需要指出的是，激光雷达与GF-2数据获取时间间隔为3—4年，研究表明研究区的优势树种兴安落叶松的年生长量平均在0.5 m左右(张璐，2014)。如果缩短激光雷达数据与GF-2数据获取时间间隔，它们之间的相关性有可能会有所提高。
|
| 图 6 GF-2的CHM与激光雷达CHM的对比 Figure 6 Comparison between GF-2 CHM and LiDAR CHM |
本研究在内蒙古自治区根河市大兴安岭林区的研究区内，以分别获取于2015年6月20日和2016年7月19日的GF-2全色数据作为立体像对，其标称的轨道侧摆角分别为0.00118°和20.4984°。通过对立体像对进行处理得到摄影测量点云，将点云进行栅格化后得到10 m分辨率的DSM，以激光雷达数据提供的DEM作为林下地形，将DSM和DEM相减得到了GF-2的CHM。
结果表明，GF-2的CHM与激光雷达的CHM存在明显的相关性，两者的空间分布格局较为一致，其中，由二者的CHM的像元得到的散点图发现，像元对像元的线性相关性(R2)可以达到0.51，均方根误差为3.6 m。因此，本研究结果表明，在林下地形已知的情况下，GF-2的异轨立体观测数据能较好地刻画森林冠层垂直结构信息，可以用于森林高度的估算。
本研究所使用的立体像对的获取时间间隔长达一年，且激光雷达数据的获取时间间隔与GF-2数据的获取时间间隔长达3—4年。在后续的研究中，将选用获取时间间隔更短的立体像对和激光雷达数据，并将针对不同的森林类型进行分析。
志 谢 本文所使用的GF-2数据由中国科学院遥感与数字地球研究所/对地观测应用技术中心提供，在此表示衷心的感谢！
| [1] | Beer C, Reichstein M, Tomelleri E, Ciais P, Jung M, Carvalhais N, Rödenbeck C, Arain M A, Baldocchi D, Bonan G B, Bondeau A, Cescatti A, Lasslop G, Lindroth A, Lomas M, Luyssaert S, Margolis H, Oleson K W, Roupsard O, Veenendaal E, Viovy N, Williams C, Woodward F I and Papale D. Terrestrial gross carbon dioxide uptake: global distribution and covariation with climate[J]. Science, 2010, 329 (5993) : 834 –838. DOI: 10.1126/science.1184984 |
| [2] | Drake J B, Dubayah R O, Clark D B, Knox R G, Blair J B, Hofton M A, Chazdon R L, Weishampel J F and Prince S. Estimation of tropical forest structural characteristics using large-footprint lidar[J]. Remote Sensing of Environment, 2002, 79 (2/3) : 305 –319. DOI: 10.1016/S0034-4257(01)00281-4 |
| [3] | Gong P, Mei X L, Biging G S and Zhang Z X. Improvement of an oak canopy model extracted from digital photogrammetry[J]. Photogrammetric Engineering and Remote Sensing, 2002, 68 (9) : 919 –924. |
| [4] | Hobi M L and Ginzler C. Accuracy assessment of digital surface models based on WorldView-2 and ADS80 stereo remote sensing data[J]. Sensors, 2012, 12 (5) : 6347 –6368. DOI: 10.3390/S120506347 |
| [5] | Kindermann G E, McCallum I, Fritz S and Obersteiner M. A global forest growing stock, biomass and carbon map based on FAO statistics[J]. Silva Fennica, 2008, 42 (3) : 387 –396. DOI: 10.14214/sf.244 |
| [6] | Lefsky M A, Harding D J, Keller M, Cohen W B, Carabajal C C, Del Bom Espirito-Santo F, Hunter M O and de Oliveira R Jr. Estimates of forest canopy height and aboveground biomass using ICESat[J]. Geophysical Research Letters, 2005, 32 (22) : L22S02 . DOI: 10.1029/2005GL023971 |
| [7] | Ni W J, Guo Z F, Zhang Z Y, Sun G Q and Huang W L. 2012. Semi-automatic extraction of digital surface model using ALOS/PRISM data with ENVI//Proceedings of 2012 IEEE International Geoscience and Remote Sensing Symposium (IGARSS). Munich, Germany: IEEE: 6557–6560 [DOI: 10.1109/IGARSS.2012.6352097] |
| [8] | Ni W J, Ranson K J, Zhang Z Y and Sun G Q. Features of point clouds synthesized from multi-view ALOS/PRISM data and comparisons with LiDAR data in forested areas[J]. Remote Sensing of Environment, 2014a, 149 : 47 –57. DOI: 10.1016/j.rse.2014.04.001 |
| [9] | Ni W J, Sun G Q and Ranson K J. Characterization of ASTER GDEM elevation data over vegetated area compared with lidar data[J]. International Journal of Digital Earth, 2015b, 8 (3) : 198 –211. DOI: 10.1080/17538947.2013.861025 |
| [10] | Ni W J, Sun G Q, Ranson K J, Pang Y, Zhang Z Y and Yao W. Extraction of ground surface elevation from ZY-3 winter stereo imagery over deciduous forested areas[J]. Remote Sensing of Environment, 2015a, 159 : 194 –202. DOI: 10.1016/j.rse.2014.12.007 |
| [11] | Ni W J, Sun G Q, Zhang Z Y, Guo Z F and He Y T. Co-registration of two DEMs: impacts on forest height estimation from SRTM and NED at mountainous areas[J]. IEEE Geoscience and Remote Sensing Letters, 2014b, 11 (1) : 273 –277. DOI: 10.1109/LGRS.2013.2255580 |
| [12] | Pan Y D, Birdsey R A, Fang J Y, Houghton R, Kauppi P E, Kurz W A, Phillips O L, Shvidenko A, Lewis S L, Canadell J G, Ciais P, Jackson R B, Pacala S W, McGuire A D, Piao S L, Rautiainen A, Sitch S and Hayes D. A large and persistent carbon sink in the world’s forests[J]. Science, 2011, 333 (6045) : 988 –993. DOI: 10.1126/science.1201609 |
| [13] | 庞勇, 李增元, 陈尔学, 车学俭, 白黎娜, 谭炳香. 干涉雷达技术用于林分高估测[J]. 遥感学报, 2003, 7 (1) : 8 –13. Pang Y, Li Z Y, Chen E X, Che X J, Bai L N and Tan B X. InSAR technology and its application to estimate stand average height[J]. Journal of Remote Sensing, 2003, 7 (1) : 8 –13. DOI: 10.11834/jrs.20030102 |
| [14] | Papathanassiou K P and Cloude S R. Single-baseline polarimetric SAR interferometry[J]. IEEE Transactions on Geoscience and Remote Sensing, 2001, 39 (11) : 2352 –2363. DOI: 10.1109/36.964971 |
| [15] | Sheng Y W, Gong P and Biging G S. Model-based conifer-crown surface reconstruction from high-resolution aerial images[J]. Photogrammetric Engineering and Remote Sensing, 2001, 67 (8) : 957 –965. |
| [16] | St-Onge B, Jumelet J, Cobello M and Véga C. Measuring individual tree height using a combination of stereophotogrammetry and lidar[J]. Canadian Journal of Forest Research, 2004, 34 (10) : 2122 –2130. DOI: 10.1139/X04-093 |
| [17] | St-Onge B, Vega C, Fournier R A and Hu Y. Mapping canopy height using a combination of digital stereo-photogrammetry and lidar[J]. International Journal of Remote Sensing, 2008, 29 (11) : 3343 –3364. DOI: 10.1080/01431160701469040 |
| [18] | Toutin T. DTM generation from Ikonos in-track stereo images using a 3D physical model[J]. Photogrammetric Engineering and Remote Sensing, 2004, 70 (6) : 695 –702. DOI: 10.14358/PERS.70.6.695 |
| [19] | 张璐. 2014. 内蒙古大兴安岭兴安落叶松林分直径结构及生长过程研究. 呼和浩特: 内蒙古农业大学 Zhang L. 2014. Research the Diametric Structure of Forest Stands and the Growth Process of Larix Gmelinii at Daxing’ Anling in Inn Mongolia. Hohhot: Inner Mongolia Agricultural University |