表观反射率(反射率、反照率)的计算

第一步、分别计算各个波段每个像元的辐射亮度L值：

L=Gain*DN+Bias

或者

LLmaxLmin*(QCALQCALmin)Lmin

QCALmaxQCALmin式中，QcaL为某一像元的DN值，即QCAL=DN。 QCALmax为像元可以取的最大值255。QCALmin为像元可以取的最小值。如果卫星数据来自LPGS(The level 1 product generation system)，则QCAL=1(Landsat-7数据属于此类型)。如果卫星数据来自美国的NLAPS ( National Landsat Archive Production System )，则QCALmin=0 (Ldsat-5的TM数据属于此类型)。

根据以上情况，对于Landsat-7来说，可以改写为(QCALmin=1)：

LLmaxLmin*(DN1)Lmin254

对于Landsat-5来说，可以改写为(QCALmin=0)：

L

LmaxLmin*DNLmin255

表1 Iandsa-7 ETM+各个反射波段的Lmax和Lmin值

Table1The values of Lmmax and Lmin for reflecting bands ETM+(W˙m-2-sr-1˙μm-1) 波段 Band 1 2 3 4 5 7

-6.2 -6.0 -4.5 -4.5 -1.0 -0.35 2000年7月1日之前 低Gain Lmin Lmax 297.5 303.4 235.5 235.5 47.7 16.6 -6.2 -6.0 -4.5 -4.5 -1.0 -0.35 高Gain Lmin Lmax 194.3 202.4 158.6 157.5 31.76 10.932 -6.2 -6.4 -5.0 -5.1 -1.0 -0.35 Lmin 低Gain Lmax 293.7 300.9 234.4 241.1 47.57 16.54 -6.2 -6.4 -5.0 -5.1 -1.0 of Landsat-7

2000年7月1日之后 高Gain Lmin Lmax 191.6 196.5 152.9 157.4 31.06 10.8 -0.35 表2 Landsat-5 TM各反射波段的Lmax和Lmin值

Table 2 The values of Lmax and Lmin for reflecting bands of Landsat-5 TM (W˙m-2-sr-1˙μm-1) 波段 1984/03/01至2003/05/04 2003/05/04之后 Band 1 2 3 4 5 7 Lmin -1.52 -2.84 -1.17 -1.51 -0.37 -0.15 Lmax 152.10 296.81 204.30 206.20 27.19 14.38 Lmin -1.52 -2.84 -1.17 -1.51 -0.37 -0.15 Lmax 193.0 365.0 264.0 221.0 30.2 16.5 为了使传感器的辐射分辨率达到最大，而又不使其达到饱和，根据地表类型(非沙漠和冰面的陆地、沙漠、冰与雪、水体、海冰、火山等6大类型)和太阳高度角状况来确定采用高增益参数或是低增益参数。一般低增益的动态范围比高增益大1.5倍，因此当地表亮度较大时，用低增益参数;其它情况用高增益参数。在非沙漠和冰面的陆地地表类型中，ETM+的1一3和5,7波段采用高增益参数，4波段在太阳高度角低于45度（天顶角>45度）时也用高增益参数，反之则用低增益参数。详见文献(NASA Landsat Project ScienceOffice , 1998b )。

第二步、计算各波段反射率（反照率、反射率）ρ：

