表觀反射率(反射率、反照率)的計算

1、表觀反射率（反射率、反照率）的計算第一步、分別計算各個波段每個像元的輻射亮度L 值：L=Gain*DN+Bias或者LmaxLmin* (QCAL QCALmin ) L minLQCALminQCALmax式中， QcaL 為某一像元的DN 值，即 QCAL=DN 。 QCALmax 為像元可以取的最大值255 。 QCALmin 為像元可以取的最小值。如果衛(wèi)星數(shù)據(jù)來自LPGS(The level 1 productgeneration system)，則 QCAL=1(Landsat-7數(shù)據(jù)屬于此類型)。如果衛(wèi)星數(shù)據(jù)來自美國的NLAPS ( National Landsat Archive

2、 Production System )，則 QCALmin=0 (Ldsat-5的 TM 數(shù)據(jù)屬于此類型)。根 據(jù) 以 上 情 況 ， 對 于Landsat-7來 說 ， 可 以 改 寫 為 (QCALmin=1)：LLmaxLmin * (DN1) Lmin254對于 Landsat-5 來說，可以改寫為 (QCALmin=0)：LLmaxLmin * DNL min255表 1Iandsa-7 ETM+ 各個反射波段的Lmax 和 Lmin 值Table1ThevaluesofLmmax andLminforreflecting bands ofLandsat-7ETM+(W m-2-s

3、r-1 m-1)波段2000年 7月 1日之前2000年 7月 1日之后Band低 Gain高 Gain低 Gain高 GainLminLmaxLminLmaxLminLmaxLminLmax1-6.2297.5-6.2194.3-6.2293.7-6.2191.62-6.0303.4-6.0202.4-6.4300.9-6.4196.53-4.5235.5-4.5158.6-5.0234.4-5.0152.94-4.5235.5-4.5157.5-5.1241.1-5.1157.45-1.047.7-1.031.76-1.047.57-1.031.067-0.3516.6-0.3510.93

4、2-0.3516.54-0.3510.8表 2 Landsat-5 TM 各反射波段的LmaxTable 2 The values of Lmax and Lmin for reflecting bands of波段1984/03/01 至 2003/05/04和 Lmin 值Landsat-5 TM (W m-2-sr-1 m-1)2003/05/04 之后BandLminLmaxLminLmax1-1.52152.10-1.52193.02-2.84296.81-2.84365.03-1.17204.30-1.17264.04-1.51206.20-1.51221.05-0.3727.19

5、-0.3730.27-0.1514.38-0.1516.5為了使傳感器的輻射分辨率達到最大，而又不使其達到飽和，根據(jù)地表類型( 非沙漠和冰面的陸地、沙漠、冰與雪、水體、海冰、火山等6 大類型 )和太陽高度角狀況來確定采用高增益參數(shù)或是低增益參數(shù)。一般低增益的動態(tài)范圍比高增益大1.5 倍，因此當(dāng)?shù)乇砹炼容^大時，用低增益參數(shù) ; 其它情況用高增益參數(shù)。在非沙漠和冰面的陸地地表類型中，ETM+ 的 1 一 3和 5,7 波段采用高增益參數(shù)，4 波段在太陽高度角低于45 度（天頂角 45 度）時也用高增益參數(shù)，反之則用低增益參數(shù)。詳見文獻(NASA Landsat Project ScienceOff

6、ice , 1998b ) 。第二步、計算各波段反射率（反照率、反射率） ：LD2iESUN(i為第 i波段）Cos( )式中， p 為人氣層頂 (TOA) 表觀反射率 (無量綱 )， 為常量 ( 球面度 str)，L 為大氣層頂進人衛(wèi)星傳感器的光譜輻射亮度 (W m-2-sr-1m-1) ，D 為日地之間距離 (天文單位 ),ESUN 為大氣層頂?shù)钠骄柟庾V輻照度 (W m-2-sr-1m-1), 為太陽的天頂角 ( =90?- ， 為太陽高度角， Cos( )也可以這樣計算： Cos( )=Sin *Sin +Cos *Cos *Cosh, 式中 甲為地理緯度， 為太陽赤緯， h 為太陽

7、的時角。太陽赤緯是太陽光與地球赤道平面的夾角)。也可以是：TE sun(L()(1 0.0167sin 2 (D 93.5) )2) cos s365其中 , s 為太陽天頂角,D 為儒略歷 (Julian)日期 , 這兩個參數(shù)可由數(shù)據(jù)頭文件讀出。 L ( ) 為入瞳輻亮度 ,Esun為外大氣層太陽輻照度。上式成立的條件是假設(shè)在大氣層頂，有一個朗勃特(Laribcitian) 反射面。太陽光以天頂角 人射到該面，該表面的輻照度為E = ESUN*Cos( )/D 2(呂斯嘩， 1981)。該表面的輻射出射度 M= L( 呂斯驟， 1981)。根據(jù) Lanbertian 反射率定義，大氣層頂?shù)谋?/p>

