CN106909750A - The computational methods and model of a kind of broad-leaved Vegetation canopy reflectivity - Google Patents

Abstract

The invention discloses the computational methods and model of a kind of broad-leaved Vegetation canopy reflectivity, computational methods are comprised the following steps：S1：The parameter of input recognizes and classified by |input paramete, is divided into blade parameter, Crown canopy parametre and soil parameters；S2：The reflectivity and transmissivity of individual blade are calculated according to blade parameter；S3：Blade parameter in Crown canopy parametre and S2 tries to achieve the extinction coefficient and scattering coefficient of canopy；S4：The associated reflections factor and reflectivity of canopy are tried to achieve according to required canopy delustring and scattering parameter；S5 tries to achieve canopy reflectance spectrum according to the associated reflections factor and reflectivity of canopy.The present invention is by PROSPECT models and SAIL Model couplings, make full use of can get parms in the case of cancel Vegetation canopy reflectivity simulation process Leaf reflectivity and transmissivity input process, it is of the invention effectively to simplify parameter acquiring problem and optimized algorithm in Vegetation canopy spectral information simulation process, speed-up computation process, while coupling model is conducive to vegetation parameter inverting.

Description

The computational methods and model of a kind of broad-leaved Vegetation canopy reflectivity

Technical field

The present invention relates to high-grade highway road domain vegetation health status evaluation, and in particular to a kind of broad-leaved Vegetation canopy reflection The computational methods and model of rate.

Background technology

Vegetation information inverting is all the time the research field of quantitative remote sensing most prospect, to carry out vegetation ecological environment Monitoring and evaluation research provides favourable basis.Domestic and foreign scholars propose numerous vegetation parameter inverse models, can be divided mainly into Statistical model and the major class of physical model two, wherein physical model enjoy favor with its outstanding broad applicability and stability.Plant Great physical basis are had by the physics inverse model of information, the change of the particular type and background environment of vegetation is not relied on Change.But physical model exist in large-scale vegetation information inverting precision it is still poor and it is computationally intensive lead to not meet The problem of production application demand, mainly due to two reasons：One is that requirement of the physical model to quantitative remote sensing is very high, is needed Reflectivity for Growing Season is finally inversed by by remote sensing image, while needing multiple to be difficult to the accurate parameters for obtaining as the defeated of physical model Enter parameter；On the other hand, current Remote sensing physical modeling is most on canopy scale, and canopy spectra is subject to Spectra of The Leaves information, The influence of the factor such as Soil Background reflectivity and canopy structure, the precision for limiting biochemical parameter inverting also increases the difficulty of inverting Degree.

Recent study and have PROSPECT models suitable for broad-leaved using most models, suitable for needle LIBERTY models and canopy SAIL models.Broad-leaved blade is regarded as and contains blade biochemical substances by N layers by PROSPECT models The compound slab model that Rough Flat Plate and N-1 layer of air are constituted is dense by being input into each pigment concentration, equivalent water thickness, white material Degree and blade construction parameter obtain the reflectivity and transmissivity of this compound slab model.LIBERTY models regard needle blade Cell is standard circular cell, it is believed that needle is formed in atmosphere by countless leaf cells stacking, by trying to achieve individual cells Optical characteristics further tries to achieve the spectral information of indefinite cell.Canopy SAIL models are to spoke on the basis of radiation transfer theory Four streams for penetrating transmission equation are approximate, model by be input into solar zenith angle, solar azimuth, view zenith angle, observed azimuth, The parameters such as the reflectivity and transmissivity of LAI, spectral reflectance, focus and single blade simulate the anti-of canopy 400nm to 2500nm Penetrate rate.Obvious can be that the output parameter of broad-leaved and needle model is realized into mould as the |input paramete of SAIL models The coupling of type, further realizes the simulation of Vegetation canopy information.Therefore, domestic and foreign scholars have successively been carried out Model coupling and its have been answered Research, current application it is more be to be coupled with SAIL canopy reflectance models using PROSPECT leaf models, the coupled mode Type using revised leaf model and canopy modeling Vegetation canopy reflectivity, with inversion accuracy higher and faster Inversion speed.However, currently used is generally the PROSAIL models of the more perfect initial stage Model coupling of development, failing To PROSPECT models and SAIL models in the case of the problems such as considering vegetation blade upper and lower surface asymmetry and hot spot-effect Coupled, although adaptability and precision have certain guarantee, vegetation parameter is still exposed in actual application The problems such as inverting is not enough so that coupling model is limited in actual application by certain procedures.

Computational methods and model it is therefore desirable to design a kind of new broad-leaved Vegetation canopy reflectivity.

The content of the invention

Technical problem solved by the invention is, for the problem that physical model is limited in actual application, also for Vegetation physical model parameter acquiring and inversion accuracy in application process are solved the problems, such as, a kind of broad-leaved Vegetation canopy reflection is proposed The computational methods and model of rate, the single blade reflectivity that single blade model is obtained by the use of acquired parameter and transmissivity as hat The |input paramete of layer model, cooperates with other |input parametes simulation Vegetation canopy reflectivity of canopy model, the vegetation hat after coupling Layer model can clipped complex parameters acquisition process, while coupling after the integrated model of model, advantageously in plant By the inverting of parameter.

The technical scheme is that：

A kind of computational methods of broad-leaved Vegetation canopy reflectivity, comprise the following steps：

S1：Parameter is recognized；

Input model parameter；And the parameter of input is tentatively divided into three major types：Blade parameter, soil parameters and canopy ginseng Number；Wherein blade parameter includes pigment content, dry matter content, equivalent water thickness and blade construction parameter；Soil parameters is soil The spectral reflectivity of earth；Crown canopy parametre is remaining other specification；

S2：Blade parameter input PROSPECT models in step S1 carry out monolithic leaf spectral simulation, calculate single The spectral reflectivity and transmissivity of blade；

S3：Soil parameters and Crown canopy parametre in step S1 calculate extinction coefficient and scattering coefficient；

S4：Extinction coefficient and scattering coefficient the input SAIL models that step S3 is calculated, the correlation for calculating canopy are anti- Penetrate the factor and reflectivity, including sun direct projection direction hemispherical reflectance, observed direction hemispherical reflectance and double direct projection reflectivity；

S5：Calculate canopy reflectance spectrum.

The step S2 specifically includes following steps：

S2.1：Calculate absorption coefficient K；Absorption coefficient is that the linear list of the various content of chemical substances of blade reaches, computing formula For：

Wherein, c (i) is the concentration of component i in blade, and component i is the chemical substances such as pigment, water, dry；K (i) is group Divide the specific absorption coefficient of i；N is blade construction parameter, represents blade hierarchy number, and N can be calculated using following empirical equation Arrive：

