Transcript Document

Retrieval of forest structural
characteristics from lidar waveforms:
new applications and methods
Ph.D. Dissertation Defense
by
Svetlana Y. Kotchenova
Dissertation Committee
Anthony B. Davis
Yuri Knyazikhin
Ranga B. Myneni
Nathan Phillips
Curtis Woodcock
October 4th, 2004
Geography Department, Boston University
Lidar remote sensing technique
Large-footprint waveform-recording laser altimeters:
- SLICER (air-borne)
8 km
- LVIS (air-borne)
- GLAS (space-borne)
Lidar waveform
Lidar = light detection and ranging
Data products:
- Canopy height
- Vertical canopy structure
(vertical distribution of nadirintercepted surfaces)
- Ground elevation
Lidar footprint
http://www.geog.umd.edu/vcl
2
SLICER (Scanning Lidar Imager of Canopies by Echo Recovery)
Data
Technical characteristics
Status
Platform
Wavelength
Pulse frequency
Pulse width
Pulse form
Footprint diameter
Transmit energy
previously-used
air-borne
1064 nm
80 Hz
4 ns
Raleigh
9-10 m
0.7 mJ
Footprint distribution pattern
Southern BOREAS study area, Saskatchewan,
Canada. July 20-30, 1996
Jack pine
Black spruce
Smithsonian Environmental Research Center.
Western shore of Chesapeake Bay, eastern
Maryland. September 7, 1995.
flight path
Mixed deciduous
forest stands with
the overstory
dominated by
tulip poplar
http://denali.gsfc.nasa.gov/research/laser/slicer/slicer.html
3
LVIS (Laser Vegetation Imaging Sensor)
(advanced version of SLICER)
Status
Platform
Wavelength
Pulse frequency
Pulse width
Pulse form
Footprint diameter
Transmit energy
currently-used
air-borne
1064 nm
500 Hz (320 Hz)
10 ns
Gaussian
1-80 m (25 m)
5 mJ
Footprint distribution pattern
flight path
. . .
. . . . . . . . . . . . . . . . . . . . . . .
. . .
27 m
Technical characteristics
. . .
9m
25 m
80 footprints
Data
La Selva Biological Research Station, Costa Rica. March,
1998. (A 1546-ha area comprised of primary (73%) and
secondary tropical rainforests and agroforestry plots.)
Harvard Forest, Massachusetts, USA. Summer 1999, summer
2003. Temperate deciduous broadleaf forest, dominant species
include white pine, hemlock, spruce, oak, and maple.
4
GLAS (Geoscience Laser Altimeter System)
GLAS will provide continuous global measurements of the Earth’s
land surface topography.
Technical characteristics
Status
Platform
Wavelength
Pulse frequency
Pulse width
Pulse form
Footprint diameter
Transmit energy
Along-track
separation
Cross-track max
Cross-track min
Repeat cycle
Life-time
launched
Jan. 2003
space-borne
1064 nm
(vegetation)
40 Hz
5 ns
Gaussian
60-70 m
5 mJ
170 m
15 km
2.5 km
183 days
3 years
http://icesat.gsfc.nasa.gov
5
Current applications of lidar measurements (1)
1. Retrieval of canopy height profiles (CHP),
or the distributions of canopy material with height
Assumptions:
1) Uniform horizontal distribution of leaves
2) Absence of multiple scattering
CHPc  ln(1  cover(h))
3) Negligible influence of non-foliated surfaces
6
Current applications of lidar measurements
2. Biomass estimation
Lidar
recorded
waveform
Canopy
height
&
CHP
(2)
Ground measurements
to relate height indices
and biomass
Calculation of
special indices
(mean, medium,
quadratic mean
height)
Regression
procedure
Biomass, basal area,
quadratic mean stem
diameter
3. Other methods
- Canopy volume method (vegetation is treated as a number of 3-dim pixels,
each of which is defined as containing canopy or not)
- Calculation of mean transmittance profiles (rate of attenuation, whole
canopy transmittance, radiation effective height)
7
Shortcomings of the existing methods
Retrieval of canopy height profiles
*
1.
Assumption of
uniform horizontal
distribution of leaves
*
2.
Assumption
of single
scattering
Foliage
clumping
is ignored
Erroneous assignment
of larger amounts of
foliage to the lower part
of the waveform
The method is not
applicable for coniferous
forests
Critical for
photosynthesis and
light transmission
studies
Biomass estimation
Requirement of extensive ground
measurements (DBH of each tree within the
footprint area, 25-30 footprint areas at least)
The method is
hardly applicable to
GLAS data
* Current applications
Numerous other applications are possible but have not been explored yet.
8
Research plan
* 1. Evaluation of the contribution of multiple scattering
* 2. Development of a new application for lidar data
* 3. Improvement of the current algorithm for retrieval of CHPs
4. Development of a new algorithm for biomass estimation
9
Task 1: Contribution of multiple scattering
a) Development of time-dependent stochastic radiative transfer theory to
describe the radiation regime in heterogeneous vegetation canopies
b) Creation of a physical model describing the propagation of lidar signals
through vegetation on the basis of this theory
c) Validation of the model with field data
d) Evaluation of multiple scattering contribution through the comparison of
reflected signals formed with and without multiple scattering
10
Time-dependent stochastic RT theory
 




 

 

 
1 I (t, r , ) 
   I (t, r , )  ( r )  ()I (t, r , ) ( r )   S (  ) I (t, r , ) d