8、觀反射率P 等于M 和 E 的比值，即iMLD 2(i 為第 i 波段）EESUNCos( )表 3 隨時間變化的日地距離（天文單位）Table 3 Earth-Sun distance at different time (Astonomical units)日數(shù)距離日數(shù)日數(shù)日數(shù)日數(shù)日數(shù)日數(shù)日數(shù)日數(shù)daydaydaydaydaydaydaydaydayday10.9832740.99451521.01402271.01283050.9925150.9836910.99931661.01582421.00923190.9892320.98531061.00331821.01672581.005

9、73350.9860460.98781211.00761961.01652741.00113490.9843600.99091351.01092131.01492880.99723650.9830表 4 Landsat-7和 Landsat-5 的大氣層頂平均太陽光譜輻照度ESUN(W m-2-sr-1m-1 )TahlP 4 Mean solar spectral iwadiance at the atmosphemic top for Landsat-7 and Landsat-5波段 Band123457Landsat-7 ESUN1969184015511044225.782.07La

10、ndsat-5 ESUN195718261554103621580.67兩步合為一步計算如下：iD 2Lma xLmin(QCALQCAL min )Lmin(i為第 i波段）ESUNCos( )QCAL maxQCAL min對于 Landsat-7 上試簡化為：iD 2LmaxL min(QCAL1)LminESUN Cos( )254對于 Landsat-5上試簡化為：iD 2LmaxLminQCALLminESUNCos( )255其中， QCAL 為圖像灰度值DN。反照率的計算：TM1 TM4波段所對應(yīng)的寬波段反照率可表示為4i(i 為 TM 第 i個波段的反射率）1Table 1.

11、 Characteristics of the Enhanced Thematic Mapper Plus (ETM+)bands.BandSpatialLowerUpperBandwidthBitsGainOffsetresolutionlimitlimit(nm)per(m)( m) ( m)pixel128.500.450.52708 0.786274521-6.1999998228.500.530.61808 0.817254878-6.0000000328.500.630.69608 0.639607867-4.5000000428.500.750.901508 0.93921568

12、6-4.5000000528.501.551.752008 0.128470589-1.0000000657.0010.4012.5021008 0.0668235330.00000000728.502.102.352508 0.044243138-0.3499999814.250.520.903808 0.786274521-6.1999998During 1G product rendering image pixels are converted to units ofabsolute radiance using 32 bit floating point calculations.

13、Pixel values are then scaled to byte values prior to media output. The followingequation is used to convert DNs in a 1G product back to radiance units:L = gain * QCAL + offsetwhich is also expressed as:L = (LMAX - LMIN ) /(QCALMAX-QCALMIN) * (QCAL-QCALMIN)+ LMINwhere:L? aperture in= Spectral Radianc

14、e at the sensorwatts/(meter squared * ster *m)gain= Rescaled gain (thedata productgaincontainedinthe Level 1 product header or ancillary data record)inwatts/(meter squared * ster *m)offset= Rescaled bias (the data product offset containedin the Level 1 product header or ancillary datarecord ) inwatt

15、s/(meter squared * ster *m)QCAL = the quantized calibrated pixel value in DNLMIN= the spectral radiance that is scaled to QCALMIN inwatts/(meter squared * ster *m)LMAX= the spectral radiance that is scaled to QCALMAX inwatts/(meter squared * ster *m)QCALMIN = the minimum quantized calibrated pixel v

16、alue(corresponding to LMIN) in DN1 (LPGS Products)0 (NLAPS Products)QCALMAX= the maximum quantized calibrated pixel value (corresponding to LMAX ) in DN= 255The LMINs and LMAXs are the spectral radiances for each band at digitalnumbers 0 or 1 and 255 (i.e QCALMIN, QCALMAX), respectively. LPGS used1

17、for QCALMIN while NLAPS used0 for QCALMIN for data products processedbefore April5, 2004. NLAPSfrom thatdate nowuses 1 for theQCALMINvalue.Other product differencesexist as well. One LMIN/LMAXsetexistsforeachgain state.These valueswill change slowly over timeas theETM+detectorsloseresponsivity.Table

18、11.2 liststwo sets ofLMINs and LMAXs.The firstsetshouldbe used for both LPGSand NLAPS1G products created beforeJuly1, 2000 and the second set for 1G products createdafterJuly 1, 2000.Please note the distinctionbetween acquisitionand processing dates.Useof the appropriate LMINs and LMAXs will ensure

19、accurate conversion toradiance units.Note for band 6:A bias was found in the pre-launchcalibrationby a team ofindependent investigatorspost launch.Thiswascorrected for in the LPGS processing system beginning Dec 20, 2000. Fordataprocessed before this, the image radiances given by the abovetransform

20、are 0.31 w/m2 ster um too high. See the official announcementfor more details.Table 11.2 ETM+ Spectral Radiance Rangewatts/(meter squared * ster *m)Before July 1, 2000After July 1, 2000BandLow GainHigh GainLow GainHigh GainNumberLMINLMAXLMINLMAXLMINLMAXLMINLMAX1-6.2297.5-6.2194.3-6.2293.7-6.2191.62-

21、6.0303.4-6.0202.4-6.4300.9-6.4196.53-4.5235.5-4.5158.6-5.0234.4-5.0152.94-4.5235.0-4.5157.5-5.1241.1-5.1157.45-1.047.70-1.031.76-1.047.57-1.031.0660.017.043.212.650.017.043.212.657-0.3516.60-0.3510.932-0.3516.54-0.3510.808-5.0244.00-5.0158.40-4.7243.1-4.7158.3For relatively clear Landsat scenes, a r

22、eduction in between-scenevariabilitycan be achieved through a normalizationforsolarirradianceby converting spectral radiance, as calculated above, to planetaryreflectanceor albedo.This combined surfaceand atmosphericreflectanceof the Earth is computed with the following formula:Where:Unitless planet

23、ary reflectanceSpectral radiance at the sensors apertureEarth-Sun distance in astronomical units from nautical handbook orinterpolated from values listed in Table11.4Mean solar exoatmospheric irradiances from Table 11.3Solar zenith angle in degreesTable 11.3ETM+ Solar Spectral IrradiancesBandwatts/(

24、meter squared *m)11969.00021840.00031551.00041044.0005225.700782.0781368.000Table 11.4 Earth-Sun Distance in Astronomical UnitsJulianDistanceJulianDistanceJulianDistanceJulianDistanceJulianDistanceDayDayDayDayDay1.983274.99451521.01402271.0128305.992515.983691.99931661.01582421.0092319.989232.98531061.0033182