N=(0.9*SLA+0.025)/(SLA-0.1)

Wherein, SLA is specific leaf area, refers to the leaf area of per unit dry weight；

S2.2：Calculate the transmissivity of interface：

The transmissivity of interface includes two：One is transmissivity of the light from air into blade, and one is light from blade The transmissivity of air inlet；

Available light is considered as non-polarized light, and perfact polarization occurs after reflection, refracted light generating unit split pole after transmission Change.Be may know that according to Snell-Descarts laws, on two media contact surface incidence angle and refraction angle and two refractive indexes it Between relation, the electromagnetic wave be given further according to Schanda through two media transmissivity by calculate light with solid angle α It is the average transmittance t of the blade of n by the air incidence refractive index that refractive index is 1_av(α, 1, n), make t₁₂=t_av(α,1,n)；

According to the achievement in research of Stern, light injects the transmissivity of the air that refractive index is 1 from refractive index for the blade of n t₂₁With t₁₂Between there is relational expression t₂₁=n^{- 2}·t₁₂, the transmissivity of interface is calculated accordingly；The meter of the transmissivity of interface It is prior art to calculate, and can be used existing software to complete to calculate；

S2.3：Calculate total transmittance τ of the light through plating media₁；

For natural light, its total transmittance τ for passing through plating media₁It is the integration in the whole spaces of the π of whole hemisphere 2, with suction Receive COEFFICIENT K and vane thickness D is related, computing formula is：

Wherein, t is intermediate variable；

S2.4：Calculate the reflectivity R of individual layer flat board_αAnd transmissivity T (1)_α(1)：

Individual layer flat board spectral simulation is carried out, the spectrum of individual layer flat board is the core of whole Spectra of The Leaves simulation, is mainly included The reflectivity and transmissivity of individual layer flat board are calculated, its computing formula is：

S2.5：Calculate monolithic foliage reflectance R_αAnd transmissivity T (N)_α(N)：

Spectra of The Leaves information simulation is carried out, Stokes is goed deep into the optical phenomena after light penetration finite layer flat board Research, has drawn light through the mass reflex rate after the flat board of N layers of identical reflectivity and transmissivity and the calculating reason of transmissivity By can now draw the simulation formula of Spectra of The Leaves information using its theory：

Wherein：

The step S3 specifically includes following steps：

S3.1 shadow compensations：

1. the geometrical factor related to delustring and scattering is calculated：First according to general Leaf angle inclination distribution probability weight and discrete Change, obtain one group of Leaf inclination of discretization；Afterwards, the corresponding delustring of each Leaf inclination and dispersion factor are calculated；For Leaf inclination θ_l, computational methods are as follows：

First, critical angle β is sought_sAnd β_o, computing formula is：

Wherein, it is horizontal plane normal direction and blade normal angular separation i.e. Leaf inclination θ_sIt is solar zenith angle, θ_oIt is observation Zenith angle；

Two critical angles are directly calculated when it is determined that denominator is not less than 1 for the result of calculation of 0 and cos β in computing formula Value；When the result of calculation of cos β is equal to 1, two are faced Jie angle and are equal to π；Wherein cos β refer to β_sAnd β_oCosine；

Then, it is θ to calculate Leaf inclination_lIndividual blade sun direct projection direction and observed direction extinction coefficient, calculate public Formula is respectively：

Wherein, L '=lai/h, lai are leaf area index, and h is canopy height；

2. auxiliary azimuthal angle beta is calculated₁、β₂、β₃：The azimuthal calculating of auxiliary depend primarily on two critical angles relatively not Etc. relation, obtaining value method is as follows：

If：
					ψ
		ψ
							ψ

Wherein, ψ is the relative bearing between solar direction and observed direction, the i.e. difference at both direction azimuth；

The single blade reflectivity multiplier related to transmissivity can be calculated after obtaining three auxiliary azimuths, for calculating pair To scattering coefficient ω；

3. two-way dispersion coefficient ω (θ are calculated_l)：

After auxiliary azimuth is obtained, the single blade reflectivity and transmittance calculation being calculated in cooperation S2 obtain double To scattering coefficient, computing formula is：

Wherein, ρ, τ, θ_lThe reflectivity of blade, transmissivity and now corresponding Leaf inclination are represented respectively；Make the ρ=R in S2_α (N)；τ=T_α(N)。

S3.2 calculates the scattering coefficient after adding leaf reflectance and transmissivity：

1. the backscattering coefficient computing formula of diffusing radiation E- and E+ is：

2. the forward scattering coefficient formulas of diffusing radiation E- and E+ are：

3. beam radia E_SBackscattering coefficient computing formula be：

4. beam radia E_SForward scattering coefficient formulas be：

5. the attenuation coefficient computing formula of diffusing radiation E- and E+ is：

6. observed direction radiates E₀Backscattering coefficient computing formula be：

7. observed direction radiates E₀Forward scattering coefficient formulas be：

S3.3 calculates general Leaf angle inclination distribution probability (blade profile profile element)：

Calculate the general Leaf angle inclination distribution under average Leaf angle inclination distribution：Leaf inclination refers to horizontal plane normal direction and blade method Line angular separation；The different type of Vegetation canopy corresponds to different vegetation blade pitch angles, while corresponding to two leaves point Cloth parameter LIDFa and LIDFb.

In the case where non-spheroid shape is distributed, general Leaf angle inclination distribution probability can be calculated according to average Leaf inclination, calculated Journey is completed by following iteration：

X=2 θ

Y=LIDFasin (x)+0.5LIDFbsin (2x)

Dx=0.5 (y-x+2 θ)

X=x+dx

Until | dx | ＜ t；

Then F (θ)=2 (y+ θ)/π.

Wherein, θ is the centrifugal pump of average Leaf inclination, and F (θ) is accumulation Leaf inclination, i.e., general Leaf angle inclination distribution probability.

For elliposoidal distribution, what parameter LIDFa was represented is distribution angle, is 30 degree, and parameter LIDFb is 0, general asking During Leaf angle inclination distribution probability, the Leaf inclination density function made using Campbell (1986) is tried to achieve；

S3.4 determines the model parameter of whole canopy；

For any canopy, blade tilt will not be fixed, be one from 0 to 90 ever-increasing continuous mistakes Journey, it is thus determined that during the model parameter of whole canopy, it by Leaf inclination is θ that computational methods are_lWhen general Leaf angle inclination distribution probability F (θ_l) be added again with after model parameter product, computing formula is：

Z=∑ F (θ_l)Z(θ_l)