c
t
4
Boundary conditions:

 



I  ( t, rxy ,0, )  f  ( t )  (  0 ) , rxy  Sf , (0 )  0
 
I  (0, r , )  0, 0  z  H







soil ()



I  ( t, rxy , H, ) 
  I  ( t, rxy , H,  )  ( ) d , ()  0

2 
0
Multiply scattered photons appear
delayed in the waveform compared
to singly scattered photons
Expected results:
H
z
Enhancement of the lower
part of the reflected signal
11
Description of canopy structure
Indicator function:

 1, if thereis a leaf at r
( r )  
0, otherwise
Integration of the RT equation
over the footprint area
Probability function p(z)
Sz
p( z )

K(z, , )
M( x, y, z)
This function defines the probability of finding
foliage elements in a horizontal plane Sz at
depth z.

Conditional probability function K(z, , )
This function defines the probability that one finds
a vegetated element at point M1 (, , ) , moving
from point M(x, y, z) along direction  , given that
there is a vegetated element at M(x, y, z) .
12
Time-dependent stochastic RT model
The integrated RT equation is solved numerically with the SOSA method (successive




orders of scattering approximations): Idn (t, r , )  J1 (t, r , )  J 2 (t, r , )  ... J n (t, r , ) ,


where J k (t, z, )  f (J k 1 (t, z, )), k  1, is the mean intensity of photons scattered k times.
Model inputs
(1) Characteristics
of incoming radiation
(2) Canopy structural parameters
(tree height, crown length, leaf area
(direction, intensity
and pulse duration)
density, leaf normal orientation distribution,
and statistical probability functions p and K)
(3) Optical
properties of
leaves and
soils
Model outputs
(1) The mean intensity of radiation over
a horizontal plane at depth z

I( t, z, ) 

I
(
t
,
x
,
y
,
z
,

)dxdy

Sz
 dxdy
Sz
(2) The mean intensity of radiation over
a vegetated area at depth z

U( t, z, ) 

 (x, y, z)I(t, x, y, z, )dxdy
Sz
 (x, y, z)dxdy
Sz
13
Model simulations for different types of forests
Averaged SLICER
waveforms
Model-simulated
waveforms
SOBS forest
Southern Old Black Spruce
 age
 height
 LAI
 location
100-155 years
10-11 m
4.0
BOREAS
Mature-aged forest
Tulip poplar association
 age
99 years
 height
32-37 m
 LAI
5.16
 location Maryland
14
Simulations with and without multiple scattering
Averaged SLICER waveforms
Model simulations:
with multiple scattering
without multiple scattering
Southern BOREAS study area
SOJP (old jack pine)
SOBS
 age, years
60-75
100-155
 height, m
16-19
10-11
2.61
4.0
 LAI
Eastern Maryland, Tulip poplar association
Mature
Intermediate
 age, years
99
41
 height, m
32-37
31-34
5.16
5.26
 LAI
15
Retrieval of canopy structural parameters
Multiply scattered photons carry information on canopy structural
parameters, namely, foliage density and gap fraction.
Single interaction: t total  t ps , where t ps is the round trip time between
the sensor and the interaction point
Several interactions: t total  t ps  t pm , where t pm is the extra time due to
multiple scattering
t pm
( N  1)d

c
1
d~ 
()

() ~ u L
d - photon mean free path

() - extinction coefficient
u L - leaf area volume density
16
Task 1: Conclusions
1. The inclusion of multiple scattering leads to a better approximation of
the return signal.
2. Multiply scattered photons magnify the signal and significantly enhance the
lower part of it.
3. In case of sparse canopies, effects of multiple scattering are insignificant and
single-scattering approximation models are expected to provide good simulations
of lidar-recorded signals.
4. Multiply scattered photons carry information on canopy structural parameters
S. Y. Kotchenova, N. V. Shabanov, Y. Knyazikhin, A. B. Davis, R. Dubayah, and R. B. Myneni,
Modeling lidar waveforms with time-dependent stochastic radiative transfer theory for remote
estimations of forests structure. Journal of Geophysical Research, 108 (D15), 2003.
17
Task 2: New application for lidar measurements
Canopy height profiles retrieved from lidar waveforms can be used
for modeling gross primary productivity of deciduous forests
18
Gross Primary Production
For several decades, monitoring and modeling the terrestrial carbon cycle have
been among the main goals for many ecological and climate change studies.
[NPP] = [GPP] – [Ra]
Autotrophic respiration
Net Primary Production
Gross Primary Production
(the carbon lost by photorespiration and by internal
plant metabolism; about half
of the total carbon uptake)
(the total amount of CO2 taken up
by land plants from the atmosphere
to participate in photosynthesis)
- Annual terrestrial GPP is estimated as about
120 Pg/yr.
- Forests cover about 26% of the total land
surface area.
Land surface: 148,300,000 km2
Forests: 38,700,000 km2
19
Overview of photosynthesis models
1) Production Efficiency Models (global and regional scales)
GPP=g(T, soil moisture, VPD).FAPAR.PAR
NPP=n(T, soil moisture, VPD).FAPAR.PAR
(g and n are the light use efficiencies, T is the ambient temperature, VPD is the vapor
pressure deficit, FAPAR is the fraction of absorbed Photosynthetically Active Radiation (PAR))
2) Models based on Farquhar’s equations (global, regional and local scales)
a) the big-leaf model
Farquhar’s model
Modeling
photosynthesis for
a unit leaf area
scaling from leaf
to canopy
b) the sunlit/shaded leaf
separation model
c) the multiple layer model
d) the combined leaf-separation
multiple-layer model
20
Advanced photosynthesis models
Unit leaf area
Farquhar’s model
+ gs(RAD, T, VPD)
Representation of the canopy
as two big leaves with different
illumination conditions
Lidar
measurements
Ph  Phsun Lsun  Phshade Lshade
Ph (RAD, T, CO2 , Ni,
VPD, O2, wind)
Two-leaf model
Use of actual foliage profiles
Ph   Lsun (i )Phsun (i )  Lshade (i )Phshade (i )
N
i
Combined two-leaf
multiple-layer model
Division of the canopy into N layers
with equal LAI + separation of sunlit
and shaded leaves in each layer
N
N
i
i
Ph  Lsun  Phsun (i )  Lshade  Phshade (i )
21
Drawbacks of the existing models
1. The two-leaf model (the leaf-separation model)
Use of the average values of diffuse and direct PAR absorbed by the canopy.
to capture the attenuation of incident PAR with height.
Fails
2. The coupled model (the combined leaf-separation multiple-layer model)
Assumption of a uniform distribution of foliage with height.
Wrong sunlit/shaded
leaf and PAR distributions except for the canopies with approximately uniform vertical
structures of foliage.
Will the accounting for the vertical distribution of
foliage lead to a better estimation of photosynthesis?
22
Distributions of PAR, sunlit/shaded leaves, and GPP
Comparison between
the uniform and actual
distributions of foliage
N
P h   P hlayer (i) 
i
L sun (i)P hsun (i) 

  


i
L shade (i)P hshade (i) 
N
23
Radiation model
Model inputs
Characteristics of the
incoming radiation



direct PAR flux
diffuse PAR flux
solar zenith angle
Canopy structural
parameters
Optical properties
of leaves and soils



tree height
vertical foliage profiles
leaf inclination
Model outputs


the distribution of direct PAR with height
the distribution of diffuse PAR with height
Reference - Shabanov et al., Stochastic modeling of radiation regime in discontinuous vegetation
canopies. Remote Sensing of Environment. 74, 125-144, 2000.
24
Comparison of photosynthesis models
Daily GPP rates calculated by
the coupled model with
uniform foliage profiles
the two-leaf model
(the sunlit/shaded leaf
separation model)
and
(the combined leaf-separation
multiple-layer model)
will be compared with the GPP rates calculated by
the coupled model with actual (lidar measured) foliage profiles
Comparison
(SLICER data collected over the mixed deciduous forest stands + environmental data)
1
100 profiles & 3 days with
different weather conditions
2
3 different profiles & 1 month
August, 1997
25
Daily GPP rates as functions of weather conditions
(1)
Measurements of incident PAR, temperature and humidity were taken at a 3-min time-step
Wye Research and Education Center, Maryland, August 1997
http://uvb.nrel.colostate.edu/UVB/uvb_climate_network.html
26
Daily GPP rates as functions of weather conditions (2)
GPP rates were calculated at a half-hour step and then integrated over a whole day-length period.
Intermediate-aged stand
- daily
GPP rates
calculated
by the twoleaf model
CHP
classification
Mature-aged stand
more than 80%
of foliage in the
first 2/3 of the
tree
more than 80%
of foliage in the
last 2/3 of the
tree
approximately
uniform
distribution of
foliage
remaining CHP
27
Daily GPP rates as functions of foliage profiles
(1)
SLICER waveforms (September 07, 1995)
PAR, temperature, humidity (August 1997)
Wye Research and Education Center, Maryland
http://uvb.nrel.colostate.edu/UVB/uvb_climate_network.html
Mean daily wind speed values (August 1997)
Andrews AFB station, Maryland
http://www.ncdc.noaa.gov/oa/climate/onlineprod/drought/
xmgr.html
28
Daily GPP rates as functions of foliage profiles
(2)
the two-leaf model
the coupled model
with uniform profiles
the coupled model
with actual profiles
(For each day, the GPP rate
was calculated at a half-hour
step and then integrated over
a whole day-length period.)
3 different canopy profiles:
(1) the majority of the foliage in
the first half of the canopy
(2) a nearly uniform distribution
(3) the majority of the foliage in
the lower part of the canopy
29
Task 2: Conclusions
1. The use of actual foliage profiles might lead to more accurate estimates of daily
GPP rates in the photosynthesis models relying on the extrapolation of Farquhar’s
equations from a unit leaf area to the whole canopy.
2. The disagreement range was almost the same for the two-leaf model and the
coupled model with uniform CHPs during the comparison with the coupled model
with actual CHPs. This means it might be useless to apply a multiple layer division
in addition to the two-leaf separation if the actual foliage profile is not known.
3. An increase of the incident diffuse PAR due to the partial cloudiness could
significantly enhance the daily carbon assimilation rate. The performed study also
demonstrates the importance of separate measurements of diffuse and direct PAR.
S. Y. Kotchenova, X. Song, N. V. Shabanov, C. S. Potter, Y. Knyazikhin, and R. B. Myneni,
Lidar remote sensing for modeling gross primary production of deciduous forests. Remote
Sensing of Environment, in Press.
30
Task 3: Improvement of the algorithm for retrieval of CHPs
Drawbacks of the current algorithm:
* 1) the assumption of a uniform horizontal distribution of needles
* 2) ignoring the effects of multiple scattering
3) the use of canopy material distribution instead of foliage distribution
31
Shoots as basic structural elements
1. Use of shoots as basic structural elements instead of needles
Description
of the canopy structure in terms of spatial and angular distribution of shoots


The shoot silhouette
area to total area ratio:


SSA ()
ST AR() 
T NA


SSA () - the shoot silhouette area in direction  ;
TNA - the total needle area of the shoot.

SSA ()


SSA ()
The mean projection of unit shoot area: G() 
T NA
STAR depends on the shoot structure, shoot orientation and position in the crown.
Reference – Stenberg, P. (1996). Simulations of the effects of shoot structure and orientation on vertical
gradients in intercepted light by coniferous canopies. Tree Physiology, 16, 99-108.
32
Accounting for multiple scattering
2a. Scattering of radiation between the needles of a shoot
(The scattering properties of needles are replaced with those of shoots)
The shoot (  sh ) and needle (  n ) scattering
 ()  p sh  n ()
coefficients are related as  sh ()  n
1  p sh  n ()
where p sh is the probability for a scattered photon to
hit the same shoot again.
Assumption: p sh = constant
Ray-tracing simulations:
psh  1  4  STAR
Reference – Smolander, S., & Stenberg, P. (2003). A method to account for shoot scale clumping in
coniferous canopy reflectance models. Remote Sensing of Environment, 88(4), 363-373.
2b. Scattering of radiation between shoots
(The RT model developed in Task 1 will be used to account for this type of scattering)
33
Modification of the algorithm
the current (M-H) algorithm
The modified algorithm
Step 1: accounting for scattering
of radiation between the
needles of a shoot (  n  sh )
Step 2: step 1 + accounting
for the shoot inclination
Step 3: steps 1-2 + accounting
for multiple scattering of
radiation between the shoots
R v (z) - the vegetation return at height z
R v (0) - the total vegetation return
R gr
- the ground return




R v (z)


CHPc (z)   ln1 

sh
R
(
0
)

R

v
gr 
gr


