Intelligent forest pest and disease damage monitoring method and system based on unmanned aerial vehicle remote sensing

Abstract

The invention relates to the technical field of remote sensing monitoring, in particular to an intelligent forest disease and pest monitoring method and system based on unmanned aerial vehicle remote sensing, and the method comprises the following steps: collecting multispectral data by an unmanned aerial vehicle, extracting reflectivity and smoothing the reflectivity, carrying out differential recognition on abnormal pixels, extracting curves and screening significant changes, segmenting scab boundaries, and classifying health states. And generating a pest and disease map layer prediction trend. According to the method, the reflectivity time sequence is constructed, differential processing is carried out, the vegetation change trend is dynamically captured, abnormal areas are identified by combining slope offset and persistence analysis, significant pixels are screened according to main peak wavelength offset, boundary information is extracted, the scab positioning precision is improved, and the recognition resolution of the lesion state is enhanced through reflectivity combined analysis; accurate description of disease spot dynamic changes is realized, static image dependence limitation is broken through, monitoring time continuity and space response capability are enhanced, and disease and insect pest change capture efficiency and state classification accuracy are effectively improved.

Description

Technical Field

[0001] The present invention relates to the field of remote sensing monitoring technology, and in particular to an intelligent forest pest monitoring method and system based on unmanned aerial vehicle (UAV) remote sensing. Background Art

[0002] The field of remote sensing monitoring technology includes the technical means of using sensors carried by aerospace platforms to obtain surface target and environmental information, and realizing the identification, classification and monitoring of surface features by perceiving and analyzing the electromagnetic radiation characteristics and their changes. The core content includes remote sensing data acquisition, remote sensing data processing and information extraction, mainly involving the application of data sources such as multispectral imaging, hyperspectral imaging, lidar, synthetic aperture radar, etc., and realizing dynamic monitoring and evaluation of natural resources, ecological environment, agriculture, forestry, urban management and other fields through image processing, feature extraction, target recognition and other methods. The technical system is based on data acquisition platform, sensor type, data processing method and application scenario, forming a complete technical framework and application system.

[0003] Among them, an intelligent monitoring method and system for forest pests and diseases based on UAV remote sensing refers to the use of UAVs equipped with multispectral imaging equipment to collect low-altitude remote sensing data of forest areas, analyze the reflection characteristics of vegetation canopies through the collected multispectral images, and spatially locate the pest and disease areas by combining image segmentation methods. At the same time, a trained convolutional neural network model is used to classify and identify the health status of vegetation, and the monitoring results are further spatially visualized and analyzed through the geographic information system. The patent subject mainly covers UAV remote sensing image acquisition, image preprocessing, multispectral feature extraction, convolutional neural network classification, spatial data display and other links, and the monitoring of forest pests and diseases is completed through a clear data processing process.

[0004] Existing technologies rely on static multispectral images to identify pest and disease areas. The processing results are easily affected by single-phase lighting conditions and changes in vegetation status. In scenes with dense forest canopies or frequent phenological changes, there are significant errors in the recognition results. Lack of continuity makes it impossible to capture the spread trajectory of lesions in a timely manner. Image segmentation results are easily interfered by areas with similar colors, resulting in unclear boundary extraction. Vegetation classification relies on single-time recognition of convolutional neural networks and lacks support from information on temporal changes. It is easy to make inaccurate judgments on health status. Under seasonal changes or frequent interference sources, it is impossible to accurately portray the formation and evolution process of lesions, which limits the timeliness and accuracy of pest and disease monitoring. Summary of the Invention

[0005] In order to solve the technical problems existing in the prior art, the present invention provides an intelligent forest pest monitoring method and system based on drone remote sensing. The technical solution is as follows:

[0006] To achieve the above objectives, the present invention adopts the following technical solution: an intelligent monitoring method for forest pests and diseases based on UAV remote sensing, comprising the following steps:

[0007] S1: Using multispectral imaging equipment carried by drones, we collected and organized forest remote sensing image data, identified the main band reflectance of the vegetation canopy area, and constructed time series data. We then processed the data using the moving average method to generate a smooth reflectance trend sequence.

[0008] S2: calling the smoothed reflection trend sequence, performing differential processing on the data of adjacent time nodes, constructing a time slope change sequence, fitting the reflectivity change trend, constructing a normal slope change interval, identifying pixels whose reflectivity changes continuously deviate from the normal interval range, analyzing the persistence of the deviation, and identifying the location of candidate abnormal areas;

[0009] S3: Based on the candidate abnormal area position image, extract the reflection curves of multiple bands, locate the main peak wavelength and calculate the abnormal offset, analyze the pixel offset degree, screen the pixels with significant change characteristics, and segment the spectral boundary image based on the spatial distribution to generate the initial boundary data of the lesion;

[0010] S4: Based on the initial boundary data of the diseased spot, the multispectral characteristic parameters of the corresponding pixel are extracted, and the reflectance combination relationship between the bands is used to classify the vegetation health status and perform classification annotation to generate classification annotation result data.

[0011] As a further solution of the present invention, the smoothed reflectance trend sequence includes the main band reflectance time series, the smoothed value of the time series data, and the moving average characteristic value; the time slope change sequence includes the reflectance difference value of adjacent nodes, the normal change interval range, and the continuous deviation index; the initial boundary data of the lesion includes the distribution of significantly changed pixels, the main peak wavelength, and the spectral boundary segmentation map; the classification annotation result data includes a multispectral feature parameter set, a vegetation health level label, and a classification area map.

[0012] As a further solution of the present invention, the step S1 is specifically as follows:

[0013] S101: Use a drone equipped with multispectral imaging equipment to collect remote sensing image data of the forest area. Number the images collected at each time point in chronological order, identify the green vegetation canopy area in the image, extract the reflectance value corresponding to the main band, and generate a vegetation canopy reflectance dataset.

[0014] S102: Based on the vegetation canopy reflectance dataset, integrating the main band reflectance values ​​corresponding to the time nodes in the order of image shooting time, establishing a time series reflectance data structure of the vegetation canopy area, and sorting the reflectance values ​​in time series to obtain a canopy reflectance time series value sequence;

[0015] S103: calling the canopy reflectance time series value sequence, setting the observation window width according to the continuously observed reflectance values ​​at each time point, smoothing each set of reflectance data sequence using a moving average method, and generating a smoothed reflectance trend sequence;

[0016] The specific formula for setting the observation window width is:

[0017]

[0018] Calculate the observation window width value;

[0019] Among them, W i Indicates the observation window width value at the i-th time point, Z c Indicates the absolute value of the canopy reflectance difference between the i+cth and icth time points when the band is λ, Q i Indicates the absolute value of the deviation of the canopy reflectance at time point i relative to the sequence mean when the band is λ, L c represents the average time interval of the cth symmetric time point relative to the central time point, U i represents the horizontal velocity of the observation platform at time point i, Θ represents the standard deviation of the canopy reflectance series with band λ, m represents the number of symmetrical time point pairs used for calculation, and c represents the index number of the time point pairs symmetrically distributed forward and backward with respect to the central time point i;

[0020] The moving average method sets a fixed time window for the vegetation canopy reflectance time series value sequence and calculates the reflectance mean value in each interval by sliding.

[0021] As a further solution of the present invention, the step S2 is specifically as follows:

[0022] S201: Calling the data values ​​of adjacent time nodes in the smooth reflection trend sequence, obtaining the reflectivity value of the pixel position at the time node, and calculating the difference of the reflectivity at the adjacent time nodes, constructing the time slope change value sequence of the pixel point, and generating the pixel-level time slope change sequence;

[0023] S202: calling the pixel-level time slope change sequence and the image main band reflectance time series data with geographic coordinate information, performing a fitting operation on the main band reflectance sequence of each pixel, extracting the trend slope value set of the pixel, and constructing the normal slope change range interval of the pixel based on the peak and valley values ​​in the trend slope value set to generate the normal range of regional slope change;

[0024] S203: Analyzing the continuity of the pixel time slope change value based on the pixel-level time slope change sequence and the normal interval of regional slope change, marking the pixel positions that exceed the normal interval of regional slope change, and counting the continuous exceeding time. Based on the exceeding time length, the pixel point positions that meet the continuity offset are screened to generate candidate abnormal area positions;

[0025] The continuity deviation refers to the situation where the slope change value of a certain pixel exceeds the normal range at multiple time nodes in the time series, and the exceeding state remains continuous and non-repeating.

[0026] As a further solution of the present invention, the step S3 is specifically as follows:

[0027] S301: Based on the image corresponding to the candidate abnormal area location, extract the pixel reflectance value sequence of the near-infrared, high green, and short-wave red edge bands in the image, construct the spectral reflectance curve of each pixel in the order of wavelength, extract and record the wavelength value corresponding to the main peak of the band spectral curve, and generate the main band peak wavelength sequence;