Wherein, F (θ_l) it is Leaf inclination θ_lCorresponding general Leaf angle inclination distribution probability, Z (θ_l) for single blade parameter, Z be through Cross Leaf inclination correction model Z=∑ F (θ_l)Z(θ_l) it is subject to the parameter of revised whole canopy；, Z (θ_l) refer to step S3.1 and K (the θ obtained in S3.2_l)、K(θ_l)、ω(θ_l)、σ(θ_l)、σ′(θ_l)、s(θ_l)、s′(θ_l)、a(θ_l)、v(θ_l)、u(θ_l) in appoint Meaning one；Z correspondingly refer to all coefficient ks in four differential equations of SAIL models for whole canopy, K, w, σ, σ ', s, Any one in s ', a, v, u.

The step S4 specifically includes following steps：

S4.1：Calculate leaf area parameter；Leaf area parameter is to evaluate the important indicator of Vegetation canopy reflectivity, its parameter meter Calculating formula is：

τ_ss=e^-klai

τ_oo=e^-Klai

ρ_dd=(e^mlai-e^-mlai)/(h₁e^mlai-h₂e^-mlai)

τ_dd=(h₁-h₂)/(h₁e^mlai-h₂e^-mlai)

Wherein, K and k represent that observed direction radiates E respectively₀Extinction coefficient and beam radia E_sExtinction coefficient； Lai represents leaf area index, is mode input parameter；

S4.2：The parameter after adding hot spot-effect is calculated, computing formula is as follows：

ρ_sd=C_S(1-τ_ssτ_dd)-D_Sρ_dd

τ_sd=D_S(τ_ss-τ_dd)-C_Sτ_ssρ_dd

ρ_do=C_o(1-τ_ooτ_dd)-D_oρ_dd

τ_do=D_o(τ_oo-τ_dd)-C_oτ_ooρ_dd

ρ_so=H_o(1-τ_ssτ_oo)-C_oτ_sdτ_oo-D_oρ_sd

Wherein, the computational methods of each intermediate variable are：

h₁=(a+m)/σ

h₂=(a-m)/σ=1/h₁

C_S=[s ' (k-a)-s σ]/(k²-m²)

C_o=[v (K-a)-u σ]/(K²-m²)

D_s=[- s (k+a)-s ' σ]/(k²-m²)

D_o=[- u (K+a)-v σ]/(K²-m²)

H_s=(uC_S+vD_s+w)/(K+k)

H_o=(sC_o+s′D_o+w)/(K+k)

Wherein, τ_ssooIt is hot spot-effect corrected parameter；Due to the problem of solar radiation, the particularly soil in sparse vegetation region The temperature difference of earth and vegetation is quite big, and this is related to its physical surface and Meteorological, it is therefore desirable to a process for temperature adjustmemt, It is likely to provide extra Vegetation canopy information simultaneously.The makeover process of hot spot-effect is：

Focus correction is carried out according to 2/ (k+K) effect, first according to given focus value calculating parameter：

The initial value of given alf is alf=10⁶,

If 1) focus value q is more than 0, then it is assumed that be there is the hot spot-effect to influence, now Wherein, θ_sIt is solar zenith angle, θ_oFor view zenith angle,It is the sun Azimuth；(its alf value is 200 if α result of calculations are more than 200)；If the result of calculation of alf is 0, calculate as follows The parameter τ gone out under the influence of the pure hot spot-effect of shadow-free_ssooAnd sumint：

τ_ssoo=τ_ss

If the result of calculation of alf is not 0, then it is assumed that be, without hot spot-effect influence, τ to be calculated by the method in 2)_ssooWith sumint；

If 2) focus value q is 0, then it is assumed that be that, without hot spot-effect influence, the method being added using circulation is calculated without focus Parameter τ under effects_ssooAnd sumint：Concretely comprise the following steps:

The first step, parameter initialization：

Second step, by following iterative calculation τ_ssooAnd sumint：

Step 1, the value for judging i：

If i is less than 20, x2=-log (1-i*fint)/alf is made；

If i=20, x2=1 is made；

If i is more than 20, terminate iteration, make τ_ssoo=f1；；

Step 2, following calculating is carried out successively：

Y2=- (K+k) * lai*x2+fhot* (- exp (- alf*x2))/alf；

F2=exp (y2)；

Sumint=sumint+ (f2-f1) * (x2-x1)/(y2-y1)；

X1=x2；

Y1=y2；

F1=f2；

Step 3, makes i=i+1, and go to step 1；

Wherein, * represents multiplication.

S4.3：Calculate bidirectional reflectance：

Bidirectional reflectance includes two parts, and scattering,single of the part with hot spot-effect influences, and another part is not with focus The Multiple Scattering influence of effect, two parts are summed to the result of calculation of the final reflectivity influence at the top of canopy：

ρ_so2=ρ_so+w*lai*sumint

S4.4：Introduce soil reflex：

Soil is located at the bottommost of whole canopy reflecting system, and electromagnetic wave can be mapped to soil ground after penetrating canopy, through anti- Canopy top is reflected back after penetrating again, a part for canopy reflectance spectrum is formed, numerical procedure is：

d_n=1-r_s*ρ_dd

S4.5：Calculate double hemispherical reflectance factors：

S4.6：Calculate sun direct projection direction hemispherical reflectance：

ρ_sd2=ρ_dd+(τ_sd+τ_ss)*r_s*τ_dd/d_n

S4.7：Calculating observation direction hemispherical reflectance：

ρ_do2=ρ_do+(τ_do+τ_oo)*r_s*τ_dd/d_n

S4.8：Calculate double direct projection reflectivity：

ρ_sod2=ρ_so+((τ_ss+τ_sd)*τ_do+(τ_sd+τ_ss*r_s*ρ_dd)*τ_oo)*r_s/d_n

ρ_so3=ρ_so2+ρ_so+τ_ssoo*r_s

Wherein, r_sIt is spectral reflectance.

The step S5 specifically includes following steps：

First, can be divided in the hope of direct sunlight line and the fundamental radiation amount of atmospheric scattering light according to database data E is not designated as it_sAnd E_d；

Then, calculate scattering in representing air and radiate the coefficient skyl for accounting for global radiation ratio：

Skyl=0.847-1.61*sin (90- θ_s)+1.04*sin(90-θ_s)²

Finally, the amount of radiation that calculating observation direction is measured to and the ratio of the amount of radiation of incidence, obtain canopy reflectance spectrum, public Formula is：

A kind of computation model of broad-leaved Vegetation canopy reflectivity, formula is as follows：

Wherein, R is canopy reflectance spectrum, E_sAnd E_dThe respectively fundamental radiation amount of direct sunlight line and atmospheric scattering light, Skyl is that scattering radiates the coefficient for accounting for global radiation ratio in representing air；ρ_do2And ρ_so3Respectively observed direction hemispherical reflectance and Direct projection reflectivity；Model parameter is solved according to the computational methods of above-mentioned broad-leaved Vegetation canopy reflectivity.

