杜伟娜， 徐爱功， 宋耀鑫， 孙华生

(1.辽宁工程技术大学,阜新 123000； 2.公安部交通管理科学研究所,无锡 214151)



新型SAR传感器一级地距产品绝对辐射定标方法

杜伟娜1， 徐爱功1， 宋耀鑫2， 孙华生1

(1.辽宁工程技术大学,阜新 123000； 2.公安部交通管理科学研究所,无锡 214151)

针对目前新型SAR传感器数据产品预处理软件缺乏的现状，全面细致地论述了ENVISAT ASAR,Radarsat2,Cosmoskymed,TerraSAR-X和Sentinel1等5种新型SAR传感器的一级地距产品(level-1 detected ground range product,L-1 DGRP)数据的绝对辐射定标方法以及软件化过程中所涉及的后向散射系数等参数获取方法； 并以Sentinel1传感器的L-1 DGRP数据为例，采用C++语言编程实现了绝对辐射定标过程； 最后将本文方法软件处理结果和欧空局S1 ToolBox软件处理结果进行对比，二者得到的后向散射系数值基本相同，证实了本文所介绍方法的正确性。

新型SAR； 绝对辐射定标； 一级地距产品(L-1 DGRP)； 后向散射系数； Sentinel1

0 引言

随着国外高分辨率合成孔径雷达(synthetic aperture Radar,SAR)数据对我国各领域的开放以及“龙计划”项目[1]的深入开展，ENVISAT ASAR,Radarsat1/2,Cosmoskymed,TerraSAR-X,ALOS PALSAR和Sentinel1等高分辨率星载SAR影像数据受到广大关注。然而，当前能够支持这些星载SAR数据预处理的商业软件还比较少； 因此，研究这些SAR传感器各级产品数据的辐射定标算法并进行软件实现，对于促进新型高分辨率合成孔径成像雷达数据的定量应用具有重要意义。与SAR数据的其他各级产品相比，一级地距产品(level-1 detected ground range product,L-1 DGRP)数据在陆地资源环境遥感中的应用更为广泛，也是卫星地面接收站和数据生产商主要提供的数据产品。L-1 DGRP数据记录的是雷达微波后向散射信号的振幅，而非地球生物物理参数定量化研究所需的地物后向散射系数[2]。在地表参数反演中，辐射定标处理是不可或缺的，该过程是将传感器接收的地物后向散射强度信息转化为地物后向散射系数的唯一途径。因此，本文根据SAR系统辐射定标的基本原理，介绍了ENVISAT ASAR,Radarsat2,Cosmoskymed,TerraSAR-X和Sentinel1等5种新型星载SAR传感器的L-1 DGRP数据的绝对辐射定标方法； 并以Sentinel1数据产品的辐射定标为例，采用C++语言编程实现了其绝对辐射定标过程。

1 SAR辐射定标基本原理

星载合成孔径雷达在探测地物目标时，雷达系统主动向目标发射无线电脉冲，其发射能量部分被目标吸收，其余大部分能量经目标散射后被雷达天线接收，形成回波信号。这一能量传输过程可用雷达方程定量表示[3-4]为

(1)

式中：Pr为雷达接收的信号功率；Pt为发射功率；Gr为总接收增益；G(θ,φ)为距离向和方位向的增益；λ为波长；S为地距；σ为目标的雷达截面积。

雷达截面积σ是在假设雷达所接收的回波信号功率密度与实际目标在同等条件下所接收的目标回波信号功率密度相同的理想条件下，定义的实际目标的等效雷达截面积，是表征地物目标散射特性的常用参数。

对于分布目标，雷达截面积σ可表示为后向散射系数的函数，即

σ=σ0Ac,

(2)

式中：σ0为单位面积的雷达截面积(即后向散射系数)，可用来表征目标对电磁波的散射能力；Ac为雷达地面分辨单元的面积。

星载SAR传感器在发射之前，通常会利用角反射器在地面试验场进行多次绝对辐射定标实验； 并且大多数SAR传感器会将辐射定标模型所需的参数简化成角度、天线增益、定标常量等[5]，并存储在数据产品的头文件中供用户使用。

2 典型SAR传感器绝对辐射定标方法

与传统SAR相比，新一代SAR传感器普遍具有高空间分辨率、高灵敏辐射强度、多极化模式和多扫描带宽等特点，表1列出5种典型SAR传感器的参数[6-11]。

表1 5种新型SAR传感器的参数

2.1 ENVISAT ASAR数据辐射定标

ENVISAT ASAR[6]是迄今为止在环境监测领域应用最为广泛的星载合成孔径雷达，为遥感应用领域的学者提供了大量高质量的数据。ASAR的L-1 DGRP数据包括： ASA_IMP_1P,ASA_IMM_1P,ASA_APP_1P, ASA_APM_1P, ASA_WSM_1P, ASA_IMG_1P和 ASA_APG_1P共7种，通常以“.N1”文件格式分发。图像中任一像元(i,j)的后向散射系数σ0为

(3)

对后向散射系数求平均，以db为单位表示，其计算公式为

(4)

式中:U和V分别为平均窗口的行和列；N=UV。

以上定标模型涉及的参数均可从“.N1”格式的文件中直接或间接获取。DN是数据产品的像元值，直接从SAR数据产品文件中的雷达数据记录中读取； 标定常数K从ASAR数据文件的元数据中读取，该参数与极化方式有关，具体读取方式见表2。

表2 ASAR数据绝对定标常数K

入射角αi,j需要通过计算得到，在注记数据集下的Geolocation Grid ADSRs中记录了SAR的入射角信息，即影像方位向的行和距离向的列均存在对应的11组入射角，各像元(i,j)的入射角通过插值得到。由于影像在方位向不超过100 km时可认为入射角随影像方位向行号变化的影响可以忽略不计，因此入射角只随列号而变化。通常选择一元二次方程作为插值函数，利用最小二乘法估计，即可计算得到列号j对应的像元入射角αi,j。

2.2 Radarsat-2数据辐射定标

2007年发射的Radarsat2[7]是Radarsat1的后续卫星，它引领了星载SAR的多项技术革命，如第一次出现了全极化模式。Radarsat2的L-1 DGRP数据以Product.xml文件组织，包含TIFF格式的影像数据和XML格式的增益偏移参数文件，其绝对辐射定标过程可以分为获取雷达截面积影像和计算入射角。

1)获取雷达截面积影像。Radarsat2影像雷达截面积σ为

σ=(DN2+A)/G ,

(5)

式中：DN为雷达影像像元值；A为增益偏移参数；G为第j像元(即第j列)的增益。其中偏移和增益参数可以在lutSigma.xml文件中获取。

2)计算像元的入射角。雷达入射角与地球半径和轨道高度有关，其几何关系为

(6)

式中：h为轨道高度；R为地球半径；s为斜距；Sj为图像中第j列的地距。对于L-1DGRP，在升轨右视和降轨左视情况下(根据Product.xml文件中的passDirection和antennaPointing值判断)，Sj的计算公式为

Sj=a+jLb+(jL)2c+(jL)3d+(jL)4e+(jL)5f ;

(7)

在降轨右视和升轨左视情况下，Sj的计算公式为

Sj=a+(J-j)Lb+[(J-j)L]2c+[(J-j)L]3d+[(J-j)L]4e+[(J-j)L]5f ,

(8)

式中：L为像元间隔；J为影像总列数；a,b,c,d,e,f分别为地距到斜距的6个转换参数，可从Product.xml文件中的groundToSlantRangeCoefficients获取。此外，Sj也可以通过最小入射角和最大入射角插值得到。

综合以上过程，可得到雷达后向散射系数图像σ0(以db为单位表示)，即

(9)

2.3 Cosmoskymed数据辐射定标

Cosmoskymed[8]是意大利研发的COSMO-skymed高分辨率雷达卫星星座的代称，该星座共有4颗卫星，重访能力高，具备全球范围观测能力。Cosmoskymed的L-1 DGRP数据(有DGM和GEC两种)通常以HDF5格式(“.H5”)分发。其绝对辐射定标过程为

(10)

式中： |imginp(i,j)|2为计算强度影像表达式；αref为参考入射角；F为尺度因子；K为定标常数。以db为单位表示，则

σ0(i,j)db=10lg[σ0(i,j)] ,

(11)

以上参数均可在HDF5文件的元数据中找到，其中αref取自Reference Slant Range；F取自Rescaling Factor；K取自Calibration Constant。

为了减少系统噪声的影响，通常可以采用一个窗口平均值替代，即

(12)

式中：r,c分别为窗口的行和列；wr,wc分别为窗口的行宽和列宽。

2.4 TerraSAR-X数据辐射定标

TerraSAR-X[9]是一颗成像分辨率和轨道精度都非常高的SAR卫星，在大范围对地观测和干涉测量领域具有明显优势。TerraSAR-X的L-1 DGRP数据也有2种： 多视地距探测产品(multi look ground range detected,MGD)和地理编码椭球校正产品(geocoding of terrain correction,GEC)，通常以“.XML”格式组织，包含TIFF影像数据、快视数据和元数据等。TerraSAR-X的L-1 DGRP数据为绝对辐射定标处理提供了定标常量、等效噪声分布和本地入射角参数，具体过程如式(13)―(16)所示，即