[0028] S302: calling the main band peak wavelength sequence, calculating the offset value between the main peak wavelength of each pixel in the three bands and the central wavelength of the reference band, comparing the offset value with a set offset threshold, screening pixel positions with offset values ​​greater than the threshold, and generating an offset abnormal pixel position group;

[0029] The reference band center wavelength is the position with the strongest spectral response extracted from the sensor response curves of the near infrared, high green and short-wave red edge bands according to the band setting of the remote sensing imaging device;

[0030] The offset threshold is a criterion set by statistically analyzing the total amount of the main peak wavelength offset of the pixels in the candidate area and according to the central tendency and discrete range of the offset degree;

[0031] S303: Based on the offset abnormal pixel position group and the two-dimensional spatial position range of the candidate abnormal area, a pixel segmentation map structure is constructed according to the spatial adjacent relationship, and the boundary contour of the connected pixel area is delineated according to the spectral boundary difference to establish the initial boundary data of the lesion.

[0032] As a further solution of the present invention, the step S4 is specifically as follows:

[0033] S401: Based on the initial boundary data of the lesion, extract the band values ​​of the multispectral image pixels corresponding to the boundary, call the visible light, near infrared and short wave infrared band values ​​of the corresponding pixel positions in the reflectance data, associate and identify the pixel data according to the spatial coordinates, and generate a spectral value set of the boundary pixels;

[0034] S402: Based on the boundary pixel spectral value set, using the reflectance combination relationship between bands, extracting multiple vegetation status assessment indicators, including the normalized vegetation index and the normalized water index, and identifying the offsets of the multiple indicators based on the typical threshold range of each indicator, classifying the vegetation health status, and generating vegetation health classification reflectance intervals;

[0035] The typical threshold range refers to the standard value range of the spectral evaluation index used to analyze the vegetation status of boundary pixels at different health levels;

[0036] The normalized vegetation index is used to quantify the difference in the ratio of red light and near-infrared reflectance of boundary pixels and to judge the photosynthetic activity and health level of vegetation;

[0037] The normalized water index is used to reflect the imbalance between the green light and near-infrared reflectance of the boundary pixels;

[0038] S403: Call the vegetation health grade reflectance interval, perform state matching on the spatial position in the initial boundary data of the diseased spot, analyze the corresponding pixel band value and classification interval value, perform spatial annotation of the vegetation state according to the classification label, and output code the vegetation pixel to generate classification annotation result data.

[0039] As a further embodiment of the present invention, the method further comprises:

[0040] S5: Calling the classification and labeling result data, superimposing the spatial location of pests and diseases onto the forest spatial coordinate system, generating a pest and disease spatial distribution layer, updating the forest area pest and disease monitoring database based on the layer data, and predicting the pest and disease development trend in combination with time information to generate a pest and disease spatial distribution trend result;

[0041] The results of the spatial distribution trend of pests and diseases include a spatial distribution layer of pests and diseases, a monitoring database update record, and a development trend forecast map.

[0042] As a further solution of the present invention, the step S5 is specifically as follows:

[0043] S501: Calling the classification and labeling result data, identifying the corresponding position of the pest and disease labeling area in the forest remote sensing image, projecting and matching the spatial coordinate information to the forest spatial coordinate system, overlaying the geometric boundaries of the classified pixels and establishing a pest and disease spatial distribution information layer, and generating pest and disease layer overlay data;

[0044] S502: Based on the overlay data of the pest and disease layer, the time label and spatial number corresponding to each pest and disease area in the layer are matched, the data unit with the same number is located in the forest area database, the pest and disease occurrence status in the database field is updated, and the spatial description field content is supplemented to obtain the pest and disease monitoring update record;

[0045] S503: Call the pest monitoring update record, analyze the pest labeling status of the spatial unit, identify the change direction and expansion rate of the pest area, predict the change trend of the pest area and mark the trend level, and generate the pest spatial distribution trend result.

[0046] As a further solution of the present invention, the specific formula for identifying the change direction and expansion rate of the pest and disease area is:

[0047]

[0048] Calculate the regional expansion rate value V;

[0049] in, Represents the expansion rate value of the area with spatial unit number j in direction type d. and represent the pest and disease labeling values ​​of the kth pixel in the spatial unit j at time t and t-1, respectively, j,k represents the pixel center distance of the kth pixel in the spatial unit j along the direction d, T j,k Represents the topological disturbance factor of the boundary of the region where the k-th pixel is located in this direction, Respectively represent the pest and disease area of ​​spatial unit j at time t and t-1, n j Represents the number of changed pixels involved in the calculation, j is the number index of the spatial unit, t is the time index, and k is the pixel index.

[0050] In another aspect, a forest pest and disease intelligent monitoring system based on drone remote sensing is provided. The system is applied to a forest pest and disease intelligent monitoring method based on drone remote sensing. The system comprises:

[0051] The multispectral time series processing module obtains vegetation canopy reflectance data through the multispectral imaging equipment carried by the UAV, uses the moving average method to perform time series smoothing on the reflectance data, generates a smoothed reflectance trend sequence, and passes it to the time series anomaly detection module;

[0052] The time series anomaly detection module calls the smoothed reflection trend sequence, performs differential processing on the data of adjacent time nodes, constructs a time slope change sequence, uses a linear regression model to fit the reflectivity change trend, calculates the slope residual standard deviation, constructs a normal slope change interval, identifies the pixel coordinates that exceed the interval for multiple consecutive cycles, generates candidate abnormal area locations, and transmits them to the lesion boundary generation module;

[0053] The lesion boundary generation module extracts the red edge band and near-infrared band reflectance curves based on the candidate abnormal area location, calculates the main peak wavelength offset, analyzes the pixel offset degree, selects pixels with significant change characteristics, and segments the spectral boundary image based on the spatial distribution to generate initial lesion boundary data, which is then passed to the health status classification module;

[0054] The health status classification module obtains the initial boundary data of the diseased spot, extracts the multispectral characteristic parameters of the corresponding pixel, classifies the vegetation health status and performs classification annotation based on the reflectance combination relationship between the bands, generates classification annotation result data, and transmits it to the spatial trend prediction module;

[0055] The spatial trend prediction module calls the classification and annotation result data, converts the coordinates of pests and diseases into the geographic coordinate system, generates a spatial distribution thermal layer, calculates the regional pest and disease density change rate, updates the forest area pest and disease monitoring database, combines time information to predict the development trend of pests and diseases, and generates the spatial distribution trend results of pests and diseases.

[0056] The beneficial effects brought about by the technical solution provided by the embodiment of the present invention include at least:

[0057] By constructing a time series and introducing a differential processing method to obtain changes in reflectivity slope, it is possible to dynamically track subtle changes in vegetation status. Combined with offset degree analysis, the early identification ability of abnormal areas is enhanced, the preliminary boundary delineation of abnormal areas is guided to improve the precise positioning of diseased areas, and the combination of reflectance features is used to enhance the classification resolution of vegetation health status. Spatially, the focus is further focused on areas of significant variation to achieve a refined presentation of distribution characteristics, effectively avoiding the problem of missing time information caused by static recognition of single-frame images, improving the continuous monitoring capability of the spread trend of pests and diseases, and improving the recognition ambiguity caused by weak image contrast in different growth periods. By constructing an analysis mechanism that combines time and space, the response sensitivity to abnormal vegetation characteristics is enhanced. BRIEF DESCRIPTION OF THE DRAWINGS

[0058] In order to more clearly illustrate the technical solutions in the embodiments of the present invention, the following briefly introduces the drawings required for use in the description of the embodiments. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without creative work.

[0059] Figure 1 It is a schematic diagram of the workflow of the present invention;

[0060] Figure 2 It is a system flow chart of the present invention. DETAILED DESCRIPTION

[0061] The technical solution of the present invention is described below in conjunction with the accompanying drawings.

[0062] In the embodiments of the present invention, words such as "exemplarily" and "for example" are used to indicate examples, illustrations, or explanations. Any embodiment or design described as an "exemplary" in the present invention should not be interpreted as being preferred or advantageous over other embodiments or designs. Rather, the use of the word "exemplary" is intended to present concepts in a concrete manner. Furthermore, in the embodiments of the present invention, "and / or" can mean both or either of the two.

[0063] In the embodiments of the present invention, the terms "image" and "picture" may be used interchangeably. It should be noted that, when the distinction between them is not emphasized, their intended meanings are the same. The terms "of," "corresponding," and "corresponding" may be used interchangeably. It should be noted that, when the distinction between them is not emphasized, their intended meanings are the same.

[0064] In the embodiments of the present invention, sometimes a subscript such as W1 may be written as a non-subscript such as W1. When the difference is not emphasized, the meanings to be expressed are the same.

[0065] In order to make the technical problems, technical solutions and advantages to be solved by the present invention clearer, a detailed description will be given below with reference to the accompanying drawings and specific embodiments.

[0066] See also Figure 1 The present invention provides a technical solution, a method for intelligent monitoring of forest pests and diseases based on UAV remote sensing, comprising the following steps:

[0067] S1: Using multispectral imaging equipment carried by drones, we collected and organized forest remote sensing image data, identified the main band reflectance of the vegetation canopy area, and constructed time series data. We then processed the data using the moving average method to generate a smooth reflectance trend sequence.

[0068] S2: Call the smoothed reflectance trend sequence, perform differential processing on the data of adjacent time nodes, construct a time slope change sequence, fit the reflectance change trend, construct a normal slope change interval, identify pixels whose reflectance changes continuously deviate from the normal interval range, analyze the persistence of the deviation, and identify the location of candidate anomaly areas;

[0069] S3: Based on the image of the candidate abnormal area, extract the reflectance curves of multiple bands, locate the main peak wavelength and calculate the abnormal offset, analyze the degree of pixel offset, screen pixels with significant change characteristics, and segment the spectral boundary image based on the spatial distribution to generate the initial boundary data of the lesion;

[0070] S4: Based on the initial boundary data of the lesion, the multispectral characteristic parameters of the corresponding pixels are extracted, and the reflectance combination relationship between the bands is used to classify the vegetation health status and perform classification annotation to generate classification annotation result data;

[0071] S5: Call the classification and labeling result data, superimpose the spatial location of pests and diseases on the forest spatial coordinate system, generate a pest and disease spatial distribution layer, update the forest area pest and disease monitoring database based on the layer data, combine the time information to predict the development trend of pests and diseases, and generate the pest and disease spatial distribution trend results.

[0072] The smoothed reflectance trend sequence includes the main band reflectance time series, the smoothed value of the time series data, and the moving average characteristic value. The time slope change sequence includes the reflectance difference value of adjacent nodes, the normal change range, and the continuous deviation index. The initial boundary data of the lesion includes the distribution of significantly changed pixels, the main peak wavelength, and the spectral boundary segmentation map. The classification and annotation result data includes the multispectral feature parameter set, the vegetation health level label, and the classification area map. The spatial distribution trend results of pests and diseases include the spatial distribution layer of pests and diseases, the monitoring database update record, and the development trend forecast map.

[0073] See also Figure 1 , the steps of S1 are as follows:

[0074] S101: Use a drone equipped with multispectral imaging equipment to collect remote sensing image data of the forest area. Number the images collected at each time point in chronological order, identify the green vegetation canopy area in the image, extract the reflectance value corresponding to the main band, and generate a vegetation canopy reflectance dataset.

[0075] To call the data values ​​of adjacent time nodes in the smoothed reflection trend sequence, we first need to clarify the sequence index of each time node and extract the main band reflectance value under the corresponding pixel number. For example, at the pixel number 1, its reflectance at time nodes t1 to t4 is 0.12, 0.15, 0.18, and 0.20. According to the time sequence, the reflectance at each adjacent node is subtracted in pairs to form a reflectance difference sequence, specifically: ΔR 12 =0.15-0.12=0.03, ΔR 23 =0.18-0.15=0.03, ΔR 34=0.20-0.18=0.02, each difference represents the change amplitude of the reflectivity per unit time within the time span, and then each difference is arranged according to the time node to form a time slope change value sequence of the pixel point, that is, {0.03, 0.03, 0.02}. During the execution process, it is necessary to traverse all pixel numbers to complete the difference sequence construction operation, and further use the pixel number as the index to construct the pixel-level time slope change sequence matrix of the entire image. In this process, it is necessary to judge whether the time difference between adjacent time nodes is uniform. If not, time standardization is required, that is, the reflectivity difference is divided by the corresponding time interval difference τ. For example, if t2-t1=10 days, t3-t2=15 days, and t4-t3=20 days, then ΔR 23 =(0.18-0.15) / 15≈0.002, rather than the aforementioned 0.03. This comparative analysis can avoid the risk of misjudgment due to inconsistent sampling intervals. In addition, if the reflectivity difference is a negative number, such as the reflectivity change of pixel 2 is: 0.35, 0.30, 0.25, 0.22, the corresponding difference is -0.05, -0.05, -0.03, then the sequence can be represented as a continuous downward trend and needs to be further classified and stored as a downward slope sequence. In actual calculations, if the reflectivity difference fluctuates in the range of -0.01 to 0.01, it can be determined as a stable change. The judgment standard is set The stable interval is [-0.01, 0.01], the rising interval is (0.01, 1], and the falling interval is [-1, -0.01]. The boundaries of each interval are set by adding or subtracting the standard deviation of the statistical mean of the maximum change amplitude of the main band in the image. For example, if the reflectivity in the historical sample is 0.75 at the maximum and 0.05 at the minimum, then its theoretical maximum fluctuation range is 0.7. Considering the standard deviation σ = 0.05, the stable interval is set within ± 2σ = 0.1. This method realizes the quantitative identification of the change trend of pixels in the entire region at different time nodes, thereby generating a pixel-level time slope change sequence of the entire image.

[0076] Table 1 Pixel time node reflectivity table

[0077]

[0078] As shown in Table 1, different pixels have different reflectivity change values ​​at different time nodes. In the subsequent steps, the time slope sequence is generated by performing pairwise difference operations on these values ​​for change analysis.

[0079] S102: Based on the vegetation canopy reflectance dataset, the main band reflectance values ​​corresponding to the time nodes are integrated in the order of image shooting time, a time series reflectance data structure of the vegetation canopy area is established, and the reflectance values ​​are sorted in time series to obtain a canopy reflectance time series value sequence;

[0080] To call the pixel-level time slope change sequence and the image main band reflectance time series data with geographic coordinate information, it is necessary to first perform a fitting operation on the main band reflectance of each pixel point. The fitting method is to use the time node and its reflectance as x and y inputs, build a time-reflectance mapping relationship, and adopt a piecewise linear regression method to fit a linear slope once in each three-node sliding window. For example, the reflectance sequence of pixel 1 is {0.12, 0.15, 0.18, 0.20}, then the linear slope calculation is performed with (t1, t2, t3) and (t2, t3, t4) as two windows respectively. The linear slope calculation formula is: slope The corresponding calculation process is If t3-t2=10 days, then k1=0.003. At the entire image level, each pixel is calculated through the window regression method to calculate its local slope value set. Then, these slope values ​​are collected and the maximum and minimum values ​​are extracted in ascending order as the slope peak and valley values ​​of the pixel. For example, the slope set of pixel 3 is {0.03, 0.02}, where the peak value is 0.03 and the valley value is 0.02. The normal variation range of the pixel is further calculated by setting the normal range as interval = [valley value-δ, peak value+δ], where δ is the adjustment coefficient, which is set to 2 times the historical slope fluctuation σ=0.005, then δ=0.01, and the normal range of the pixel is [0. 01,0.04]. In actual application, this interval can be used for subsequent change anomaly detection scenarios. By comparing whether the local slope in the current time period exceeds this range, pixel-level screening and judgment are performed. If the current slope of a pixel is 0.045 and exceeds [0.01,0.04], it is marked as an abnormal change area. This judgment method compares the current slope value of the participating item with the upper and lower bounds of the stable slope obtained by historical change fitting. If the current value is greater than the upper limit or less than the lower limit, the abnormal flag 1 is output, otherwise it is 0. The normal slope change range interval of the pixel constructed in this way can eventually generate a normal interval map of geographic slope fluctuation on a regional scale, realizing the spatial expression of the change interval of each pixel point.

[0081] S103: calling the canopy reflectance time series value sequence, setting the observation window width according to the continuously observed reflectance values ​​at each time point, smoothing each set of reflectance data sequence using the moving average method, and generating a smoothed reflectance trend sequence;

[0082] The specific formula for setting the observation window width is:

[0083]

[0084] Calculate the observation window width value;

[0085] Among them, W i Indicates the observation window width value at the i-th time point, Zc Indicates the absolute value of the canopy reflectance difference between the i+cth and icth time points when the band is λ, Q i Indicates the absolute value of the deviation of the canopy reflectance at time point i relative to the sequence mean when the band is λ, L c represents the average time interval of the cth symmetric time point relative to the central time point, U i represents the horizontal velocity of the observation platform at time point i, Θ represents the standard deviation of the canopy reflectance series with band λ, m represents the number of symmetrical time point pairs used for calculation, and c represents the index number of the time point pairs symmetrically distributed forward and backward with respect to the central time point i;

[0086] formula:

[0087]

[0088] Detailed explanation of the formula and the process of formula calculation and derivation:

[0089] The formula is used to calculate the observation window width value W at the i-th time point i , which is used to determine the data window span around time point i when smoothing canopy reflectance time series data;

[0090] Parameter meaning and setting value:

[0091] Z c : When the band is λ, the absolute value of the canopy reflectance difference between the i+cth time point and the icth time point. The reflectance at the i+cth time point is set to 0.48 and the reflectance at the icth time point is set to 0.42. Then Z c =|0.48-0.42|=0.06;

[0092] Q i : When the band is λ, the absolute value of the deviation of the canopy reflectance at the i-th time point relative to the sequence mean. If the reflectance at the i-th time point is 0.45 and the sequence mean reflectance is 0.40, then Q i =|0.45-0.40|=0.05;