Beneficial effect：

The present invention will simulate broad-leaved based on newest PROSPECT broad-leaveds model and SAIL canopy scale-model investigation achievements The PROSPECT models of Spectra of The Leaves information are coupled with the SAIL models of simulation canopy spectrum information, and making full use of to obtain Parameter save Spectra of The Leaves and simulate this process, cancel Vegetation canopy reflectivity simulation process Leaf reflectivity and transmissivity The process of input, adds canopy SAIL models to complete coupling process leaf model simulation algorithm, directly by vegetation physical and chemical parameter Simulation canopy spectra.Carry out coupling the essence for being greatly improved physical model simulated spectra using physical model newest research results Degree, improves model business Motor ability, while the parameter acquiring effectively simplified in Vegetation canopy spectral information simulation process is asked Inscribe and optimized algorithm, speed-up computation process, improve the inverting ability of vegetation information.

Brief description of the drawings

Fig. 1 is principle of the invention figure；

Fig. 2 is the model parameter of present invention input；

Fig. 3 is experimental result of the present invention

Specific embodiment

The present invention is further described below in conjunction with accompanying drawing.

Based on the present invention is using the SAILH canopy models after hot spot-effect improvement, given parameter is divided first Class, is divided into blade parameter and Crown canopy parametre two parts, and afterwards using the physical and chemical parameter of blade, the base absorption for calculating blade is made With in conjunction with reflectivity and transmissivity that individual layer flat board is derived without the transmittance calculation under the absorption of interface, finally utilization point Shelf theory derives the spectral information of given structural parameters lower blade；Simultaneously for Crown canopy parametre first with given Leaf inclination Distributed model parameter calculates blade profile profile element, calculates the extinction coefficient of canopy model in conjunction with the geometric element for calculating and dissipates Penetrate coefficient and combined with blade profile distributed model, the complete Crown canopy parametre under intact leaf is distributed is obtained by cyclic process, finally, The reflectivity for determining canopy is sought using two radiation parameters and sky scattering coefficient tried to achieve.

The implementation case is by taking the tree that highway side growing state is good and canopy spectra is easily determined as an example to this The technical scheme of invention is further described.As shown in drawings, improved physical model PROSAILH involved in the present invention Computational methods it is as follows：

By determining, the parameter such as following table of the tree：

Chlorophyll content	45.80	Leaf angle inclination distribution type	Apsacline
				Equivalent water thickness	0.02	Leaf area index	4.30
Dry matter content	0.02	Solar zenith angle	30.00
				Blade construction parameter	1.30	View zenith angle	10.00
Focus	8.00	Relative bearing	0.00

Choice experiment wavelength is 400nm-2399nm

S1 parametric classifications：It is divided into two parts after the relevant parameter for obtaining vegetation, wherein blade parameter includes the color of blade Cellulose content, equivalent water thickness, dry matter content and maximum incident angle (generally 59 degree)；Crown canopy parametre is joined including Leaf angle inclination distribution The background reflectivity of number, leaf area index, solar zenith angle, view zenith angle, relative bearing and focus and soil.

S2：Spectra of The Leaves information is resolved：After relevant parameter is obtained, it is necessary first to calculate the spectral reflectivity of individual blade And transmissivity, numerical procedure is as follows：

(1) absorption coefficient K is calculated：Absorption coefficient is that the linear list of the various content of chemical substances of blade reaches, and computing formula is：

(2) without transmittance calculation under the absorption of interface：Available light is considered as non-polarized light, occurs after reflection complete Polarization, after transmission there is partial polarization in refracted light.Be may know that on two media contact surface according to Snell-Descarts laws Relation between incidence angle and refraction angle and two refractive indexes, two media is passed through further according to the electromagnetic wave that Schanda is given Transmissivity can calculate the average transmittance of blade of the light with α solid angles by air incidence refractive index that refractive index is 1 as n t_av(α, 1, n), now assign its code t₁₂, according to the achievement in research of Stern, light is from refractive index for the blade of n injects refraction Rate is 1 air transmissivity t₂₁With t₁₂Between there is relational expression t₂₁=n^{- 2}t₁₂, the transmissivity of interface is calculated accordingly.

(3) transmission of the light in isotropic medium：For natural light, its total transmittance τ for passing through plating media is Integration in the whole π spaces of hemisphere 2, related to absorption coefficient K and vane thickness D, computing formula is：

(4) individual layer flat board spectral simulation：The spectrum of individual layer flat board is the core of whole Spectra of The Leaves simulation, main to include list The reflectivity and transmissivity of layer flat board, its computing formula is：

(5) monolithic foliage reflectance and transmissivity are simulated：Stokes enters to the optical phenomena after light penetration finite layer flat board Go further investigation, draw light through the mass reflex rate and transmissivity after the flat board of N layers of identical reflectivity and transmissivity, The simulation formula of Spectra of The Leaves information can be now drawn using its theory：

Wherein：

S3：Calculate extinction coefficient and scattering coefficient

(6) geometric element is calculated：

(7) blade profile profile element is calculated：

General Leaf angle inclination distribution under average Leaf angle inclination distribution：The different type of Vegetation canopy corresponds to different vegetation leaves Piece angle of inclination, while corresponding to two leaf distributed constants LIDFa and LIDFb.In the case where non-spheroid shape is distributed, can basis Average Leaf inclination calculates general inclination angle, and calculating process can be completed with a simple iteration, and false code is：

X=2 θ

Y=LIDFasin (x)+0.5LIDFbsin (2x)

Dx=0.5 (y-x+2 θ)

X=x+dx

Until | dx | ＜ t；

Then F (θ)=2 (y+ θ)/π

Wherein θ is the centrifugal pump of average Leaf inclination, and F (θ) is accumulation Leaf inclination.

For elliposoidal distribution, what parameter LIDFa was represented is distribution angle, is 30 degree, and parameter LIDFb is 0, general asking During angular distribution, the Leaf inclination density function made using Campbell (1986) is tried to achieve.

(8) delustring and scattering coefficient are calculated：

Shadow compensation：

1. geometrical factor related to delustring and scattering is calculated：First to Leaf angle inclination distribution weighting and discretization, one is obtained The central angle angle value Leaf angle inclination distribution of group discretization.Afterwards, each Leaf inclination correspondence calculates one group of delustring and dispersion factor, single Numerical procedure is as follows：

Critical angle β is sought first_sAnd β_o, computing formula is：