iLD2ESUNCos()(i为第i波段）

式中，p为人气层顶(TOA)表观反射率(无量纲)，π为常量(球面度str)，L为大气层顶进人卫星传感器的光谱辐射亮度(W˙m-2-sr-1˙μm-1)，D为日地之间距离(天文单位),ESUN为大气层顶的平均太阳光谱辐照度(W˙m-2-sr-1˙μm-1),θ为太阳的天顶角(θ=90˚-β，β为太阳高度角， Cos(θ)也可以这样计算：Cos(θ)=Sinφ*Sinδ+Cosφ*Cosδ*Cosh,式中φ甲为地理纬度，φ为太阳赤纬，h为太阳的时角。太阳赤纬是太阳光与地球赤道平面的夹角)。

也可以是：

TL()2(D93.5)2(10.0167sin)Esun()coss365

其中,θs 为太阳天顶角, D 为儒略历(Julian) 日期,这两个参数可由数据头文件读

出。L (λ) 为入瞳辐亮度, Esun为外大气层太阳辐照度。

上式成立的条件是假设在大气层顶，有一个朗勃特(Laribcitian)反射面。太阳光以天顶角θ人射到该面，该表面的辐照度为E = ESUN*Cos(θ)/D2(吕斯哗，1981)。该表面的辐射出射度M=πL(吕斯骤，1981)。根据Lanbertian反射率定义，大气层顶的表观反射率P等于M和E的比值，即

MLD2iEESUNCos()(i为第i波段）

表 3 随时间变化的日地距离（天文单位）

Table 3 Earth-Sun distance at different time (Astonomical units) 日数 day 1 15 32 46 距离 day 0.9832 0.9836 0.9853 0.9878 日数 day 74 91 106 121 日数 day 0.9945 0.9993 1.0033 1.0076 日数 day 152 166 182 196 日数 day 1.0140 1.0158 1.0167 1.0165 日数 day 227 242 258 274 日数 day 1.0128 1.0092 1.0057 1.0011 日数 day 305 319 335 349 日数 day 0.9925 0.9892 0.9860 0.9843 60 0.9909 135 1.0109 213 1.0149 288 0.9972 365 0.9830 表 4 Landsat-7 和Landsat-5的大气层顶平均太阳光谱辐照度ESUN(W˙m-2-sr-1˙μm-1) TahlP 4 Mean solar spectral iwadiance at the atmosphemic top for Landsat-7 and Landsat-5 波段Band Landsat-7 ESUN Landsat-5 ESUN 1 1969 1957 2 1840 1826 3 1551 1554 4 1044 1036 5 225.7 215 7 82.07 80.67 两步合为一步计算如下：

iLmaxLmin(QCALQCAL)LminminESUNCos()QCALmaxQCALmin对于Landsat-7上试简化为：

D2(i为第i波段）iLmaxLmin(QCAL1)Lmin

ESUNCos()254LmaxLmin QCALLminESUNCos()255D2对于Landsat-5上试简化为：

iD2 其中，QCAL为图像灰度值DN。

反照率的计算：

TM1～TM4波段所对应的宽波段反照率可表示为

i14 (i为TM第i个波段的反射率）

Table 1. Characteristics of the Enhanced Thematic Mapper Plus (ETM+)

bands. Spatial Lower Upper Bits

Bandwidth

Band resolution limit limit per

(nm)

(m) (µm) (µm) pixel 1 2 3 4 5 6 7 8

28.50 28.50 28.50 28.50 28.50 28.50 14.25

0.45 0.53 0.63 0.75 1.55 2.10 0.52

0.52 0.61 0.69 0.90 1.75 2.35 0.90

Gain Offset

70 80 60 150 200 2100 250 380

8 0.786274521 -6.1999998 8 0.817254878 -6.0000000 8 0.639607867 -4.5000000 8 0.939215686 -4.5000000 8 0.128470589 -1.0000000 8 0.066823533 0.00000000 8 0.044243138 -0.3499999 8 0.786274521 -6.1999998

57.00 10.40 12.50

11.3.1 Conversion to Radiance

During 1G product rendering image pixels are converted to units of absolute radiance using 32 bit floating point calculations. Pixel values are then scaled to byte values prior to media output. The following equation is used to convert DN's in a 1G product back to radiance units:

Lλ = \"gain\" * QCAL + \"offset\"

which is also expressed as:

Lλ = ((LMAXλ - LMINλ)/(QCALMAX-QCALMIN)) * (QCAL-QCALMIN) + LMINλ where:

Lλ

= Spectral Radiance at the sensorճ aperture in

watts/(meter squared * ster * μm)

\"gain\" = Rescaled gain (the data product \"gain\" contained in

the Level 1 product header or ancillary data record)

in watts/(meter squared * ster * μm)

\"offset\" = Rescaled bias (the data product \"offset\" contained

in the Level 1 product header or ancillary data

record ) in watts/(meter squared * ster * μm)

QCAL = the quantized calibrated pixel value in DN LMINλ = the spectral radiance that is scaled to QCALMIN in

watts/(meter squared * ster * μm)

LMAXλ = the spectral radiance that is scaled to QCALMAX in

watts/(meter squared * ster * μm)

QCALMIN = the minimum quantized calibrated pixel value

(corresponding to LMINλ) in DN = 1 (LPGS Products) = 0 (NLAPS Products)

QCALMAX = the maximum quantized calibrated pixel value

(corresponding to LMAXλ) in DN = 255

The LMINs and LMAXs are the spectral radiances for each band at digital numbers 0 or 1 and 255 (i.e QCALMIN, QCALMAX), respectively. LPGS used 1 for QCALMIN while NLAPS used 0 for QCALMIN for data products processed before April 5, 2004. NLAPS from that date now uses 1 for the QCALMIN value. Other product differences exist as well. One LMIN/LMAX set exists for each gain state. These values will change slowly over time as the ETM+ detectors lose responsivity. Table 11.2 lists two sets of LMINs and LMAXs. The first set should be used for both LPGS and NLAPS 1G products created before July 1, 2000 and the second set for 1G products created after July 1, 2000. Please note the distinction between acquisition and processing dates. Use of the appropriate LMINs and LMAXs will ensure accurate conversion to radiance units. Note for band 6: A bias was found in the pre-launch calibration by a team of independent investigators post launch. This was corrected for in the LPGS processing system beginning Dec 20, 2000. For data processed before this, the image radiances given by the above transform are 0.31 w/m2 ster um too high. See the official announcement for more details.

Table 11.2 ETM+ Spectral Radiance Range watts/(meter squared * ster * μm) Before July 1, 2000 Low Gain Band Number LMIN LMAX 1 2 3 4 High Gain LMIN LMAX After July 1, 2000 Low Gain High Gain LMIN LMAX LMIN LMAX -6.2 297.5 -6.2 194.3 -6.2 293.7 -6.2 191.6 -6.0 303.4 -6.0 202.4 -6.4 300.9 -6.4 196.5 -4.5 235.5 -4.5 158.6 -5.0 234.4 -5.0 152.9 -4.5 235.0 -4.5 157.5 -5.1 241.1 -5.1 157.4 5 6 7 8

-1.0 47.70 -1.0 31.76 -1.0 47.57 -1.0 31.06 0.0 17.04 3.2 12.65 0.0 17.04 3.2 12.65 -0.35 16.60 -0.35 10.932 -0.35 16.54 -0.35 10.80 -5.0 244.00 -5.0 158.40 -4.7 243.1 -4.7 158.3 11.3.2 Radiance to Reflectance

For relatively clear Landsat scenes, a reduction in between-scene

variability can be achieved through a normalization for solar irradiance by converting spectral radiance, as calculated above, to planetary reflectance or albedo. This combined surface and atmospheric reflectance of the Earth is computed with the following formula:

Where:

= Unitless planetary reflectance

= Spectral radiance at the sensor's aperture = Earth-Sun distance in astronomical units from nautical handbook or

interpolated from values listed in Table 11.4

= Mean solar exoatmospheric irradiances from Table 11.3

= Solar zenith angle in degrees

Table 11.3 ETM+ Solar Spectral Irradiances Band 1 2 watts/(meter squared * μm) 1969.000 1840.000 3 4 5 7 8 1551.000 1044.000 225.700 82.07 1368.000

Table 11.4 Earth-Sun Distance in Astronomical Units Julian Julian Julian Julian Julian Distance Distance Distance Distance Distance Day Day Day Day Day 1 15 32 46 60 .9832 .9836 .9853 .9878 .9909 74 91 106 121 135 .9945 .9993 1.0033 1.0076 1.0109 152 166 182 196 213 1.0140 1.0158 1.0167 1.0165 1.0149

11.3.3 Band 6 Conversion to Temperature

ETM+ Band 6 imagery can also be converted from spectral radiance (as

described above) to a more physically useful variable. This is the effective at-satellite temperatures of the viewed Earth-atmosphere system under an assumption of unity emmissivity and using pre-launch calibration constants listed in Table 11.5. The conversion formula is:

227 242 258 274 288 1.0128 1.0092 1.0057 1.0011 .9972 305 319 335 349 365 .9925 .9892 .9860 .9843 .9833

Where: T

= Effective at-satellite temperature in Kelvin

K2 = Calibration constant 2 from Table 11.5

K1 = Calibration constant 1 from Table 11.5 L

= Spectral radiance in watts/(meter squared * ster * ?m)

Table 11.5 ETM+ and TM Thermal Band Calibration Constants Landsat 7 Landsat 5 Constant 1- K1 watts/(meter squared * ster * μm) 666.09 607.76

Constant 2 - K2 Kelvin 1282.71 1260.56