CN105910556A - Leaf area vertical distribution information extraction method - Google Patents

Abstract

The invention discloses a leaf area vertical distribution information extraction method for a corn plant. The method comprises the steps of utilizing a laser radar as a data source, identifying a corn leaf echo point cloud by means of a normal-based differential method; voxelizing a corn plant echo point cloud, and extracting leaf area vertical distribution information by means of a canopy analysis method. The leaf area vertical distribution information extraction method has advantages of quickly and accurately extracting leaf area vertical distribution information of the corn plant, simplifying a corn plant structure parameter extracting process, improving automatic degree and furthermore realizing no damage of the corn leaf.

Description

A kind of leaf area vertical distribution information extracting method

Technical field

The present invention relates to agricultural remote sensing technical field, be more particularly to a kind of Estimating Leaf Area In Maize vertical Distributed intelligence extracting method.

Background technology

Corps leaf surface amasss vertical distribution situation, is agricultural remote sensing theorem inverting research, crops One of the important parameter in mechanism model research, crops heredity and breeding field.Blade is plant Carrying out the important supporting body of the physiological process such as photosynthesis, quantitative description blade vertically divides accurately Cloth situation coerces monitoring for crop growth monitoring, Potential hazards, the yield by estimation, hereditary effect are divided Analysis etc. has very important significance.

At present, it is big that the corps leaf surface commonly used in agricultural application amasss vertical distribution information extracting method Cause can be divided into two classes: the direct method of measurement and the indirect method of measurement.Wherein, although the direct method of measurement Certainty of measurement is higher, but needs manually gather in the crops and measure vaned area, so can only be Measuring on finite point, time-consuming, laborious and plant itself is had destructiveness, this is for for a long time The measurement of the large-scale structure characteristic in sequence dynamic monitoring and big region is unpractical, thus Correction and checking frequently as the indirect method of measurement.Indirect measurement method includes utilizing between optical instrument Connecing contact type measurement and utilize the non-contact measurement of remote sensing technology, wherein instrument contacts formula is measured Mode still measure for point is upper, though non-cpntact measurement based on remote sensing technology is metering system on face, But being difficult to penetrate corps canopy, the most current inverse model is to regard corps canopy as one Individual homosphere, the difference ignored in canopy vertical direction, it is impossible to express canopy inside and inferior leads Sheet distribution situation.

Summary of the invention

(1) to solve the technical problem that

As a example by milpa, the technical problem to be solved in the present invention is the most accurate, large area Measure corps leaf surface and amass vertical distribution information, obtain corps leaf surface and amass index, the most right Crops blade causes damage.

(2) technical scheme

In order to solve above-mentioned technical problem, the invention provides a kind of leaf area vertical distribution information Extracting method, comprises the following steps:

S1, utilize three-dimensional laser scanner to obtain milpa laser point cloud, and carry out multi-site Cloud registrates；

S2, the normal of calculating milpa echo point cloud, and utilize normal difference method, pass through Cluster operation identifies blade echo point cloud；

S3, maize leaf echo point cloud is carried out voxelization, and the transmission path according to laser is true Determine the attribute of voxel；

S4, according to grid height by plant be layered, according to each point of the property calculation of described voxel Sharp Light interception-rate on layer, laser penetration rate；Laser is calculated according to site location and some cloud coordinate Angle of incidence, sets up the G-function expressing Leaf positional distribution；

S5, based on laser light incident angle, swash Light interception-rate, laser penetration rate, G-function, level Layer height, builds Estimating Leaf Area In Maize vertical distribution model, calculates Semen Maydis blade face on each level course Long-pending size, thus obtain Estimating Leaf Area In Maize vertical distribution result.

Preferably, Multiple Scattering is first removed before described step S1 carrying out multistation point cloud registering Point cloud and rough error point cloud.