Two critical angles are directly calculated when it is determined that denominator is not less than 1 for the result of calculation of 0 and cos β in computing formula Value；When the result of calculation of cos β is more than 1, two are faced Jie angle and are equal to π；Wherein cos β refer to β_sAnd β_oCosine；

The extinction coefficient that sun direct projection direction and observed direction can be just calculated behind two interim boundaries angle is calculated, calculates public Formula is respectively：

Wherein, L '=lai/h, lai are leaf area index, and h is canopy height；

2. it is used to calculate the auxiliary azimuthal angle beta of two-way dispersion coefficient W₁、β₂、β₃Calculate：The azimuthal calculating of auxiliary mainly takes Certainly in the relation such as not relatively of two critical angles, obtaining value method is as follows：

If：
					ψ
		ψ
							ψ

Wherein, ψ is the relative bearing between solar direction and observed direction, the i.e. difference at both direction azimuth；

The single blade reflectivity multiplier related to transmissivity can be calculated after obtaining three auxiliary azimuths, for calculating pair To scattering coefficient W.

3. two-way dispersion coefficient is calculated：After auxiliary azimuth is obtained, coordinate the single blade reflection being calculated in S2 Rate and transmissivity can be calculated two-way dispersion coefficient, and computing formula is：

Wherein, ρ, τ, θ_lThe reflectivity of blade, transmissivity and now corresponding Leaf inclination are represented respectively；Make the ρ=R in S2_α (N)；τ=T_α(N)。

Model parameter is resolved

(9) calculating (addition of leaf reflectance and transmissivity) of scattering coefficient：

1. the backscattering coefficient computing formula of diffusing radiation E- and E+ is：

2. the forward scattering coefficient formulas of diffusing radiation E- and E+ are：

3. beam radia E_SBackscattering coefficient computing formula be：

4. beam radia E_SForward scattering coefficient formulas be：

5. the attenuation coefficient computing formula of diffusing radiation E- and E+ is：

6. observed direction E₀Backscattering coefficient computing formula be：

7. observed direction E₀Forward scattering coefficient formulas be：

For any canopy, blade tilt will not be fixed, be one from 0 to 90 ever-increasing continuous mistakes Journey, it is thus determined that during the model parameter of whole canopy, computational methods are Leaf angle inclination distribution function as probability function F (θ) and mould It is added again after shape parameter product, computing formula is：

Z=∑ F (θ_l)Z(θ_l)

Wherein, F (θ_l) it is Leaf inclination θ_lCorresponding general Leaf angle inclination distribution probability, Z (θ_l) for single blade parameter, Z be through Cross Leaf inclination correction model Z=∑ F (θ_l)Z(θ_l) it is subject to the parameter of revised whole canopy；, Z (θ_l) refer to step 1) and step It is rapid 2) in the k (θ that obtain_l)、K(θ_l)、ω(θ_l)、σ(θ_l)、σ′(θ_l)、s(θ_l)、s′(θ_l)、a(θ_l)、v(θ_l)、u(θ_l) in appoint Meaning one；Z correspondingly refer to all coefficient ks in four differential equations of SAIL models for whole canopy, K, w, σ, σ ', s, Any one in s ', a, v, u.

S4：Calculate the associated reflections factor and reflectivity of canopy, including sun direct projection direction hemispherical reflectance, observed direction Hemispherical reflectance and double direct projection reflectivity；

(1) addition of leaf-area coefficient：Leaf area index is to evaluate the important indicator of Vegetation canopy reflectivity, its parameter meter Calculating formula is：

τ_ss=e^-klai

τ_oo=e^-Klai

ρ_dd=(e^mlai-e^-mlai)/(h₁e^mlai-h₂e^-mlai)

τ_dd=(h₁-h₂)/(h₁e^mlai-h₂e^-mlai)

(2) parameter of solution is calculated：Parameter after addition hot spot-effect is more complicated, and computing formula is as follows：

ρ_sd=C_S(1-τ_ssτ_dd)-D_Sρ_dd

τ_sd=D_S(τ_ss-τ_dd)-C_Sτ_ssρ_dd

ρ_do=C_o(1-τ_ooτ_dd)-D_oρ_dd

τ_do=D_o(τ_oo-τ_dd)-C_oτ_ooρ_dd

ρ_so=H_o(1-τ_ssτ_oo)-C_oτ_sdτ_oo-D_oρ_sd

The computational methods of wherein each variable are：

h₁=(a+m)/σ

h₂=(a-m)/σ=1/h₁

C_S=[s ' (k-a)-s σ]/(k²-m²)

C_o=[v (K-a)-u σ]/(K²-m²)

D_s=[- s (k+a)-s ' σ]/(k²-m²)

D_o=[- u (K+a)-v σ]/(K²-m²)

H_s=(uC_S+vD_s+w)/(K+k)

H_o=(sC_o+s′D_o+w)/(K+k)

(3) amendment of hot spot-effect：Due to the problem of solar radiation, the particularly soil in sparse vegetation region and vegetation The temperature difference is quite big, and this is related to its physical surface and Meteorological, it is therefore desirable to a process for temperature adjustmemt, while being likely to Extra Vegetation canopy information is provided.The makeover process of hot spot-effect is：

Focus correction is carried out according to 2/ (k+K) effect, first according to given focus value calculating parameter：

The initial value of given alf is alf=10⁶,

If 1) focus value q is more than 0, then it is assumed that be there is the hot spot-effect to influence, now Wherein, θ_sIt is solar zenith angle, θ_oFor view zenith angle,It is the sun Azimuth；(its alf value is 200 if α result of calculations are more than 200)；If the result of calculation of alf is for 0, shadow-free is calculated pure Parameter τ under the influence of hot spot-effect_ssooAnd sumint, computing formula is：

τ_ssoo=τ_ss

If the result of calculation of alf is not 0, then it is assumed that be, without hot spot-effect influence, τ to be calculated by the method in 2)_ssooWith sumint；

If 2) focus value q is 0, then it is assumed that be that, without hot spot-effect influence, the method being added using circulation is calculated without focus Parameter τ under effects_ssooAnd sumint：Concretely comprise the following steps:

The first step provides the initial parameter in circulation：

Fhot=lai*sqrt (K*k)；

X1=0；

Y1=0；

F1=1；

Fint=(- exp (- alf)) * 0.05；

Sumint=0；

Second step calculating parameter：

Finally make τ_ssoo=f1；

* in above-mentioned code represents multiplication；

In program calculating, intermediate parameters x2 is made according to given formula result of calculation in preceding 19 circulations and is substituted into below Sumint is calculated in formula, intermediate parameters x2 values are made in the 20th circulation to be calculated in 1 and substitution formula below Sumint, makes τ_ssooEqual to the f1 that the 20th cycle calculations go out.

