CN107330898B - Quantitative marking calculation method and system for vegetation vertical zone - Google Patents

Quantitative marking calculation method and system for vegetation vertical zone Download PDF

Info

Publication number
CN107330898B
CN107330898B CN201710467468.5A CN201710467468A CN107330898B CN 107330898 B CN107330898 B CN 107330898B CN 201710467468 A CN201710467468 A CN 201710467468A CN 107330898 B CN107330898 B CN 107330898B
Authority
CN
China
Prior art keywords
vegetation
ndvi
vertical zone
boundary
dem
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.)
Active
Application number
CN201710467468.5A
Other languages
Chinese (zh)
Other versions
CN107330898A (en
Inventor
王心源
项波
甄静
廖颖
骆磊
刘传胜
杨瑞霞
朱岚巍
常纯
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Institute of Remote Sensing and Digital Earth of CAS
Original Assignee
Institute of Remote Sensing and Digital Earth of CAS
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Institute of Remote Sensing and Digital Earth of CAS filed Critical Institute of Remote Sensing and Digital Earth of CAS
Priority to CN201710467468.5A priority Critical patent/CN107330898B/en
Publication of CN107330898A publication Critical patent/CN107330898A/en
Application granted granted Critical
Publication of CN107330898B publication Critical patent/CN107330898B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/11Region-based segmentation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00Image enhancement or restoration
    • G06T5/50Image enhancement or restoration using two or more images, e.g. averaging or subtraction
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00Image enhancement or restoration
    • G06T5/70Denoising; Smoothing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/40Analysis of texture
    • G06T7/41Analysis of texture based on statistical description of texture
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10032Satellite or aerial image; Remote sensing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20212Image combination
    • G06T2207/20221Image fusion; Image merging

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Probability & Statistics with Applications (AREA)
  • Image Processing (AREA)

Abstract

A vegetation vertical zone quantitative scoring calculation method and system. The method comprises the following steps: (1) acquiring remote sensing image data of a vegetation vertical zone area to be extracted through a satellite, and extracting corresponding slope and sloping data through a digital elevation model DEM of the area; (2) performing image processing on the image data, and performing sample area screening; (3) constructing a sample area scatter diagram, and extracting a vertical zone boundary value of initial vegetation; (4) extracting a boundary of a vegetation vertical zone based on neighborhood statistical analysis; (5) and outputting the result image and data. The system comprises a data reading module, a region selection module, a scatter diagram module, a neighborhood analysis module and an output module. The vegetation vertical zone quantitative scribing method and system have objective and accurate results and high application value.

Description

Quantitative marking calculation method and system for vegetation vertical zone
Technical Field
The application relates to the technical field of data processing, in particular to a vegetation vertical zone quantitative scribing calculation method and system.
Background
Under the influence of factors such as air temperature, precipitation, illumination, soil and the like, the vegetation type changes along with the rise of the altitude and presents a gradual change characteristic. Meanwhile, the vegetation types vary with the elevation and are different due to the difference of regions. These present significant challenges given the vertical banding pattern of the delineated vegetation.
At present, methods for extracting the vertical zone of the alpine vegetation at home and abroad are mainly divided into three types: the method comprises the traditional ground investigation method, the remote sensing interpretation method and the ecological model method.
The ground investigation method utilizes local ecological characteristics to represent regional ecological characteristics, has certain errors and has the defects of difficult sampling point acquisition, long time consumption, high cost and the like. For high mountain vegetation, these deficiencies are more evident.
The method for extracting vegetation vertical zones by using remote sensing interpretation is mainly divided into the following six methods:
1) visual interpretation method
The visual interpretation method is mainly based on high-resolution remote sensing image data and carries out visual interpretation on distribution positions of vegetation vertical belts. The method is simple and easy to implement, but has large artificial subjective influence.
2) Land utilization classification method
The land utilization classification method is mainly based on land utilization classification of remote sensing images and combined with visual interpretation, morphology, geostatistics, landscape principles and other methods to extract vegetation vertical zones. The method depends on the precision of image classification, and the quantitative scoring error of the vegetation vertical band is large.
3) Number classification method
The quantity classification method is to take the pixels of the remote sensing image data as the sample, take the spectral attributes of the pixels as the quantity characteristics of the sample, and classify the pixels on the sample line by adopting the vegetation quantity classification method. The method is complex in calculation, is only suitable for extracting the position of the boundary line of the vegetation vertical zone on a specific sample line, and cannot extract the distribution boundary line of the vegetation vertical zone which is continuously distributed in a large range.
4) Neighborhood statistical method
The neighborhood statistical method mainly divides the remote sensing image into two types; then selecting a neighborhood window with a proper size, and solving the average value of all grids in the neighborhood of the remote sensing image data after bisection; finally, defining a distribution threshold value, and quantitatively marking out a boundary of the vegetation vertical zone. The method can extract various types of vegetation vertical zone boundary lines, but has strong subjectivity on the definition of the distribution threshold value.
5) Edge detection method
The edge detection method mainly comprises the steps of classifying remote sensing image data according to vegetation types, then performing classification post-processing, morphological processing and the like on classification results, and finally automatically identifying boundaries among different vegetation types by using the edge detection method. The method is suitable for extracting the rapidly-changed vegetation vertical zone boundary line with a simple vegetation boundary line shape, has a good extraction effect on image data with large scale and low spatial resolution, and is not suitable for the gradually-changed vegetation vertical zone boundary line with a complex vegetation boundary line shape.
6) Correlation analysis method the correlation analysis method is mainly to analyze the position of the boundary of a vertical zone of vegetation by analyzing the correlation between a normalized vegetation index (NDVI) and environmental factors such as altitude. The method comprises the steps of firstly obtaining geographic attributes by utilizing ground sampling, then obtaining an NDVI value in a sampling point remote sensing image according to a geographic position, and finally analyzing ecological environment characteristics by utilizing the sampling points, wherein the method is not separated from ground sampling point investigation.
The ecological model method is not mature enough at present, and the most key factor to be considered in modeling is how to improve the application range and precision of the model. When the dynamic change of the boundary of the vegetation vertical zone is predicted, no ideal prediction model exists at present.
Disclosure of Invention
The purpose of the present application is to overcome the above disadvantages, and provide a fine calculation method for quantitatively ruling a vertical zone of vegetation and a system for quantitatively ruling a vertical zone of vegetation, which can more comprehensively reflect the change characteristics of NDVI with altitude.
The following technical scheme is adopted in the application:
the technical scheme of the fine vegetation vertical zone quantitative ruling calculation method is as follows:
a vegetation vertical zone quantitative scoring calculation method is used for quantitatively scoring vegetation vertical zone distribution characteristics, and is characterized by comprising the following steps:
(1) acquiring remote sensing image data of a vegetation vertical zone area to be extracted through a satellite, and extracting corresponding slope and sloping data through a digital elevation model DEM of the area;
(2) carrying out image processing on the remote sensing image data, and screening sample areas;
(3) constructing a sample area scatter diagram, and extracting a vertical zone boundary value of initial vegetation;
(4) extracting a boundary of a vegetation vertical zone based on neighborhood statistical analysis;
(5) and outputting the result image and data.
The vegetation vertical zone quantitative scoring calculation method further comprises the following preferred scheme:
in the step (1), the remote sensing image, the digital elevation model DEM, the slope and the slope direction need to be the same vegetation vertical zone area to be extracted.
In the step (2), the image processing of the remote sensing image data includes the following steps:
2.1.1, carrying out radiometric calibration, atmospheric correction and terrain correction preprocessing operations on the remote sensing image;
2.1.2 extracting a vegetation coverage index NDVI based on the preprocessed remote sensing image;
and 2.1.3, carrying out image fusion on the NDVI, the DEM, the gradient and the slope data.
Wherein, the calculation formula of the vegetation coverage index NDVI in the step 2.1.2 is as follows:
Figure BDA0001326259420000031
in the formula, NIR and R are reflectance values of the vegetation at a near infrared band and an infrared band, respectively.
In the step (2), the sample region screening includes the following steps:
2.2.1, carrying out slope screening on the NDVI, the DEM and the fusion data of the slope and the slope direction, and screening an area with the slope not less than 5 degrees as a vertical zone area of vegetation to be extracted;
2.2.2 carrying out slope screening on the mountain areas of the vertical zone of the vegetation to be extracted, and ensuring that the areas of the vertical zone of the vegetation to be extracted are in the same slope.
The step (3) comprises the following steps:
3.1 constructing a DEM-NDVI scatter distribution chart by utilizing a digital elevation model DEM and a vegetation coverage index NDVI extracted from a remote sensing image based on the fusion data of the vegetation vertical zone area to be extracted;
3.2, performing density segmentation on the DEM-NDVI scatter distribution diagram;
3.3 based on the DEM-NDVI scatter distribution diagram after density segmentation, obtaining the NDVI average value and the sliding average value curve of the specific altitude by using a sliding average method;
3.4 determining a proper fitting curve function;
and 3.5, analyzing the change trend of the fitted curve, finding out the demarcation points of different sections of the curve, and extracting the initial vegetation vertical zone demarcation value.
In step 3.2, the density segmentation of the DEM-NDVI scatter distribution map includes the following steps:
3.2.1 in a DEM-NDVI scatter diagram, counting the number of scatter points in an interval range by taking 0.01 as an interval for a certain specific altitude;
3.2.2 calculate the NDVI density value in the range of every 0.01 interval according to the formula
Figure BDA0001326259420000032
In the formula, Mi is an NDVI density value in the ith NDVI interval range, Ni is the number of scattering points in the ith NDVI interval range, and Nmax is the maximum number of scattering points in all the NDVI interval ranges;
3.2.3, performing the operation on all the altitudes in the DEM-NDVI scattergram to obtain a DEM-NDVI density map;
3.2.4 selecting the area with the NDVI density value more than 0.95 as the core area of the DEM-NDVI scatter distribution diagram.
In step 3.3, the moving average method includes the following steps:
3.3.1 defining a window of size 100m to find the average value of NDVI in the window;
3.3.2 sliding the window along the direction of the elevation according to the sliding distance of 5m so as to obtain the NDVI average value of the equal-spacing elevation;
3.3.3 connecting all the obtained NDVI mean points to obtain a sliding mean curve.
In step 3.4, the moving average curve is fitted with an N-th order polynomial regression equation; correlation coefficient r for fitting degree of fitting curve and DEM-NDVI sliding average value curve2To measure, N can be determined by multiple comparison experiments, and is selected such that r2The minimum N value at which the maximum is reached is taken as the N value of the fitted curve.
In step 3.5, the dividing points of the different sections of the curve are determined by solving the inflection point and the second derivative extreme point of the fitted curve function.
In the step (4), based on neighborhood statistical analysis, the extraction of the boundary of the vegetation vertical zone comprises the following contents:
4.1, acquiring an initial threshold value and a corresponding experimental threshold value range of a vertical zone boundary NDVI of the vegetation;
4.2, respectively carrying out iterative optimization within each experimental threshold range to obtain an optimal threshold of the vegetation vertical zone boundary NDVI;
4.3 performing optimal thresholding (i.e., performing image segmentation on the NDVI image based on the NDVI optimal threshold) on the sample region NDVI image (i.e., representing different types of NDVI image segmentation results by a single value) and assigning;
4.4 selecting a neighborhood window with a specific size, and obtaining a vegetation type probability map by a neighborhood analysis method;
and 4.5, defining probability distribution thresholds of different vegetation types, and extracting a final vegetation vertical zone boundary.
In the step 4.1, the initial threshold value of the boundary NDVI of the vertical vegetation zone is an NDVI value corresponding to the boundary point of different sections of the fitting curve of the vertical vegetation zone, and the experimental threshold range is an interval range in which the corresponding initial threshold value of the NDVI fluctuates by 0.1 up and down;
in the step 4.2, the maximum inter-class variance method is applied to each experimental threshold range to respectively obtain the corresponding optimum threshold of the boundary NDVI of the vegetation vertical zone.
In the step 4.3, the assigning of the NDVI image after the optimal threshold segmentation includes the following steps:
4.3.1 the assignment range of the NDVI image after the optimal threshold segmentation is [0,1 ];
4.3.2 the different categories of NDVI images are assigned a value of 1, …, 0, respectively, at equal intervals from low to high altitude. (for example, in a specific research area, the original remote sensing image, the altitude distribution characteristics and the segmented NDVI image are combined, the segmented NDVI image can be reclassified into a needle-broad mixed forest, a coniferous forest, a mountain bush meadow and a mountain bare rock bare soil, and is respectively assigned with 1, 0.5 and 0.) generally, the vegetation dense area can be assigned with the value of 1, the vegetation non-vegetation area can be assigned with the value of 0, and the intermediate transition vegetation area can be assigned with a certain value within the range of [0,1] at equal intervals.
The neighborhood window of a specific size in said step 4.4 is a circular window with a radius of 200 m.
The extraction of the boundary of the vegetation vertical zone in the step 4.5 comprises the following contents:
4.5.1 generating isolines with intervals of 0.1 based on the vegetation type probability map;
4.5.2 acquiring an initial boundary of the vegetation vertical zone in the research area based on the initial threshold of the NDVI;
4.5.3 overlapping the contour line, the initial boundary of the vegetation vertical zone and the original remote sensing image of the vegetation vertical zone area to be extracted, and selecting the contour line which is optimally matched with the initial boundary of the vegetation vertical zone and the type of vegetation visually interpreted by the remote sensing image as the optimized boundary of the vegetation vertical zone;
4.5.4, the optimized boundary line of the vegetation vertical zone is deleted by combining factors such as the theoretical distribution altitude and the slope direction of the vegetation, and the final boundary line of the vegetation vertical zone is obtained.
The application also discloses a vegetation vertical belt quantitative scoring system, and the technical scheme of the invention is as follows:
the vegetation vertical zone quantitative marking system is used for quantitatively marking vegetation vertical zone distribution characteristics and comprises a data reading module, a zone selection module, a scatter diagram module, a neighborhood analysis module and an output module; the method is characterized in that:
the data reading module is used for reading a remote sensing image and a Digital Elevation Model (DEM) of a vertical zone area of vegetation to be extracted, carrying out image processing on the read remote sensing image, extracting vegetation coverage index, gradient and slope data of the preprocessed remote sensing image, and providing NDVI, the DEM and fusion data of the gradient and the slope to the area selection module;
the region selection module performs slope screening and slope zoning on the remote sensing image processed by the data reading module so as to obtain a vegetation vertical zone region to be extracted as a research region, and provides screened vegetation vertical zone mountain region fusion data to be extracted to the scatter diagram module;
the scatter diagram module constructs a DEM-NDVI scatter diagram according to the mountain area fusion data of the vegetation vertical zone to be extracted after being screened by the selection module, obtains the NDVI of the vegetation vertical zone and a critical value of the DEM, and provides an initial threshold value of the NDVI of the boundary of the vegetation vertical zone for the neighborhood analysis module;
the neighborhood analysis module is used for acquiring a vegetation vertical zone boundary and providing final vegetation vertical zone boundary vector data and NDVI and DEM critical average values to the output module;
and the output module is used for outputting the vertical belt boundary vector data of the vegetation and the NDVI and DEM critical average values.
The vegetation vertical zone quantitative scoring system further comprises the following preferred scheme:
in the data reading module, the remote sensing image, the DEM, the gradient and the slope need to be in the same research area range.
The data reading module comprises an image preprocessing unit, an NDVI extraction unit and an image fusion unit:
the image preprocessing unit is used for carrying out radiometric calibration, atmospheric correction and terrain correction preprocessing operations on the remote sensing image;
the NDVI extraction unit is used for extracting the vegetation coverage index NDVI of the preprocessed remote sensing image;
the image fusion unit is used for carrying out image fusion on the NDVI, the DEM, the gradient and the slope data.
The NDVI is a normalized difference vegetation index, and the calculation formula is as follows:
Figure BDA0001326259420000061
in the formula, NIR and R are reflectance values of the vegetation at a near infrared band and an infrared band, respectively.
The selection area module comprises a slope screening unit and a slope screening unit:
the slope screening unit is used for carrying out slope screening on NDVI, DEM, slope and slope fusion data, and screening an area with the slope not less than 5 degrees as a vertical mountainous area of vegetation to be extracted;
the slope screening unit is used for conducting slope screening on the mountain areas of the vertical zones of the vegetation to be extracted, and ensuring that the mountain areas of the vertical zones of the vegetation to be extracted are in the same slope.
The scatter diagram module comprises a scatter distribution diagram generating unit, a density dividing unit, a moving average unit, a fitting curve unit and a demarcation point extracting unit:
the scattered point distribution diagram generating unit is used for reading the mountain area fusion data of the vegetation vertical zone to be extracted after being screened by the area selection module and constructing a DEM-NDVI scattered point distribution diagram by utilizing the DEM and the NDVI;
the density segmentation unit is used for performing density segmentation on the DEM-NDVI scatter distribution diagram;
the moving average unit is used for moving average the DEM-NDVI scatter distribution diagram after the density division, and acquiring an NDVI average value and a moving average value curve of a specific altitude;
the fitting curve unit is used for determining a proper fitting curve function;
and the demarcation point extraction unit is used for analyzing the variation trend of the fitted curve and finding out the demarcation points of different sections of the curve.
The density segmentation unit comprises a scattered point statistics subunit, a density calculation subunit, a density map subunit and a core area generation subunit;
the scatter counting subunit is used for reading the scatter at a certain specific altitude in the DEM-NDVI scatter diagram, and counting the number of the scatters within a certain altitude interval range by taking an NDVI value of 0.01 as an interval;
the density calculation subunit is used for calculating the NDVI density value within the 0.01 interval range, and the calculation formula is
Figure BDA0001326259420000062
In the formula, Mi is an NDVI density value in the ith NDVI interval range, Ni is the number of scattering points in the ith NDVI interval range, and Nmax is the maximum number of scattering points in all the NDVI interval ranges;
the density map subunit is used for performing the above operation on all the altitudes in the DEM-NDVI scattergram to obtain a DEM-NDVI density map;
and the core area generating subunit is used for selecting an area with an NDVI density value larger than 0.95 as a DEM-NDVI scatter distribution diagram core area.
The moving average unit comprises a window average value operator unit, an equidistant average value calculation subunit and a curve generation subunit;
the window average value operator unit is used for defining a window with the size of 100m to obtain the NDVI average value in the window range;
the equidistant average value operator unit is used for sliding the window along the altitude direction according to the sliding distance of 5m to obtain the NDVI average value of the equidistant altitude;
and the curve generation subunit is used for connecting all the obtained NDVI mean value points to obtain a sliding mean value curve.
In the fitting curve unit, the moving average value is curved by using an N-th order polynomial regression equationFitting the curve, and using the correlation coefficient r to fit the fitting curve with the DEM-NDVI sliding average value curve2Performing a measurement, wherein N is determined by multiple comparison experiments, and r is selected2The minimum N value at which the maximum is reached is taken as the N value of the fitted curve.
In the boundary point extraction unit, the boundary points of different sections of the fitting curve are determined by solving inflection points and second-order derivative extreme points of the fitting curve function.
The neighborhood analysis module comprises a threshold range acquisition unit, an optimal threshold extraction unit, a reclassification unit, a probability map generation unit and a boundary extraction unit;
the threshold range obtaining unit is used for obtaining an initial threshold and a corresponding experimental threshold range of a vertical belt boundary NDVI of the vegetation;
the optimal threshold extraction unit is used for respectively carrying out iterative optimization within each experimental threshold range to obtain an optimal threshold of a vegetation vertical zone boundary NDVI;
the reclassification unit is used for performing optimal threshold segmentation on the NDVI data in the sample area range, and performing reclassification and assignment on the NDVI data;
the probability map generating unit is used for selecting a neighborhood window with a specific size and acquiring a vegetation type probability map by a neighborhood analysis method;
and the boundary extraction unit is used for defining probability distribution thresholds of different vegetation types and extracting a final boundary of the vegetation vertical zone.
In the threshold range obtaining unit, the NDVI values corresponding to the demarcation points of different sections of the fitting curve determined by the demarcation point extracting unit in the scatter diagram module are used as the NDVI initial threshold of the vertical boundary of the vegetation zone, and the range of the corresponding NDVI initial threshold which fluctuates up and down by 0.1 is used as the experimental threshold range.
The optimal threshold extraction unit comprises an inter-class variance subunit, and the inter-class variance subunit is used for respectively solving the corresponding vegetation vertical zone dividing line NDVI optimal threshold in each experiment threshold range of the NDVI image in the sample area range by applying the idea of the maximum inter-class variance method.
And the reclassification unit is used for determining an assignment range and assigning the reclassified images.
In the probability map generating unit, the neighborhood window with a specific size is a circular window with a radius of 200 m.
The boundary extraction unit comprises an isoline generation subunit, an initial boundary acquisition subunit, a boundary initial optimization subunit and a boundary fine extraction subunit;
the contour line generation subunit is used for generating contour lines with the interval of 0.1 based on the vegetation type probability map acquired by the probability map generation unit;
the initial boundary acquisition subunit is used for acquiring a vegetation vertical zone initial boundary of the research area based on an NDVI initial threshold of the vegetation vertical zone boundary;
the boundary initial optimization subunit is used for superposing the isoline, the vegetation vertical zone initial boundary and the original remote sensing image of the research area, and selecting the isoline which is optimally matched with the vegetation vertical zone initial boundary and the remote sensing image visual interpretation vegetation type as the optimized vegetation vertical zone boundary;
and the boundary fine-extraction subunit is used for deleting the optimized boundary of the vegetation vertical zone by combining factors such as theoretical distribution altitude and slope direction of the vegetation to obtain the final boundary of the vegetation vertical zone.
The vegetation vertical zone quantitative ruling algorithm and system disclosed by the application have the following advantages:
1) compared with the traditional ground investigation method, the method has the advantages that all pixels in the research area are used as sample points, so that on one hand, the number of the sample points is increased, and the sample can better represent the overall distribution condition of the vegetation vertical zone; on the other hand, the time and labor input of ground investigation are reduced; in addition, the method and the device carry out slope and slope screening on the research area, ensure that the experimental sample area can more truly represent the distribution characteristics of the vegetation vertical zone in the research area, and ensure that the quantitative scoring result is more accurate and reliable.
2) Compared with the existing remote sensing interpretation method, the extracted boundary of the vegetation vertical zone is based on trend fitting and rule statistics, so that the distribution structure of the vegetation vertical zone is more regular, and the influence of artificial subjectivity is reduced.
3) Before curve fitting is carried out, the density segmentation is carried out on the scatter diagram, and random fluctuation in sample values is effectively eliminated.
4) In the invention, the idea of the maximum inter-class variance method is applied in each experimental threshold range to respectively obtain the corresponding optimum threshold of the dividing line NDVI of the vegetation vertical zone, so that the calculation amount of iterative optimization is reduced on one hand; on the other hand, the influence of noise is reduced, and the limitation of the maximum inter-class variance method is weakened; in addition, the boundary of the extracted vegetation vertical belt is more accurate.
Drawings
FIG. 1 is a schematic flow chart of a method for calculating a quantitative score of a vertical zone of vegetation according to the present disclosure;
fig. 2 is a schematic structural view of a vertical vegetation zone quantitative scoring system disclosed in the present application.
Detailed Description
As shown in the attached figure 1, the application discloses a vegetation vertical zone quantitative scoring calculation method for quantitatively scoring vegetation vertical zone distribution characteristics, which comprises the following steps:
step 1: the method comprises the steps of obtaining remote sensing image data of a vegetation vertical zone area to be extracted through a satellite, and extracting corresponding slope and sloping data through a digital elevation model DEM of the area.
In this embodiment, the remote sensing image, the DEM, the gradient and the slope direction need to be in the same research area range. In order to ensure that the experimental data are in the same research area range, the original remote sensing image and the DEM need to be cut by using the vector data of the research area. And the gradient and slope data are both cut DEM extraction results.
Step 2: and carrying out image processing on the image data, and screening sample areas.
In this embodiment, the step 2 includes the following specific steps:
step 2.1: carrying out image processing on the remote sensing image data;
step 2.2: and carrying out sample area screening on the remote sensing image after image processing.
Wherein, the step 2.1 specifically comprises the following procedures:
step 2.1.1: carrying out radiometric calibration, atmospheric correction and terrain correction preprocessing operations on the remote sensing image;
due to the influences of non-uniform sensing media, sensor parameter change, topographic relief, earth rotation, atmospheric refraction and the like, deformation errors occur to remote sensing images, and meanwhile, spectral data are distorted. Therefore, the remote sensing image data is subjected to preprocessing operations such as radiometric calibration, atmospheric correction, terrain correction, image fusion, image enhancement and the like before being used.
Step 2.1.2: and extracting the vegetation coverage index NDVI based on the preprocessed remote sensing image. The NDVI is calculated by the formula:
Figure BDA0001326259420000091
in the formula, NIR and R are reflectance values of the vegetation at a near infrared band and an infrared band, respectively.
Step 2.1.3: and carrying out image fusion on the NDVI, the DEM, the gradient and the slope data.
Step 2.2 comprises the following:
step 2.2.1: gradient screening is carried out on the fusion data, and an area with the gradient not less than 5 degrees is screened as a vertical mountainous area of vegetation to be extracted;
step 2.2.2: and (4) carrying out slope screening on the mountain areas of the vertical zones of the vegetation to be extracted, and ensuring that the mountain areas of the vertical zones of the vegetation to be extracted are in the same slope.
The experimental area is required to be a single hillside in the embodiment, the selected area is naturally transited from the bottom of the mountain to the top of the mountain, and all vegetation vertical zones are uniformly distributed. The slope of the mountain land can be divided into the following categories: flat slope (below 5 degrees), gentle slope (5-15 degrees), slope (16-25 degrees), steep slope (26-35 degrees), abrupt slope (36-45 degrees), dangerous slope (above 45 degrees) and the like. Therefore, when the sample region is screened, the gradient is required to be not less than 5 degrees on the one hand, and the same gradient direction is required on the other hand.
And step 3: and constructing a sample area scatter diagram, and extracting a vertical zone boundary value of the initial vegetation.
In this embodiment, the step 3 includes the following specific steps:
step 3.1: constructing a DEM-NDVI scatter distribution diagram by utilizing the DEM and the NDVI based on the vegetation vertical zone mountain area fusion data to be extracted;
step 3.2: performing density segmentation on the DEM-NDVI scatter distribution diagram;
at different altitudes, DEM and NDVI have different distribution characteristics. The boundary of a vegetation vertical zone does not exist in a region with lower altitude, and the NDVI value is relatively higher; the part of the mountain top with higher altitude is completely free of vegetation, and the NDVI value is relatively low; at the waists, NDVI values decrease with increasing altitude. In general, the fitted curve of NDVI as a function of altitude follows the trend of a parabola, but there is inevitably some random fluctuation. In this embodiment, the DEM-NDVI scattergram is visually observed, and it is found that the NDVI span of the core area of the visual scattergram is 0.2, so that the scattergram is obtained by counting the number of scattergrams within the interval range at intervals of 0.01, and the core area of the DEM-NDVI scattergram can be obtained by screening the density values. Thus, random fluctuations in the sample values can be effectively eliminated.
In this embodiment, the density division includes the following contents:
3.2.1 in a DEM-NDVI scatter diagram, counting the number of scatter points in an interval range by taking an NDVI value of 0.01 as an interval for a certain specific altitude;
3.2.2 calculate the NDVI density values over a 0.01 interval, using the formula
Figure BDA0001326259420000101
In the formula, Mi is an NDVI density value in the ith NDVI interval range, Ni is the number of scattering points in the ith NDVI interval range, and Nmax is the maximum number of scattering points in all the NDVI interval ranges;
3.2.3, performing the operation on all the altitudes in the DEM-NDVI scattergram to obtain a DEM-NDVI density map;
3.2.4 selecting the area with the NDVI density value more than 0.95 as the core area of the DEM-NDVI scatter distribution diagram.
Step 3.3: based on the DEM-NDVI scatter distribution diagram after density segmentation, obtaining an NDVI average value and a sliding average value curve of a specific altitude by using a sliding average method;
in this embodiment, the moving average method includes the following steps:
3.3.1 defining a window of size 100m to find the average value of NDVI in the window;
3.3.2 sliding the window along the direction of the elevation according to the sliding distance of 5m so as to obtain the NDVI average value of the equal-spacing elevation;
3.3.3 connecting all the obtained NDVI mean points to obtain a sliding mean curve.
In the embodiment, the vegetation type shows vertical terrain characteristic with rising altitude, the vertical terrain has certain fluctuation, and the moving average method can effectively eliminate random fluctuation in the observed value.
In order to determine the influence of sliding windows with different sizes on the experimental precision, experiments with different window sizes are respectively carried out on DEM-NDVI scattered points, and the fact that when the window size is within 100m, the curve variation trends of the sliding average values with different window sizes are basically consistent, only differences exist in the aspect of local curve volatility, and the smaller the window is, the larger the curve volatility is. When the window size is 100m, the curve is the most smooth, and the curve is completely overlapped with the curve of the window size within 100m, and the change trend is completely consistent. When the sliding window is larger than 100m (e.g. 200m), the curve is still smooth, but the trend of the curve has deviated from the sliding mean curve with the window size within 100 m. Thus, a window of 100m is defined to average the NDVI over the window.
In order to determine the influence of different sliding distances on the experimental precision, experiments of different sliding distances are respectively carried out on the DEM-NDVI scattering point, and the sliding curve starts to translate to a high-altitude area when the sliding distance is larger than 20 m. Theoretically, the smaller the sliding distance, the more the obtained sliding average curve is consistent with the actual situation. To ensure accuracy, the sliding window is set to 100m, and the sliding distance is set to one fourth of the vertical accuracy (20m) of the DEM, namely the sliding window for calculating the average value of the NDVI is in an elevation range of 100m, and the sliding distance is 5m each time. If the average value of the NDVI within the altitude of 3000-3100 m is taken as the average NDVI of 3050m, after sliding for 5m, the average value of the NDVI within the altitude of 3005-3105 m is taken as the average NDVI of 3055 m.
Step 3.4: determining a suitable fitted curve function;
in this embodiment, the moving average curve may be fitted using an nth order polynomial regression equation. Available correlation coefficient r of fitting degree of fitting curve and DEM-NDVI sliding average value curve2In this regard, N can be determined by multiple comparative experiments. Is selected so that r2The minimum N value at which the maximum is reached is taken as the N value of the fitted curve.
Step 3.5: and analyzing the change trend of the fitted curve, and finding out the demarcation points of different sections of the curve.
The NDVI of the same vegetation type has similarity (consistency) with the change trend of the altitude, and the NDVI shows similar slope and slope change rate on a fitting curve; in a transition region between different vegetations, the structure of the vegetations is relatively complex, so that the NDVI of the vegetations has relatively unstable variation trend along with the altitude, and the NDVI shows a slope with larger fluctuation and a slope variation rate on a fitting curve. In a low altitude area, the fluctuation range of the NDVI is small, and the NDVI is sharply reduced after a certain altitude is reached. Thus, the demarcation points for the different segments of the curve are determined by inflection points and second derivative extreme points of the fitted curve function. For the same type of vegetation, the NDVI decline rates are closer. And the demarcation point for different vegetation types should theoretically be the point at which the NDVI rate of decline is fastest.
In this embodiment, the fitted curve of the lower altitude region is similar to a straight line, and the first derivative (i.e. slope) of the region has only a little fluctuation, so it can be regarded as the approximate region of NDVI descent rate, i.e. the region of the same kind of vegetation. When the elevation rises above a certain elevation A, the NDVI descending rate begins to become fast, that is, the interior vegetation type changes, so the NDVI value corresponding to the position A can be judged as a boundary NDVI value of the vegetation vertical zone. When the elevation rises from A to B, the NDVI still becomes smaller, but the trend of the descending rate of the NDVI begins to slow down, because the dividing point of different vegetation types should be the point with the fastest descending rate of the NDVI theoretically, the corresponding NDVI value at B can be judged as another dividing line NDVI value of the vegetation vertical zone. As the altitude increases from B to C, the second derivative reaches a maximum at C, i.e. the NDVI rate of descent reaches the fastest, after which the NDVI rate of descent begins to slow. In the altitude interval above C, the NDVI descending rate tends to be flat, the fluctuation amplitude of the first derivative is similar to that of the lower altitude interval, and the fluctuation amplitude belongs to small-amplitude fluctuation, so that the NDVI value corresponding to the position where C is judged can be regarded as another boundary NDVI value of the vegetation vertical zone. Thus, all the demarcation points of different sections of the fitting curve can be found out.
And 4, step 4: and extracting the boundary of the vegetation vertical zone based on neighborhood statistical analysis.
In this embodiment, the step 4 includes the following specific steps:
step 4.1: acquiring an initial threshold value and a corresponding experimental threshold value range of a boundary NDVI of a vegetation vertical zone;
in this embodiment, the initial threshold value of the vertical belt boundary NDVI of the vegetation is the NDVI value corresponding to the demarcation point of different sections of the curve in the step 3.5, and the experimental threshold range is the range of the corresponding NDVI initial threshold value fluctuating by 0.1 up and down;
step 4.2: respectively carrying out iterative optimization within each experimental threshold range to obtain an optimal threshold of a vegetation vertical zone boundary NDVI;
in this embodiment, based on the NDVI images of the research area, the maximum inter-class variance method is applied to each experimental threshold range to obtain the corresponding optimum threshold of the NDVI of the vertical zone boundary of the vegetation.
The basic idea of the maximum inter-class variance method is to divide data into two classes by using a certain threshold, and select the threshold that maximizes the variance between the two classes as the optimal threshold. When the area difference between the two types is not large, the image can be effectively segmented; however, when the two types of areas are different greatly, the segmentation effect is poor.
In this embodiment, in order to avoid the above defect of the maximum inter-class variance, the idea of applying the maximum inter-class variance method within each experimental threshold range is selected to respectively obtain the corresponding optimal thresholds, which not only can reduce the calculation amount, but also can ensure that the area difference between two adjacent classes is not large, thereby more effectively segmenting the image. In addition, the optimal threshold extraction is carried out on the basis of the initial threshold of the vegetation vertical zone dividing line NDVI, so that the final quantitative scoring result is more accurate.
Step 4.3: performing optimal threshold segmentation on the NDVI image of the research area, and performing reclassification and assignment on the NDVI image;
in this embodiment, the assignment range is [0,1], and the assignment ranges are 1, 0.6, 0.3, and 0 for different categories of images at equal intervals from low to high according to the altitude.
Step 4.4: selecting a neighborhood window with a specific size, and obtaining a vegetation type probability map by a neighborhood analysis method;
in this embodiment, the neighborhood window of a specific size is a circular window with a radius of 200 m.
The neighborhood statistics is to take the grid to be calculated as the center, to perform statistical analysis on the grid values in a certain range (field) around the grid, and then to output the result value to the grid position to form a new grid layer. The replacement of the vertical zone of vegetation is mainly caused by the change of the elevation gradient, the grids in the vertical zone are relatively homogeneous, and the heterogeneous grids gradually increase along with the rise of the elevation until the transition is another vertical zone type of vegetation. And carrying out neighborhood statistics on the vegetation type probability map to obtain an average value in a neighborhood, so that the grid attribute values in the same vegetation vertical zone have small difference, and the grid attribute values at different vegetation vertical zones are close to the average value of two vegetation vertical zones. When the radius of the adopted field is too small compared with the transition zone, the difference between different vegetation types cannot be increased to protrude the boundary; when the radius of the adopted domain is too large, some larger plaques of the present implanted type are smoothed, and the resolution of the neighborhood statistical map is reduced.
In this embodiment, through experiments on neighborhood radii of different sizes, it is found that when a neighborhood window is a circular window with a radius of 200m, the pixel attribute values of the vegetation of the same type do not change much overall, but patches of smaller non-native vegetation zones are smoothed out, and the contrast between vertical zones of vegetation is more obvious.
Step 4.5: and defining probability distribution thresholds of different vegetation types, and extracting a final vegetation vertical zone boundary.
In this embodiment, the extraction of the vegetation vertical zone boundary includes the following contents:
4.5.1 generating isolines with the interval of 0.1 based on the vegetation type probability map obtained in the step 4.4;
4.5.2 acquiring an initial boundary of the vegetation vertical zone in the research area based on the initial threshold of the NDVI;
4.5.3 overlapping the contour line, the initial boundary of the vegetation vertical zone and the original remote sensing image of the research area, and selecting the contour line which is optimally matched with the initial boundary of the vegetation vertical zone and the type of vegetation visually interpreted by the remote sensing image as the optimized boundary of the vegetation vertical zone;
4.5.4, the optimized boundary line of the vegetation vertical zone is deleted by combining factors such as the theoretical distribution altitude and the slope direction of the vegetation, and the final boundary line of the vegetation vertical zone is obtained.
In the embodiment, the boundary of the vegetation vertical zone is extracted only according to the distribution density of the vegetation by neighborhood statistics and contour line segmentation, and according to field ground investigation, the evaporation capacity of a north slope (a shade slope) in a research area is small, the water condition is good, the vegetation is vigorous in growth, and a south slope is a leeward slope with less rainfall, so that a 'north forest south grass' landscape is formed. Therefore, there is a large misjudgment in extracting the vegetation vertical zone boundary only according to the vegetation distribution density. When the vegetation vertical zone is quantitatively scribed in the research area, the isoline obtained in the step 4.5.3 is judged and deleted by combining factors such as theoretical distribution altitude and slope direction of the vegetation, and the final boundary of the vegetation vertical zone is obtained.
And 5: and outputting the result image and data.
As shown in fig. 2, the application also discloses a vegetation vertical zone quantitative scribing system for quantitatively scribing the vegetation vertical zone distribution characteristics, which comprises a data reading module, a region selection module, a scatter diagram module, a neighborhood analysis module and an output module.
The data reading module is mainly used for reading data and primarily processing images, and provides NDVI, DEM, gradient and slope fusion data for the region selection module; the region selection module is mainly used for performing slope screening and slope direction region division on the fusion data and providing screened fusion data of the vegetation vertical zone mountainous region to be extracted to the scatter diagram module; the scatter diagram module is mainly used for constructing a DEM-NDVI scatter diagram and providing an NDVI initial threshold of a vertical zone boundary of vegetation for the neighborhood analysis module; the neighborhood analysis module is mainly used for acquiring a vegetation vertical zone boundary and providing final vegetation vertical zone boundary vector data and NDVI and DEM critical average values to the output module; and the output module is mainly used for outputting the vertical belt boundary vector data of the vegetation and the NDVI and DEM critical average values.
Module 1: and the data reading module is used for reading the remote sensing image, the digital elevation model DEM, the slope and the slope data and carrying out image processing.
In this embodiment, the remote sensing image, the DEM, the gradient and the slope direction need to be in the same research area range. In order to ensure that the experimental data are in the same research area range, the original remote sensing image and the DEM need to be cut by using the vector data of the research area. And the gradient and slope data are both cut DEM extraction results.
The module (1) comprises:
1.1 an image preprocessing unit, which is used for carrying out radiometric calibration, atmospheric correction and terrain correction preprocessing operations on the remote sensing image;
1.2NDVI extraction unit for extracting vegetation coverage index NDVI of the preprocessed remote sensing image,
and 1.3, an image fusion unit for carrying out image fusion on the NDVI, the DEM, the gradient and the slope data.
In the above unit 1.2, NDVI is a normalized differential vegetation index, and the calculation formula is:
Figure BDA0001326259420000141
in the formula, NIR and R are reflectance values of the vegetation at a near infrared band and an infrared band, respectively.
In this embodiment, due to the influences of non-uniformity of the sensing medium, variation of sensor parameters, topographic relief, earth rotation, atmospheric refraction and the like, a deformation error occurs in the remote sensing image, and the spectral data is also distorted. Therefore, the remote sensing image data is subjected to preprocessing operations such as radiometric calibration, atmospheric correction, terrain correction, image fusion, image enhancement and the like before being used.
And (3) module 2: and the region selection module is used for performing slope screening and slope zoning on the research region.
In this embodiment, the selection module includes a slope screening unit and a slope screening unit. The slope screening unit is used for carrying out slope screening on NDVI, DEM, slope and slope fusion data, and screening an area with the slope not less than 5 degrees as a vertical mountainous area of vegetation to be extracted; and the slope screening unit is used for carrying out slope screening on the mountain areas of the vertical zones of the vegetation to be extracted, and ensuring that the mountain areas of the vertical zones of the vegetation to be extracted are in the same slope.
In the region selection module, the slope screening is system automatic screening, and the slope partition is semi-automatic partition. By manually inputting N core slopes (in this embodiment, southeast and northwest slopes), the system automatically divides the data into N parts, performs subsequent operations, and outputs corresponding results.
And a module 3: and the scatter diagram module is used for constructing a DEM-NDVI scatter diagram and acquiring the critical values of the NDVI and the DEM of the vegetation vertical zone.
In this embodiment, the scattergram module includes the following units:
unit 3.1: the scattered point distribution diagram generating unit is used for reading the screened vegetation vertical zone mountain area fusion data to be extracted and constructing a DEM-NDVI scattered point distribution diagram by utilizing the DEM and the NDVI;
unit 3.2: the density segmentation unit is used for performing density segmentation on the DEM-NDVI scatter distribution diagram;
in this embodiment, the unit 3.2 includes the following sub-units:
subunit 3.2.1: the scatter counting subunit is used for reading the scatter of a certain specific altitude in the DEM-NDVI scatter diagram, and counting the number of the scatters within a certain altitude interval range by taking 0.01 as an interval;
subunit 3.2.2: a density calculating subunit for calculating the NDVI density value within the interval range of NDVI value 0.01, the calculation formula is
Figure BDA0001326259420000151
In the formula, Mi is an NDVI density value in the ith NDVI interval range, Ni is the number of scattering points in the ith NDVI interval range, and Nmax is the maximum number of scattering points in all the NDVI interval ranges;
subunit 3.2.3: the density map subunit is used for performing the operation on all the altitudes in the DEM-NDVI scattergram to obtain a DEM-NDVI density map;
subunit 3.2.4: and the core area generating subunit is used for selecting an area with the NDVI density value larger than 0.95 as the core area of the DEM-NDVI scatter distribution diagram.
In this embodiment, the DEM and NDVI have different distribution characteristics at different altitudes. The boundary of a vegetation vertical zone does not exist in a region with lower altitude, and the NDVI value is relatively higher; the part of the mountain top with higher altitude is completely free of vegetation, and the NDVI value is relatively low; at the waists, NDVI values decrease with increasing altitude. In general, the fitted curve of NDVI as a function of altitude follows the trend of a parabola, but there is inevitably some random fluctuation. In this embodiment, the DEM-NDVI scattergram is visually observed, and it is found that the NDVI span of the core area of the visual scattergram is 0.2, so that the number of scatters in the interval range is counted at intervals of NDVI values of 0.01 to obtain a scattergram, and the core area of the DEM-NDVI scattergram can be obtained by screening the density values. Thus, random fluctuations in the sample values can be effectively eliminated.
Unit 3.3: the moving average unit is used for moving average the DEM-NDVI scatter distribution diagram after the density division is carried out, and obtaining an NDVI average value and a moving average value curve of a specific altitude;
in this embodiment, the unit 3.3 includes the following sub-units:
subunit 3.3.1: the window mean value calculating subunit is used for defining a window with the size of 100m to obtain the NDVI mean value in the window range;
subunit 3.3.2: the equidistant mean value calculating subunit is used for sliding the window along the altitude direction according to the sliding distance of 5m to obtain the equidistant altitude NDVI mean value;
subunit 3.3.3: and the curve generation subunit is used for connecting all the obtained NDVI mean values to obtain a sliding mean value curve.
In the embodiment, the vegetation type shows vertical terrain characteristic with rising altitude, the vertical terrain has certain fluctuation, and the moving average method can effectively eliminate random fluctuation in the observed value.
In order to determine the influence of sliding windows with different sizes on the experimental precision, experiments with different window sizes are respectively carried out on DEM-NDVI scattered points, and the fact that when the window size is within 100m, the curve variation trends of the sliding average values with different window sizes are basically consistent, only differences exist in the aspect of local curve volatility, and the smaller the window is, the larger the curve volatility is. When the window size is 100m, the curve is the most smooth, and the curve is completely overlapped with the curve of the window size within 100m, and the change trend is completely consistent. When the sliding window is larger than 100m (e.g. 200m), the curve is still smooth, but the trend of the curve has deviated from the sliding mean curve with the window size within 100 m. Thus, a window of 100m is defined to average the NDVI over the window.
In order to determine the influence of different sliding distances on the experimental precision, experiments of different sliding distances are respectively carried out on the DEM-NDVI scattering point, and the sliding curve starts to translate to a high-altitude area when the sliding distance is larger than 20 m. Theoretically, the smaller the sliding distance, the more the obtained sliding average curve is consistent with the actual situation. To ensure accuracy, the sliding window is set to 100m, and the sliding distance is set to one fourth of the vertical accuracy (20m) of the DEM, namely the sliding window for calculating the average value of the NDVI is in an elevation range of 100m, and the sliding distance is 5m each time. If the average value of the NDVI within the altitude of 3000-3100 m is taken as the average NDVI of 3050m, after sliding for 5m, the average value of the NDVI within the altitude of 3005-3105 m is taken as the average NDVI of 3055 m.
Unit 3.4: a fitting curve unit for determining a suitable fitting curve function;
in this embodiment, the moving average curve may be fitted using an nth order polynomial regression equation. Available correlation coefficient r of fitting degree of fitting curve and DEM-NDVI sliding average value curve2In this regard, N can be determined by multiple comparative experiments. Is selected so that r2The minimum N value at which the maximum is reached is taken as the N value of the fitted curve.
Unit 3.5: and the demarcation point extraction unit is used for analyzing the change trend of the fitted curve and finding out the demarcation points of different sections of the curve.
The NDVI of the same vegetation type has similarity (consistency) with the change trend of the altitude, and the NDVI shows similar slope and slope change rate on a fitting curve; in a transition region between different vegetations, the structure of the vegetations is relatively complex, so that the NDVI of the vegetations has relatively unstable variation trend along with the altitude, and the NDVI shows a slope with larger fluctuation and a slope variation rate on a fitting curve. In a low altitude area, the fluctuation range of the NDVI is small, and the NDVI is sharply reduced after a certain altitude is reached. Thus, the demarcation points for the different segments of the curve are determined by inflection points and second derivative extreme points of the fitted curve function. For the same type of vegetation, the NDVI decline rates are closer. And the demarcation point for different vegetation types should theoretically be the point at which the NDVI rate of decline is fastest.
In this embodiment, the fitted curve of the lower altitude region is similar to a straight line, and the first derivative (i.e. slope) of the region has only a little fluctuation, so it can be regarded as the approximate region of NDVI descent rate, i.e. the region of the same kind of vegetation. When the elevation rises above a certain elevation A, the NDVI descending rate begins to become fast, that is, the interior vegetation type changes, so the NDVI value corresponding to the position A can be judged as a boundary NDVI value of the vegetation vertical zone. When the elevation rises from A to B, the NDVI still becomes smaller, but the trend of the descending rate of the NDVI begins to slow down, because the dividing point of different vegetation types should be the point with the fastest descending rate of the NDVI theoretically, the corresponding NDVI value at B can be judged as another dividing line NDVI value of the vegetation vertical zone. As the altitude increases from B to C, the second derivative reaches a maximum at C, i.e. the NDVI rate of descent reaches the fastest, after which the NDVI rate of descent begins to slow. In the altitude interval above C, the NDVI descending rate tends to be flat, the fluctuation amplitude of the first derivative is similar to that of the lower altitude interval, and the fluctuation amplitude belongs to small-amplitude fluctuation, so that the NDVI value corresponding to the position where C is judged can be regarded as another boundary NDVI value of the vegetation vertical zone. Thus, all the demarcation points of different sections of the fitting curve can be found out.
And (4) module: a neighborhood analysis module for obtaining boundary of vegetation vertical zone
In this embodiment, the neighborhood analysis module includes the following units:
unit 4.1: the threshold range acquisition unit is used for acquiring an initial threshold and a corresponding experimental threshold range of the vertical belt boundary NDVI of the vegetation;
in this embodiment, the initial threshold of the vertical belt boundary NDVI of the vegetation is the NDVI value corresponding to the boundary point of different sections of the curve in the 22.5 boundary point extraction unit, and the experimental threshold range is the range of the corresponding NDVI initial threshold fluctuating by 0.1 up and down;
unit 4.2: the optimal threshold extraction unit is used for respectively carrying out iterative optimization within each experimental threshold range to obtain an optimal threshold of the vegetation vertical zone dividing line NDVI;
in this embodiment, the optimal threshold extraction unit mainly uses the idea of the maximum inter-class variance method to respectively obtain the corresponding optimal thresholds of the vertical zone boundary NDVI of the vegetation within each experimental threshold range of the NDVI image in the research area.
The basic idea of the maximum inter-class variance method is to divide data into two classes by using a certain threshold, and select the threshold that maximizes the variance between the two classes as the optimal threshold. When the area difference between the two types is not large, the image can be effectively segmented; however, when the two types of areas are different greatly, the segmentation effect is poor.
In this embodiment, in order to avoid the above defect of the maximum inter-class variance, the idea of applying the maximum inter-class variance method within each experimental threshold range is selected to respectively obtain the corresponding optimal thresholds, which not only can reduce the calculation amount, but also can ensure that the area difference between two adjacent classes is not large, thereby more effectively segmenting the image. In addition, the optimal threshold extraction is carried out on the basis of the initial threshold of the vegetation vertical zone dividing line NDVI, so that the final quantitative scoring result is more accurate.
Unit 4.3: the reclassification unit is used for performing optimal threshold segmentation on the NDVI images of the research area and reclassifying and assigning the NDVI images;
in this embodiment, the assignment range is [0,1], and the assignment ranges are 1, 0.6, 0.3, and 0 for different categories of images at equal intervals from low to high according to the altitude.
Unit 4.4: the probability map generating unit is used for selecting a neighborhood window with a specific size and acquiring a vegetation type probability map by a neighborhood analysis method;
in this embodiment, the neighborhood window of a specific size is a circular window with a radius of 200 m.
The neighborhood statistics is to take the grid to be calculated as the center, to perform statistical analysis on the grid values in a certain range (field) around the grid, and then to output the result value to the grid position to form a new grid layer. The replacement of the vertical zone of vegetation is mainly caused by the change of the elevation gradient, the grids in the vertical zone are relatively homogeneous, and the heterogeneous grids gradually increase along with the rise of the elevation until the transition is another vertical zone type of vegetation. And carrying out neighborhood statistics on the vegetation type probability map to obtain an average value in a neighborhood, so that the grid attribute values in the same vegetation vertical zone have small difference, and the grid attribute values at different vegetation vertical zones are close to the average value of two vegetation vertical zones. When the radius of the adopted field is too small compared with the transition zone, the difference between different vegetation types cannot be increased to protrude the boundary; when the radius of the adopted domain is too large, some larger plaques of the present implanted type are smoothed, and the resolution of the neighborhood statistical map is reduced.
In this embodiment, through experiments on neighborhood radii of different sizes, it is found that when a neighborhood window is a circular window with a radius of 200m, the pixel attribute values of the vegetation of the same type do not change much overall, but patches of smaller non-native vegetation zones are smoothed out, and the contrast between vertical zones of vegetation is more obvious.
Unit 4.5: and the boundary extraction unit is used for defining probability distribution thresholds of different vegetation types and extracting a final boundary of the vegetation vertical zone.
In this embodiment, the unit 4.5 includes the following sub-units:
subunit 4.5.1: the contour line generation subunit is used for generating contour lines with the interval of 0.1 based on the vegetation type probability map acquired by the 27.4 probability map generation unit;
subunit 4.5.2: the initial boundary acquisition subunit is used for acquiring a vegetation vertical zone initial boundary of the research area based on an NDVI initial threshold of the vegetation vertical zone boundary;
subunit 4.5.3: the boundary initial optimization subunit is used for superposing the isoline, the vegetation vertical zone initial boundary and the original remote sensing image of the research area, and selecting the isoline which is optimally matched with the vegetation vertical zone initial boundary and the remote sensing image visual interpretation vegetation type as the optimized vegetation vertical zone boundary;
subunit 4.5.4: and the boundary fine extraction subunit is used for deleting the optimized boundary of the vegetation vertical zone by combining factors such as theoretical distribution altitude and slope direction of the vegetation and the like to obtain the final boundary of the vegetation vertical zone.
In the embodiment, the boundary of the vegetation vertical zone is extracted only according to the distribution density of the vegetation by neighborhood statistics and contour line segmentation, and according to field ground investigation, the evaporation capacity of a north slope (a shade slope) in a research area is small, the water condition is good, the vegetation is vigorous in growth, and a south slope is a leeward slope with less rainfall, so that a 'north forest south grass' landscape is formed. Therefore, there is a large misjudgment in extracting the vegetation vertical zone boundary only according to the vegetation distribution density. When the vegetation vertical zone is quantitatively scribed in the research area, the isoline obtained in the step 4.5.3 is judged and deleted by combining factors such as theoretical distribution altitude and slope direction of the vegetation, and the final boundary of the vegetation vertical zone is obtained.
And a module 5: and the output module is used for outputting the vector data of the vertical zone boundary of the vegetation and the NDVI and DEM critical average values.
The embodiments of the present application have been described and illustrated in detail in the drawings of the specification by the applicant, but it should be understood by those skilled in the art that the above embodiments are only the preferred embodiments of the present application, and the detailed description is only for the purpose of helping the reader to better understand the spirit of the present invention, and not for limiting the scope of the present application, but rather, any improvement or modification made based on the spirit of the present application should fall within the scope of the present application.