[0093] L c : The average time interval between the cth symmetrical time point and the central time point. Let the i+cth time point be the 10th day and the icth time point be the 6th day. Then L c = (10-6) / 2 = 2 days;

[0094] U i : The horizontal velocity of the observation platform at the i-th time point, set U i =5 m / s;

[0095] Θ: standard deviation of the canopy reflectance sequence for the λ band, obtained by calculating the square root of the average of the squares of the deviations of the reflectance from the mean value at each time point in the entire sequence. If the sequence reflectance is set to [0.38, 0.42, 0.45, 0.40, 0.43], then Θ≈0.039;

[0096] m: the number of symmetric time point pairs used for calculation, set m=2;

[0097] Substitute the parameters into the formula for calculation:

[0098]

[0099] Result 2 shows that when smoothing canopy reflectance time series data at time point i, it is recommended to use a data window centered at that time point and containing two time points before and after it, that is, a data window with a total span of five time points, to ensure the effectiveness and accuracy of the smoothing process.

[0100] The moving average method sets a fixed time window for the vegetation canopy reflectance time series value sequence and calculates the mean reflectance value in each interval by sliding.

[0101] See also Figure 1 , the steps of S2 are as follows:

[0102] S201: Calling the data values ​​of adjacent time nodes in the smooth reflection trend sequence, obtaining the reflectivity value of the pixel position at the time node, and calculating the difference of the reflectivity at adjacent time nodes, constructing the time slope change value sequence of the pixel point, and generating the pixel-level time slope change sequence;

[0103] Call the data values ​​of adjacent time nodes in the smooth reflection trend sequence. First, based on the multi-phase main band reflectance image sequence recorded in the image, locate the data of each pixel at each time node. For each pixel number, extract its reflectance at each time point from t1 to t4. For example, the time series of pixel number 2 is 0.35, 0.30, 0.25, and 0.22. Calculate the reflectance difference between adjacent time nodes, i.e., ΔR 12 =-0.05, ΔR 23 =-0.05, ΔR 34 =-0.03, and a time slope change value sequence {-0.05, -0.05, -0.03} is formed in sequence. If the time intervals are equally spaced (for example, all 10 days), the difference can be used directly to represent the slope change; if the time intervals are uneven, the difference needs to be divided by the corresponding time span to reflect the reflectivity change rate after standardization, such as ΔR 23=(0.25-0.30) / 15=-0.0033. This operation is performed pixel by pixel, and on this basis, a pixel-level time slope change matrix of the entire image is generated. Furthermore, the main peak wavelength information of the main band is introduced to record the main peak wavelength change sequence of each pixel at each time node. For example, the main peak wavelength of pixel 2 is 682nm, 685nm, 688nm, and 691nm, and the corresponding change Δλ sequence is {+3nm, +3nm, +3nm}. This information will serve as an auxiliary criterion for subsequent combination with the time slope to form a two-dimensional index system. This sequence will be collaboratively judged with the time slope sequence. If the two show an offset trend at multiple nodes at the same time, the ability to recognize the dynamic change of the pixel point is enhanced, and finally the pixel-level time slope change sequence and the main peak wavelength offset combined data are formed, providing a basis for subsequent fitting and screening.

[0104] Table 2 Pixel time node reflectivity and main peak wavelength

[0105]

[0106]

[0107] As shown in Table 2, each pixel not only has a temporal reflectivity sequence, but also comes with main peak wavelength information, thus supporting the coordinated analysis of the two-dimensional dynamic changes of slope and wavelength.

[0108] S202: Calling the pixel-level time slope change sequence and the image main band reflectance time series data with geographic coordinate information, performing a fitting operation on the main band reflectance sequence of each pixel, extracting the trend slope value set of the pixel, and constructing the normal slope change range interval of the pixel based on the peak and valley values ​​in the trend slope value set to generate the normal range of regional slope change;

[0109] Call the pixel-level time slope change sequence and the image main band reflectance time series data with geographic coordinate information, perform sliding window linear fitting on the reflectance sequence of each pixel, extract the local change slope value, and form a trend slope set. In the fitting operation, the window length is selected as three time nodes, and the difference slope of two points is constructed in sequence. Taking pixel 1 as an example, if its main band reflectance time series is {0.12, 0.15, 0.18, 0.20}, and the node times are t1 = 0, t2 = 10, t3 = 20, t4 = 30 (unit: day), then the three groups of slopes are calculated as Construct a slope set {0.003, 0.003, 0.002}, then sort the set and determine the minimum value as the valley value k min =0.002, the maximum value is the peak k max =0.003, set the normal interval to Ik =[k min -δ,k max +δ], where δ is the deviation tolerance parameter, which is set to the historical sample standard deviation σ k For example, if the standard deviation of the pixel sample slope is σ k =0.0005, then set δ=2σ k =0.001, and the normal interval I is obtained k =[0.002-0.001,0.003+0.001]=[0.001,0.004]. This interval is used as the normal pixel change slope range for the next stage of offset judgment. In addition, in order to enhance the discrimination dimension, the main band center wavelength change judgment is introduced, and the main band center wavelength λ of each time node is recorded. i , construct the change sequence Δλ i =λ i+1 -λ i , assuming that the main band drift judgment threshold is λ0 = 2nm, if the main band wavelength of the pixel continuously changes to Δλ = {3, 2.5, 2}nm, then the main band offset condition Δλ is satisfied i >λ0, combined with the slope change results, it can assist in the subsequent determination of continuous abnormal areas and realize the joint reference of the main wavelength offset and the time slope sequence.

[0110] S203: Analyze the continuity of the pixel time slope change value based on the pixel-level time slope change sequence and the normal interval of regional slope change, mark the pixel positions that exceed the normal interval of regional slope change, and count the continuous exceeding time. Based on the exceeding time length, select the pixel points that meet the continuity offset and generate the candidate abnormal area position;

[0111] According to the pixel-level time slope change sequence and the normal interval of regional slope change, the slope value k of each time node of all pixels is calculated in turn. i With its normal variation range I k =[k min -δ,k max +δ] to perform a comparison operation, if the slope of a node That is k i <k min -δ or k i >k max +δ, it is marked as an offset point, and further sliding window scanning is performed to count the number of consecutive offset nodes L, that is, when multiple consecutive time nodes meet And when there is no interruption, it is recorded as a continuous offset sequence, and the continuous offset length judgment threshold T is set. When L≥T, the current pixel position is judged as an abnormal candidate pixel point. For example, the slope sequence of pixel 2 is {-0.05, -0.05, -0.03}, and the normal interval is [-0.06, -0.02]. Its three points are not in this interval. The continuous offset time is L=3, and the threshold T=3 nodes is set. Because L≥T, it is marked as an abnormal point. If the slope of a node returns to the normal interval, the continuity is broken, and the accumulated value is reset and recalculated. By executing this judgment logic within the entire pixel set, all pixel number sets that meet the continuous offset conditions can be obtained as candidate areas. At the same time, if the central wavelength change of the main band continuously exceeds the set threshold λ0, that is, Δλ i >λ0 holds true and coincides with the time node of the slope offset segment, the pixel is further marked as a collaborative offset pixel. This dual-index collaborative mechanism improves the discrimination accuracy and spatial consistency of candidate region screening;

[0112] Table 2 Pixel slope continuous offset determination table

[0113]

[0114] Table 2 lists the continuous deviation judgment results after comparing the slope change with the normal range. Combined with the main wavelength continuous deviation exceeding the threshold, it will further constitute the main slope-wavelength coordinated anomaly candidate area;

[0115] Continuous offset means that the slope change value of a certain pixel exceeds the normal range at multiple time nodes in the time series, and the exceeding state remains continuous and non-repeating.

[0116] See also Figure 1 , the specific steps of S3 are:

[0117] S301: Based on the image corresponding to the candidate abnormal area location, extract the pixel reflectance value sequence of the near-infrared, high green, and short-wave red edge bands in the image, construct the spectral reflectance curve of each pixel in the order of wavelength, extract and record the wavelength value corresponding to the main peak of the band spectral curve, and generate the main band peak wavelength sequence;