(4) bidirectional reflectance is calculated：Bidirectional reflectance includes two parts, part band hot spot-effect, the non-band of another part Hot spot-effect, two parts are summed to final result of calculation:

ρ_so2=ρ_so+w*lai*sum

(5) introducing of soil reflex：Soil is located at the bottommost of whole canopy reflecting system, and electromagnetic wave penetrates canopy After can be mapped to soil ground, be reflected back again after reflection canopy top, formed canopy reflectance spectrum a part, numerical procedure For：1-r_s*ρ_dd

(6) double hemispherical reflectance factors are calculated：

(7) sun direct projection direction hemispherical reflectance is calculated：

ρ_sd2=ρ_dd+(τ_sd+τ_ss)*r_s*τ_dd/ (1-rs* ρ_dd)

(8) observed direction hemispherical reflectance is calculated：

ρ_do2=ρ_do+(τ_do+τ_oo)*r_s*τ_dd/ (1-rs* ρ_dd)

(9) double direct projection reflectivity are calculated：

ρ_sod2=ρ_so+((τ_ss+τ_sd)*τ_do+(τ_sd+τ_ss*r_s*ρ_dd)*τ_oo)*r_s/ (1-rs* ρ_dd)

ρ_so3=ρ_so2+ρ_so+τ_ssoo*r_s

S5：Calculate canopy reflectance spectrum：

Direct sunlight line can be tried to achieve with the fundamental radiation amount of atmospheric scattering light according to database data, respectively E_SAnd E_d, secondly calculate atmospheric scattering coefficient skyl parameters：

Sskyl=0.847-1.61*sin (90- θ_s)+1.04*sin(90-θ_s)²

The ratio of the amount of radiation that the computational methods of canopy reflectance spectrum are measured to for observed direction and the amount of radiation of incidence, formula For：

The vegetation reflectivity such as Fig. 3 institutes for calculating.

This High-grade highway area vegetation physical and chemical parameter chosen and contrast spectrum used are measured value, modeling knot Fruit determines availability and precision problem using RMSE.

After testing, the spectral information root-mean-square error (RMSE) of model this example simulation is 0.15706, simulation algorithm simulation Result is available and precision is higher.

Described example is a part of example of the invention, rather than whole examples, it is impossible to be interpreted as to of the invention Limitation.Based on the example in the present invention, those of ordinary skill in the art are obtained on the premise of innovative labor is not made Every other implementation method belong to protection scope of the present invention.

Claims (6)

1. a kind of computational methods of broad-leaved Vegetation canopy reflectivity, it is characterised in that comprise the following steps：

S1：Parameter is recognized；

Input model parameter；And the parameter of input is tentatively divided into three major types：Blade parameter, soil parameters and Crown canopy parametre；

S2：Blade parameter input PROSPECT models in step S1 carry out monolithic leaf spectral simulation, calculate individual blade Spectral reflectivity and transmissivity；

S3：The spectral reflectivity of the individual blade that soil parameters and Crown canopy parametre, step S2 in step S1 are obtained and thoroughly Penetrate rate and calculate extinction coefficient and scattering coefficient；

S4：Step S3 is calculated extinction coefficient and scattering coefficient input SAIL models, calculate canopy associated reflections because Son and reflectivity；

S5：Calculate canopy reflectance spectrum.

2. computational methods of broad-leaved Vegetation canopy reflectivity according to claim 1, it is characterised in that the step S2 tools Body is comprised the following steps：

S2.1：Calculate absorption coefficient K：

$K=\Σ\frac{c\left(i\right)*k\left(i\right)}{N}$

Wherein, c (i) is the concentration of component i in blade；K (i) is the specific absorption coefficient of component i；N is blade construction parameter, table Show blade hierarchy number；N is calculated using following empirical equation：

N=(0.9*SLA+0.025)/(SLA-0.1)

Wherein, SLA is specific leaf area, refers to the leaf area of per unit dry weight；

S2.2：Calculate the transmissivity of interface：

Calculate the average transmittance t of blade of the light with solid angle α by air incidence refractive index that refractive index is 1 as n_av(α, 1, n), make t₁₂=t_av(α,1,n)；

S2.3：Calculate total transmittance τ of the light through plating media₁；

${\mathrm{\τ}}_1=(1-KD)\exp (-KD)+{\left(KD\right)}^2{\∫}_{KD}^{\mathrm{\∞}}{t}^{-1}{\exp }^{-t}dt$

Wherein, t is intermediate variable, and K is absorption coefficient, and D is vane thickness；

S2.4：Calculate the reflectivity R of individual layer flat board_αAnd transmissivity T (1)_α(1)：

$R_{\mathrm{\α}}\left(1\right)=1-t_{12}+\frac{{{\mathrm{\τ}}_1}^2{t}_{12}^2(n^2-t_{12})}{n^4-{{\mathrm{\τ}}_1}^2{(n^2-t_{12})}^2}$

$T_{\mathrm{\α}}\left(1\right)=\frac{n^2{{\mathrm{\τ}}_1}^2{t}_{12}^2}{n^4-{{\mathrm{\τ}}_1}^2{(n^2-t_{12})}^2}$

Wherein, n is the refractive index of blade；

S2.5：Calculate monolithic foliage reflectance R_αAnd transmissivity T (N)_α(N), i.e. light is through N layers of identical reflectivity and transmissivity Mass reflex rate and transmissivity after flat board：

$\begin{array}{c}{R}_{\mathrm{\α}}\left(N\right)\\ =\frac{R_{\mathrm{\α}}\left(1\right)({\mathrm{ab}}^{N-1}-a^{-1}{b}^{1-N})+\left(T_{\mathrm{\α}}\right(1\left)T_{90}\right(1)-R_{\mathrm{\α}}(1\left)R_{90}\right(1\left)\right)(b^{N-1}-b^{1-N})}{{\mathrm{ab}}^{N-1}-a^{-1}{b}^{1-N}-R_{90}\left(1\right)(b^{N-1}-b^{1-N})}\end{array}$

$T_{\mathrm{\α}}\left(N\right)=\frac{T_{\mathrm{\α}}\left(1\right)(a-a^{-1})}{{\mathrm{ab}}^{N-1}-a^{-1}{b}^{1-N}-R_{90}\left(1\right)(b^{N-1}-b^{1-N})}$

Wherein：

$\mathrm{\δ}=\sqrt{(R_{90}{\left(1\right)}^2-T_{90}{\left(1\right)}^2-1)-4T_{90}{\left(1\right)}^2}.$

3. computational methods of broad-leaved Vegetation canopy reflectivity according to claim 2, it is characterised in that the step S3 tools Body is comprised the following steps：