Claims (30)

1. A vegetation vertical zone quantitative scoring calculation method is used for quantitatively scoring vegetation vertical zone distribution characteristics, and is characterized by comprising the following steps:
(1) acquiring remote sensing image data of a vegetation vertical zone area to be extracted through a satellite, and extracting corresponding slope and sloping data through a digital elevation model DEM of the area;
(2) carrying out image processing on the remote sensing image data, extracting a vegetation coverage index NDVI from the preprocessed remote sensing image, then carrying out image fusion on the NDVI, the DEM, the slope and the slope direction data, and carrying out sample area screening;
(3) constructing a DEM-NDVI scatter point distribution diagram of a sample area, performing density segmentation on the DEM-NDVI scatter point distribution diagram, then obtaining an NDVI average value and a sliding average value curve of a specific altitude, performing curve fitting on the sliding average value through an N-th polynomial regression equation to obtain a fitting curve function, and extracting a vertical zone boundary value of initial vegetation according to the fitting curve function;
(4) extracting a boundary of a vegetation vertical zone based on neighborhood statistical analysis;
4.1, acquiring an initial threshold value and a corresponding experimental threshold value range of a vertical zone boundary NDVI of the vegetation;
4.2, respectively carrying out iterative optimization within each experimental threshold range to obtain an optimal threshold of the vegetation vertical zone boundary NDVI;
4.3 performing optimal threshold segmentation on the NDVI image in the sample area, namely NDVI data in the sample area range, and performing reclassification and assignment on the NDVI image, wherein the optimal threshold segmentation refers to performing image segmentation on the NDVI image based on an NDVI optimal threshold, and the reclassification refers to representing different types of NDVI image segmentation results by using a single value;
4.4 selecting a neighborhood window with a specific size, and obtaining a vegetation type probability map by a neighborhood analysis method;
4.5, defining probability distribution thresholds of different vegetation types, and extracting a final boundary of a vegetation vertical zone; (5) and outputting the result image and data.
2. The vegetation vertical zone quantitative scoring computation method of claim 1, wherein:
in the step (1), the remote sensing image, the digital elevation model DEM, the slope and the slope direction need to be within the same area of the vegetation vertical zone to be extracted.
3. The vegetation vertical zone quantitative scoring computation method of claim 1, wherein:
in the step (2), the image processing of the remote sensing image data includes the following steps:
2.1.1, carrying out radiometric calibration, atmospheric correction and terrain correction preprocessing operations on the remote sensing image;
2.1.2 extracting a vegetation coverage index NDVI based on the preprocessed remote sensing image;
and 2.1.3, carrying out image fusion on the NDVI, the DEM, the gradient and the slope data.
4. The vegetation vertical zone quantitative scoring computation method of claim 3, wherein:
in step 2.1.2, the vegetation coverage index NDVI is calculated as:
Figure FDA0002243877180000021
in the formula, NIR and R are reflectance values of the vegetation at a near infrared band and an infrared band, respectively.
5. The vegetation vertical zone quantitative scoring computation method of claim 3, wherein:
in the step (2), the sample region screening includes the following steps:
2.2.1, carrying out slope screening on the NDVI, the DEM and the fusion data of the slope and the slope direction, and screening an area with the slope not less than 5 degrees as a vertical zone area of vegetation to be extracted;
2.2.2 carrying out slope screening on the mountain areas of the vertical zone of the vegetation to be extracted, and ensuring that the areas of the vertical zone of the vegetation to be extracted are in the same slope.
6. The vegetation vertical zone quantitative scoring computation method of claim 1, wherein:
the step (3) comprises the following steps:
3.1 constructing a DEM-NDVI scatter distribution chart by utilizing a digital elevation model DEM and a vegetation coverage index NDVI extracted from a remote sensing image based on the fusion data of the vegetation vertical zone area to be extracted;
3.2, performing density segmentation on the DEM-NDVI scatter distribution diagram;
3.3 based on the DEM-NDVI scatter distribution diagram after density segmentation, obtaining the NDVI average value and the sliding average value curve of the specific altitude by using a sliding average method;
3.4 determining a proper fitting curve function;
and 3.5, analyzing the change trend of the fitted curve, finding out the demarcation points of different sections of the curve, and extracting the initial vegetation vertical zone demarcation value.
7. The vegetation vertical zone quantitative scoring computation method of claim 6, wherein:
in the step 3.2, the density segmentation of the DEM-NDVI scatter distribution diagram comprises the following contents:
3.2.1 in a DEM-NDVI scatter diagram, counting the number of scatter points in an interval range by taking 0.01 as an interval for a certain specific altitude;
3.2.2 calculate the NDVI density value in the range of every 0.01 interval according to the formula
Figure FDA0002243877180000022
In the formula, Mi is an NDVI density value in the ith NDVI interval range, Ni is the number of scattering points in the ith NDVI interval range, and Nmax is the maximum number of scattering points in all the NDVI interval ranges;
3.2.3, performing the operation on all the altitudes in the DEM-NDVI scattergram to obtain a DEM-NDVI density map;
3.2.4 selecting the area with the NDVI density value more than 0.95 as the core area of the DEM-NDVI scatter distribution diagram.
8. The vegetation vertical zone quantitative scoring computation method of claim 6, wherein:
in step 3.3, the moving average method includes the following steps:
3.3.1 defining a window of size 100m to find the average value of NDVI in the window;
3.3.2 sliding the window along the direction of the elevation according to the sliding distance of 5m so as to obtain the NDVI average value of the equal-spacing elevation;
3.3.3 connecting all the obtained NDVI mean points to obtain a sliding mean curve.
9. The vegetation vertical zone quantitative scoring computation method of claim 6, wherein:
in the step 3.4, fitting the sliding average curve by using an N-th-order polynomial regression equation; the fitting degree of the fitting curve and the DEM-NDVI sliding average value curve is measured by a correlation coefficient r2, N can be determined through multiple comparison experiments, and the minimum N value which enables r2 to reach the maximum value is selected as the proper fitting curve N value.
10. The vegetation vertical zone quantitative scoring computation method of claim 6, wherein:
in step 3.5, the dividing points of the different sections of the curve are determined by solving the inflection point and the second derivative extreme point of the fitted curve function.
11. The vegetation vertical zone quantitative scoring computation method of claim 1, wherein:
in the step 4.1, the initial threshold value of the boundary NDVI of the vertical vegetation zone is the NDVI value corresponding to the boundary point of different sections of the fitting curve of the vertical vegetation zone, and the experimental threshold range is the range of the corresponding NDVI initial threshold value fluctuating by 0.1.
12. The method of claim 11, wherein the method comprises:
in the step 4.2, the maximum inter-class variance method is applied to each experimental threshold range to respectively obtain the corresponding optimum threshold of the boundary NDVI of the vegetation vertical zone.
13. The method of claim 12, wherein the method comprises:
in the step 4.3, the assigning of the NDVI image after the optimal threshold segmentation includes the following steps:
4.3.1 the assignment range of the NDVI image after the optimal threshold segmentation is [0,1 ];
4.3.2 assigning the different categories of the NDVI images to be 1, …, 0 at equal intervals according to the altitude from low to high; namely, the vegetation dense area is assigned with the value of 1, the vegetation non-dense area is assigned with the value of 0, and the intermediate transition vegetation area is assigned with a certain value within the range of [0,1] at equal intervals.
14. The method of claim 13, wherein the method comprises:
the neighborhood window of a specific size in said step 4.4 is a circular window with a radius of 200 m.
15. The method of claim 14, wherein the method comprises:
the extraction of the boundary of the vegetation vertical zone in the step 4.5 comprises the following contents:
4.5.1 generating isolines with intervals of 0.1 based on the vegetation type probability map;
4.5.2 acquiring an initial boundary of the vegetation vertical zone in the research area based on the initial threshold of the NDVI;
4.5.3 overlapping the contour line, the initial boundary of the vegetation vertical zone and the original remote sensing image of the vegetation vertical zone area to be extracted, and selecting the contour line which is optimally matched with the initial boundary of the vegetation vertical zone and the type of vegetation visually interpreted by the remote sensing image as the optimized boundary of the vegetation vertical zone;
4.5.4, the optimized boundary line of the vegetation vertical zone is deleted by combining factors such as the theoretical distribution altitude and the slope direction of the vegetation, and the final boundary line of the vegetation vertical zone is obtained.
16. The utility model provides a vegetation vertical zone ration ruling system for the distribution characteristic is taken to ration ruling vegetation vertical zone, vegetation vertical zone ration ruling system includes data reading module, district selection module, scatter plot module, neighborhood analysis module and output module, its characterized in that:
the data reading module is used for reading a remote sensing image and a digital elevation model DEM of a vertical zone of vegetation to be extracted, carrying out image processing on the read remote sensing image, extracting vegetation coverage index, gradient and slope data of the preprocessed remote sensing image, and providing fusion data of the vegetation coverage index NDVI, the digital elevation model DEM and the gradient and slope to the area selection module;
the region selection module performs slope screening and slope zoning on the remote sensing image processed by the data reading module so as to obtain a vegetation vertical zone region to be extracted as a research region, and provides screened vegetation vertical zone mountain region fusion data to be extracted to the scatter diagram module;
the scatter diagram module constructs a DEM-NDVI scatter point distribution diagram according to mountain area fusion data of a vegetation vertical zone to be extracted after being screened by the selection module, density segmentation is carried out on the DEM-NDVI scatter point distribution diagram, then an NDVI average value and a sliding average value curve of a specific altitude are obtained, a fitting curve function is obtained by fitting the sliding average value curve through an N-th-order polynomial regression equation, a critical value of the vegetation vertical zone NDVI and the DEM is obtained according to the fitting curve function, and an initial threshold value of the vegetation vertical zone dividing line NDVI is provided for the neighborhood analysis module;
the neighborhood analysis module is used for acquiring a vegetation vertical zone boundary and providing final vegetation vertical zone boundary vector data and NDVI and DEM critical average values to the output module;
the output module is used for outputting the vector data of the boundary of the vegetation vertical zone and the NDVI and DEM critical average values;
the neighborhood analysis module comprises a threshold range acquisition unit, an optimal threshold extraction unit, a reclassification unit, a probability map generation unit and a boundary extraction unit;
the threshold range obtaining unit is used for obtaining an initial threshold value of a vegetation vertical zone dividing line NDVI and a corresponding experiment threshold range;
the optimal threshold extraction unit is used for respectively carrying out iterative optimization within each experimental threshold range to obtain an optimal threshold of a vegetation vertical zone boundary NDVI;
the reclassification unit is used for performing optimal threshold segmentation on the NDVI data in the sample area range, and performing reclassification and assignment on the NDVI data;
the probability map generating unit is used for selecting a neighborhood window with a specific size and acquiring a vegetation type probability map by a neighborhood analysis method;
and the boundary extraction unit is used for defining probability distribution thresholds of different vegetation types and extracting a final boundary of the vegetation vertical zone.
17. The vertical strip of vegetation quantitative scoring system of claim 16, wherein:
in the data reading module, the remote sensing image, the DEM, the gradient and the slope need to be in the same research area range.
18. The vertical strip of vegetation quantitative scoring system of claim 16, wherein:
the data reading module comprises an image preprocessing unit, an NDVI extraction unit and an image fusion unit:
the image preprocessing unit is used for carrying out radiometric calibration, atmospheric correction and terrain correction preprocessing operations on the remote sensing image;
the NDVI extraction unit is used for extracting the vegetation coverage index NDVI of the preprocessed remote sensing image;
the image fusion unit is used for carrying out image fusion on the NDVI, the DEM, the gradient and the slope data.
19. The vertical strip of vegetation quantitative scoring system of claim 18, wherein:
the NDVI is a normalized difference vegetation index, namely a vegetation coverage index, and the calculation formula is as follows:
Figure FDA0002243877180000051
in the formula, NIR and R are reflectance values of the vegetation at a near infrared band and an infrared band, respectively.
20. The vertical strip of vegetation quantitative scoring system of claim 16, wherein:
the selection area module comprises a slope screening unit and a slope screening unit:
the slope screening unit is used for carrying out slope screening on NDVI, DEM, slope and slope fusion data, and screening an area with the slope not less than 5 degrees as a vertical mountainous area of vegetation to be extracted;
the slope screening unit is used for conducting slope screening on the mountain areas of the vertical zones of the vegetation to be extracted, and ensuring that the mountain areas of the vertical zones of the vegetation to be extracted are in the same slope.
21. The vertical strip of vegetation quantitative scoring system of claim 16, wherein:
the scatter diagram module comprises a scatter distribution diagram generating unit, a density dividing unit, a moving average unit, a fitting curve unit and a demarcation point extracting unit:
the scattered point distribution diagram generating unit is used for reading the mountain area fusion data of the vegetation vertical zone to be extracted after being screened by the area selection module and constructing a DEM-NDVI scattered point distribution diagram by utilizing the DEM and the NDVI;
the density segmentation unit is used for performing density segmentation on the DEM-NDVI scatter distribution diagram;
the moving average unit is used for moving average the DEM-NDVI scatter distribution diagram after the density division, and acquiring an NDVI average value and a moving average value curve of a specific altitude;
the fitting curve unit is used for determining a proper fitting curve function;
and the demarcation point extraction unit is used for analyzing the variation trend of the fitted curve and finding out the demarcation points of different sections of the curve.
22. The vertical strip of vegetation quantitative scoring system of claim 21, wherein:
the density segmentation unit comprises a scattered point statistics subunit, a density calculation subunit, a density map subunit and a core area generation subunit;
the scatter counting subunit is used for reading scatter points at a certain specific altitude in the DEM-NDVI scatter diagram, and counting the number of scatter points in a certain altitude interval range by taking 0.01 as an interval;
the density calculation subunit is used for calculating the NDVI density value within the 0.01 interval range, and the calculation formula is
Figure FDA0002243877180000061
In the formula, Mi is an NDVI density value in the ith NDVI interval range, Ni is the number of scattering points in the ith NDVI interval range, and Nmax is the maximum number of scattering points in all the NDVI interval ranges;
the density map subunit is used for performing the above operation on all the altitudes in the DEM-NDVI scattergram to obtain a DEM-NDVI density map;
and the core area generating subunit is used for selecting an area with an NDVI density value larger than 0.95 as a DEM-NDVI scatter distribution diagram core area.
23. The vertical strip of vegetation quantitative scoring system of claim 22, wherein:
the moving average unit comprises a window average value operator unit, an equidistant average value calculation subunit and a curve generation subunit;
the window average value operator unit is used for defining a window with the size of 100m to obtain the NDVI average value in the window range;
the equidistant average value operator unit is used for sliding the window along the altitude direction according to the sliding distance of 5m to obtain the NDVI average value of the equidistant altitude;
and the curve generation subunit is used for connecting all the obtained NDVI mean value points to obtain a sliding mean value curve.
24. The vertical strip of vegetation quantitative scoring system of claim 21, wherein:
in a fitting curve unit, an N-th-order polynomial regression equation is adopted to perform curve fitting on the sliding average value curve, and the fitting degree of the fitting curve and the DEM-NDVI sliding average value curve is related to a coefficient r2Performing a measurement, wherein N is determined by multiple comparison experiments, and r is selected2The minimum N value at which the maximum is reached is taken as the N value of the fitted curve.
25. The vertical strip of vegetation quantitative scoring system of claim 22 or 24, wherein:
in the boundary point extraction unit, the boundary points of different sections of the fitting curve are determined by solving inflection points and second-order derivative extreme points of the fitting curve function.
26. The vertical strip of vegetation quantitative scoring system of claim 25, wherein:
in the threshold range obtaining unit, the NDVI values corresponding to the demarcation points of different sections of the fitting curve determined by the demarcation point extracting unit in the scatter diagram module are used as the NDVI initial threshold of the vertical boundary of the vegetation zone, and the range of the corresponding NDVI initial threshold which fluctuates up and down by 0.1 is used as the experimental threshold range.
27. The vertical strip of vegetation quantitative scoring system of claim 16 or 26, wherein:
the optimal threshold extraction unit comprises an inter-class variance subunit, and the inter-class variance subunit is used for respectively solving the corresponding vegetation vertical zone dividing line NDVI optimal threshold in each experiment threshold range of the NDVI image in the sample area range by applying the idea of the maximum inter-class variance method.
28. The vertical strip of vegetation quantitative scoring system of claim 16, wherein:
and the reclassification unit is used for determining an assignment range and assigning the reclassified images.
29. The vertical strip of vegetation quantitative scoring system of claim 16, wherein:
in the probability map generating unit, the neighborhood window with a specific size is a circular window with a radius of 200 m.
30. The vertical strip of vegetation quantitative scoring system of claim 16, wherein:
the boundary extraction unit comprises an isoline generation subunit, an initial boundary acquisition subunit, a boundary initial optimization subunit and a boundary fine extraction subunit;
the contour line generation subunit is used for generating contour lines with the interval of 0.1 based on the vegetation type probability map acquired by the probability map generation unit;
the initial boundary acquisition subunit is used for acquiring a vegetation vertical zone initial boundary of the research area based on an NDVI initial threshold of the vegetation vertical zone boundary;
the boundary initial optimization subunit is used for superposing the isoline, the vegetation vertical zone initial boundary and the original remote sensing image of the research area, and selecting the isoline which is optimally matched with the vegetation vertical zone initial boundary and the remote sensing image visual interpretation vegetation type as the optimized vegetation vertical zone boundary;
and the boundary fine-extraction subunit is used for deleting the optimized boundary of the vegetation vertical zone by combining factors such as theoretical distribution altitude and slope direction of the vegetation to obtain the final boundary of the vegetation vertical zone.
CN201710467468.5A 2017-06-20 2017-06-20 Quantitative marking calculation method and system for vegetation vertical zone Active CN107330898B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710467468.5A CN107330898B (en) 2017-06-20 2017-06-20 Quantitative marking calculation method and system for vegetation vertical zone

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710467468.5A CN107330898B (en) 2017-06-20 2017-06-20 Quantitative marking calculation method and system for vegetation vertical zone