[0118] Based on the image corresponding to the candidate anomaly area, first extract the near-infrared, high green, and short-wave red edge bands from the image data, and read the pixel reflectance values ​​of the corresponding bands respectively. For each pixel in the image, the reflectance value of each band is extracted in sequence according to the wavelength order of high green (550nm), short-wave red edge (705nm), and near-infrared (850nm), and the spectral reflectance curve of the pixel is constructed. In the process of constructing the spectral curve, the wavelength sequence is set to 550, 705, 850, nm, and the reflectance value sequence is the extracted pixel value. For example, for a certain pixel, its reflectance values ​​are 0.32, 0.45, and 0.51 respectively. When constructing the curve, the wavelength is used as the horizontal axis and the reflectance value is used as the horizontal axis. The ordinate is plotted, and then in the spectral curve, the maximum reflectance value position is identified by traversing the reflectance value sequence, and the corresponding main peak wavelength value is extracted. That is, for the above reflectance sequence, the maximum value is 0.51, and the corresponding wavelength is 850nm. This wavelength value is recorded in the main band peak wavelength sequence, and the above operation is performed to traverse all pixels to finally generate the main band peak wavelength sequence. During the extraction process, special attention should be paid to the fact that the wavelength sequence is fixed at 550nm, 705nm, and 850nm, and the main peak is identified by direct maximum value comparison without introducing fitting or interpolation methods. All wavelength-reflectance data are directly read from the original pixels of the image data without the need for interpolation or resampling. Taking a candidate abnormal area in a specific image area as an example, the area size is 100×100 pixels. After extracting the three-band data from the image data, three matrix data are formed, namely the near-infrared reflectivity matrix, the high green reflectivity matrix, and the short-wave red edge reflectivity matrix. The matrix dimensions are all 100×100. Then, for the pixel at the (25,30) position, its high green reflectivity is read as 0.28, the short-wave red edge reflectivity is 0.42, and the near-infrared reflectivity is 0.39. The reflectivity curve is drawn in the order of wavelength, and the sequence is 0.28, 0.42, and 0.39. The maximum reflectivity is 0.42, corresponding to a wavelength of 705nm. Therefore, 705nm is recorded as the peak wavelength of the main band of the pixel. For ease of understanding, the pixel data of some candidate areas are shown in the following table:

[0119] Table 3 Main band peak wavelength extraction table

[0120] Line number Column number High green reflectivity Shortwave red edge reflectivity Near-infrared reflectivity Main peak wavelength (nm) 25 30 0.28 0.42 0.39 705 50 45 0.35 0.47 0.44 705 75 80 0.32 0.38 0.55 850

[0121] As shown in Table 3, the peak wavelength of the main band is extracted and determined according to the wavelength corresponding to the maximum reflectivity.

[0122] S302: Calling the main band peak wavelength sequence, calculating the offset value between the main peak wavelength of each pixel in the three bands and the central wavelength of the reference band, comparing the offset value with the set offset threshold, screening the pixel positions with offset values ​​greater than the threshold, and generating an offset abnormal pixel position group;

[0123] Call the generated main band peak wavelength sequence, read the main peak wavelength value of each pixel in turn, and calculate the offset value between the main peak wavelength and the reference center wavelength of the band to which it belongs according to the preset reference band center wavelength, which are high green 550nm, short wave red edge 705nm, and near infrared 850nm for each pixel. The specific operation is to subtract the reference center wavelength of the corresponding band from the main peak wavelength. The offset value is the absolute value of the difference between the main peak wavelength and the reference center wavelength. After the offset value set is generated, it is traversed and screened to remove pixels whose offset value is greater than the offset degree threshold. The positions are screened out to form an offset abnormal pixel position group. The offset degree threshold is set based on the statistical result of the total offset value distribution of all pixels in the candidate area. For example, there are 10,000 pixels in the candidate area. After calculating the offset values ​​of all pixels, a frequency histogram of the offset value is drawn. The offset threshold is set to the mean plus twice the standard deviation by statistically analyzing the mean and standard deviation. For example, if the mean is 10nm and the standard deviation is 5nm, the offset threshold is set to 20nm. Pixels with offset values ​​greater than 20nm are considered abnormal pixels and added to the offset abnormal pixel position group. Taking the pixels in Table 3 as an example, the main peak wavelength of the (25, 30) pixel is 705 nm. Compared with the short-wave red edge center wavelength of 705 nm, the offset value is |705-705| = 0 nm. The main peak wavelength of the (75, 80) pixel is 850 nm. Compared with the near-infrared center wavelength of 850 nm, the offset value is |850-850| = 0 nm. However, in other pixels, if the main peak wavelength is 780 nm, compared with the near-infrared center wavelength of 850 nm, the offset value is |780-850| = 70 nm, which exceeds the set threshold of 20 nm and is determined to be an abnormal pixel.

[0124] The reference band center wavelength is the position with the strongest spectral response extracted from the sensor response curve of the near-infrared, high green and short-wave red edge bands according to the band setting of the remote sensing imaging equipment;

[0125] The offset threshold is a criterion set by statistically analyzing the total amount of the main peak wavelength offset of the pixels in the candidate area and based on the central trend and discrete range of the offset degree.

[0126] S303: Based on the offset abnormal pixel position group and the two-dimensional spatial position range of the candidate abnormal area, a pixel segmentation map structure is constructed according to the spatial adjacent relationship, and the boundary contour of the connected pixel area is delineated according to the spectral boundary difference to establish the initial boundary data of the lesion;

[0127] Based on the position group of offset abnormal pixels and the two-dimensional spatial distribution of candidate abnormal areas, a two-dimensional matrix structure was established according to the row and column coordinate system. Clustering was performed according to the eight-neighborhood connectivity of the pixels, and adjacent offset abnormal pixels were merged into connected areas. By traversing the boundaries of each connected area, the reflectance difference between the boundary pixels and their outer neighboring pixels was calculated based on the difference in spectral reflectance curves of each pixel. Specifically, the absolute difference in reflectance between the boundary pixels and the outer neighboring pixels was calculated at three wavelengths of 550nm, 705nm, and 850nm. The difference judgment threshold was set to 0.15. If the reflectance difference in any band was greater than 0.15, it was judged as a boundary pixel, and the boundary contour of the connected area was drawn to finally generate the initial boundary data of the lesion. Taking a cluster connected area as an example, which contains 50 offset abnormal pixels, the boundary pixels are traversed and the reflectance difference between each boundary pixel and its neighboring pixels is calculated in three bands. Assuming that the reflectance of the boundary pixel in the near-infrared band is 0.50 and that of the neighboring pixel is 0.32, the difference is |0.50-0.32|=0.18, which is greater than 0.15 and is determined to be the boundary. The boundary determination of all boundary pixels is completed step by step, and finally the boundary contour is closed to obtain the initial boundary data of the lesion.

[0128] See also Figure 1 , the steps of S4 are as follows:

[0129] S401: Based on the initial boundary data of the lesion, the band values ​​of the multispectral image pixels corresponding to the boundary are extracted, the visible light, near infrared and short wave infrared band values ​​of the corresponding pixel positions in the reflectance data are called, the pixel data are associated and identified according to the spatial coordinates, and a spectral value set of the boundary pixel is generated;

[0130] Based on the initial boundary data of the lesion, the boundary vector coordinate set of the lesion area needs to be extracted first. For each boundary point (x, y), the multispectral image data after georeferencing is called, and the image pixel number that completely corresponds to its spatial position is read in sequence. The reflectance value of the pixel in the visible light, near infrared and short wave infrared bands is obtained to form a spectral value record with a geographic location index. The extraction process needs to traverse all boundary points and perform spatial precise matching at the pixel level. For example, the boundary of a lesion area contains 3 boundary pixels, numbered 101, 102, and 103, and the corresponding band values ​​on the image are: red light reflectance, near infrared reflectance, and short wave infrared reflectance. The emissivity is 0.12, 0.18, and 0.25, the green reflectivity is 0.20, 0.22, and 0.24, the near infrared is 0.45, 0.35, and 0.30, and the short-wave infrared is 0.28, 0.33, and 0.30. The boundary pixel spectral value set is generated by combining these pixel spatial coordinates and band values. The fields in the data structure include "pixel number", "geographic coordinates (x, y)", "red light reflectivity", "green light reflectivity", "near infrared reflectivity", "short-wave infrared reflectivity", etc. Each pixel point is used as a recording unit and corresponds one-to-one to the geographic data of the lesion boundary.

[0131] S402: Based on the spectral value set of the boundary pixels and the reflectance combination relationship between the bands, multiple vegetation status assessment indicators are extracted, including the normalized vegetation index and the normalized water index. Based on the typical threshold range of each indicator, the offset of multiple indicators is identified, the vegetation health status is classified, and the vegetation health classification reflectance interval is generated;

[0132] According to the spectral value set of the boundary pixels, the reflectance value of each pixel is sequentially calculated. At this stage, the normalized vegetation index (NDVI) and normalized water index (NDWI) need to be extracted to reflect the photosynthetic activity and water status of the vegetation, respectively. The calculation of NDVI uses red light and near-infrared reflectance as input, and the calculation formula is: Taking pixel 101 as an example, its red light is 0.12 and near infrared is 0.45, then The calculation formula of NDWI is: The green light of pixel 101 is 0.20 and the near infrared is 0.45, then After completing the calculation of NDVI and NDWI for all boundary pixels, the vegetation health level is determined according to the standard division interval. The NDVI health level interval is set as follows: healthy is [0.5, 1], sub-healthy is [0.3, 0.5), unhealthy is [0, 0.3), and the NDWI health level interval is set as follows: healthy is [-0.4, -0.2], sub-healthy is [-0.2, -0.1], unhealthy is [-0.1, 0.0]. Through the double-index cross-matching rule, if NDVI and NDWI are If both WI fall into the healthy range, the vegetation health level is marked as "healthy". For example, if the NDVI is 0.58 and the NDWI is -0.38, it matches the healthy level. If the NDVI is 0.32 and falls within the range [0.3, 0.5), and the NDWI is -0.23 and falls within the range [-0.4, -0.2], it is marked as "sub-healthy". If the NDVI is 0.09 and falls within the range [0, 0.3), and the NDWI is -0.11, it is "unhealthy". The above levels correspond to the reflectance intervals of vegetation health classification.

[0133] Table 4 Boundary pixel multispectral and evaluation index table

[0134]

[0135]

[0136] As shown in Table 4 , there are differences in the NDVI and NDWI indicators of different boundary pixels, which are divided into different vegetation health levels after being judged by the standard threshold interval;

[0137] The typical threshold range refers to the standard value range of the spectral evaluation index used to analyze the vegetation status of boundary pixels at different health levels;

[0138] The Normalized Difference Vegetation Index is used to quantify the difference in the ratio of red and near-infrared reflectance of boundary pixels and to judge the photosynthetic activity and health level of vegetation;

[0139] The normalized water index is used to reflect the imbalance between the reflectance of green light and near-infrared bands in boundary pixels.

[0140] S403: Calling the vegetation health classification reflectance interval, performing state matching on the spatial position in the initial boundary data of the diseased spot, analyzing the corresponding pixel band value and the classification interval value, spatially annotating the vegetation state according to the classification label, and outputting the vegetation pixel encoding to generate classification annotation result data;

[0141] The vegetation health classification reflectance interval is called to match the classification status of each pixel spatial position in the initial boundary data of the diseased spot. In the specific execution process, the current red, green, and near-infrared reflectance values ​​of each boundary pixel are located with the index values ​​shown in Table 3, and the NDVI and NDWI results are compared to see if they fall into the standard grade interval. If they match, the classification result is assigned according to the corresponding label of the interval. In this process, it is necessary to traverse in sequence according to the pixel number. For example, for pixel number 101, it is known that its NDVI is 0.58 and NDWI is -0.38, both of which are within the health grade standard interval, then it is assigned a value of "healthy" and the generated code value is 1; Pixel 102, NDVI is 0.32, NDWI is -0.23, matches sub-health, coded as 2; pixel 103, NDVI is 0.09, NDWI is -0.11, matches unhealthy, coded as 3. Finally, the spatial positions, evaluation indicators and health codes of all pixels are output as vector or raster result datasets, where the code field is used to identify different categories of pixels in the subsequent processing module. For example, healthy area 1 can be used as a retained object, and unhealthy area 3 is a removal target. Through this step, the health status of the lesion boundary pixels is spatially graded and labeled, and the data foundation is laid for subsequent disease spread modeling or intervention evaluation.

[0142] See also Figure 1 , the specific steps of S5 are:

[0143] S501: Calling the classification and annotation result data, identifying the corresponding position of the pest and disease annotated area in the forest remote sensing image, projecting and matching the spatial coordinate information to the forest spatial coordinate system, overlaying the geometric boundaries of the classified pixels and establishing a pest and disease spatial distribution information layer, and generating pest and disease layer overlay data;

[0144] To call the classification annotation result data, we first need to extract the spatial location information and health grade classification fields of the classified pixels, and use them as the key geometry and attribute content in the vector dataset. Then, we perform a spatial projection transformation operation to convert the image coordinate system (such as WGS 84) of the classification result to the geographic projection coordinate system (such as CGCS2000 projection) used by the forest resource database. When performing the transformation, we need to read the longitude and latitude of the pixel center point and apply the projection transformation parameters to convert the geographic coordinates (X_lon, Y_lat) of each classified pixel into projection coordinates (X_proj, Y_proj). For example, the center point of the pixel numbered A001 is (122.12°, 42.31°). After the projection transformation, we get (432100.0, 3121000.0). Then, based on the merging operation of adjacent pixels of the same category, we construct their spatial boundaries with surface vectors to form pest and disease classification blocks. It is then overlaid on the existing forest remote sensing image base map. The current classification result is embedded in the GIS environment as a new layer using a layer management structure and named "Pest and Disease Layer 2025-05-01". Each pest and disease patch in the layer is identified by an area number, such as A001, A002, and A003, corresponding to the unhealthy, sub-healthy, and healthy classification results, respectively. All surface graphics and attribute fields are combined to form an independent data layer. The fields include "area number", "projection X", "projection Y", "classification label", etc., and are exported in shapefile or GeoJSON format as pest and disease layer overlay data.

[0145] S502: Based on the overlay data of the pest and disease layer, the time label and spatial number corresponding to each pest and disease area in the layer are matched, the data unit with the same number is located in the forest area database, the pest and disease occurrence status in the database field is updated, and the spatial description field content is supplemented to obtain the pest and disease monitoring update record;

[0146] Based on the overlay data of the pest and disease layer, the contents of the "layer time" and "area number" fields in the layer attributes are first extracted, and the corresponding data units in the forest database are located through the spatial attribute matching mechanism. The specific process is to match the spatial index field in the database, match the number A001 in the layer with the field "spatial number" in the database, and at the same time perform spatial overlap judgment on the layer coordinate point (X_proj, Y_proj) and the geometric boundary of the forest patch recorded in the database. After confirming that it falls into the same management unit, the attribute update operation is performed on the record, marking the "pest and disease status" field as "pest and disease status". If the layer is marked as "healthy", this field will not be updated. In addition, the "Remarks" field is added, and the pest type description is added according to the accompanying description in the classification annotation record. For example, the A001 area is added with the note "spread of red pine caterpillars", and the A002 area is added with "spread of larch leaf blight". A complete pest and disease monitoring update record is formed in the attribute table. The fields include "spatial number", "layer time", "classification label", "database status update", "spatial position (X, Y)" and "remarks". This information can eventually be exported as an update log record document for updating the pest and disease database of the forestry data platform;

[0147] Table 5 Spatial update record of pest and disease layer

[0148]

[0149] As shown in Table 5, each area number in the pest and disease layer overlay data has completed the location matching and status update operations in the corresponding forest database, and is accompanied by explanatory information to form a structured monitoring record.

[0150] S503: Calling pest monitoring update records, analyzing pest labeling of spatial units, identifying the direction of change and expansion rate of pest areas, predicting the change trend of pest areas and marking the trend level, and generating pest spatial distribution trend results;

[0151] The specific formula for identifying the direction of change and expansion rate of pest and disease areas is:

[0152]

[0153] Calculate the regional expansion rate value V;

[0154] in, Represents the expansion rate value of the area with spatial unit number j in direction type d. and represent the pest and disease labeling values ​​of the kth pixel in the spatial unit j at time t and t-1, respectively, j,k represents the pixel center distance of the kth pixel in the spatial unit j along the direction d, Tj,k Represents the topological disturbance factor of the boundary of the region where the k-th pixel is located in this direction, Respectively represent the pest and disease area of ​​spatial unit j at time t and t-1, n j Represents the number of pixels involved in the calculation, j is the number index of the spatial unit, t is the time index, and k is the pixel index;

[0155] formula:

[0156]

[0157] Detailed explanation of the formula and the process of formula calculation and derivation:

[0158] The formula is used to calculate the pest area expansion rate value of spatial unit j in direction d The results are used to assess the expansion trend and speed of pests and diseases in a specific direction;

[0159] Parameter meaning and setting value:

[0160] T j,k : The topological disturbance factor of the boundary of the region where the k-th pixel is located in direction d, which is calculated by analyzing the difference in pest and disease status between the pixel and its neighboring pixels, and the value range is between 0 and 1;

[0161] The area of ​​pests and diseases in spatial unit j at time t is obtained by counting the number of damaged pixels in the remote sensing image and multiplying it by the area of ​​a single pixel;

[0162] The area of ​​pests and diseases in spatial unit j at time t-1 is obtained by counting the number of damaged pixels in the remote sensing image and multiplying it by the area of ​​a single pixel;

[0163] Assume that in spatial unit j, there are 5 pixels along direction d where changes in pest and disease status occur. The parameters are as follows:

[0164] Pixel 1: D j,1 =10 m, T j,1 =0.3;

[0165] Pixel 2: D j,2 =10 m, T j,2 =0.4;

[0166] Pixel 3: D j,3 =10 m, T j,3 =0.2;

[0167] Pixel 4: D j,4=10 m, T j,4 =0.5;

[0168] Pixel 5: D j,5 =10 m, T j,5 =0.3;

[0169] Set the pest and disease area:

[0170] square meters;

[0171] square meters;

[0172] Substitute the parameters into the formula for calculation:

[0173]

[0174] The result of 1.723 indicates that the pest area expansion rate in spatial unit j in direction d is 1.723. This result is used to assess the expansion trend and speed of pests and diseases in a specific direction and helps to formulate targeted prevention and control measures.