S3.1 shadow compensations：

1. the geometrical factor related to delustring and scattering is calculated：

First according to general Leaf angle inclination distribution probability weight and discretization, one group of Leaf inclination of discretization (i.e. horizontal plane method is obtained Line direction and blade normal angular separation)；Afterwards, the corresponding delustring of each Leaf inclination and dispersion factor are calculated；For Leaf inclination θ_l, computational methods are as follows：

First, critical angle β is sought_sAnd β_o, computing formula is：

${\mathrm{\β}}_s=arccos(-\frac{{\mathrm{cos\θ}}_l{\mathrm{cos\θ}}_s}{{\mathrm{sin\θ}}_l{\mathrm{sin\θ}}_s})$

${\mathrm{\β}}_o=arccos(-\frac{{\mathrm{cos\θ}}_l{\mathrm{cos\θ}}_o}{{\mathrm{sin\θ}}_l{\mathrm{sin\θ}}_o})$

Wherein, θ_sIt is solar zenith angle, θ_oIt is view zenith angle；

Two values of critical angle are directly calculated when it is determined that denominator is not less than 1 for the result of calculation of 0 and cos β in computing formula； When the result of calculation of cos β is equal to 1, two are faced Jie angle and are equal to π；Wherein cos β refer to β_sAnd β_oCosine；

Then, it is θ to calculate Leaf inclination_lIndividual blade sun direct projection direction and observed direction extinction coefficient, computing formula difference For：

$k\left({\mathrm{\θ}}_l\right)=\frac{2}{\mathrm{\π}}{L}^{\′}\[({\mathrm{\β}}_s-\frac{\mathrm{\π}}{2}){\mathrm{cos\θ}}_l+{\mathrm{sin\β}}_s{\mathrm{tan\θ}}_s{\mathrm{sin\θ}}_l\]$

$K\left({\mathrm{\θ}}_l\right)=\frac{2}{\mathrm{\π}}{L}^{\′}\[({\mathrm{\β}}_o-\frac{\mathrm{\π}}{2}){\mathrm{cos\θ}}_l+{\mathrm{sin\β}}_o{\mathrm{tan\θ}}_o{\mathrm{sin\θ}}_l\]$

Wherein, L '=lai/h, lai are leaf area index, and h is canopy height；

2. auxiliary azimuthal angle beta is calculated₁、β₂、β₃：Obtaining value method is as follows：

If： ψ ψ ψ

Wherein, ψ ψ are the relative bearing between solar direction and observed direction, the i.e. difference at both direction azimuth；

3. two-way dispersion coefficient ω (θ are calculated_l)：

After auxiliary azimuth is obtained, coordinate the single blade reflectivity R being calculated in S2_αAnd transmissivity T (N)_α(N) and thoroughly The rate of penetrating is calculated two-way dispersion coefficient, and computing formula is：

$\begin{array}{c}\mathrm{\ω}\left({\mathrm{\θ}}_l\right)=\frac{L^{\′}}{2\mathrm{\π}}\{\[\mathrm{\π}\mathrm{\ρ}-{\mathrm{\β}}_2(\mathrm{\τ}+\mathrm{\ρ})\](2{\cos }^2{\mathrm{\θ}}_l+{\sin }^2{\mathrm{\θ}}_l{\mathrm{tan\θ}}_s{\mathrm{tan\θ}}_o\cos \mathrm{\Ψ})+(\mathrm{\ρ}\\ +\mathrm{\τ}){\mathrm{sin\β}}_2\[\frac{2{\cos }^2{\mathrm{\θ}}_l}{{\mathrm{cos\β}}_s{\mathrm{cos\β}}_o}+{\mathrm{cos\β}}_1{\mathrm{cos\β}}_3{\sin }^2{\mathrm{\θ}}_l{\mathrm{tan\θ}}_s{\mathrm{tan\θ}}_o\]\}\end{array}$

Wherein, ρ, τ, θ_lThe reflectivity of blade, transmissivity and now corresponding Leaf inclination are represented respectively；ρ=R_α(N)；τ=T_α (N)；

S3.2 calculates the scattering coefficient after adding leaf reflectance and transmissivity：

1. the backscattering coefficient computing formula of diffusing radiation E- and E+ is：

$\mathrm{\σ}\left({\mathrm{\θ}}_l\right)=L^{\′}(\frac{\mathrm{\τ}+\mathrm{\ρ}}{2}+\frac{\mathrm{\ρ}-\mathrm{\τ}}{2}{\cos }^2{\mathrm{\θ}}_l)$

2. the forward scattering coefficient formulas of diffusing radiation E- and E+ are：

${\mathrm{\σ}}^{\′}\left({\mathrm{\θ}}_l\right)=L^{\′}(\frac{\mathrm{\τ}+\mathrm{\ρ}}{2}-\frac{\mathrm{\ρ}-\mathrm{\τ}}{2}{\cos }^2{\mathrm{\θ}}_l)$

3. beam radia E_SBackscattering coefficient computing formula be：

$s\left({\mathrm{\θ}}_l\right)=\frac{\mathrm{\ρ}+\mathrm{\τ}}{2}k\left({\mathrm{\θ}}_l\right)-\frac{\mathrm{\ρ}-\mathrm{\τ}}{2}{L}^{\′}{\cos }^2{\mathrm{\θ}}_l$

4. beam radia E_SForward scattering coefficient formulas be：

$s^{\′}\left({\mathrm{\θ}}_l\right)=\frac{\mathrm{\ρ}+\mathrm{\τ}}{2}k\left({\mathrm{\θ}}_l\right)+\frac{\mathrm{\ρ}-\mathrm{\τ}}{2}{L}^{\′}{\cos }^2{\mathrm{\θ}}_l$

5. the attenuation coefficient computing formula of diffusing radiation E- and E+ is：

$a\left({\mathrm{\θ}}_l\right)=L^{\′}(1-\frac{\mathrm{\ρ}+\mathrm{\τ}}{2}+\frac{\mathrm{\ρ}-\mathrm{\τ}}{2}{\cos }^2{\mathrm{\θ}}_l)$

6. observed direction radiates E₀Backscattering coefficient computing formula be：

$v\left({\mathrm{\θ}}_l\right)=\frac{\mathrm{\ρ}+\mathrm{\τ}}{2}K\left({\mathrm{\θ}}_l\right)+\frac{\mathrm{\ρ}-\mathrm{\τ}}{2}{L}^{\′}{\cos }^2{\mathrm{\θ}}_l$

7. observed direction radiates E₀Forward scattering coefficient formulas be：