By the point outside solid space shared by whole plant before removal Multiple Scattering point cloud and rough error point cloud Cloud is removed, and then utilizes some cloud strength information to remove soil echo point cloud；Point cloud registering is to two The point cloud coordinate scanned of standing sets up transition matrix, carries out position and apart from upper conversion.

Preferably, equidistant to plant to be scanned in described step S1, angularly set up and swash Optical radar, scanning field of view angle is 360 ° × 317 °, places 6 marks around plant to be scanned Target ball；

Preferably, described step S2 specifically includes following steps:

S21, calculate the normal of pretreated milpa echo point cloud.Search maize leaf All echo point clouds in 1 p contiguous range r, and utilize method of least square plane fitting p The curved surface formed a little in vertex neighborhood scope r, calculates the phase tangent plane normal of curved surface；

S22, the normal difference calculated in two contiguous range.Two kinds of yardsticks calculate some p's Normal, setting search radius r₁<r₂<r₃<…<r_n, then the normal difference under two kinds of described yardsticks is:

$\mathrm{\&Delta;}n(p,r_1,r_2)=\frac{n(p,r_1)-n(p,r_2)}{2}$

S23, normal ambiguity eliminate.Choose z value is maximum in contiguous range point as starting point, Using the normal direction of described starting point as standard normal direction, seek the method calculating remaining described some cloud If the inner product (l of line direction and described standard normal direction_i,l_j), if (l_i,l_j) ≈ 1, the most do not change described point The normal direction of cloud, otherwise by the normal direction reverse process of described some cloud；

S24, utilize normal difference identification maize leaf echo point cloud.Judge the method under two kinds of yardsticks Line difference delta n (p, r₁,r₂) whether in the threshold range set, if in threshold range, belong to phase With the some cloud of attribute, otherwise belong to the some cloud of different attribute, thus identify blade echo point cloud and Non-photosynthetic affected tissues echo point cloud.

Preferably, described step S3 is particularly as follows: carry out voxelization by described blade echo point cloud Formula is:

$i=\mathrm{int}\left(\frac{x-x_{\min }}{\mathrm{\&Delta;}i}\right)+1,j=\mathrm{int}\left(\frac{y-y_{\min }}{\mathrm{\&Delta;}j}\right)+1,k=\mathrm{int}\left(\frac{z-z_{\min }}{\mathrm{\&Delta;}k}\right)+1$

Wherein, (i, j, k) be a cloud voxel coordinate, and int is bracket function, and (x, y, after z) being registration Some cloud coordinate, (x_min,y_min,z_min) it is that (x, y, z) minima of coordinate, (Δ i, Δ j, Δ k) they are voxels Size.

Preferably, described step S4 specifically includes following steps:

S41, milpa echo point cloud carry out horizontal slice.Set level course height Δ H, root According to milpa height H, whole strain is divided into H/ Δ H layer；

S42, the sharp Light interception-rate calculated on each level course.According to laser in plant voxel Transmission route, layering judges the attribute of plant voxel, and calculates Laser Interception in every layer of plant Number of voxel n_IK number of voxel n that (), laser pass through_P(k), and utilize following formula calculating laser to cut Obtain rate:

n_I(k)/(n_I(k)+n_P(k))

Wherein, K is the number of plies.

S43, the G-function of foundation expression Leaf positional distribution.G (θ) is to be perpendicular to beam direction (θ Angle) averaging projection's area of upper unit leaf area.

Preferably, building Estimating Leaf Area In Maize vertical distribution information extraction model in described step S5 is:

$LAD(h,\mathrm{\&Delta;}H)=\frac{cos\mathrm{\&theta;}}{G\left(\mathrm{\&theta;}\right)}\&CenterDot;\frac{1}{\mathrm{\&Delta;}H}\underset{k=m_h}{\overset{m_{h+\mathrm{\&Delta;}H}}{\&Sigma;}}\frac{n_I\left(k\right)}{n_I\left(k\right)+n_P\left(k\right)}$