σ0=(ks|DN|2-B)sinθloc,

(13)

(14)

(15)

(16)

式(13)―(16)中：ks为辐射定标常量，对应XML格式元数据文件中calibrationConstant节点下calFactor的参数值；DN为像元的后向散射强度值，影像数据位于IMAGEDATA文件夹；B为强度等效噪声，反映不同噪声分布模型对信号的影响，所涉及的参数均可在元数据文件中获取(由于B的影响很小，在有些情况下该参数可以忽略)；θloc为本地入射角，其计算所需的参数GGIM可以在georef.xml文件中获取；deg和coffi分别为多项式拟合的阶数和系数；τ为快时间；τref为参考快时间；τmin和τmax分别为快时间的最小值和最大值。

2.5 Sentinel1数据辐射定标

欧空局于2014年发射的新一代SAR卫星Sentinel1[11](亦称“哨兵1号”)是ENVISAT卫星的后继卫星，继续担负着对地环境观测的任务，其数据可以在欧空局官方网站免费申请获得。Sentinel1的L-1 DGRP数据代码为GRD，以“manifest.safe”索引文件组织，包含影像数据(measurement文件夹)、快视数据(preview文件夹)、轨道及标定参数数据(annotation文件夹)以及说明数据(support文件夹)。对于Sentinel1的L-1 DGRP数据的绝对辐射定标，数据产品本身提供了一个定标矢量，将影像的强度值转化为后向散射系数，具体转换过程为

(17)

式中：DN为TIFF格式地距影像中对应像元的像素值；Ssigma为定标参数，该参数可通过“annotation”文件夹下的“Calibration”文件夹中的XML元数据文件sigmaNought域得到。元数据中Ssigma是一个查找表(LUT)，为行方向和列方向具有一定间隔的矢量数组； 对于任意像元的Ssigma，需要通过插值得到，插值方法往往选择双线性插值方法。

3 实验与分析

本文以Sentinel1的L-1 DGRP数据为例，采用VC++编程语言和GDAL栅格数据读写库，验证了绝对辐射定标方法。实验数据为2015年3月15日获取的覆盖环渤海地区Sentinel1的L-1 DGRP，影像大小为19 342行×25 267列，空间分辨率为5 m×20 m，干涉宽幅扫描(IW)成像模式，扫描带宽为250 km。根据2.5节介绍的方法对Sentinel1的L-1 DGRP数据进行绝对辐射定标。首先读取影像的像元值，Sentinel1影像采用16 bit记录强度信息，因此像元值在0～65 535之间； 然后获取每个像元定标参数值S，由于数据产品提供了30行×633列的定标参数矩阵以及对应的行列坐标矩阵，因此采用双线性内插方法计算其他行列坐标像元的定标参数值； 最后根据式(17)计算得到每个像元的后向散射系数值，即得到辐射定标处理结果。

图1(a)为用本文方法的软件实现的Sentinel1的L-1 DGRP数据辐射定标结果，图1(b)为用欧空局的S1-ToolBox软件对相同数据进行辐射定标得到的结果，两者在视觉上差别甚小。

(a) 本文方法 (b) 欧空局S1-ToolBox软件

图1 2种不同方法软件辐射标定结果

Fig.1 Calibration results by using two kinds of software from different methods

为了进一步分析本文方法绝对辐射定标处理结果的精度，从以上2个辐射定标结果中随机选取了12组样本点的后向散射系数值进行对比，误差分析结果如表3所示。

表3 2种方法标定精度对比

可以看出，2种方法处理结果相对误差的最大值为0.003，最小值为0。考虑到计算机插值处理精度的问题，以及忽略有效数字的差别，本文方法的软件处理得到的结果与欧空局S1-ToolBox软件处理的结果影像是基本相同的，从而证实了本文所介绍方法的正确性。

4 结论

辐射定标是SAR影像地表参数定量反演中的关键步骤，本文在明确几种新型SAR一级地距产品(L-1 DGRP)数据组织结构特点的基础上，结合辐射定标基本原理，阐述了几种新型SAR的L-1 DGRP辐射定标方法，结论如下：

1) 针对ENVISAT ASAR,Radarsat2,Cosmoskymed,TerraSAR-X和Sentinel1的L-1 DGRP数据，发展了利用产品的定标参数和入射角插值方法、进而利用辐射定标公式进行后向散射系数计算。

2) 以Sentinel1的L-1 DGRP数据为例，通过对2种方法的软件处理结果对比研究，验证了本文所述算法的正确性。

3)下一步的研究将充分考虑地形因子对新型SAR的L-1 DGRP数据辐射定标的影响，并发展更为严密的辐射定标模型。

志谢： 感谢欧洲空间局提供ESA S1 ToolBox软件。

[1] Desnos Y L,Li Z Y,Gao Z H,et al.The dragon programme-status and achievements[C]//Proceedings of the Envisat Symposium 2007.Montreux,Switzerland:ESA,2007.

[2] Laur H,Bally P,Meadows P,et al.ERS SAR Calibration.Derivation of the Backscattering Coefficientσ0in ESA ERS SAR PRI Products[R].ESA/ESRIN ES-TN-RS-PM-HL09.2009:3-6.

[3] Cumming I G,Wong F H.Digital Processing of Synthetic Aperture Radar Data:Algorithms and Implementation[M].Boston,MA:Artech House Remote Sensing Library,2005:56-123.

[4] 袁孝康.星载合成孔径雷达导论[M].北京:国防工业出版社,2003. Yuan X K.Introduction to Synthetic Aperture Radar[M].Beijing:National Defense Industry Press,2003.

[5] Schwerdt M,Brautigam B,Bachmann M,et al.TerraSAR-X calibration-first results[C]//Proceedings of the 2007 IEEE International Geoscience and Remote Sensing Symposium.Barcelona:IEEE,2007:3932-3935.

[6] Rosich B,Meadows P.Absolute Calibration of ASAR Level 1 Products Generated with PF-ASAR[R].ESA-ESRIN,ENVI-CLVL-EOPG-TN-03-0010,2004:5-8.

[7] Slade B.RADARSAT-2 Product Description[R].RN-SP-52-1238,2009:2-25.

[8] Agenzia Spaziale Italiana.COSMO-SkyMed SAR Products Handbook[M].Agenzia:ASI,2007:399-405.

[9] Fritz T,Werninghaus R.TerraSAR-X Ground Segment Level 1b Product Format Specification[R].Berlin:Clustert Applied Remote Sensing(CAF),German Aerospace Center(DLR),2007:156-170.

[10]Rosenqvist A,Shimada M,Suzuki S,et al.Operational performance of the ALOS global systematic acquisition strategy and observation plans for ALOS-2 PALSAR-2[J].Remote Sensing of Environment,2014,155:3-12.

[11]Geudtner D,Torres R,Snoeij P,et al.Sentinel1 system capabilities and applications[C]//Proceedings of the 2014 IEEE Geoscience and Remote Sensing Symposium.Quebec City,QC:IEEE,2014:1457-1460.

(责任编辑： 邢宇)

Absolute radiometric calibration of level-1 detected ground range products of new SAR sensors

DU Weina1， XU Aigong1， SONG Yaoxin2， SUN Huasheng1

(1.LiaoningTechnicalUniversity，Fuxin123000，China； 2.TrafficManagementResearchInstituteofMinistryofPublicSecurity，Wuxi214151,China)

For the current situation of the lack of the new SAR sensor data preprocessing software, this paper introduced in detail the methods of absolute radiometric calibration and the parameter acquisition for several new SAR sensor level-1 detected products, such as ENVISAT ASAR，Radarsat2，Cosmoskymed，TerraSAR-X and Sentinel1. In addition, the absolute radiometric calibration process was achieved by programming with the level-1 detected ground range products(L-1 DGRP) data of Sentinel1 sensor, and C++ programming language was used to achieve the absolute radiation of the calibration process. At last, the radiometric calibration results produced by the method developed in this paper and implemented in the authors’ software were compared with those by ESA S1 ToolBox, the freely distributed SAR data processing tool by European Space Agency, and it is shown that the two numerical back scattering systems are basically the same. The radiometric calibration method developed in this paper is proved to be correct by the program implementation.

new SAR sensors； absolute radiometric calibration； level-1 detected ground range products(L-1 DGRP)； back scattering coefficient； Sentinel1

10.6046/gtzyyg.2016.04.05

杜伟娜,徐爱功,宋耀鑫,等.新型SAR传感器一级地距产品绝对辐射定标方法[J].国土资源遥感,2016,28(4)：30-34.(Du W N,Xu A G,Song Y X,et al.Absolute radiometric calibration of level-1 detected ground range products of new SAR sensors[J].Remote Sensing for Land and Resources,2016,28(4)：30-34.)

2015-05-13；

2015-07-08

国家自然科学基金项目“遥感数据的空间分辨率和波段数对土地覆盖制图的影响研究”(编号： 41201454)资助。

TP 751.1

A

1001-070X(2016)04-0030-05

杜伟娜(1988-),女,硕士研究生,主要研究方向为SAR图像处理与应用。Email： RSwendu@163.com。