Interactive comment on “ Accurate LAI retrieval method based on PROBA / CHRIS data ”

Abstract. Leaf area index (LAI) is one of the key structural variables in terrestrial vegetation ecosystems. Remote sensing offers an opportunity to accurately derive LAI at regional scales. The anisotropy of canopy reflectance, variations in background characteristics, and variability in atmospheric conditions constitute three factors that can strongly constrain the accuracy of retrieved LAI. Based on a hybrid canopy reflectance model, a new hyperspectral directional second derivative method (DSD) is proposed in this paper. This method can estimate LAI accurately through analyzing the canopy anisotropy. The effect of the background can also be effectively removed. With the aid of a widely-accepted atmospheric model, the influence of atmospheric conditions can be minimized as well. Thus the inversion precision and the dynamic range can be markedly improved, which has been proved by numerical simulations. As the derivative method is very sensitive to random noise, we put forward an innovative filtering approach, by which the data can be de-noised in spectral and spatial dimensions synchronously. It shows that the filtering method can remove random noise effectively; therefore, the method can be applied to hyperspectral images. The study region was situated in Zhangye, Gansu Province, China; hyperspectral and multi-angular images of the study region were acquired via the Compact High-Resolution Imaging Spectrometer/Project for On-Board Autonomy (CHRIS/PROBA), on 4 June 2008. After the pre-processing procedures, the DSD method was applied, and the retrieved LAI was validated by ground reference data at 11 sites. Results show that the new LAI inversion method is accurate and effective with the aid of the innovative filtering method.


Introduction
Leaf area index (LAI) is one of the important geometricstructure parameters in terrestrial vegetation ecosystems.Accurate LAI retrieval is always a key task of such remotesensing applications, and it is also a precondition of cropyield estimation, vegetative evapotranspiration calculation and the energy-exchange evaluation of terrestrial vegetation (Badhwar et al., 1986;Townshend et al., 1988;Baret et al., 1991;Wood et al., 1993;Bonan, 1995;Bicheron et al., 1999;Buermann et al., 2001;Kamalesh et al., 2008).In general, there are three categories of LAI retrieval methods.The first category is experiential inversion method based on Vegetation Index (VI), such as retrieving LAI according to the correlation between the Normalized Difference Vegetation Index (NDVI) and LAI (Turner, 1999;Brown et al., 2000;Haboudane et al., 2004;Wang et al., 2005).Although such an approach is simple and widely-used, its accuracy is very low, and some VI can easily get saturated.The second category can be called the reflectance model inversion method, by which the anisotropic effect caused by sun-object-sensor geometry can be removed, but the effect caused by background variations cannot be avoided, especially within mixed pixels.The complex background inhibits accurate LAI retrieval (Li et al., 1986;Nilson et al., 1989;Chen et al., 1997;Xu, 2005).In contrast, the hyperspectral second derivative method, which was founded at the end of the 1980s, can remove the effect of the soil background effectively, and this hyperspectral derivative method was established to retrieve LAI (Pu et al., 2000).Since the method is very sensitive to noise even if the canopy anisotropy is neglected, it is difficult to be used operationally.In this research, we have found that through the use of the proper filtering method which can remove noises brought about by factors other than the variation of LAI or Chlorophyll concentration, the directional second derivative (DSD) method is more reliable and accurate.
On 22 October 2001, PROBA (Project On-Board Autonomy) was launched as a technology demonstrator within ESA's General Support Technology Programme.CHRIS (Compact High-Resolution Imaging Spectrometer), the prime instrument of the PROBA mission, is hyperspectral, and its objective is data collection of BRDF (Bidirectional Reflectance Distribution Function) for a better understanding of the spectral reflectance.CHRIS acquires a set of images with five view angles (−55 • , −36 • , 0 • , 36 • and 55 • ) of each target during each acquisition sequence in 2.5 min (Verrelst et al., 2007;Kamalesh et al., 2008; CHRIS/PROBA Website), and the 62-band images are collected in the visible and near-infrared regions from 400 nm to 1050 nm.This new hyperspectral sensor offers a new opportunity to effectively improve the inversion accuracy of LAI.The main objective of this paper is to establish a useful method that can best estimate LAI by multi-angular hyperspectral data, as such an approach can constrain the effect of background variations and the canopy anisotropic reflectance.

The BRDF model
Supposing the height of canopy is H , there are four components within the sensor field-of-view, which are illuminated ground and canopy, as well as shadowed ground and canopy; K g ,K c ,K z ,K t represent the area ratios of the four parts respectively.Reflectance of targets can be described as ρ=ρ 1 +ρ m , where ρ 1 stands for the reflectance contributed by the single scattering and ρ m multiple scattering (Li et al., 1986;Nilson et al., 1989;Nilson, 1991;Chen, 1997).ρ 1 describes the anisotropic characteristics of the object, and according to the geometric model, ρ can be expressed as (Xu, 2009): Here, ρ g and ρ v are the hemispherical albedos of the soil background and the leaf, G v and G s mean the G function of the view direction and the solar direction, respectively (Xu, 2005); µ v =cosθ v , µ 0 =cosθ 0 , θ v and θ 0 are the viewing zenith angle and the solar zenith angle, respectively.Nilson parameter (λ 0 ) is used to describe the clumping effect of foliage (Nilson et al., 1989).(φ) is an empirical function used to describe the hot-spot effect.φ is the angle between light and the viewing directions.(0)=1, (π)=0, when 0<φ≤π, and (φ)=e − φ π −φ .µ 0 F 0 is the direct irradiance of the sun, and E d is the diffuse irradiance of atmosphere; the contribution of multi-scattering can be expressed by the Hapke model (Hapke, 1981(Hapke, , 1986)).
In order to evaluate the hybrid model, we used the field test data measured in Shunyi District, Beijing 2001, in an area covered with winter wheat (Wang et al., 2006).Measurements were carried out using an Analytical Spectral Devices, Inc. (ASD) spectroradiometer with a spectral range from 400 nm to 2500 nm.For the same points, several measurements were taken in different viewing directions along the principle plane, varying from −40 • (backward) to 40 • (forward) at an interval of 5 • .The bi-directional reflectance spectra simulated by the model are almost consonant with the measured spectra (Fig. 1a-f).The angular distribution of the wheat leaves is taken as Erectophile.The input leaf reflectance, transmittance spectrum and the soil reflectance spectrum were measured using ASD (Fig. 1g).The ratio of sky diffuse light data was also obtained in the field campaign (Fig. 1h).

The Directional Second Derivative (DSD)
According to Eq. ( 1), the second derivative of wavelength for the object was the sum of second derivative of wavelength for the leaf and the soil spectrum, when ρ v ρ g , Eq. ( 1) can be written as: Where ρ , ρ v and ρ g are the second derivative of wavelength for the object, the leaf, and the soil spectrum, respectively.According to Eq. ( 2), if ρ v ρ g , ρ is independent of the soil spectrum.

Band choice
Comparing the second derivative spectra simulated by the new hybrid canopy reflectance model with the second derivative spectra measured with the different backgrounds (Fig. 2), ρ v ρ g can be satisfied within 0.68 µm-0.71µm and 0.73 µm-0.75 µm bands, not only for the soil background, but also for the water body, concrete, and so on.So the    hyperspectral second derivative method can successfully remove the background effect.At the same time, Fig. 3 shows that multiple scattering can be ignored within the 0.68 µm-0.71µm-range, and thus it is the best band area choice.

Inversion formula
As mentioned earlier, within the 0.68 µm-0.71µm-range, multiple scattering can be neglected completely, so Eq. ( 2) can be transformed.Suppose X= ρ ρ v , then Equation ( 3) is the DSD inversion formula for LAI.Since ρ v ρ g , within the spectral range, Therefore, the influence of backgrounds can be removed by using the DSD method.When φ=0, (0)=1, that is the hot-spot direction; if b=λ 0 ,we have Equation ( 4) shows the relationship between the second derivative and LAI, as the target is supposed to be lambertian.
In fact, the vegetation-soil system is non-isotropic, and the relationship between the second derivative and LAI changes with the solar and view angles.It is just the difference between Eqs. ( 3) and ( 4).The range of G v is 0.3∼0.8, and λ 0 can be gained as a priori knowledge.0<λ 0 <1, θ v is known, and θ v <60 • .So (φ) is the only parameter that can dramatically affect the relationship between X and LAI.Impacts of multi-scattering on reflectance spectra (LAI=4.87; the leaf angular distribution is spherical; the leaf and soil reflectance spectrum was measured using ASD; the sky diffuse light ratio was measured using CE318).Figure 4 shows the relationship between X and LAI as (φ) changes.It can also be found that sensitivity of DSD method is very high, and that the dynamic range of the method is also large.The optimal view angle for LAI inversion is near the hot spot.When the view angle is near the cold spot, X can quickly saturate when LAI equals 3, so the direction is not suitable for LAI retrieval.Numerical simulations also show that the inversion error can be less than 5% in the hot-spot direction when LAI is larger than 0.5, even when 15% random errors are added (Fig. 5).The effective range of the method is LAI>0.5.When the LAI is very low, the vegetation information is also very low.Therefore, it is impossible to retrieve LAI accurately.