Wherein, θ is laser light incident angle；Δ H is level course height；m_hAnd m_h+ΔHExist for voxel Coordinate in vertical direction, equal to h and h+ Δ H (the h=Δ k × m in rectangular coordinate system_h)；n_I(k) And n_PK () is the number of voxel that Laser Interception and laser pass through on kth layer respectively；G (θ) be It is perpendicular to averaging projection's area of the upper unit leaf area of beam direction (θ angle).

Preferably, the method for the present invention also includes step S7, utilizes leaf area vertical distribution information Extraction model calculates the leaf area on each level course of plant, and with the actual measurement of field measurement Mean absolute percentage error between leaf area carries out plant leaf area vertical distribution information retrieval Precision evaluation.

(3) beneficial effect

The invention provides a kind of leaf area vertical distribution information extracting method, the method for the present invention With ground laser radar as data source, method identification maize leaf based on normal difference is utilized to return Wave point cloud, then by all maize leaf echo point cloud voxelizations, sets up leaf area vertical distribution Information extraction model, extract leaf area vertical distribution information, can quickly, extracted with high accuracy jade Rice plant leaf area vertical distribution information, simplifies plant structural parameters and extracts process, improve Automaticity, will not cause damage to maize leaf simultaneously.

Accompanying drawing explanation

In order to be illustrated more clearly that the embodiment of the present invention or technical scheme of the prior art, below The accompanying drawing used required in embodiment or description of the prior art will be briefly described, below Accompanying drawing in description is only some embodiments of the present invention, for those of ordinary skill in the art From the point of view of, on the premise of not paying creative work, it is also possible to obtain other according to these accompanying drawings Accompanying drawing.

Fig. 1 is the flow chart of a kind of leaf area vertical distribution information extracting method of the present invention；

Fig. 2 a is that ground laser radar website of the present invention arranges schematic diagram；

Fig. 2 b be the present invention step S1 in scanning field of view angle schematic diagram；

Fig. 3 is the cloud voxelization of maize leaf echo point and leaf area vertical distribution parameter in the present invention Calculate schematic diagram；

Fig. 4 is the leaf area vertical distribution information extracting method utilizing the present invention, extracts Semen Maydis and plants The result schematic diagram of strain leaf area body vertical distribution information.

Detailed description of the invention

With embodiment, the present invention is described in further detail below in conjunction with the accompanying drawings.Following example For the present invention is described, but can not be used for limiting the scope of the present invention.

Fig. 1 is the proposed by the invention flow chart extracting Estimating Leaf Area In Maize vertical distribution information, Described content comprises the following steps:

S1, according to laser radar, obtain milpa three-dimensional point cloud, carry out multistation point cloud registering；

S2, the normal of calculating milpa echo point cloud, calculate normal difference method, utilizes poly- Generic operation obtains maize leaf echo point cloud；

S3, described maize leaf echo point cloud is carried out voxelization, and according to laser in voxel Transmission path determines the attribute of voxel；

S4, according to grid height by milpa be layered, according to each point of the property calculation of described voxel Sharp Light interception-rate on layer, laser penetration rate；Calculate laser light incident angle, set up expression blade and divide The G-function (expression Leaf positional distribution) of cloth；

S5, based on laser light incident angle, swash Light interception-rate, laser penetration rate, G-function, level Layer height, builds Estimating Leaf Area In Maize vertical distribution model, calculates Semen Maydis blade face on each level course Long-pending size, thus obtain Estimating Leaf Area In Maize vertical distribution result.

By the method for the present invention can quickly, extracted with high accuracy angles of corn plant leaves area normal divides Cloth information, simplifies plant structural parameters and extracts process, improve automaticity, the most not Maize leaf can be caused damage.

The present invention obtains Semen Maydis using RIEGL VZ-1000 laser scanner as laser radar, scanning The individual plant Semen Maydis in greenhouse, corn growing season is tasseling stage, builds leaf area vertical distribution information and carries Delivery type, utilizes said method to measure leaf area vertical distribution information.