[0175] See also Figure 2 A forest pest and disease intelligent monitoring system based on UAV remote sensing is provided. The forest pest and disease intelligent monitoring system based on UAV remote sensing is used to implement the above-mentioned forest pest and disease intelligent monitoring method based on UAV remote sensing. The system includes:

[0176] The multispectral time series processing module obtains vegetation canopy reflectance data through the multispectral imaging equipment carried by the UAV, uses the moving average method to perform time series smoothing on the reflectance data, generates a smoothed reflectance trend sequence, and passes it to the time series anomaly detection module;

[0177] The time series anomaly detection module calls the smoothed reflectance trend sequence, performs differential processing on the data of adjacent time nodes, constructs a time slope change sequence, uses a linear regression model to fit the reflectance change trend, calculates the standard deviation of the slope residual, constructs the normal slope change interval, identifies the pixel coordinates that exceed the interval for multiple consecutive cycles, generates candidate anomaly area locations, and passes them to the lesion boundary generation module;

[0178] The lesion boundary generation module extracts the red edge band and near-infrared band reflectance curves based on the location of the candidate abnormal area, calculates the main peak wavelength offset, analyzes the degree of pixel offset, screens pixels with significant change characteristics, and segments the spectral boundary image based on spatial distribution to generate initial lesion boundary data, which is then passed to the health status classification module.

[0179] The health status classification module obtains the initial boundary data of the diseased spots, extracts the multispectral characteristic parameters of the corresponding pixels, and uses the reflectance combination relationship between the bands to classify the vegetation health status and perform classification annotation. It generates classification and annotation result data and passes it to the spatial trend prediction module;

[0180] The spatial trend prediction module calls the classification and labeling result data, converts the coordinates of pests and diseases into the geographic coordinate system, generates a spatial distribution thermal layer, calculates the regional pest and disease density change rate, updates the forest area pest and disease monitoring database, combines time information to predict the development trend of pests and diseases, and generates the spatial distribution trend results of pests and diseases.

[0181] The above embodiments can be implemented in whole or in part through software, hardware (such as circuits), firmware, or any other combination. When implemented using software, the above embodiments can be implemented in whole or in part in the form of a computer program product. The computer program product includes one or more computer instructions or computer programs. When the computer instructions or computer program are loaded or executed on a computer, the processes or functions described in accordance with the embodiments of the present invention are generated in whole or in part. The computer can be a general-purpose computer, a special-purpose computer, a computer network, or other programmable device. The computer instructions can be stored in a computer-readable storage medium or transmitted from one computer-readable storage medium to another. For example, the computer instructions can be transmitted from one website, computer, server, or data center to another website, computer, server, or data center via a wired method (such as infrared, wireless, microwave, etc.). The computer-readable storage medium can be any available medium that can be accessed by a computer or a data storage device such as a server or data center that contains a collection of one or more available media. The available medium can be a magnetic medium (such as a floppy disk, hard disk, or magnetic tape), an optical medium (such as a DVD), or a semiconductor medium. The semiconductor medium can be a solid-state drive.

[0182] It should be understood that the term "and / or" as used herein simply describes a relationship between associated objects, indicating that three possible relationships exist. For example, "A and / or B" can represent: A alone, A and B together, or B alone. A and B can be singular or plural. Furthermore, the character " / " as used herein generally indicates an "or" relationship between the associated objects, but it may also indicate an "and / or" relationship. For specific understanding, please refer to the context.

[0183] In this disclosure, "at least one" means one or more, and "plurality" means two or more. "At least one of the following" or similar expressions refers to any combination of these items, including any combination of single or plural items. For example, "at least one of a, b, or c" can mean: a, b, c, ab, ac, bc, or abc, where a, b, and c can be single or plural.

[0184] It should be understood that in various embodiments of the present invention, the size of the serial numbers of the above-mentioned processes does not mean the order of execution. The execution order of each process should be determined by its function and internal logic, and should not constitute any limitation on the implementation process of the embodiments of the present invention.

[0185] Those skilled in the art will appreciate that the units and algorithm steps of each example described in conjunction with the embodiments disclosed herein can be implemented in electronic hardware, or a combination of computer software and electronic hardware. Whether these functions are performed in hardware or software depends on the specific application and design constraints of the technical solution. Professionals and technicians can use different methods to implement the described functions for each specific application, but such implementation should not be considered beyond the scope of the present invention.

[0186] Those skilled in the art will clearly understand that, for the convenience and brevity of description, the specific working processes of the above-described equipment, devices and units can refer to the corresponding processes in the aforementioned method embodiments and will not be repeated here.

[0187] In the several embodiments provided by the present invention, it should be understood that the disclosed devices, apparatuses and methods can be implemented in other ways. For example, the device embodiments described above are merely illustrative. For example, the division of the units is merely a logical function division. In actual implementation, there may be other division methods, such as multiple units or components can be combined or integrated into another device, or some features can be ignored or not executed. Another point is that the mutual coupling or direct coupling or communication connection shown or discussed can be through some interface, indirect coupling or communication connection of the device or unit, which can be electrical, mechanical or other forms.

[0188] The units described as separate components may or may not be physically separate, and the components shown as units may or may not be physical units, that is, they may be located in one place or distributed across multiple network units. Some or all of these units may be selected to achieve the purpose of this embodiment according to actual needs.

[0189] In addition, each functional unit in each embodiment of the present invention may be integrated into one processing unit, or each unit may exist physically separately, or two or more units may be integrated into one unit.

[0190] If the functions are implemented in the form of software functional units and sold or used as independent products, they can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of the present invention, or the part that contributes to the prior art, or the part of the technical solution, can be embodied in the form of a software product. The computer software product is stored in a storage medium and includes several instructions for enabling a computer device (which can be a personal computer, server, or network device, etc.) to execute all or part of the steps of the method described in each embodiment of the present invention. The aforementioned storage medium includes various media that can store program codes, such as a USB flash drive, a mobile hard disk, a read-only memory (ROM), a random access memory (RAM), a magnetic disk, or an optical disk.

[0191] The above description is merely a specific embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any modifications or substitutions that can be easily conceived by a person skilled in the art within the technical scope disclosed in the present invention should be included within the scope of protection of the present invention. Therefore, the scope of protection of the present invention should be based on the scope of protection of the claims.

Claims

1. An intelligent monitoring method for forest pests and diseases based on UAV remote sensing, characterized in that: The method comprises: S1: Using multispectral imaging equipment carried by drones, we collected and organized forest remote sensing image data, identified the main band reflectance of the vegetation canopy area, and constructed time series data. We then processed the data using the moving average method to generate a smooth reflectance trend sequence. S2: calling the smoothed reflection trend sequence, performing differential processing on the data of adjacent time nodes, constructing a time slope change sequence, fitting the reflectivity change trend, constructing a normal slope change interval, identifying pixels whose reflectivity changes continuously deviate from the normal interval range, analyzing the persistence of the deviation, and identifying the location of candidate abnormal areas; S3: Based on the candidate abnormal area position image, extract the reflection curves of multiple bands, locate the main peak wavelength and calculate the abnormal offset, analyze the pixel offset degree, screen the pixels with significant change characteristics, and segment the spectral boundary image based on the spatial distribution to generate the initial boundary data of the lesion; S4: Based on the initial boundary data of the diseased spot, the multispectral characteristic parameters of the corresponding pixel are extracted, and the reflectance combination relationship between the bands is used to classify the vegetation health status and perform classification annotation to generate classification annotation result data.

2. The method for intelligent monitoring of forest pests and diseases based on UAV remote sensing according to claim 1 is characterized in that: The smoothed reflectance trend sequence includes the main band reflectance time series, the smoothed value of the time series data, and the moving average characteristic value; the time slope change sequence includes the reflectance difference value of adjacent nodes, the normal change interval range, and the continuous deviation index; the initial boundary data of the lesion includes the distribution of significantly changed pixels, the main peak wavelength, and the spectral boundary segmentation map; the classification annotation result data includes the multispectral feature parameter set, the vegetation health level label, and the classification area map.

3. The method for intelligent monitoring of forest pests and diseases based on UAV remote sensing according to claim 1 is characterized in that: The steps of S1 are as follows: S101: Use a drone equipped with multispectral imaging equipment to collect remote sensing image data of the forest area. Number the images collected at each time point in chronological order, identify the green vegetation canopy area in the image, extract the reflectance value corresponding to the main band, and generate a vegetation canopy reflectance dataset. S102: Based on the vegetation canopy reflectance dataset, integrating the main band reflectance values ​​corresponding to the time nodes in the order of image shooting time, establishing a time series reflectance data structure of the vegetation canopy area, and sorting the reflectance values ​​in time series to obtain a canopy reflectance time series value sequence; S103: calling the canopy reflectance time series value sequence, setting the observation window width according to the continuously observed reflectance values ​​at each time point, smoothing each set of reflectance data sequence using a moving average method, and generating a smoothed reflectance trend sequence; The specific formula for setting the observation window width is: Calculate the observation window width value; Among them, W i Indicates the observation window width value at the i-th time point, Z c Indicates the absolute value of the canopy reflectance difference between the i+cth and icth time points when the band is λ, Q i Indicates the absolute value of the deviation of the canopy reflectance at time point i relative to the sequence mean when the band is λ, L c represents the average time interval of the cth symmetric time point relative to the central time point, U i represents the horizontal velocity of the observation platform at time point i, Θ represents the standard deviation of the canopy reflectance series in band λ, m represents the number of symmetrical time point pairs used for calculation, and c represents the index number of the time point pairs symmetrically distributed forward and backward with respect to the central time point i; The moving average method sets a fixed time window for the vegetation canopy reflectance time series value sequence and calculates the reflectance mean value in each interval by sliding.