The innovative filtering approach
As mentioned before, the DSD method is very sensitive to noise, so it is necessary to conduct pre-processing to remove the effect of noise.An innovation de-noising method is proposed that can filter the noise in both the spectral dimension and the spatial dimension.First, the Minimum Noise Fraction (MNF) transform is used in the spatial dimension, and then the Fast Fourier Transform (FFT) in spectral dimension is adopted.
As far as the maize and wheat canopy spectra in the study area are concerned, the spectra variation in the visible-near infrared wavelength region (400-900 nm) are mainly caused by LAI and chlorophyll concentration of the leaf according to related studies (Liu et al., 2000;Pu et al., 1997Pu et al., , 2000)).The variation of canopy spectra through the Fourier transform caused by LAI and the chlorophyll concentration can be found only within the low frequency domain, much lower than the noise frequency, so a threshold can be found in the frequency domain to distinguish between the useful signal and noise.In this case, noise in the hyperspectral image can be removed by a low-pass filter.The Butterworth lowpass filter is chosen, but the threshold of the low-pass filter function should be determined first.In order to determine frequency ranges caused by the variation of LAI and chlorophyll concentration, the spectra simulated by PRO-SAIL model are adopted as the benchmark without the noise.The Fourier Transform of the canopy spectrum is analyzed to choose the threshold for the Butterworth low-pass filter.Because CHRIS mode 3 only acquired 18 bands, the spectrum is linearly interpolated to 180 bands first.By comparing the simulated reflectance with image samples in the frequency domain, the minimal value of Butterworth low-pass filter is given as 0.01283; frequencies higher than this can be recognized and wiped out.
The result using the FFT smoothing is shown in Fig. 6.By using the two-step de-noising method, the result is more reasonable.Since CHRIS images have only 18 bands, the effect of filtering is small.In order to demonstrate the ability of the innovative filtering method to remove noise brought about by factors other than the variation of LAI or chlorophyll concentration, the method was applied to a Hyperion image, which has 70 bands in the VNIR spectral range and 172 bands in the SWIR spectral range.From the results (Fig. 7), it can be seen that the innovative filtering method can remove the noises effectively.3 Study region and dataset

Study region
The study region is an oasis, situated in Zhangye City, in the middle reach of the Heihe River Basin, the second largest inland river basin ( 37 It is one of the most extremely arid areas of the country.The type and distribution of landscapes and the plant communities in the region can be found at website http://heihe.westgis.ac.cn/.The crops sowed in the oasis are mainly maize, wheat, barley and benne, and the growth of these crops always requires irrigation.The "Watershed Allied Telemetry Experiment Research (WATER)" remote sensing experiment was carried out in the Heihe River Basin, 2008 (Li et al., 2009).Our field experiment constituted one part of the WATER experiment.

Pre-processing of spaceborne CHRIS image
The CHRIS images were acquired on 4 June 2008, which are the 1A products.Mode 3 for land studies results in a spatial resolution of approximately 17 m at nadir and a band set consisting of 18 bands.At first, the HDFCLEAN software was used to de-stripe the image.The radiation of the image was corrected using the gain and bias provided by the CHRIS program.Geometric correction was performed using 1:10000 scale topographic maps.The geocorrection root mean squared error (RMSE) of 5 angles' images is all less than one pixel.As the pixel sizes in the 5 view-angle images are different, the resolutions of 5 images are over an average of 51 m to remove the effect of the varying view angles.Four images of the study region were intercepted to retrieve LAI (Fig. 8).
Atmospheric correction of the CHRIS radiance data products was performed to obtain surface reflectance using atmospheric sounding data and ACORN (1.5 mode), which is a specialized atmospheric-correction software for the highspectral-resolution images.ACORN is capable of processing data from tilted sensors by accounting for varying path lengths through the atmosphere and varying transmittance.Atmospheric-correction results in the retrieval of HDRF (Hemispherical Directional Reflectance Factor) data sets for the various CHRIS view angles.Figure 9 shows the spectra before and after the atmospheric processing.

Field measurements and validation
Ground reference data were collected in a maize and wheat field, the same as prior knowledge, G v =0.6, λ 0 =0.6 in the maize canopy, and G v =0.1, λ 0 =0.97 in the wheat canopy.Ground data collection included spectral radiometric measurements using an ASD FieldSpecPro FR. and canopy structure manual measurement.Two methods were used to measure LAI: the canopy was measured using a LI-COR LAI-2000 Plant Canopy Analyzer, or the plant was manually harvested and measured using LI-COR LAI 3100 area meter.Leaf water and chlorophyll content were measured in the laboratory.The reflectance of the 706 nm-band was used to calculate the directional second derivative, θ s =25 • , ϕ s =137.00 • when CHRIS transited the study area on 4 June and 24 June 2008.First the sowing area of maize and wheat was classified, and LAI was retrieved from CHRIS data at 5 view angles, based on the least-square method.The retrieved LAI is shown in Fig. 10.In order to evaluate the LAI retrieval performance, the estimated LAI (Fig. 10) was compared with the LAI retrieved by NDVI method according to the field measurements (Fig. 11).Since there are only 11 simultaneous ground LAI in the study area on the day when the remote sensing image was acquired, they are used only as validating dataset.The relationship between NDVI and LAI was established by using other ground measurement data measured during different growth periods, LAI changed from 0.5 to 6 in other 4 sites.
From the two figures, obvious changes in gradient are shown in Fig. 10, which represents a larger dynamic LAI range.In order to compare the two retrieval methods   quantitatively, ground measured LAI at 11 sites (Table 1) were taken as the validation dataset.Scatter of measured LAI and retrieved results of the two methods are shown in Fig. 12. LAI retrieved by using the DSD method shows a better agreement with the true values.The error (RMSE) has been analyzed to quantify the agreement between measured LAI at 11 sites and estimated values (Table 2).The two figures show that the inversion-derived LAI using DSD method is more accurate.

Conclusion and discussion
In this paper, the relationship between the second derivative of different view angles and LAI was established based on a new hybrid canopy reflectance model.Numerical simulation showed that this method can not only eliminate the effect of background and the anisotropy of canopy reflectance, but it also improved the dynamic range and inversion accuracy.As the DSD method is sensitive to noise, an innovative denoising method was proposed, which can reduce noise in the spatial and spectral dimensions.Additionally, the threshold of the low-pass filter function was determined according to the spectra simulated by PROSAIL model as a function of changes in LAI and chlorophyll concentration.The study demonstrated that the new filtering method can eliminate the random noise of the CHRIS image.
In order to validate the DSD method, Zhangye, Gansu Province was chosen as the study region.LAI was retrieved using multi-angular and hyperspectral CHRIS images on 4 June 2008.In order to evaluate the LAI retrieval performance, the estimated LAI using DSD method was compared with the LAI retrieved using NDVI method according to the field measurements.It shows that the DSD method can retrieve LAI accurately when applying the innovative denoising method.
If the sensor entails more than three view angles, the G function, λ 0 , and LAI can be retrieved simultaneously, and the requirement of a priori knowledge can be reduced, but the inversion accuracy will be influenced.Addressing this issue will constitute the next step of our research effort.

Figure 2 .
Figure 2. Comparison of second derivative of leaf spectrum and different background spectra,all the spectra are cited from the Johns Hopkins University Spectral Library, it can be found at http://asterweb.jpl.nasa.gov/speclib/.

19Figure 3 .
Figure3.Impacts of multi-scattering on reflectance spectra (LAI=4.87; the leaf angular distribution is spherical; the leaf and soil reflectance spectrum was measured using ASD; the sky diffuse light ratio was measured using CE318)

Figure 3 .Fig. 4 .
Figure 3. Impacts of multi-scattering on reflectance spectra (LAI=4.87; the leaf angular distribution is spherical; the leaf and soil reflectance spectru was measured using ASD; the sky diffuse light ratio was measured using CE318)

20Figure 5 .Fig. 5 .
Figure 5.The inversion error of the method in numerical simulation

20Figure 5 .Fig. 6 .
Figure 5.The inversion error of the method in numerical simulation

Figure. 7 Fig. 7 .
Figure.7 Hyperion spectrum before and after filtering by using the innovative fi method

22Figure 9 .Fig. 10 .
Figure 9. Soil and canopy spectra on CHRIS image before and after atmospheric correction

Table 1 .
List of in situ-measured LAI.

Table 2 .
Error analysis of retrieval results using CHRIS images.