Preferably, equidistant to plant to be scanned in described step S1, angularly set up and swash Optical radar, scanning field of view angle is 360 ° × 317 °, as shown in Fig. 2 a, 2b；For ensureing difference station Point-clouds Registration precision between point, places 6 target balls around plant to be scanned；Enter First Multiple Scattering point cloud and rough error point cloud is removed before row multistation point cloud registering.Remove described many When secondary scattering point cloud and rough error point cloud, first the some cloud outside the solid space of milpa place is gone Remove, then utilize some cloud strength information to remove soil surface echo point cloud；Described point cloud registering is The point cloud coordinate of two websites is rotated and rigid body translation, i.e. sets up transition matrix and carry out position Put and apart from upper conversion.

RiScan software is utilized to carry out multi-site Registration of Measuring Data, splicing, and by writing C++ journey Sequence removes Multiple Scattering point cloud and rough error point cloud.

Described step S2 specifically includes following steps:

S21, calculate the normal of pretreated milpa echo point cloud.The solution of some cloud normal Scheme is i.e. to analyze the characteristic vector of a covariance matrix and eigenvalue (or PCA main constituent divides Analysis), its covariance matrix creates from the point of proximity cloud of query point, utilizes PCA Solving unit normal vector, solving of the problems referred to above is converted to solve one and half by PCA The eigen vector of positive definite covariance matrix；Described some cloud expresses object three dimensional structure Coordinate is respectively x, y, z, calculates the spy of 3 × 3 covariance matrix M according to PCA Value indicative and characteristic vector, wherein minimal eigenvalue characteristic of correspondence vector is the normal of described some cloud Direction, place；

$M=\frac{1}{k}\underset{i=1}{\overset{k}{\&Sigma;}}(p_i-\stackrel{\&OverBar;}{p}){(p_i-\stackrel{\&OverBar;}{p})}^T$

Wherein, p_iFor point to be analyzed,For a p_iNeighborhood point cloud barycenter,For row Vector,Transposition for described column vector.

S22, the concept of introducing metric space, calculate the normal difference in two contiguous range.? The normal of some p, setting search radius r is calculated on two kinds of yardsticks₁<r₂<r₃<…<r_n, then described in two kinds Normal difference under yardstick is:

$\mathrm{\&Delta;}n(p,r_1,r_2)=\frac{n(p,r_1)-n(p,r_2)}{2}$

Wherein, Δ n (p, r₁,r₂) representing the difference of normal direction under two kinds of yardsticks, this difference is one Three-dimensional vector, including x, y, z coordinate, the x, y, z coordinate of some cloud normal and the method for a cloud Line curvature.

S23, normal ambiguity eliminate.For the incisal plane of certain point in a cloud, it is corresponding Normal points to the both direction that angle is 180 ° that is expressed as that may be random, and i.e. there is some cloud normal The direction ambiguity of vector.It is therefore desirable to unified normal direction a little so that it is point to same side To.Due to topographical surface relative smooth, so the normal direction between adjacent two some clouds should approximate In parallel, theoretical based on this, choose z value in contiguous range (height value, it is simply that above (x, y, z) In z value) maximum some p_iAs starting point, for a p_iWith its point of proximity cloud p_j, 2 right The normal answered is respectively l_i,l_j, take inner product (l_i,l_j), if (l_i,l_j) ≈ 1 then thinks that both normal directions refer to To same direction, otherwise it is assumed that both normal orientation are contrary, need one of them normal is negated.

S24, utilize normal difference identification maize leaf echo point cloud.Judge the method under two kinds of yardsticks Line difference delta n (p, r₁,r₂) whether in the threshold range set, if in threshold range, belong to phase With the some cloud of attribute, otherwise belong to the some cloud of different attribute, thus identify blade echo point cloud and Non-photosynthetic affected tissues echo point cloud.