4. The method for intelligent monitoring of forest pests and diseases based on UAV remote sensing according to claim 3 is characterized in that: The specific steps of S2 are: S201: Calling the data values ​​of adjacent time nodes in the smooth reflection trend sequence, obtaining the reflectivity value of the pixel position at the time node, and calculating the difference of the reflectivity at the adjacent time nodes, constructing the time slope change value sequence of the pixel point, and generating the pixel-level time slope change sequence; S202: calling the pixel-level time slope change sequence and the image main band reflectance time series data with geographic coordinate information, performing a fitting operation on the main band reflectance sequence of each pixel, extracting the trend slope value set of the pixel, and constructing the normal slope change range interval of the pixel based on the peak and valley values ​​in the trend slope value set to generate the normal range of regional slope change; S203: Analyzing the continuity of the pixel time slope change value based on the pixel-level time slope change sequence and the normal interval of regional slope change, marking the pixel positions that exceed the normal interval of regional slope change, and counting the continuous exceeding time. Based on the exceeding time length, the pixel point positions that meet the continuity offset are screened to generate candidate abnormal area positions; The continuity deviation refers to the situation where the slope change value of a certain pixel exceeds the normal range at multiple time nodes in the time series, and the exceeding state remains continuous and non-repeating.

5. The method for intelligent monitoring of forest pests and diseases based on UAV remote sensing according to claim 4 is characterized in that: The specific steps of S3 are: S301: Based on the image corresponding to the candidate abnormal area location, extract the pixel reflectance value sequence of the near-infrared, high green, and short-wave red edge bands in the image, construct the spectral reflectance curve of each pixel in the order of wavelength, extract and record the wavelength value corresponding to the main peak of the band spectral curve, and generate the main band peak wavelength sequence; S302: calling the main band peak wavelength sequence, calculating the offset value between the main peak wavelength of each pixel in the three bands and the central wavelength of the reference band, comparing the offset value with a set offset threshold, screening pixel positions with offset values ​​greater than the threshold, and generating an offset abnormal pixel position group; The reference band center wavelength is the position with the strongest spectral response extracted from the sensor response curves of the near infrared, high green and short-wave red edge bands according to the band setting of the remote sensing imaging device; The offset threshold is a criterion set by statistically analyzing the total amount of the main peak wavelength offset of the pixels in the candidate area and according to the central tendency and discrete range of the offset degree; S303: Based on the offset abnormal pixel position group and the two-dimensional spatial position range of the candidate abnormal area, a pixel segmentation map structure is constructed according to the spatial adjacent relationship, and the boundary contour of the connected pixel area is delineated according to the spectral boundary difference to establish the initial boundary data of the lesion.

6. The method for intelligent monitoring of forest pests and diseases based on UAV remote sensing according to claim 5 is characterized in that: The specific steps of S4 are: S401: Based on the initial boundary data of the lesion, extract the band values ​​of the multispectral image pixels corresponding to the boundary, call the visible light, near infrared and short wave infrared band values ​​of the corresponding pixel positions in the reflectance data, associate and identify the pixel data according to the spatial coordinates, and generate a spectral value set of the boundary pixels; S402: Based on the boundary pixel spectral value set, using the reflectance combination relationship between bands, extracting multiple vegetation status assessment indicators, including the normalized vegetation index and the normalized water index, and identifying the offsets of the multiple indicators based on the typical threshold range of each indicator, classifying the vegetation health status, and generating vegetation health classification reflectance intervals; The typical threshold range refers to the standard value range of the spectral evaluation index used to analyze the vegetation status of boundary pixels at different health levels; The normalized vegetation index is used to quantify the difference in the ratio of red light and near-infrared reflectance of boundary pixels and to judge the photosynthetic activity and health level of vegetation; The normalized water index is used to reflect the imbalance between the green light and near-infrared reflectance of the boundary pixels; S403: Call the vegetation health grade reflectance interval, perform state matching on the spatial position in the initial boundary data of the diseased spot, analyze the corresponding pixel band value and classification interval value, perform spatial annotation of the vegetation state according to the classification label, and output code the vegetation pixel to generate classification annotation result data.

7. The method for intelligent monitoring of forest pests and diseases based on UAV remote sensing according to claim 1, characterized in that: The method further comprises: S5: Calling the classification and labeling result data, superimposing the spatial location of pests and diseases onto the forest spatial coordinate system, generating a pest and disease spatial distribution layer, updating the forest area pest and disease monitoring database based on the layer data, and predicting the pest and disease development trend in combination with time information to generate a pest and disease spatial distribution trend result; The results of the spatial distribution trend of pests and diseases include a spatial distribution layer of pests and diseases, a monitoring database update record, and a development trend forecast map.

8. The method for intelligent monitoring of forest pests and diseases based on UAV remote sensing according to claim 7 is characterized in that: The specific steps of S5 are: S501: Calling the classification and labeling result data, identifying the corresponding position of the pest and disease labeling area in the forest remote sensing image, projecting and matching the spatial coordinate information to the forest spatial coordinate system, overlaying the geometric boundaries of the classified pixels and establishing a pest and disease spatial distribution information layer, and generating pest and disease layer overlay data; S502: Based on the overlay data of the pest and disease layer, the time label and spatial number corresponding to each pest and disease area in the layer are matched, the data unit with the same number is located in the forest area database, the pest and disease occurrence status in the database field is updated, and the spatial description field content is supplemented to obtain the pest and disease monitoring update record; S503: Call the pest monitoring update record, analyze the pest labeling status of the spatial unit, identify the change direction and expansion rate of the pest area, predict the change trend of the pest area and mark the trend level, and generate the pest spatial distribution trend result.

9. The method for intelligent monitoring of forest pests and diseases based on UAV remote sensing according to claim 8, characterized in that: The specific formula for identifying the change direction and expansion rate of the pest and disease area is: Calculate the regional expansion rate value V; in, Represents the expansion rate value of the area with spatial unit number j in direction type d. and represent the pest and disease labeling values ​​of the kth pixel in the spatial unit j at time t and t-1, respectively, j,k represents the pixel center distance of the kth pixel in the spatial unit j along the direction d, T j,k Represents the topological disturbance factor of the boundary of the region where the k-th pixel is located in this direction, Respectively represent the pest and disease area of ​​spatial unit j at time t and t-1, n j Represents the number of changed pixels involved in the calculation, j is the number index of the spatial unit, t is the time index, and k is the pixel index.

10. An intelligent forest pest monitoring system based on UAV remote sensing, characterized by: The system is used to implement the intelligent forest pest monitoring method based on drone remote sensing according to any one of claims 1 to 9, and the system comprises: The multispectral time series processing module obtains vegetation canopy reflectance data through the multispectral imaging equipment carried by the UAV, uses the moving average method to perform time series smoothing on the reflectance data, generates a smoothed reflectance trend sequence, and passes it to the time series anomaly detection module; The time series anomaly detection module calls the smoothed reflection trend sequence, performs differential processing on the data of adjacent time nodes, constructs a time slope change sequence, uses a linear regression model to fit the reflectivity change trend, calculates the slope residual standard deviation, constructs a normal slope change interval, identifies the pixel coordinates that exceed the interval for multiple consecutive cycles, generates candidate abnormal area locations, and transmits them to the lesion boundary generation module; The lesion boundary generation module extracts the red edge band and near-infrared band reflectance curves based on the candidate abnormal area location, calculates the main peak wavelength offset, analyzes the pixel offset degree, selects pixels with significant change characteristics, and segments the spectral boundary image based on the spatial distribution to generate initial lesion boundary data, which is then passed to the health status classification module; The health status classification module obtains the initial boundary data of the diseased spot, extracts the multispectral characteristic parameters of the corresponding pixel, classifies the vegetation health status and performs classification annotation based on the reflectance combination relationship between the bands, generates classification annotation result data, and transmits it to the spatial trend prediction module; The spatial trend prediction module calls the classification and annotation result data, converts the coordinates of pests and diseases into the geographic coordinate system, generates a spatial distribution thermal layer, calculates the regional pest and disease density change rate, updates the forest area pest and disease monitoring database, combines time information to predict the development trend of pests and diseases, and generates the spatial distribution trend results of pests and diseases.

Images