$u\left({\mathrm{\θ}}_l\right)=\frac{\mathrm{\ρ}+\mathrm{\τ}}{2}K\left({\mathrm{\θ}}_l\right)-\frac{\mathrm{\ρ}-\mathrm{\τ}}{2}{L}^{\′}{\cos }^2{\mathrm{\θ}}_l$

Each Leaf inclination θ that S3.3 is obtained first against 1. middle discretization, calculates its corresponding general Leaf angle inclination distribution probability F (θ)：

In the case where non-spheroid shape is distributed, general Leaf angle inclination distribution probability is calculated according to average Leaf inclination, calculating process passes through Following iteration is completed：

X=2 θ

Y=LIDFasin (x)+0.5LIDFbsin (2x)

Dx=0.5 (y-x+2 θ)

X=x+dx

Until | dx | ＜ t；

Then F (θ)=2 (y+ θ)/π；

Wherein, LIDFa and LIDFb is leaf distributed constant；θ is the centrifugal pump of average Leaf inclination, and F (θ) is accumulation Leaf inclination, i.e., one As Leaf angle inclination distribution probability；

For elliposoidal distribution, general Leaf angle inclination distribution probability F (θ) is asked using Campbell Leaf inclination density functions；

S3.4 determines the model parameter of whole canopy；

1. the corresponding Leaf angle inclination distribution probability of each Leaf inclination that middle discretization is obtained be multiplied with model parameter respectively after phase again Plus, new model parameter is obtained, computing formula is：

Z=∑ F (θ_l)Z(θ_l)

Wherein, F (θ_l) it is Leaf inclination θ_lCorresponding general Leaf angle inclination distribution probability, Z (θ_l) refer to and obtain in step S3.1 and S3.2 K (θ_l)、K(θ_l)、ω(θ_l)、σ(θ_l)、σ′(θ_l)、s(θ_l)、s′(θ_l)、a(θ_l)、v(θ_l)、u(θ_l) in any one；Z phases In all coefficient ks, K, w, σ, σ ', s, s ', a, v, u in answering ground to refer to four differential equations of SAIL models for whole canopy Any one.

4. computational methods of broad-leaved Vegetation canopy reflectivity according to claim 3, it is characterised in that the step S4 tools Body is comprised the following steps：

S4.1：Calculate leaf area parameter：

τ_ss=e^-klai

τ_oo=e^-Klai

ρ_dd=(e^mlai-e^-mlai)/(h₁e^mlai-h₂e^-mlai)

τ_dd=(h₁-h₂)/(h₁e^mlai-h₂e^-mlai)

Wherein, lai represents leaf area index, is mode input parameter；

S4.2：The parameter after adding hot spot-effect is calculated, computing formula is as follows：

ρ_sd=C_S(1-τ_ssτ_dd)-D_Sρ_dd

τ_sd=D_S(τ_ss-τ_dd)-C_Sτ_ssρ_dd

ρ_do=C_o(1-τ_ooτ_dd)-D_oρ_dd

τ_do=D_o(τ_oo-τ_dd)-C_oτ_ooρ_dd

ρ_so=H_o(1-τ_ssτ_oo)-C_oτ_sdτ_oo-D_oρ_sd

Wherein, the computational methods of each intermediate variable are：

$m=\sqrt{a^2-{\mathrm{\σ}}^2}$

h₁=(a+m)/σ

h₂=(a-m)/σ=1/h₁

C_S=[s ' (k-a)-s σ]/(k²-m²)

C_o=[v (K-a)-u σ]/(K²-m²)

D_s=[- s (k+a)-s ' σ]/(k²-m²)

D_o=[- u (K+a)-v σ]/(K²-m²)

H_s=(uC_S+vD_s+w)/(K+k)

H_o=(sC_o+s′D_o+w)/(K+k)

Wherein, τ_ssooIt is hot spot-effect corrected parameter with sumint；

S4.3：Calculate bidirectional reflectance：

Bidirectional reflectance includes two parts, and scattering,single of the part with hot spot-effect influences, and another part is not with hot spot-effect Multiple Scattering influence, two parts are summed to the result of calculation of the final reflectivity influence at the top of canopy：

ρ_so2=ρ_so+w*lai*sumint

S4.4：Introduce soil reflex：

d_n=1-r_s*ρ_dd

S4.5：Calculate double hemispherical reflectance factors：

${\mathrm{\ρ}}_{dd2}={\mathrm{\ρ}}_{dd}+(x_{dd}^2*r_s)/d_n$

S4.6：Calculate sun direct projection direction hemispherical reflectance：

ρ_sd2=ρ_dd+(τ_sd+τ_ss)*r_s*τ_dd/d_n

S4.7：Calculating observation direction hemispherical reflectance：

ρ_do2=ρ_do+(τ_do+τ_oo)*r_s*τ_dd/d_n

S4.8：Calculate double direct projection reflectivity：

ρ_sod2=ρ_so+((τ_ss+τ_sd)*τ_do+(τ_sd+τ_ss*r_s*ρ_dd)*τ_oo)*r_s/d_n

ρ_so3=ρ_so2+ρ_so+τ_ssoo*r_s

Wherein, r_sIt is spectral reflectance.

5. computational methods of broad-leaved Vegetation canopy reflectivity according to claim 4, it is characterised in that the step S5 tools Body is comprised the following steps：

First, the fundamental radiation amount of direct sunlight line and atmospheric scattering light is tried to achieve according to database data, E is designated as respectively_sWith E_d；

Then, calculate scattering in representing air and radiate the coefficient skyl for accounting for global radiation ratio：

Skyl=0.847-1.61*sin (90- θ_s)+1.04*sin(90-θ_s)²

Finally, the amount of radiation that calculating observation direction is measured to and the ratio of the amount of radiation of incidence, obtain canopy reflectance spectrum, formula For：

$R=\frac{{\mathrm{\ρ}}_{do2}*skyl*E_d+{\mathrm{\ρ}}_{so3}*(1-skyl)*E_s}{skyl*E_d+(1-skyl)*E_s}.$

6. a kind of computation model of broad-leaved Vegetation canopy reflectivity, it is characterised in that the computation model is：

$R=\frac{{\mathrm{\ρ}}_{do2}*skyl*E_d+{\mathrm{\ρ}}_{so3}*(1-skyl)*E_s}{skyl*E_d+(1-skyl)*E_s};$

Wherein, R is canopy reflectance spectrum, E_sAnd E_dRespectively the fundamental radiation amount of direct sunlight line and atmospheric scattering light, skyl For scattering radiates the coefficient for accounting for global radiation ratio in representing air；ρ_do2And ρ_so3Respectively observed direction hemispherical reflectance and direct projection Reflectivity；The computational methods of model parameter broad-leaved Vegetation canopy reflectivity according to claim 5 are solved.