As it is shown on figure 3, described step S3 is particularly as follows: carry out voxel by described blade echo point cloud The formula changed is

$i=\mathrm{int}\left(\frac{x-x_{min}}{\mathrm{\&Delta;}i}\right)+1,j=\mathrm{int}\left(\frac{y-y_{min}}{\mathrm{\&Delta;}j}\right)+1,k=\mathrm{int}\left(\frac{z-z_{min}}{\mathrm{\&Delta;}k}\right)+1$

Wherein, (i, j, k) be a cloud voxel coordinate, and int is bracket function, and (x, y z) are a cloud number According to the coordinate after registration, (x_min,y_min,z_min) it is that (x, y, z) minima of coordinate, (Δ i, Δ j, Δ k) are The size of voxel.

As shown in Figure 3, it is preferable that described step S4 specifically includes following steps:

S41, milpa echo point cloud carry out horizontal slice.Set level course height Δ H, root According to milpa height H, whole strain is divided into H/ Δ H layer；

S42, swash Light interception-rate be according to laser in plant voxel transmission path determine, really During Ding, all voxels in plant are involved in calculating；One layer of the parallel correspondence of laser scanner, It is defined as k=0 layer；The part of the part that laser beam is upwards launched, i.e. k=0+；Laser beam is downward The part launched, i.e. the part of k=0-；According to laser transmission route in plant voxel, point Layer judges the attribute of maize leaf laser point cloud voxel, and the attribute of described voxel includes that laser beam blocks Cut, laser beam passes through, laser beam arrives, and marks the attribute of each voxel, and 1 is Laser Interception, 2 pass through for laser, and 3 arrive this voxel for not having laser；Calculate Laser Interception in every layer of plant Number of voxel n_IK number of voxel n that (), laser pass through_P(k), and utilize following formula to calculate each leaf The sharp Light interception-rate of sheet level course:

n_I(k)/(n_I(k)+n_P(k))

Wherein, k is the number of plies；n_IK () is the summation of the number of voxel of Laser Interception in kth layer, i.e. Attribute is set to the summation of the voxel of 1；n_PK () is the number of voxel that in kth layer, laser passes through Summation, all of laser is launched from laser scanner station, it is considered to k=m_bLayer is all to be had been marked as All of voxel outside the voxel of attribute 1, if an at least laser beam in certain voxel Arrive, then the attribute of this voxel is labeled as 2, i.e. laser beam passes through；If in certain voxel Do not have any laser beam to arrive, then the attribute of this voxel is labeled as 3, represent there is no laser Bundle arrives this voxel；K=m_bIn Ceng, all of voxel is all assigned to the attribute of 1,2,3, due to Attribute 1 represents that laser beam is intercepted, and all laser beams through voxel that attribute is 1 are the most no longer Transmitting to next layer, transmission path stops；Then k=m is judged_bThe voxel properties of+1 layer, for The laser that laser is not intercepted by last layer, passes k=m according to it_bThe situation of+1 layer of voxel judges Attribute is the voxel of 2,3, and every layer of later judgement is by that analogy, the most all of until having traveled through Voxel；

S43, the G-function of calculation expression Leaf positional distribution.G (θ) is to be perpendicular to beam direction (θ Angle) averaging projection's area of upper unit leaf area.

Preferably, in described step S5, Estimating Leaf Area In Maize vertical distribution information extraction model is:

$LAD(h,\mathrm{\&Delta;}H)=\frac{cos\mathrm{\&theta;}}{G\left(\mathrm{\&theta;}\right)}\&CenterDot;\frac{1}{\mathrm{\&Delta;}H}\underset{k=m_h}{\overset{m_{h+\mathrm{\&Delta;}H}}{\&Sigma;}}\frac{n_I\left(k\right)}{n_1\left(k\right)+n_P\left(k\right)}$