34
Conifer canopy models
1. The horizontal distribution of shoots is uniform; all layers are characterized, on average,
by the same ST AR value – the standard simulation.
mod 1
2. Shoots are uniformly distributed, but the mean ST AR increases from 0.1 to 0.25 with the
degree of shading. The increase is modeled in two ways:
a) as a linear function of a degree of shading:
mod 2a
STAR  0.1  0.15(DS)
b) as a quadratic function of the degree of shading exceeding 50%
ST AR  0.1 for DS  0.5
mod 2b
ST AR  0.1  0.15((DS  0.5) / 0.5) for DS  0.5
2
3. ...
DS is calculated from the field measurements of PAR transmittance. Reference - W. Ni, X. Li, C. Woodcock,
J.-L. Roujean, & R. E. Davis (1997). Transmission of solar radiation in boreal conifer forests: measurements
and models. Journal of Geophysical Research, 102 (D24), 29555-29566.
35
Retrieval of Canopy Height Profiles
(1)
SLICER data:
the southern old black
spruce site, BOREAS,
Canada. July 1996.
the M-H algorithm
The modified
algorithm
mod 1
mod 2a
mod 2b
36
Retrieval of Canopy Height Profiles
(2)
Average canopy area indices (CAIs) of 50 SOBS CHPs retrieved by
the modified algorithm
the M-H algorithm
mod 1
mod 2a
mod 2b
0.45
Step 1
Step 2
Step 3
final CAI
0.56 (23.0% )
2.09
(14.1% )
1.80
0.50 (11.0% )
1.56
(12.9% )
1.36
0.56 (24.9% )
2.61
(14.8% )
2.22
Field measurements of CAI: 1.87
Reference - J. M. Chen, P. M. Rich, S. T. Gower, J. M. Norman, & S. Plummer (1997). Leaf area index of
boreal forests: Theory, techniques, and measurements. Journal of Geophysical Research, 102 (D24),
29429-29443.
37
Task 3: Conclusions
1. Clumping of needles into a shoot, shoots inclination and multiple scattering of
radiation between the shoots are the only effects that can be taken into account
without contradicting the basis of the M-H algorithm.
2. The main advantage of the modified algorithm is the use of the mathematically
corrected approach for transferring from the distribution of needles to the distribution
of shoots.
S. Y. Kotchenova, N. V. Shabanov, Y. Knyazikhin, and R. B. Myneni (2004). Retrieval of canopy
height profiles from lidar measurements over coniferous forests. IEEE Transactions on Geoscience
and Remote Sensing, in Review.
38
Case study 4. New algorithm for biomass estimation *
Drawbacks of the current algorithm:
Requirement of extensive ground measurements to relate canopy
height profiles and biomass volumes (in particular, each tree diameter
within a footprint should be measured).
Research plan for this study:
A new algorithm is planned to be developed with the help of allometric
scaling theory of plant energetics and geometry.
33
Allometric scaling theory
A general model for the structure and allometry of plant vascular systems
Assumptions: (1) the branch structure of a tree is considered
a space-filling fractal system; (2) energy required to distribute
resources is minimized.
Results: the rate of fluid transport through this system scales
as ¾ of the tree mass
Reference – West et al. (1998). A general model for the origin
of allometric scaling laws in biology. Science, 276, 122-126.
Relationships between density and mass in resource-limited plants
N max - the max number of individuals that can be supported by unit
R area
Q - the rate of resource supply per unit area
M - the average rate of resource use per individual
Mtot - the average plant mass
- the total plant mass
3 / 4
3/ 4
R  NmaxQ  NmaxM
Nmax  M
Mtot  NmaxM  ( Nmax )1/ 3
Reference – Enquist et al. (1998). Allometric scaling of plant energetics
and population density. Nature, 305 (10), 163-165.
34
Application of the scaling model
Problems:
1. Lidar instruments do not distinguish between foliage and woody material.
How to extract foliage distribution from the retrieved CHP?
2. The general model for the structure and allometry of plant vascular systems is
developed for an individual tree, while the lidar footprint seizes more than one tree.
How to determine the number of trees within a footprint?
Solutions:
1. Evaluation of the percentage of branch surface in the CHP using the structure
and allometry general model.
2. Combination of the lidar technique with other remote sensing tools capable to
capture the horizontal structure of vegetation.
35
Mapping of the horizontal structure of vegetation
1. Aerial photography
It is one of the most widely used remote sensing tools
in forestry. A large number of air photos are available.
Reference – http://www.eoc.csiro.au
2. Polarimetric SAR Interferometry
Reference - S. R. Cloude, D. Corr. Tree height
retrieval using single baseline polarimetric
interferometry. Proceedings of ESA Workshop,
POLInSAR 2003, Frascati, Italy, 14-16 January.
Scots Pine forest model
Simulated SAR image
36
Conclusions
- The proposed research will potentially make an important contribution
into further development of lidar remote sensing technique.
- The lack of lidar data is a significant obstacle of the developed program.
37
Thank you!
Questions?..
38