CN102914501A - Method for calculating extinction coefficients of three-dimensional forest canopy by using laser-point cloud - Google Patents
Method for calculating extinction coefficients of three-dimensional forest canopy by using laser-point cloud Download PDFInfo
- Publication number
- CN102914501A CN102914501A CN2012102606356A CN201210260635A CN102914501A CN 102914501 A CN102914501 A CN 102914501A CN 2012102606356 A CN2012102606356 A CN 2012102606356A CN 201210260635 A CN201210260635 A CN 201210260635A CN 102914501 A CN102914501 A CN 102914501A
- Authority
- CN
- China
- Prior art keywords
- line
- projection
- dimensional
- coefficient
- canopy
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 230000008033 biological extinction Effects 0.000 title claims abstract description 49
- 238000000034 method Methods 0.000 title claims abstract description 43
- 238000009826 distribution Methods 0.000 claims abstract description 42
- 238000005070 sampling Methods 0.000 claims abstract description 34
- 238000004458 analytical method Methods 0.000 claims abstract description 14
- 230000005855 radiation Effects 0.000 claims abstract description 12
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 11
- 230000008859 change Effects 0.000 claims abstract description 4
- 238000004364 calculation method Methods 0.000 claims description 18
- 238000012935 Averaging Methods 0.000 claims description 10
- 230000008569 process Effects 0.000 claims description 4
- 238000004088 simulation Methods 0.000 claims description 3
- 238000006467 substitution reaction Methods 0.000 claims description 3
- 230000004907 flux Effects 0.000 claims description 2
- 239000011159 matrix material Substances 0.000 claims description 2
- 230000035772 mutation Effects 0.000 claims 1
- 238000011160 research Methods 0.000 abstract description 12
- 230000006378 damage Effects 0.000 abstract description 2
- 238000007781 pre-processing Methods 0.000 abstract 1
- 238000010586 diagram Methods 0.000 description 19
- 239000010410 layer Substances 0.000 description 17
- 241000196324 Embryophyta Species 0.000 description 6
- 238000005516 engineering process Methods 0.000 description 6
- 238000000605 extraction Methods 0.000 description 5
- 241001269238 Data Species 0.000 description 3
- 230000007246 mechanism Effects 0.000 description 3
- 239000000203 mixture Substances 0.000 description 3
- 241000208140 Acer Species 0.000 description 2
- 238000005315 distribution function Methods 0.000 description 2
- 235000013399 edible fruits Nutrition 0.000 description 2
- 238000011156 evaluation Methods 0.000 description 2
- 239000000284 extract Substances 0.000 description 2
- 239000000463 material Substances 0.000 description 2
- 238000013178 mathematical model Methods 0.000 description 2
- 239000002356 single layer Substances 0.000 description 2
- NCGICGYLBXGBGN-UHFFFAOYSA-N 3-morpholin-4-yl-1-oxa-3-azonia-2-azanidacyclopent-3-en-5-imine;hydrochloride Chemical compound Cl.[N-]1OC(=N)C=[N+]1N1CCOCC1 NCGICGYLBXGBGN-UHFFFAOYSA-N 0.000 description 1
- 239000002028 Biomass Substances 0.000 description 1
- 235000014466 Douglas bleu Nutrition 0.000 description 1
- 240000001416 Pseudotsuga menziesii Species 0.000 description 1
- 235000005386 Pseudotsuga menziesii var menziesii Nutrition 0.000 description 1
- 241000209140 Triticum Species 0.000 description 1
- 235000021307 Triticum Nutrition 0.000 description 1
- 239000000654 additive Substances 0.000 description 1
- 230000000996 additive effect Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008901 benefit Effects 0.000 description 1
- 238000004173 biogeochemical cycle Methods 0.000 description 1
- 238000012512 characterization method Methods 0.000 description 1
- 239000003086 colorant Substances 0.000 description 1
- 239000012141 concentrate Substances 0.000 description 1
- 230000002950 deficient Effects 0.000 description 1
- 238000007598 dipping method Methods 0.000 description 1
- 230000009931 harmful effect Effects 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 238000011835 investigation Methods 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 108700002783 roundabout Proteins 0.000 description 1
Images
Landscapes
- Image Processing (AREA)
Abstract
The invention provides a method for calculating the extinction coefficients of a three-dimensional forest canopy under any given incident ray by using a geometric algorithm, belonging to the research field of a method of acquiring the structural parameters of the forest canopy. The method for calculating the extinction coefficients of the three-dimensional forest canopy by using the geometric algorithm comprises the following steps: acquiring and preprocessing the three-dimensional laser-point cloud data of a plant canopy; carrying out three-dimensional gridding on the point cloud data; carrying out point cloud slice arithmetic based on a voxel data structure; carrying out linear sampling analysis on the point cloud slices; calculating an average projection coefficient of a kth layer of slices; and calculating the extinction coefficient of any given incident ray in each layer of slices and the extinction coefficient of whole canopy. Compared with the traditional observing way, the invention has the characteristics of being small in working amount, and free from the contact observation and the damage of the canopy structure and radiation, and has objectiveness, high efficiency and precision; a method of extracting the three-dimensional structures and biological and physical diversity information from the laser radar data is developed to characterize the horizontal and vertical distribution change rules of the leaves.
Description
One, technical field
The present invention relates to the method that a kind of cloud data that utilizes the Three Dimensional Ground laser scanner to obtain calculates the Forest Canopy extinction coefficient, specifically, refer to a kind of improved method (flow process as shown in Figure 1) of utilizing the computational geometry algorithm to come the extinction coefficient under any given incident ray of Calculation of Three Dimensional Forest Canopy.
Two, background technology
The Forest Canopy structure occupies an important position in the interaction of soil-canopy-atmosphere, and closely related with the Exchange of material and energy in the biogeochemical cycle.Extinction coefficient is defined as unit leaf area perpendicular to the averaging projection's area on the plane of radiation direction, is the quantitative important factor of understanding Forest Canopy structure and radiation delivery.He is based on the directional distribution function (comprising inclination angle and azimuthal distribution function) at leaf inclination angle, calculates for the incident ray of any given angle.
The Leaf angle inclination distribution density function that obtains indirectly, exactly the vegetation canopy is a challenging job all the time, and traditional method mainly is to utilize the method for directly observing by blade to obtain the Leaf angle inclination distribution density function of given isolated tree.Such as: the equipment making such as people's utilization hornwork such as Norman and compass a kind of instrument (accompanying drawing 2 left sides) of measuring comparatively simply and easily blade angle and space distribution; The people such as Lang utilize the dirigibility of mechanical arm and cleverly how much transitive relations made a kind of instrument that is used for measuring blade space and angular distribution that can move freely and stretch (accompanying drawing 2 right sides).
But because its workload is large, subjectivity is strong, and restricted to the height of measured target trees, so can not be widely used in the practical study.Even more important a bit is, directly observation tends to affect even the structure of vegetation destruction canopy, and then destroys the inner and following radiation profiles situation of canopy.And with regard to forest, be difficult to utilize direct observation instrument to pursue blade in the reality and measure.Therefore, need a kind of new fast, accurately, round-about way obtains angle and the space distribution situation of canopy leaves.In recent years, also someone begins to attempt indirectly obtaining the blade tilt distribution.For example, Shilbayama in 2007 and Watanabe are published in " Estimating the mean leaf inclination angle of wheat canopies using reflected polarized light " literary composition of periodical " Plant Production Science " the 10th volume, propose to utilize the reflection polarized light to average the evaluation method at leaf inclination angle..
At present, also having a kind of comparatively popular method is to utilize the mathematics geometric model to come the approx blade actual distribution of approaching to reality trees, and wherein widely accepted is what model of spheroid that Campbell proposed in the eighties in 20th century.Its core concept is that the different distributions at each tangent plane inclination angle, spheroid surface is come the angular distribution of the actual blade of approximate expression, and regulates the variation of different angle distribution density function by the long and short axle of adjusting spheroid.The method can be used following equation expression:
Wherein
When χ=1, this distribution will become sphere and distribute, this moment Λ=2, a, b are respectively half short, semi-major axis of spheroid.But the actual distribution of approaching blade that this method just is similar in theory is not enough to describe protean truth.
The appearance of ground laser radar scanner is so that we can directly obtain the Forest Canopy structural information of high precision (grade) from three-dimensional perspective, and then the extinction coefficient of calculating canopy.At present domestic about the applied research of laser radar in forest, particularly the research for three-dimensional canopy structure also is in starting and exploratory stage.Existing research major part concentrates on utilizes ground laser radar to carry out parameter (comprising the diameter of a cross-section of a tree trunk 1.3 meters above the ground, the height of tree and biomass etc.) the extraction aspect of substantially every wooden dipping of higher or forest sample ground level.Such as: the people such as Ma Liqun had introduced laser radar in the application aspect the extraction Forest Vertical structural parameters in " application of laser radar in the estimation of Forest Vertical structural parameters " literary composition of " World Forestry research " the 1st volume in 2011; Forest application facet at laser radar, generating high-precision digital terrain model is a wherein crucial step, and the people such as Zhou Shufang have inquired into method and the exact evaluation that utilizes airborne laser radar system generating digital ground model in " DEM based on airborne laser radar data obtains and uses " literary composition of " remote sensing technology and application " the 22nd volume in 2007.In addition, the people such as Pang Yong has carried out a large amount of more deep researchs to the height of tree estimation that utilizes laser radar system to extract forest Dan Mu and sample prescription level between 2006-2009.But at present domestic that the three dimensional point cloud that utilizes laser radar to obtain is extracted research theory and the technology of canopy extinction coefficient is deep not enough, needs further to strengthen.
Three, summary of the invention
The objective of the invention is:
The three dimensional point cloud that provides a cover directly to utilize the ground laser radar system to generate in conjunction with the method for computational geometry, calculates the individual plant trees of any given incident light or the algorithm of forest sample prescription extinction coefficient from three-dimensional perspective.And provide the contactless method of quick and precisely observing the Forest Canopy structural characteristic parameter.
Principle of the present invention is as follows:
Utilize newer remote sensing technology means (Three Dimensional Ground Laser Radar Scanning system), technological means in conjunction with computational geometry, directly obtain individual plant trees under the assigned direction incident light and the extinction coefficient of forest sample size blade from the angle of three-dimensional, and then calculate the extinction coefficient of vegetation canopy.The three-dimensional mathematical model of model is accurately located the position of individual blade, adopts on this basis the computational geometry method, and three dimensional point cloud is carried out gridding, utilizes the point cloud slicing algorithm that the some cloud of differing heights is classified.Then carry out the line sampling analysis, obtain the relevant information that blade is intercepted and captured light, and then calculate the projection coefficient of each slicing layer, i.e. extinction coefficient, thus obtain the extinction coefficient of whole canopy.Finish the calculating of averaging projection's coefficient of extinction coefficient under the given incident ray and canopy integral body.
Technical scheme of the present invention mainly may further comprise the steps:
(1) at first utilizes the ground laser radar scanning system, obtain the three dimensional point cloud of vegetation canopy.Wherein comprised the space geometry of scanning impact point and the energy information that laser beam is rebounded, three-dimensional coordinate directly provides the locus coordinate information of any point, and this also is the basis of the mathematical model in the traditional optical theory.The cloud data that obtains at first carries out Image Mosaics, and manually removes the ground point cloud.
(2) three dimensional network is formatted.At first define an X in a territory, cloud sector, Y, the cartesian coordinate system of three coordinate axis of Z comprises all cloud datas wherein.And cloud data is divided into limited zonule, and to set up take volume elements (voxel) as basic data structure, this process is called three dimensional network and formats.Each volume elements determines its size by length (l) wide (w), high (h) three parameters, and all cloud datas are divided into m * n * p volume elements, wherein, and m=(X
Max-X
Min)/w, n=(Y
Max-Y
Min)/l, P=(Z
Max-Z
Min)/h.If, l=w=h, then volume elements is a cube.With the three-dimensional center point of these cloud datas, i.e. [(X
Max-X
Min)/2, (Y
Max-Y
Min)/2, (Z
Max-Z
Min)/2] be new initial point, set up cartesian coordinate system.Coordinate axis Z is the direction of growth of trunk, and perpendicular to surface level, X, Y-axis are positioned on the plane vertical with Z axis, and Z axis is axle longitudinally, and X, Y-axis are horizontal axle.Shown in accompanying drawing 3 (a, b, c).
(3) point cloud slicing algorithm.After the three dimensional network process of formatting finished, cloud data was divided in length and breadth numerous data Layers of direction, i.e. section, and cloud data can be considered the in length and breadth stack of direction section.Dropping cut slice shown in accompanying drawing 3 (d).All volume elements be made as (i, j, k) (i=1,2 ..., m; J=1,2 ... n; K=1,2 ..., p).Work as k=1, i and j are the arbitrary value in the given area, then can be expressed as ground floor or first cross section all volume elements.
For the relative position relation of approximate simulation incident ray and Forest Canopy, developed the slicing mode of another point of rotation cloud.At first, the cloud data rotation to the relative position that we will simulate, is then carried out three dimensional network to it and formatted.The such plane parallel at dropping cut slice and XY axle place, the plane parallel at vertical section and Z axis place.We suppose that direct sunlight injects from zenith direction, and are parallel with Z-direction.When the incident light direction changes, can realize by the rotation cloud data.
Specifically, the central point of three dimensional point cloud is the initial point of cartesian coordinate system, and cloud data then rotates around X, Y, Z axis according to the needs of a cloud with the section relative position.For example, cloud data around 30 ° of Y-axis rotations, is equivalent to that a cloud is cut into slices or sun incident light is cut into slices with 30 ° pitch angle.Therefore, some cloud at an arbitrary position (0 °-360 ° of levels, vertical 0 °-90 °) is cut into slices, and is called the directivity section.Any cloud data all can generate thin plate to extract three-dimensional canopy information by the directivity section.
(4) line sampling analysis.The space distribution of blade in section analyzed by the line method of sampling.The section of each position of cloud data is comprised of many line-transects, and these line-transects are perpendicular or parallel in sun incident light r (θ, β).When research during along the distribution of r direction blade, be to consider the line-transect parallel with r.If the three-dimensional parameter of each volume elements is set to 1.5 times of the laser radar sampled distance, we can guarantee that each volume elements at most only comprises a sampling point.So, volume elements (N) number of non-NULL represents the number of times that line-transect and blade intersect in every line-transect, and N obeys stochastic distribution.If greater than 1.5 times of sampled distance, when comprising a plurality of sampling point in the volume elements, we need to calculate different N values with spot projections all on the line-transect to the face parallel with the line-transect direction, as the number of times that intersects with blade.
Definition P
nWhen passing whole cloud data for line-transect, intersect n time possibility with blade.P
0What the expression line-transect passed through all is the space, and summation exponent is 0.By finding the solution the mean value of all line-transects N when passing vegetation, can intersect in the hope of all line-transects and blade in the whole cloud data mean value (m) and the variance (σ of number of times
2).This number of times that intersects not is actual blade point cloud density, but projects to the projection on the plane with the line-transect perpendicular direction.Variance (the σ of every line-transect
2The relative variation of)/average (m) can reflect that blade is in the distribution situation of real space.The Leaf positional distribution of isolated tree or standing forest can change (σ relatively according to it
2/ m=1,＞1, and＜1) being divided into accordingly three kinds of situations: regular distribution, stochastic distribution and gathering distribute.
The line sampling analysis is the key of calculating whole canopy and differing heights section extinction coefficient, is conducive to simultaneously a cloud density, radiation flux and the radiation profiles situation is isoparametric obtains.
(5) sun incident light pitch angle (θ) calculates averaging projection's coefficient of k layer section fixedly the time.When only studying the sun incident light of a direction of r (θ, β), that is to say the projection coefficient G of k layer section
k(r) refer to the projection coefficient of r (θ, the β) direction of fixing.When not considering the incident light position angle, and the k layer averaging projection coefficient G of pitch angle when being θ
k(θ), be by the integral and calculating of azimuthal angle beta on [0,2 π]:
(6) projection coefficient under any given incident ray of calculating, i.e. extinction coefficient.At first consider to cut into slices with fixed-direction r (θ, β), finish after three dimensional network formats at cloud data, all three-dimensional point clouds all are enclosed in the three-dimensional border, generate section and the line-transect of horizontal and vertical.If all sections are parallel with direction r (θ, β), so all line-transects also are parallel to direction r (θ, β), and the direction of the longitudinal axis namely represents this direction.Each volume elements can make up one perpendicular to the plane of this main shaft, and all sampling points can project on this plane.In the 3D grid framework, this projecting plane be two perpendicular to the volume elements border of line-transect direction it.If the sampling point number that comprises in each non-NULL volume elements is n:
If n=1 a., the square of each two dimension is take user's sampling interval as the length of side, and then sampling point represents the center of this square leaf area.If sampling interval is s, then the leaf area area is s * s.Projection coefficient (is s * s), with the ratio of the area of the length (l) of volume elements * width (w) for area take Sampling Distance as the length of side.
If b. n=2 can make up a line-transect and represent leaf area, foursquare unit area be multiply by for the projected length of this line-transect in projected area of blade.Projection coefficient is after the projection and the front line-transect Length Ratio of projection.
If point is at first detected all on a face or line in n 〉=3 c. whether all.If on same line-transect, then these points are decomposed (n=2) on the less yardstick.If on same, then make up a triangle and represent leaf area, projected area of blade is calculated by the projection coordinate on three summits.Three projection coordinates are made as (x
1, y
1), (x
2, y
2), (x
3, y
3), leg-of-mutton projected area A
pBe projected area of blade, by following matrix computations:
If the institute d. in the volume elements has a few not coplanar or on a line, then make up the convex closure of a three-dimensional, half of convex closure surface area namely represents the blade area (A) of reality.Then, all spot projections to the plane parallel with line-transect, are set up the convex closure of two dimension with these points, its area namely represents the blade area (A of projection
p)
In sum, projection coefficient A
p/ A is calculated by the non-NULL volume elements.The k layer is perpendicular to the projection coefficient G of incident light
k(r) can be calculated by following formula:
Wherein N is total number of k layer section non-NULL volume elements, and k is the number of plies of section.
Obtain G
k(r) after, substitution formula (2) is obtained averaging projection's coefficient of each slicing layer, and the projection coefficient of whole canopy (being extinction coefficient) is the summation of each layer section projection coefficient.
Compared with prior art, advantage of the present invention is embodied in: laser radar technique provides effective technological means for obtain the Forest Canopy structural parameters from three-dimensional perspective, carry out indirect observation in conjunction with the computational geometry method, compare with traditional observation method, workload is little, need not contact observation, do not destroy canopy structure and radiation characteristic, have objective efficient accurate characteristics; Developed the method for from laser radar data, extracting three-dimensional structure and biophysics diversity information, with the Concentration Changing Pattern characterization of blade.
Concrete beneficial effect is as follows:
The three dimensional point cloud that the present invention utilizes laser scanner to obtain take the extinction coefficient of volume elements as the basic calculation Forest Canopy, provides a kind of method that indirectly need not the observation Forest Canopy biophysical parameters that contacts.The present invention keeps sunshine from the tree crown top along the incident of longitudinal axis parallel direction by the rotation cloud data, come sunshine in the simulating reality situation at the sight of different angles irradiation tree crown, not only easy but also efficient, and can not cause any harmful effect to forest structure and radiation characteristic.The present invention not only can calculate the extinction coefficient of canopy integral body, can also calculate the extinction coefficient of tree crown differing heights section, and this is conducive to further to study light radiation the holding back and distribution mechanism of canopy inside, and the research of other biological physical parameter.The present invention observes Forest Canopy mechanism from three dimensions, its key is to utilize the geometric projection technology that three dimensional point cloud is converted to the two-dimensional grid image, calculate from different directions the extinction coefficient on canopy differing heights plane, not only greatly improve estimation precision, and provide bulk information for the research of the radiation delivery mechanism of canopy inside.
Practical application shows, the invention provides the effective ways of indirect observation Forest Canopy extinction coefficient, no matter is individual plant trees or the extinction coefficient of forest sample prescription canopy, all can directly use this method to try to achieve.The method has overcome classic method need spend plenty of time and a large amount of manpower and materials, affects the Forest Canopy architectural characteristic, and the large defective of error.Improve the efficient of vegetation biophysical parameters estimation, strengthened popularization and validity that three-dimensional laser scanning technique is used.Can serve better the resource environment research projects such as forest inventory investigation, vegetation ecological remote sensing.
Four, description of drawings
Fig. 1 is schematic flow sheet of the present invention;
Fig. 2 is traditional blade angle and space distribution observation instrument;
Fig. 3 is voxel data structural representation and dropping cut slice schematic diagram
A. two independent volume elements, and determine three parameters of its big or small length;
B. one 2 * 2 * 23 d grid space;
C. the three-dimensional space grid schematic diagram of volume elements structure;
D. based on the some cloud dropping cut slice algorithm schematic diagram of volume elements structure;
Fig. 4 is the three dimensional point cloud of a strain North America pesudotsuga taxifolia;
Fig. 5 is that three dimensional network is formatted and the point cloud slicing schematic diagram
A. for the some cloud dropping cut slice algorithm schematic diagram based on the volume elements structure of individual plant North America pesudotsuga taxifolia;
B. the volume elements three-dimensional space grid of monolayer slices (5 * 5 * 5 volume elements) schematic diagram;
Fig. 6 is the line sampling schematic diagram of three dimensional point cloud on the basis of three-dimensional network of higher
A. extraction and the analysis of " line sampling " (red rectangular parallelepiped among the figure);
B. line is sampled and is extracted the transversal section of sample;
C. the histogram analysis of sample is extracted in the line sampling.
Fig. 7 is single volume elements inner vanes extinction coefficient computing method schematic diagram
A. single volume elements inner vanes three-dimensional point cloud schematic diagram;
B. blade 3-D out inclusion fruit and face are expressed;
C. blade 3-D out inclusion fruit and point, line are expressed;
D. the two-dimentional convex hull computation result of blade projection point set on surface level and face are expressed.
Fig. 8 is the point cloud slicing layer extinction coefficient distribution characteristics schematic diagram of 14m-16m height
A. this slicing layer point cloud distribution schematic diagram;
B. the three-dimensional extinction coefficient spatial distribution characteristic of this slicing layer schematic diagram;
C. this slicing layer two-dimensional points cloud density spatial distribution feature schematic diagram;
D. this slicing layer two dimension extinction coefficient spatial distribution characteristic schematic diagram.
Five, embodiment
Below by example the present invention is further explained:
Take a strain North America pesudotsuga taxifolia (Douglas-fir) as research object (height of tree is 25m approximately), use Three Dimensional Ground laser scanner Leica ScanStation 2 (its parameter is as shown in table 1) and high-precision GPSs, the collection of three dimensional point cloud is carried out in a side in tree, and sampling interval is 2cm.Manually remove ground point cloud and other noise spot clouds, obtain the three dimensional point cloud of individual plant North America pesudotsuga taxifolia, as shown in Figure 4.
Table 1 three-dimensional laser scanner Leica ScanStation 2 parameters
After the cloud data of obtaining pesudotsuga taxifolia and pre-service, as shown in Figure 1, adopt the computational geometry technology that three dimensional point cloud is processed, the first step is that three dimensional network is formatted, and second step is the point cloud slicing based on volume elements, and the 3rd step was the line sampling analysis.
According to technical scheme steps (2), cloud data is carried out gridding, set up the voxel data structure, each voxel size is 0.5m * 0.5m * 1m.As shown in accompanying drawing 5 (b), be volume elements three-dimensional space grid (5 * 5 * 5 volume elements) schematic diagram of monolayer slices, the some Yun Yanse among the figure represents differing heights, shown in the some cloud be highly to be some cloud between 17 to 18 meters in the accompanying drawing 5 (a).
Cloud data after three dimensional network formatted is cut into slices, be the level point cloud Slicing Algorithm schematic diagram based on the volume elements structure for a strain North America pesudotsuga taxifolia such as accompanying drawing 5 (b), grey is 7 slice plane of level, the different height of different colours representative of some cloud.Shine from different directions the situation of tree crown for simulated solar irradiation, we with sun incident light be fixed on parallel with Z axis directly over, cloud data by the rotation pesudotsuga taxifolia comes the variation of simulated solar incident light and tree crown relative angle, to calculate the projection coefficient under the arbitrarily angled incident light.
According to technical scheme steps (4), the cloud data that three dimensional network is formatted carries out line sampling extraction and analysis.We carry out the line sampling extraction and analysis of vertical direction to the pesudotsuga taxifolia cloud data, shown in accompanying drawing 6 (a), red rectangular parallelepiped is " line sampling ", the transversal section that accompanying drawing 6 (b) extracts sample for the line sampling, accompanying drawing 6 (c) extracts the histogram analysis of sample for the line sampling, the cloud data between visible 14m-18m is maximum.
For describing better the extinction coefficient computing method of single volume elements inner vanes, the below describes as an example of the blade of Da Ye maple example, and as shown in Figure 7,7 (a) are the three-dimensional point cloud of maple blade.According to d in the technical scheme steps (6), blade point cloud is made up a three-dimensional convex closure, shown in accompanying drawing 7 (b), for result of calculation and the face of convex closure are expressed.Accompanying drawing 7 (c) is that the result of calculation of convex closure and point, line are expressed.The surface area of Calculation of Three Dimensional convex closure, half of its surface area namely represent actual blade area (A).Then three-dimensional point cloud is projected to the vertical plane of sun incident light on (i.e. vertical with Z axis, parallel with XY axle plane), set up the convex closure of two dimension with these points, its area namely represents the blade area (A of projection
p), shown in accompanying drawing 7 (d).A
p/ A i.e. the extinction coefficient of this blade.
The algorithm that proposes according to the present invention is analyzed pesudotsuga taxifolia point cloud, and is described according to technical scheme steps (5) and (6), obtains some cloud density spatial distribution and the extinction coefficient of each volume elements in the differing heights section.Take highly as the section at 14m-16m place as example, as shown in Figure 8: Fig. 8 (a) is this slicing layer point cloud distribution schematic diagram, Fig. 8 (b) is the three-dimensional extinction coefficient spatial distribution characteristic of this slicing layer schematic diagram, and Fig. 8 (c) and 8 (d) are respectively this slicing layer two-dimensional points cloud density and extinction coefficient spatial distribution characteristic schematic diagram.
As can be seen from the above results, the present invention can also provide the three-dimensional information of tree crown leaves density vertical change and certain height horizontal distribution thereof, is that additive method hardly matches.
Claims (10)
1. method of utilizing laser point cloud Calculation of Three Dimensional Forest Canopy extinction coefficient, it mainly may further comprise the steps:
(1) three-dimensional laser point cloud data of vegetation canopy obtains and pre-service;
(2) three dimensional network of cloud data is formatted: define one with X, Y, Z is the cartesian coordinate system of axle, cloud data is divided into limited zonule, foundation is take the data structure of volume elements (voxel) as the basis, and each volume elements determines its size by length (l), wide (w), high (h) three parameters;
(3) based on the point cloud slicing algorithm of volume elements structure: after the three dimensional network process of formatting finished, cloud data was divided in length and breadth numerous data Layers of direction, i.e. section, and cloud data can be considered the in length and breadth stack of direction section; All volume elements be made as (i, j, k) (i=1,2 ..., m; J=1,2 ... n; K=1,2 ..., p); Work as k=1, i and j are the arbitrary value in the given area, then can be expressed as ground floor or first cross section all volume elements; Suppose that direct sunlight injects from zenith direction, parallel with Z-direction; By point of rotation cloud, but the variation of the relative position of approximate simulation incident ray and Forest Canopy, and some cloud at an arbitrary position (0 °-360 ° of levels, vertical 0 °-90 °) is cut into slices, and is called the directivity section;
(4) line sampling analysis: the space distribution of blade in section analyzed by the line method of sampling, and volume elements (N) number of non-NULL represents the number of times that line-transect and blade intersect in every line-transect, and N meets stochastic distribution; Definition P
nFor when the parallel line-transect of sun incident light direction passes whole cloud data, intersect n time possibility with blade; P
0What the expression line-transect passed through all is the space, and summation exponent is 0; By finding the solution the mean value of all line-transects N when passing vegetation, can intersect in the hope of all line-transects and blade in the whole cloud data mean value (m) and the variance (σ of number of times
2); Variance (the σ of every line-transect
2The relative variation of)/average (m) can reflect that blade is in the distribution situation of real space; The Leaf positional distribution of isolated tree or standing forest can change (σ relatively according to it
2/ m=1,＞1, and＜1) being divided into accordingly three kinds of situations: regular distribution, stochastic distribution and gathering distribute;
(5) sun incident light pitch angle (θ) calculates averaging projection's coefficient of k layer section fixedly the time: the projection coefficient G of k layer section
k(r) refer to the projection coefficient of r (θ, the β) direction of fixing; When not considering the incident light position angle, and the k layer averaging projection coefficient G of pitch angle when being θ
k(θ), be by the integral and calculating of azimuthal angle beta on [0,2 π]:
(6) projection coefficient under any given incident ray of calculating, be extinction coefficient: with direction r (θ, sun incident light β) is fixed as parallel with Z axis, the cloud data that relative position according to incident light and canopy is rotated is cut into slices, generate section and the line-transect of horizontal and vertical; If the sampling point number that comprises in each non-NULL volume elements is n
If n=1 a., the square of each two dimension is take user's sampling interval as the length of side, and then sampling point represents the center of this square leaf area; If sampling interval is s, then the leaf area area is s * s; Projection coefficient is: the area take Sampling Distance as the length of side (is s * s) and the ratio of the area of the length (l) of volume elements * width (w);
If b. n=2 can make up a line-transect and represent leaf area, foursquare unit area be multiply by for the projected length of this line-transect in projected area of blade; Projection coefficient is after the projection and the front line-transect Length Ratio of projection;
If point is at first detected all on a face or line in n 〉=3 c. whether all; If on same line-transect, then these points are decomposed (n=2) on the less yardstick; If on same, then make up a triangle and represent leaf area, projected area of blade is calculated by the projection coordinate on three summits; Three projection coordinates are made as (x
1, y
1), (x
2, y
2), (x
3, y
3), leg-of-mutton projected area A
pBe projected area of blade, by following matrix computations:
If the institute d. in the volume elements has a few not coplanar or on a line, then make up the convex closure of a three-dimensional, half of convex closure surface area namely represents the blade area (A) of reality; Then, all spot projections to the plane parallel with line-transect, are set up the convex closure of two dimension with these points, its area namely represents the blade area (A of projection
p)
In sum, projection coefficient A
p/ A is calculated by the non-NULL volume elements; The k layer is perpendicular to the projection coefficient G of incident light
k(r) can be calculated by following formula:
Wherein N is total number of k layer section non-NULL volume elements, and k is the number of plies of section;
Obtain G
k(r) after, substitution formula (1) is obtained averaging projection's coefficient of each slicing layer, and the projection coefficient of whole canopy (being extinction coefficient) is the summation of each layer section projection coefficient.
2. a kind of method of utilizing laser point cloud Calculation of Three Dimensional Forest Canopy extinction coefficient according to claim 1, its feature is in step (1), described three-dimensional laser point cloud data is the forest cover canopy point cloud that is obtained by the Three Dimensional Ground laser scanner, the space geometry of scanning impact point and the energy information that laser beam is rebounded have wherein been comprised, and the locus coordinate information of each point, the point cloud that obtains is carried out Image Mosaics, and manually remove the ground point cloud, as the data source of extracting canopy structure information.
3. a kind of method of utilizing laser point cloud Calculation of Three Dimensional Forest Canopy extinction coefficient according to claim 1 and 2, it is characterized in that in the step (2), define one with X in a territory, cloud sector, Y, Z is the cartesian coordinate system of coordinate axis, the cloud data that obtains is carried out three dimensional network format, cloud data is divided into limited zonule, set up the data structure based on volume elements.
4. according to claim 1 or 3 described a kind of methods of utilizing laser point cloud Calculation of Three Dimensional Forest Canopy extinction coefficient, it is characterized in that in the step (3) volume elements being divided into numerous sections; For the relative position relation of approximate simulation incident ray and Forest Canopy, direct sunlight is set as parallel with Z-direction, simulate the variation of its relative position by point of rotation cloud.
5. according to claim 1,3 or 4 described a kind of methods of utilizing laser point cloud Calculation of Three Dimensional Forest Canopy extinction coefficient, it is characterized in that in the step (3), point cloud at an arbitrary position (0 °-360 ° of levels, vertical 0 °-90 °) is cut into slices, therefore be called the directivity section.
6. according to claim 1,3,4 or 5 described a kind of methods of utilizing laser point cloud Calculation of Three Dimensional Forest Canopy extinction coefficient, it is characterized in that in the step (4), utilize the line-transect parallel with sun incident light direction that cloud data is carried out the line sampling analysis, volume elements (N) number of non-NULL represents the number of times that line-transect and blade intersect in every line-transect, and N meets stochastic distribution; By finding the solution the mean value of all line-transects N when passing vegetation, can intersect in the hope of all line-transects and blade in the whole cloud data mean value (m) and the variance (σ of number of times
2); Utilize the variance (σ of every line-transect
2The relative mutation analysis blade of)/average (m) is at the distribution situation of real space, σ
2/ m=1,＞1, distribute corresponding with regular distribution, stochastic distribution and the gathering of blade respectively with＜1.
7. a kind of method of utilizing laser point cloud Calculation of Three Dimensional Forest Canopy extinction coefficient described in according to claim 6, it is characterized in that in the step (4), the line sampling analysis is the key of calculating whole canopy and differing heights section extinction coefficient, is conducive to simultaneously a cloud density, radiation flux and the radiation profiles situation is isoparametric obtains.
8. according to claim 1, a kind of method of utilizing laser point cloud Calculation of Three Dimensional Forest Canopy extinction coefficient described in 3,4,5 or 6, it is characterized in that in the step (5), sun incident light pitch angle (θ) is fixedly the time, the projection coefficient G of k layer section
k(r) refer to the projection coefficient of r (θ, the β) direction of fixing, when not considering the incident light position angle, the coefficient G of k layer averaging projection the when pitch angle is θ
k(θ), by the integral and calculating of azimuthal angle beta on [0,2 π].
9. according to claim 1, a kind of method of utilizing laser point cloud Calculation of Three Dimensional Forest Canopy extinction coefficient described in 3,4,5,6,7 or 8, it is characterized in that in the step (6), with direction r (θ, sun incident light β) is fixed as parallel with Z axis, the cloud data that relative position according to incident light and canopy is rotated is cut into slices, if the sampling point number that comprises in each non-NULL volume elements is n
(1) during n=1, projection coefficient is the area s take Sampling Distance as the length of side * s, with the ratio of the area of the length (l) of volume elements * width (w);
(2) during n=2, projection coefficient be after the projection and projection before the line-transect Length Ratio;
(3) n 〉=3 o'clock, if all points all on same line-transect, the situation when then these points being decomposed into n=2 is calculated projection coefficient; If on same, then make up a triangle and represent leaf area, calculate projected area of blade by the projection coordinate on three summits; If not on same or the same line, then make up a three-dimensional convex closure, half of convex closure surface area is actual blade area, then with all spot projections to the plane parallel with line-transect, set up the convex closure of two dimension with these points, its area is the projected area of blade; Projection coefficient be after the projection and projection before the ratio of blade area.
10. a kind of method of utilizing laser point cloud Calculation of Three Dimensional Forest Canopy extinction coefficient described in according to claim 9, it is characterized in that in step (5) and (6), with the r (θ that calculates in (6), β) averaging projection's coefficient of each slicing layer of calculating in the projection coefficient substitution (5) of direction, addition obtains the projection coefficient of whole canopy, i.e. extinction coefficient.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201210260635.6A CN102914501B (en) | 2012-07-26 | 2012-07-26 | Method for calculating extinction coefficients of three-dimensional forest canopy by using laser-point cloud |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201210260635.6A CN102914501B (en) | 2012-07-26 | 2012-07-26 | Method for calculating extinction coefficients of three-dimensional forest canopy by using laser-point cloud |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN102914501A true CN102914501A (en) | 2013-02-06 |
| CN102914501B CN102914501B (en) | 2015-01-14 |
Family
ID=47612963
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201210260635.6A Expired - Fee Related CN102914501B (en) | 2012-07-26 | 2012-07-26 | Method for calculating extinction coefficients of three-dimensional forest canopy by using laser-point cloud |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN102914501B (en) |
Cited By (17)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105103195A (en) * | 2013-02-27 | 2015-11-25 | 阿瑞斯梅迪卡有限责任公司 | Image processing |
| CN105371789A (en) * | 2015-10-09 | 2016-03-02 | 南京大学 | Method for utilizing aviation laser point cloud to calculate effective leaf area index |
| CN106383998A (en) * | 2016-09-09 | 2017-02-08 | 厦门大学 | Ground laser radar scanning-based automatic tree breast-height diameter calculation method |
| CN107643048A (en) * | 2017-07-28 | 2018-01-30 | 北京农学院 | Factors of enumeration extraction method based on cloud data |
| CN108171761A (en) * | 2017-12-13 | 2018-06-15 | 北京大学 | A kind of point cloud inner frame coding method and device that transformation is schemed based on Fourier |
| CN109472850A (en) * | 2018-11-26 | 2019-03-15 | 广东精鹰传媒股份有限公司 | A kind of implementation method of 3 D stereo electric current light effect |
| CN110702166A (en) * | 2019-09-29 | 2020-01-17 | 北京农业信息技术研究中心 | Device for measuring plant canopy parameters |
| CN111006645A (en) * | 2019-12-23 | 2020-04-14 | 青岛黄海学院 | Unmanned aerial vehicle surveying and mapping method based on motion and structure reconstruction |
| CN111551530A (en) * | 2020-04-23 | 2020-08-18 | 江苏大学 | Method and device for acquiring three-dimensional distribution information of chlorophyll fluorescence of canopy of crop group |
| CN111581710A (en) * | 2020-05-19 | 2020-08-25 | 北京数字绿土科技有限公司 | Automatic acquiring method and device for deflection of overhead transmission line tower |
| CN111580131A (en) * | 2020-04-08 | 2020-08-25 | 西安邮电大学 | Method for identifying vehicles on expressway by three-dimensional laser radar intelligent vehicle |
| CN112254820A (en) * | 2020-10-14 | 2021-01-22 | 中国科学院空天信息创新研究院 | Discrete forest scene thermal infrared radiation transmission simulation method |
| CN112287859A (en) * | 2020-11-03 | 2021-01-29 | 北京京东乾石科技有限公司 | Object recognition method, device and system, computer readable storage medium |
| CN112946669A (en) * | 2021-02-05 | 2021-06-11 | 国际竹藤中心 | Moso bamboo forest quantity identification method and device based on foundation laser radar |
| CN112964172A (en) * | 2020-12-08 | 2021-06-15 | 聚时科技(上海)有限公司 | Aviation blade surface measuring method and measuring equipment based on structured light camera |
| CN114663786A (en) * | 2022-03-21 | 2022-06-24 | 南京林业大学 | Forest stand radiation flux calculation method based on point cloud data and computer graphics |
| CN117434517A (en) * | 2023-12-21 | 2024-01-23 | 西南林业大学 | Forest canopy height estimation method based on extinction coefficient optimization |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101520307A (en) * | 2008-02-26 | 2009-09-02 | 中国计量科学研究院 | Method for measuring tree-crown volume fractal dimension by applying three-dimensional laser image-scanning system |
| CN102305622A (en) * | 2011-06-14 | 2012-01-04 | 北京林业大学 | Arbor three-dimensional green quantity measuring method based on three-dimensional laser scanner |
| CN102608620A (en) * | 2012-03-12 | 2012-07-25 | 北京北科安地科技发展有限公司 | Laser scanning point cloud vegetation filtering method on basis of reflection strength and terrain |
-
2012
- 2012-07-26 CN CN201210260635.6A patent/CN102914501B/en not_active Expired - Fee Related
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101520307A (en) * | 2008-02-26 | 2009-09-02 | 中国计量科学研究院 | Method for measuring tree-crown volume fractal dimension by applying three-dimensional laser image-scanning system |
| CN102305622A (en) * | 2011-06-14 | 2012-01-04 | 北京林业大学 | Arbor three-dimensional green quantity measuring method based on three-dimensional laser scanner |
| CN102608620A (en) * | 2012-03-12 | 2012-07-25 | 北京北科安地科技发展有限公司 | Laser scanning point cloud vegetation filtering method on basis of reflection strength and terrain |
Non-Patent Citations (4)
| Title |
|---|
| SOPHIE MOREAU ET AL.: "Biomass quantification of Andean wetland forages using ERS satellite SAR data for optimizing livestock management", 《REMOTE SENSING OF ENVIRONMENT》 * |
| 冯永康等: "多光谱数据波段选择方法试验研究——以湖北神农架林区为例", 《遥感应用》 * |
| 罗旭: "基于三维激光扫描测绘系统的森林计测学研究", 《中国博士学位论文全文数据库》 * |
| 郑光等: "结合树龄信息的遥感森林生态系统生物量制图", 《遥感学报》 * |
Cited By (26)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105103195A (en) * | 2013-02-27 | 2015-11-25 | 阿瑞斯梅迪卡有限责任公司 | Image processing |
| CN105103195B (en) * | 2013-02-27 | 2019-06-07 | 阿瑞斯梅迪卡有限责任公司 | Image processing |
| CN105371789A (en) * | 2015-10-09 | 2016-03-02 | 南京大学 | Method for utilizing aviation laser point cloud to calculate effective leaf area index |
| CN105371789B (en) * | 2015-10-09 | 2019-02-19 | 南京大学 | A method of utilizing airborne laser point cloud computing effective leaf area index |
| CN106383998B (en) * | 2016-09-09 | 2019-03-12 | 厦门大学 | A kind of tree breast-height diameter automatic calculating method based on ground laser radar scanning |
| CN106383998A (en) * | 2016-09-09 | 2017-02-08 | 厦门大学 | Ground laser radar scanning-based automatic tree breast-height diameter calculation method |
| CN107643048A (en) * | 2017-07-28 | 2018-01-30 | 北京农学院 | Factors of enumeration extraction method based on cloud data |
| CN108171761A (en) * | 2017-12-13 | 2018-06-15 | 北京大学 | A kind of point cloud inner frame coding method and device that transformation is schemed based on Fourier |
| CN108171761B (en) * | 2017-12-13 | 2020-10-16 | 北京大学 | Point cloud intra-frame coding method and device based on Fourier image transformation |
| CN109472850A (en) * | 2018-11-26 | 2019-03-15 | 广东精鹰传媒股份有限公司 | A kind of implementation method of 3 D stereo electric current light effect |
| CN110702166A (en) * | 2019-09-29 | 2020-01-17 | 北京农业信息技术研究中心 | Device for measuring plant canopy parameters |
| CN111006645A (en) * | 2019-12-23 | 2020-04-14 | 青岛黄海学院 | Unmanned aerial vehicle surveying and mapping method based on motion and structure reconstruction |
| CN111580131A (en) * | 2020-04-08 | 2020-08-25 | 西安邮电大学 | Method for identifying vehicles on expressway by three-dimensional laser radar intelligent vehicle |
| CN111551530A (en) * | 2020-04-23 | 2020-08-18 | 江苏大学 | Method and device for acquiring three-dimensional distribution information of chlorophyll fluorescence of canopy of crop group |
| CN111581710A (en) * | 2020-05-19 | 2020-08-25 | 北京数字绿土科技有限公司 | Automatic acquiring method and device for deflection of overhead transmission line tower |
| CN111581710B (en) * | 2020-05-19 | 2021-04-13 | 北京数字绿土科技有限公司 | Automatic acquiring method and device for deflection of overhead transmission line tower |
| CN112254820A (en) * | 2020-10-14 | 2021-01-22 | 中国科学院空天信息创新研究院 | Discrete forest scene thermal infrared radiation transmission simulation method |
| CN112254820B (en) * | 2020-10-14 | 2021-04-27 | 中国科学院空天信息创新研究院 | Discrete forest scene thermal infrared radiation transmission simulation method |
| CN112287859A (en) * | 2020-11-03 | 2021-01-29 | 北京京东乾石科技有限公司 | Object recognition method, device and system, computer readable storage medium |
| CN112964172A (en) * | 2020-12-08 | 2021-06-15 | 聚时科技(上海)有限公司 | Aviation blade surface measuring method and measuring equipment based on structured light camera |
| CN112964172B (en) * | 2020-12-08 | 2022-08-26 | 聚时科技(上海)有限公司 | Aviation blade surface measuring method and measuring equipment based on structured light camera |
| CN112946669A (en) * | 2021-02-05 | 2021-06-11 | 国际竹藤中心 | Moso bamboo forest quantity identification method and device based on foundation laser radar |
| CN114663786A (en) * | 2022-03-21 | 2022-06-24 | 南京林业大学 | Forest stand radiation flux calculation method based on point cloud data and computer graphics |
| CN114663786B (en) * | 2022-03-21 | 2024-03-22 | 南京林业大学 | Stand radiation flux calculating method based on point cloud data and computer graphics |
| CN117434517A (en) * | 2023-12-21 | 2024-01-23 | 西南林业大学 | Forest canopy height estimation method based on extinction coefficient optimization |
| CN117434517B (en) * | 2023-12-21 | 2024-02-27 | 西南林业大学 | Forest canopy height estimation method based on extinction coefficient optimization |
Also Published As
| Publication number | Publication date |
|---|---|
| CN102914501B (en) | 2015-01-14 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN102914501B (en) | Method for calculating extinction coefficients of three-dimensional forest canopy by using laser-point cloud | |
| CN103413133B (en) | Automatically-extracting power line method in random laser point cloud data | |
| CN106248003B (en) | A kind of method of three-dimensional laser point cloud extraction Vegetation canopy concentration class index | |
| Zheng et al. | Computational-geometry-based retrieval of effective leaf area index using terrestrial laser scanning | |
| CN106127857B (en) | The on-board LiDAR data modeling method of integrated data driving and model-driven | |
| CN102314546B (en) | Method for estimating plant growth biomass liveweight variation based on virtual plants | |
| CN100485662C (en) | Characteristic analytical method for product point clouds surface based on dynamic data access model | |
| CN107748871A (en) | A kind of three-dimensional face identification method based on multiple dimensioned covariance description with the sparse classification of local sensitivity Riemann's core | |
| CN104699952B (en) | A kind of wetland aquatic vegetation canopy BRDF Monte Carlo models | |
| Ma et al. | Retrieving forest canopy extinction coefficient from terrestrial and airborne lidar | |
| CN109446986B (en) | Effective feature extraction and tree species identification method for tree laser point cloud | |
| CN102401898A (en) | Quantified simulation method for forest remote sensing data of synthetic aperture radar | |
| CN105371789A (en) | Method for utilizing aviation laser point cloud to calculate effective leaf area index | |
| CN105389538A (en) | Method for estimating forest leaf-area index based on point cloud hemisphere slice | |
| CN109446691A (en) | Based on laser point cloud and aerodynamic live standing tree wind resistance analysis method | |
| CN105910556A (en) | Leaf area vertical distribution information extraction method | |
| CN106408608A (en) | Method for extracting trunk diameter from ground laser radar point cloud data | |
| CN105678236B (en) | A kind of land vehicles canopy reflection of polarization modeling method | |
| CN107479065A (en) | A kind of three-dimensional structure of forest gap method for measurement based on laser radar | |
| CN103106632A (en) | Fusion method of different-accuracy three-dimension point cloud data based on mean shift | |
| CN106909750B (en) | A kind of computational methods of broad-leaved Vegetation canopy reflectivity | |
| CN112068153A (en) | Crown clearance rate estimation method based on foundation laser radar point cloud | |
| CN105844057A (en) | Fast laser scanning imaging simulation method based on light beam and triangular patch intersection | |
| Li et al. | An iterative-mode scan design of terrestrial laser scanning in forests for minimizing occlusion effects | |
| CN104915982A (en) | Canopy layer illumination distribution prediction model construction method and illumination distribution detection method |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| C14 | Grant of patent or utility model | ||
| GR01 | Patent grant | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20150114 Termination date: 20210726 |
|
| CF01 | Termination of patent right due to non-payment of annual fee |