Wherein, θ is laser light incident angle；Δ H is level course height；m_hAnd m_h+ΔHExist for voxel Coordinate in vertical direction, equal to h and h+ Δ H (the h=Δ k × m in rectangular coordinate system_h)；n_I(k) And n_PK () is the number of voxel that Laser Interception and laser pass through on kth layer respectively；G (θ) be It is perpendicular to averaging projection's area of the upper unit leaf area of beam direction (θ angle).

Preferably, the method for the present invention can comprise the further steps of:

S6, by utilizing above-mentioned model to calculate the leaf area on each level course of plant, and with open country Mean absolute percentage error between the actual measurement leaf area of outer field survey carries out plant leaf area Vertical distribution information retrieval precision evaluation, result is as shown in Figure 4.

Laser radar (Light Detection And Ranging, LiDAR) is a kind of active remote sensing Technology, by the distance between the laser pulse determination sensor and the object that are sent by laser instrument, Can be used for obtaining high accuracy crops three dimensional structure information.Laser radar has millimetre-sized measurement Precision, it is provided that a kind of nondestructive corps canopy high precision three-dimensional measurement means so that Understand large area leaf area vertical stratification in depth to be possibly realized.

The present invention with laser radar as data source, can quickly, extracted with high accuracy individuality plant Leaf area vertical distribution information.The present invention utilizes some cloud recognition methods based on normal difference, essence Really identify there is in curling, blade out-of-flatness, blade arteries and veins to the blade of the features such as leaf front side recess； Based on plant voxel, extract leaf area vertical distribution information, thus improve leaf area and vertically divide The precision of cloth parameter measurement and speed, simplify plant structural parameters and extract process, and improve Automaticity.

Embodiment of above is merely to illustrate the present invention, rather than limitation of the present invention.Although ginseng According to embodiment, the present invention is described in detail, it will be understood by those within the art that, Technical scheme is carried out various combination, amendment or equivalent, without departure from this The spirit and scope of inventive technique scheme, all should contain in the middle of scope of the presently claimed invention.

Claims (7)

1. an angles of corn plant leaves area normal distributed intelligence extracting method, it is characterised in that bag Include following steps:

S1, according to laser radar, obtain milpa three-dimensional point cloud, carry out multistation point cloud registering；

S2, the normal of calculating milpa three-dimensional point cloud, calculate normal difference method, utilizes poly- Generic operation obtains maize leaf echo point cloud；

S3, described maize leaf echo point cloud is carried out voxelization, and according to laser in voxel Transmission path determines the attribute of voxel；

S4, according to grid height by milpa be layered, form several level courses, according to institute State the sharp Light interception-rate in each layering of property calculation of voxel, laser penetration rate；Calculate laser Angle of incidence, sets up the G-function expressing Leaf positional distribution；

S5, based on laser light incident angle, swash Light interception-rate, laser penetration rate, G-function, level Layer height, builds Estimating Leaf Area In Maize vertical distribution model, calculates Semen Maydis blade face on each level course Long-pending size, thus obtain Estimating Leaf Area In Maize vertical distribution result.

Method the most according to claim 1, it is characterised in that right in described step S1 Plant to be scanned is equidistant, angularly set up laser radar, and scanning field of view angle is 360 ° × 317 °, around plant to be scanned, place 6 target balls；Carry out multi-site cloud to join Before standard, need to remove Multiple Scattering point cloud and rough error point cloud.

Method the most according to claim 1, it is characterised in that described step S2 is concrete Comprise the following steps:

S21, the normal of calculating milpa echo point cloud: search maize leaf certain point p neighborhood All echo point clouds in scope r, and utilize in method of least square plane fitting contiguous range all The curved surface of some composition, calculates the phase tangent plane normal of curved surface；

S22, the normal difference calculated in two contiguous range: on two kinds of yardsticks, calculate some p's Normal, setting search radius r₁<r₂<r₃<…<r_n, then the normal difference under two kinds of described yardsticks is: n