Publications (2)

Publication Number Publication Date
CN107330898A CN107330898A (en) 2017-11-07
CN107330898B true CN107330898B (en) 2020-04-28

Family

ID=60194425

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710467468.5A Active CN107330898B (en) 2017-06-20 2017-06-20 Quantitative marking calculation method and system for vegetation vertical zone

Country Status (1)

Country Link
CN (1) CN107330898B (en)

Families Citing this family (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111861836B (en) * 2020-07-20 2022-10-18 云南财经大学 Three-dimensional mountain land planning method and device, storage medium and computer equipment
CN112347978B (en) * 2020-11-25 2022-06-03 中国科学院东北地理与农业生态研究所 Ecological system attribute component composition structure description method based on remote sensing vegetation index
CN112907587B (en) * 2021-04-01 2022-03-01 西南石油大学 High mountain forest line extraction method based on Otsu and edge detection algorithm of GEE
CN113433553B (en) * 2021-06-23 2022-08-02 哈尔滨工程大学 Precise navigation method for multi-source acoustic information fusion of underwater robot
CN114061591B (en) 2021-11-18 2022-07-12 东南大学 Contour line matching method based on sliding window data backtracking
CN115019184B (en) * 2022-07-28 2023-02-07 北京卫星信息工程研究所 Remote-sensing-image-based stony desertification degree automatic grading method and device
CN116468321A (en) * 2023-04-23 2023-07-21 河北省科学院地理科学研究所 Topography comprehensive index method for representing vegetation growth distribution rule

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105957115A (en) * 2016-05-18 2016-09-21 河北工程大学 Crop planting structure remote sensing extraction method under generalized DEM idea
CN106384116A (en) * 2016-08-29 2017-02-08 北京农业信息技术研究中心 Terahertz imaging based plant vein recognition method and device

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105957115A (en) * 2016-05-18 2016-09-21 河北工程大学 Crop planting structure remote sensing extraction method under generalized DEM idea
CN106384116A (en) * 2016-08-29 2017-02-08 北京农业信息技术研究中心 Terahertz imaging based plant vein recognition method and device

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
"基于DEM-NDVI的高山植被带定量刻划";常纯,王心源等;《地理研究》;20151130;第34卷(第11期);第2113-2123页 *
"长白山植被垂直带地形控制机制研究";郭聃;《中国博士学位论文全文数据库 基础科学辑》;20150215;第1、15、20-21、33、40、50页 *

Also Published As

Publication number Publication date
CN107330898A (en) 2017-11-07

Similar Documents

Publication Publication Date Title
CN107330898B (en) Quantitative marking calculation method and system for vegetation vertical zone
CN101840581B (en) Method for extracting profile of building from satellite remote sensing image
US9064151B2 (en) Device and method for detecting plantation rows
CN105243367A (en) Method and device for monitoring scope of water body based on satellite remote sensing data
Pirotti et al. A comparison of tree segmentation methods using very high density airborne laser scanner data
CN111091079A (en) TLS-based method for measuring dominant single plant structural parameters of vegetation in alpine and fragile regions
CN102855485A (en) Automatic wheat earing detection method
Hofmeister et al. Intercomparison of Sentinel-2 and modelled snow cover maps in a high-elevation Alpine catchment
CN113935366A (en) Automatic classification method for point cloud single wood segmentation
Lei et al. A novel algorithm of individual tree crowns segmentation considering three-dimensional canopy attributes using UAV oblique photos
CN115205528A (en) Feature selection method for geographic object-oriented image analysis
CN107103295A (en) Optical remote sensing image cloud detection method of optic
CN113705441A (en) High-spatial-temporal-resolution surface water body extraction method cooperating with multispectral and SAR images
Pan et al. Remote sensing inversion of soil organic matter by using the subregion method at the field scale
CN116541658A (en) Urban near-ground kilometer height wind profile measurement analysis method and device
CN111476434B (en) GIS-based soil heavy metal fractal dimension spatial variation analysis method
Volkova et al. Cloud Cover and Precipitation Monitoring Based on Data from Polar Orbiting and Geostationary Satellites
Kanjir et al. Application of object based approach to heterogeneous land cover/use
Overton et al. Flood modelling and vegetation mapping in large river systems
Bouchayer Synthesis of distributed snowpack simulations relevant for avalanche hazard forecasting
CN113900097B (en) Glacier quantity detection method based on satellite remote sensing data
Liu et al. Assessment of tree attributes extraction algorithms
Cui et al. A novel quantitative analysis for diurnal dynamics of Ulva prolifera patch in the Yellow Sea from Geostationary Ocean Color Imager observation
Kearsley Characterizing spatial variability of tropical rainforest structure using hemispherical photography, in the reserves of Yangambi and Yoko (Democratic Republic of Congo)
Liu Aboveground biomass estimation of individual trees with airborne Lidar data

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant