(5) 
(6) 
> remotely 
| + coss)]} e. 
7.1) (7) 
jn method 
ecular and 
al. 1996). 
S, such as 
by firstly 
988) with 
xponential 
target is 
n as: 
1+00s5)} e 
(u,)- 8) 
liance L is 
led reflief 
etermined 
M). Si 
and solar 
cal depths 
th known 
vell as the 
bedo ®, 
which are 
previously 
al surface 
the Junge 
, 1983) in 
it can be 
JV by the 
Liu e: al. 
ve aerosol 
aw image 
with the 
j| type is 
ges from 
[0.01 » 10]um. Then @, and P,(®) can be determined from the 
Mie theory. 
Here we propose the robust algorithm to retrieve aerosol optical 
depth from the image itself. One can select pixels of a given 
canopy with different altitude z, slope s and incidence angle / 
from the ancillary data, such as airphoto, base map and DTM. 
The merit function is defined as the ratio of standard deviation 
to average of derived surface reflectances of those selected 
pixels. By minimizing the merit function such as Powell 
minimization method (Press et al. 1992), optimized aerosol 
optical depth can be determined. The algorithm is so-called 
“robust” because of its robustness. Certainly, the algorithm is 
limited by the a-priori knowledge of the geographic coordinates 
of the selected canopy. The dynamic range of the selected 
pixels also limits the application of the robust algorithm. 
The altitude-dependent aerosol optical depth is modeled as 
747 Tao-Co Z, Where Tao, Co are the sea-level optical depth and the 
linear decay rate (km). In this study, the variance of surface 
reflectances of selected pixels in robust algorithm is too large, if 
one use exponential decay model of aerosol. It is also noted 
that the altitude-dependent aerosol optical depth shows linear 
decrease from simulation study of LOWTRAN 7 (Liu 1995). 
After knowing «, and P4(O), py and optical depth z, are the 
left unknowns in equation (8). The optimized z,, and c, can 
be retrieved by using the proposed robust algorithm. Penalty 
function suggested by Liang and Strahler (1994) is added in the 
merit function in order to force the retrieved z,, and c, to fall in 
the reasonable ranges which can be determined from 
LOWTRAN 7 with reasonable visibility ranges 5-336 km (Liu 
1995). After c, , P,(O) are assumed and the initial guesses of 
Tuo and c, are given, the four dimensional look-up table L —/( p, 
z, S, Sauj) is built up. Then, every p, of the selected pixels in 
the robust algorithm is determined from the interpolation of the 
above look-up table with known z, s, S,u; and L. The 
optimized z,, and c, are finally iteratively retrieved by 
minimizing the variation of reflectances of the selected pixels. 
The retrieved ones are so called "optimized" because those 
deviated values from the optimized ones would result in larger 
variation of the reflectances of the selected canopy (see section 
4) And then p, image of the topographic-effect corrected can 
be generated by the optimzed aerosol optical depth and the 
look-up table. 
By comparing equation (7) and equation (8), one can obtains 
the remotely sensed ground reflectance p which can be written 
as: 
| [1-053 1-axs]e*^ +04) xor, 
id [^ [1+05<o>(1—00oss)] pr, 
SN 
>) ©) 
where <p> is simply determined by averaging p, of the local 
window (Richter ,1990). The window size is selected to be 
about I km which is the same as the previous study (Liu et al. 
1996). If the terrrain is flat, i.e. s=0 » equation (9) will be 
certainly reduced to be equation (15) of Liu et al. (1996). 
3. DATA DESCRIPTION 
The test site is located in Taoyuan county, northern Taiwan and 
is about 7.5km X 7.5km. The altitudes range from 100~400m. 
The slopes range from 0°~49° with mean slope 12°. The area 
is covered mainly by hardwoods, acacia, tea, grassland, bare 
soil, high-reflective land and urban. 
A SPOT image taken on 30 December 1993 is used to testify 
the algorithm. Its apparent reflectance image is shown in Fig. 
la. The solar and viewing zenith angles are 51.9° and 15.8° 
respectively. The relative azimuth angle is 53°. 
DTM provided by Aerial Survey Institute with 40m spatial 
resolution is resampled to correspond with that of satellite 
image by using cubic convolution method (Kawata et al. 1988). 
DTM is firstly gaussian-filtered to remove much of noise 
(Seemuller 1989). Slope s and aspect a are then determined 
by the following equations with two partial derivatives 
computed by the third-order finite difference method (Skidmore 
1989): 
[02/0] (Zi 1 jr1)+2(Zis1 j)H (Zi 1j1)] 
  
-"[Giij 2G) j-1))/8AX, (10) 
[G2/0y].j*[(i1j1)*2(2ij)*(2134)] 
-[Gi 3) *2(2;j-1)* Gi j:)))/8Ay. (11) 
tans = [(@ / 2)? + (@ 1 @)*] (12) 
N (Oz ^x). (13) 
(Oz Gy) 
The incidence angle / of beam irradiance is computed as follows 
(Duguay and LeDrew 1992): 
diz cos '(cos9, coss + sin@, sin scos(a —@ E (14) 
where Q, ¢, are solar zenith and azimuth angles respectively. 
Ray tracing method is applied to determine cast-shadow 
function S; (Proy et al. 1989). Fig. 2 shows the shaded relief 
Sa; image determined from DTM data and corresponded to the 
solar geometry at the time of the SPOT overpass. 
4. RESULTS AND DISCUSSIONS 
To account for the altitude dependence of the molecular 
scattering and absorption, radiosonde data of the Panchiao city 
located from about 40km is used as inputs of LOWTRAN 7. 
The robust algorithm is applied to retrieve the optimized aerosol 
optical depth. In this study, pixels of hardwood are selected 
carefully by the base maps. To avoid mixed pixel. central 
pixel of homogeneous slope is selected. — The dynamic ranges 
107 
International Archives of Photogrammetry and Remote Sensing. Vol. XXXI, Part B7. Vienna 1996