$\mathrm{\&Delta;}n(p,r_1,r_2)=\frac{n(p,r_1)-n(p,r_2)}{2}$

S23, normal ambiguity eliminate: choose in contiguous range the point of z value maximum as starting point, Using the normal direction of described starting point as standard normal direction, seek the method calculating remaining described some cloud If the inner product (l of line direction and described standard normal direction_i,l_j), if (l_i,l_j) ≈ 1, the most do not change described point The normal direction of cloud, otherwise by the normal direction reverse process of described some cloud；

S24, utilize normal difference identification maize leaf echo point cloud: judge the method under two kinds of yardsticks Line difference delta n (p, r₁,r₂) whether in the threshold range set, if in threshold range, belong to phase With the some cloud of attribute, otherwise belong to the some cloud of different attribute, thus identify blade echo point cloud and Non-photosynthetic affected tissues echo point cloud.

Method the most according to claim 1, it is characterised in that described step S3 is concrete For: by the formula that described blade echo point cloud carries out voxelization be

$i=\mathrm{int}\left(\frac{x-x_{\min }}{\mathrm{\&Delta;}i}\right)+1,j=\mathrm{int}\left(\frac{y-y_{\min }}{\mathrm{\&Delta;}j}\right)+1,k=\mathrm{int}\left(\frac{z-z_{\min }}{\mathrm{\&Delta;}k}\right)+1$

Wherein, (i, j, k) be a cloud voxel coordinate, and int is bracket function, and (x, y z) are a cloud number According to the coordinate after registration, (x_min,y_min,z_min) it is that (x, y, z) minima of coordinate, (Δ i, Δ j, Δ k) are The size of voxel.

Method the most according to claim 1, it is characterised in that described step S4 is concrete Comprise the following steps:

S41, milpa echo point cloud carry out horizontal slice: set level course height Δ H, root According to milpa height H, whole strain is divided into H/ Δ H layer；

S42, the sharp Light interception-rate calculated on each level course: according to laser in plant voxel Transmission route, layering judges the attribute of plant voxel, and calculates Laser Interception in every layer of plant Number of voxel n_IK number of voxel n that (), laser pass through_P(k), and utilize following formula calculating laser to cut Obtain rate:

n_I(k)/(n_I(k)+n_P(k))

Wherein, k is the number of plies；

S43, the G-function of calculation expression Leaf positional distribution: G (θ) are to be perpendicular in beam direction Averaging projection's area of unit leaf area, θ is laser light incident angle.

Method the most according to claim 1, it is characterised in that beautiful in described step S5 Rice leaf area vertical distribution model is:

$LAD(h,\mathrm{\&Delta;}H)=\frac{cos\mathrm{\&theta;}}{G\left(\mathrm{\&theta;}\right)}\&CenterDot;\frac{1}{\mathrm{\&Delta;}H}\underset{k=m_h}{\overset{m_{h+\mathrm{\&Delta;}H}}{\mathrm{\&Sigma;}}}\frac{n_I\left(k\right)}{n_I\left(k\right)+n_P\left(k\right)}$

Wherein, θ is laser light incident angle；Δ H is level course height；m_hAnd m_h+ΔHExist for voxel Coordinate in vertical direction, equal to h and h+ Δ H (the h=Δ k × m in rectangular coordinate system_h)；n_I(k) And n_PK () is the number of voxel that Laser Interception and laser pass through on kth layer respectively；G (θ) be It is perpendicular to averaging projection's area of unit leaf area in beam direction.

7. according to the method described in any one of claim 1 to 6, it is characterised in that also include Precision evaluation step:

Utilize leaf area vertical distribution information model, calculate the leaf on each level course of milpa Area, and with the mean absolute percentage error between the actual measurement leaf area of field measurement, Carry out angles of corn plant leaves area normal distributed intelligence extraction accuracy evaluation.