A crop monitoring method integrating lidar and multispectral imaging

Abstract

The present invention provides a crop monitoring method integrating laser radar and multispectral imagery, extracting laser radar point cloud data and multispectral satellite imagery data of target crops and performing data preprocessing, unifying the data into the same coordinate system for data fusion, extracting the structural parameters of the target crop canopy, estimating the physiological parameters of the target crop canopy, comprehensively analyzing the structural parameters and the physiological parameters, constructing quantitative crop growth assessment indicators, simulating the ideal growth state of the target crop through a crop growth model, judging whether the growth of the target crop is abnormal, generating abnormal diagnosis results, and finally obtaining growth warning information and management decision recommendations of the target crop. The present invention performs a coordinated characterization of the three-dimensional structure and spectral characteristics of crops, constructs a quantitative comprehensive evaluation index, and dynamically evaluates the growth of crops, which can significantly improve the accuracy, timeliness, and quantification level of crop growth monitoring.

Description

Technical Field

[0001] The present invention belongs to the technical field of crop monitoring, and in particular relates to a crop monitoring method integrating laser radar and multispectral imaging. Background Art

[0002] In the existing agricultural science field, crop growth monitoring is often carried out to detect crop growth and evaluate crop growth conditions. However, traditional crop monitoring solutions are often not comprehensive and the evaluation of crop conditions is not accurate enough. The specific problems are as follows:

[0003] (1) Existing crop growth monitoring technologies mainly use a single remote sensing data source, such as optical images or radar data, which can only obtain limited crop canopy information and is difficult to fully reflect the structural characteristics and physiological conditions of crops;

[0004] (2) Traditional crop growth analysis mostly uses a single indicator or a simple combination of indicators, lacking comprehensive quantitative indicators and making it difficult to comprehensively evaluate crop growth status;

[0005] (3) Previous crop growth monitoring was mostly based on instantaneous observations in a specific or single period, lacking dynamic analysis of the crop growth process and having limitations in identifying abnormal crop growth.

[0006] Therefore, we need to develop a crop monitoring method that integrates lidar and multispectral imaging, which can monitor the comprehensive growth of crops and dynamically evaluate the growth status of crops based on multiple indicators, which can significantly improve the accuracy, timeliness and quantification level of crop growth monitoring. Summary of the invention

[0007] The purpose of the present invention is to provide a crop monitoring method integrating lidar and multispectral imaging to solve the problems of incomplete monitoring, inaccurate evaluation and inability to dynamically evaluate in the existing crop growth monitoring schemes mentioned in the above background technology.

[0008] To achieve the above object, the present invention provides a crop monitoring method integrating laser radar and multispectral imagery, the method is specifically as follows:

[0009] Extracting laser radar point cloud data and multispectral satellite image data of the target crop, and performing data preprocessing on the laser radar point cloud data and the multispectral satellite image data to obtain first point cloud data and first image data of the target crop respectively;

[0010] Unifying the first point cloud data and the first image data into the same coordinate system through data registration, and performing data fusion to obtain first fused data of the target crop;

[0011] Based on the first fused data, extracting structural parameters of the canopy of the target crop for describing the three-dimensional growth state of the target crop, wherein the structural parameters include but are not limited to canopy height, canopy coverage and leaf area index;

[0012] Based on the first fused data, the physiological parameters of the canopy of the target crop are estimated by a radiation transfer model inversion method to characterize the biochemical condition of the target crop, wherein the physiological parameters include but are not limited to chlorophyll content and leaf nitrogen content;

[0013] Comprehensively analyzing the structural parameters and the physiological parameters, constructing a quantitative crop growth assessment index, simulating the ideal growth state of the target crop through a crop growth model, determining whether the growth of the target crop is abnormal, and generating an abnormal diagnosis result;

[0014] Generating growth warning information and management decision suggestions for the target crop according to the abnormal diagnosis results;

[0015] The step of judging whether the growth of the target crop is abnormal includes: assuming that the crop growth evaluation index of the target crop in the i-th monitoring period is , calculated by the principal component analysis model, represents the actual growth state of the target crop; the simulation result of the crop growth model of the target crop in the same period is , represents the ideal growth state of the target crop, and the relative difference between the actual growth state and the ideal growth state It is expressed as: ,in, It is used to indicate the abnormal degree of growth of the target crop in the i-th monitoring period. A positive value indicates that the growth is lower than the ideal state, and a negative value indicates that the growth is better than the ideal state.

[0016] Based on the above scheme, the data fusion adopts the Gram-Schmidt functional transformation method to fuse the point cloud digital surface model after the data registration and the first image data in the feature space to obtain the first fused data, which specifically includes:

[0017] Resampling the point cloud digital surface model to a coordinate system consistent with the first image data, and performing standardization processing on the point cloud digital surface model and the first image data respectively;

[0018] Orthogonalizing the feature vectors of the point cloud digital surface model and the first image data by using a Gram-Schmidt functional transformation method, wherein the first principal component direction is selected as the first orthogonal direction during the orthogonal transformation;

[0019] The orthogonally transformed feature vectors are weighted fused, and the feature vectors are converted into the original data space through the inverse Gram-Schmidt functional transformation method to obtain the fused point cloud digital surface model and the first image data, that is, the first fused data.

[0020] Based on the above scheme, the extraction of the structural parameters of the target crop canopy specifically includes:

[0021] Using a point cloud digital surface model construction algorithm, a high-resolution point cloud digital surface model is generated from the first fused data, and a difference method is used to extract crop canopy height information from the high-resolution point cloud digital surface model;

[0022] extracting vegetation coverage information from the first fused data using a normalized difference vegetation index threshold method;

[0023] spatially superimposing the crop canopy height information and the vegetation coverage information, calculating the canopy coverage of each pixel using the height information and coverage information of the crop canopy height information, and calculating the canopy coverage of each pixel in the first fused data;

[0024] Based on the crop canopy height information and the normalized difference vegetation index information, an empirical regression model is used to estimate the leaf area index.

[0025] Based on the above scheme, the estimating of the physiological parameters of the target crop canopy specifically includes:

[0026] Based on the PROSAIL canopy reflectance model, the radiation transfer process of the target crop canopy under different parameter conditions is simulated, the input parameter range and step size of the PROSAIL canopy reflectance model are set according to prior knowledge, a series of parameter-canopy reflectance pairs are obtained through forward simulation, and a lookup table of the parameter-canopy reflectance is constructed;

[0027] The physiological parameter inversion is realized through the lookup table, the canopy reflectance in the first fused data is extracted, and the matching spectrum is searched in the lookup table to obtain the parameter inversion value corresponding to the canopy reflectance, and the single pixel parameter inversion is completed. The physiological parameter inversion is performed on the full image element in the first fused data to obtain a physiological parameter inversion map.

[0028] Based on the above scheme, the comprehensive analysis is specifically as follows: preprocessing and standardizing the structural parameters and the physiological parameters, wherein the standardization is implemented using the maximum and minimum values, and the conversion is specifically performed according to the following formula: ,in, is the original parameter value, and are the minimum and maximum values ​​of each of the structural parameter and the physiological parameter of the target crop, respectively, is the standardized parameter value, and its value range is .

[0029] Based on the above scheme, the quantitative crop growth assessment index is constructed as follows:

[0030] The standardized structural parameters and physiological parameters are subjected to dimensionality reduction and integration through a principal component analysis model to extract principal component variables reflecting crop growth, wherein the principal component analysis model converts the original variables into linearly independent new variables, namely, principal component variables, through orthogonal transformation;

[0031] The principal component score matrix F is calculated by the principal component analysis model, and the first k columns of the principal component score matrix F are used as the crop growth assessment index.

[0032] Based on the above scheme, the process of extracting the principal component variables by the principal component analysis model includes:

[0033] Features Crop parameter indicators, each indicator has observations, forming The original data matrix X = [ x 1 , x 2 , ..., x n ] T ;

[0034] For the original data matrix Centering is done by subtracting the mean of each row to obtain a centralized matrix ; Calculate the centralization matrix The covariance matrix of : ;

[0035] For the covariance matrix Perform eigenvalue decomposition: ,in, For the eigenvalues, is the corresponding eigenvector; sort the eigenvalues ​​from large to small: ; Take the front The eigenvector corresponding to the largest eigenvalue , forming the eigenvector matrix P =[ v 1 , v 2 , ..., v k ] ;

[0036] The centralization matrix Projection onto the eigenvector matrix On the above, we get the principal component score matrix : ; The principal component score matrix Each row of is the principal component score of an observed sample, and each column is a principal component. The principal component score reflects the position of the original variable in the new coordinate system;

[0037] Calculate the variance contribution of each principal component : ; The variance contribution rate Indicates The principal components can explain the covariance matrix proportion;

[0038] Before The principal component with a cumulative variance contribution rate greater than 80% is used as the principal component variable.

[0039] Based on the above scheme, based on the principal component analysis model, the crop growth model is introduced to dynamically evaluate the growth of the target crop, and the crop growth model is a mathematical model based on crop physiological and ecological processes that quantitatively describes crop growth and development and yield formation, including but not limited to WOFOST;

[0040] The crop growth model includes the following modules: a phenological development module, a photosynthetic production module, a material distribution module, a soil moisture module and a soil nutrient module.

[0041] Based on the above solution, the data preprocessing specifically includes:

[0042] For the laser radar point cloud data, a statistical filtering algorithm is used to remove noise points and outliers, and an octree-based spatial indexing method is used to perform spatial segmentation into point cloud units, the point cloud units are classified into ground points and non-ground points, ground points are extracted to generate a digital terrain model, and the non-ground point data is extracted as the canopy information of the target crop, wherein the first point cloud data includes the digital terrain model and the canopy information;

[0043] For the multispectral satellite image data, the DN value in the multispectral satellite image data is converted into apparent reflectance through satellite radiation calibration parameters to eliminate the influence of sensor response differences, and based on the 6S atmospheric radiation transfer model and the atmospheric parameters when the satellite passes, the multispectral satellite image data is atmospherically corrected to eliminate the influence of atmospheric scattering and absorption, and the multispectral satellite image data after processing is used as the first image data.

[0044] Based on the above scheme, the data registration adopts an image registration algorithm based on mutual information, specifically: converting the first point cloud data into a point cloud digital surface model in a raster format, and having the same resolution and data range as the first image data; extracting scale-invariant feature points from the point cloud digital surface model and the first image data respectively through a scale-invariant feature transformation algorithm, and calculating the initial transformation matrix between the point cloud digital surface model and the first image data through matching the scale-invariant feature points; based on the initial transformation matrix, using the Powell optimization algorithm to minimize the mutual information measure between the point cloud digital surface model and the first image data, to obtain a registration transformation matrix, and based on the registration transformation matrix, resampling the point cloud digital surface model to a coordinate system consistent with the first image data to complete the data registration.

[0045] Compared with the prior art, the present invention has the following advantages and effects:

[0046] (1) Based on the complementary characteristics of LiDAR and multispectral imaging, the target crop canopy is characterized from two dimensions: three-dimensional structure and spectral characteristics. LiDAR can obtain the height, volume and other structural parameters of the target crop canopy with high precision, while multispectral imaging can invert the crop's leaf area index, chlorophyll content and other physiological parameters, thus achieving a more comprehensive and accurate data characterization of the target crop's growth.

[0047] (2) Based on the canopy structure parameters and physiological parameters of the target crops obtained from multiple sources, and through data mining methods such as principal component analysis, a quantitative comprehensive evaluation index for crop growth was constructed, which can quantitatively reflect the growth quality of the target crops from multiple dimensions and can more accurately evaluate the growth of the target crops;

[0048] (3) Design a crop growth model to simulate the dynamic changes of crops throughout their growth period. Compare the actual monitored values ​​of the target crop growth with the ideal values ​​simulated by the model to quantitatively determine the degree of deviation of the target crop growth. This can timely detect abnormal seedling conditions and analyze the causes, realize continuous monitoring and early warning of crop growth, and provide technical support for efficient and sustainable crop management. BRIEF DESCRIPTION OF THE DRAWINGS

[0049] In order to more clearly illustrate the specific embodiments of the present invention or the technical solutions in the prior art, the following is a brief introduction to the drawings required for the specific embodiments or the description of the prior art. In all the drawings, similar elements or parts are generally identified by similar reference numerals. In the drawings, the elements or parts are not necessarily drawn according to the actual scale.

[0050] Figure 1This is a flow chart of a crop monitoring method integrating laser radar and multispectral imaging provided by an embodiment of the present invention. DETAILED DESCRIPTION

[0051] In order to more clearly explain the purpose, technical solutions and advantages of the present invention, the technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, rather than all of the embodiments. The example implementation methods can be implemented in various forms and should not be understood as being limited to the examples described herein. On the contrary, these implementation methods are provided to make the present invention more comprehensive and complete, and to fully convey the concepts of the example implementation methods to those skilled in the art.

[0052] In addition, the described features, structures or characteristics may be combined in one or more embodiments in any suitable manner. In the following description, many specific details are provided to provide a full understanding of the embodiments of the present invention. However, those skilled in the art will appreciate that the technical solution of the present invention can be practiced without one or more of the specific details, or other methods, components, devices, steps, etc. may be adopted. In other cases, known methods, devices, implementations or operations are not shown or described in detail to avoid blurring various aspects of the present invention.

[0053] The block diagrams shown in the accompanying drawings are merely functional entities and do not necessarily correspond to physically independent entities. That is, these functional entities may be implemented in software form, or in one or more hardware modules or integrated circuits, or in different networks and / or processor devices and / or microcontroller devices.

[0054] The flowcharts shown in the accompanying drawings are only exemplary and do not necessarily include all the contents and operations / steps, nor must they be executed in the order described. For example, some operations / steps can be decomposed, and some operations / steps can be combined or partially combined, so the actual execution order may change according to actual conditions.

[0055] The present invention is described in detail below in conjunction with specific embodiments:

[0056] As attached Figure 1 As shown, an embodiment of the present invention provides a crop monitoring method integrating laser radar and multispectral imagery, and the specific steps of the method are as follows:

[0057] Step S1, extracting laser radar point cloud data and multispectral satellite image data of a target crop, and performing data preprocessing on the laser radar point cloud data and the multispectral satellite image data to obtain first point cloud data and first image data of the target crop respectively;

[0058] Step S2: unifying the first point cloud data and the first image data into the same coordinate system through data registration, and performing data fusion to obtain first fused data of the target crop;

[0059] Step S3, based on the first fusion data, extracting structural parameters of the canopy of the target crop, for describing the three-dimensional growth state of the target crop, wherein the structural parameters include but are not limited to canopy height, canopy coverage and leaf area index;

[0060] Step S4: based on the first fused data, estimating physiological parameters of the canopy of the target crop by using a radiation transfer model inversion method, so as to describe the biochemical condition of the target crop, wherein the physiological parameters include but are not limited to chlorophyll content and leaf nitrogen content;

[0061] Step S5, comprehensively analyzing the structural parameters and the physiological parameters, constructing a quantitative crop growth assessment index, simulating the ideal growth state of the target crop through a crop growth model, determining whether the growth of the target crop is abnormal, and generating an abnormal diagnosis result;

[0062] Step S6: generating growth warning information and management decision suggestions for the target crop based on the abnormal diagnosis results.

[0063] Preferably, in step S1, the data preprocessing specifically includes:

[0064] Step S101: for the laser radar point cloud data, use a statistical filtering algorithm to remove noise points and outliers, and use an octree-based spatial indexing method to perform spatial segmentation into point cloud units, classify the point cloud units into ground points and non-ground points, extract ground points to generate a digital terrain model, extract the non-ground point data as the canopy information of the target crop, and the first point cloud data includes the digital terrain model and the canopy information;

[0065] Specifically, the statistical filtering algorithm includes but is not limited to Gaussian filtering and median filtering, and appropriate filtering parameters are selected according to the characteristics of the laser radar point cloud data; the segmentation parameters of the octree need to be optimized according to the point cloud density and computer hardware performance. In this embodiment, the point cloud is divided into a number of cube units with a side length of 0.5 to 1 meter;

[0066] Specifically, a progressive TIN densification filtering algorithm is adopted to iteratively construct a triangulated network to approximate the ground and gradually eliminate non-ground points to extract the ground points. The key parameters of the progressive TIN densification filtering algorithm include but are not limited to an iteration distance threshold and an iteration angle threshold, which are set specifically according to different terrain conditions.

[0067] Step S102: For the multispectral satellite image data, the DN value in the multispectral satellite image data is converted into an apparent reflectance through satellite radiation calibration parameters to eliminate the influence of sensor response differences, and based on the 6S atmospheric radiation transfer model and the atmospheric parameters when the satellite passes, the multispectral satellite image data is atmospherically corrected to eliminate the influence of atmospheric scattering and absorption, and the multispectral satellite image data after processing is used as the first image data.

[0068] Specifically, the multispectral satellite image data is accompanied by calibration coefficients such as gain, offset and solar radiation intensity table. These parameters are used to convert DN value into radiation brightness. The conversion formula is: , where: L is the radiance; G is the gain; B is the offset or intercept; DN is the digital value of the image; after the radiance is calculated, it is converted into the apparent reflectance, and the conversion formula is: ,in To show reflectivity; is the radiance; is the distance between the Earth and the Sun, calculated based on the date of observation; is the solar radiation intensity, is the solar zenith angle.

[0069] Preferably, in step S2, the data registration adopts an image registration algorithm based on mutual information, specifically: converting the first point cloud data into a point cloud digital surface model in a raster format, and having the same resolution and data range as the first image data; extracting scale-invariant feature points from the point cloud digital surface model and the first image data respectively through a scale-invariant feature transformation algorithm, and calculating an initial transformation matrix between the point cloud digital surface model and the first image data through matching the scale-invariant feature points; based on the initial transformation matrix, using a Powell optimization algorithm to minimize the mutual information measure between the point cloud digital surface model and the first image data, obtaining a registration transformation matrix, and based on the registration transformation matrix, resampling the point cloud digital surface model to a coordinate system consistent with the first image data, thereby completing the data registration;

[0070] Specifically, commonly used feature point matching algorithms include Euclidean distance matching, RANSAC matching, etc. The mutual information metric is calculated using the Parzen window method, and the window size is generally set to 11x11 pixels.

[0071] Preferably, in step S2, the data fusion adopts the Gram-Schmidt functional transformation method to fuse the point cloud digital surface model after the data registration and the first image data in the feature space to obtain the first fused data, which is specifically as follows:

[0072] Step S201, resampling the point cloud digital surface model to a coordinate system consistent with the first image data, and performing standardization processing on the point cloud digital surface model and the first image data respectively to eliminate the influence of dimensional differences;

[0073] Step S202, orthogonalizing the feature vectors of the point cloud digital surface model and the first image data by using the Gram-Schmidt functional transformation method to eliminate the correlation between the features, wherein the first principal component direction is selected as the first orthogonal direction during the orthogonal transformation;

[0074] Step S203, performing weighted fusion on the feature vectors after the orthogonal transformation, and converting the feature vectors into the original data space by using the inverse Gram-Schmidt functional transformation method, to obtain the fused point cloud digital surface model and the first image data, that is, the first fused data;

[0075] Specifically, the fusion weight of the weighted fusion may be set according to the importance of the feature. For example, in the present embodiment, the weight of the height feature is set to be greater than the spectral feature.

[0076] Preferably, in step S3, extracting the structural parameters of the target crop canopy specifically includes:

[0077] Step S301: using a point cloud digital surface model construction algorithm to generate a high-resolution point cloud digital surface model in the first fused data, and using a difference method to extract crop canopy height information from the high-resolution point cloud digital surface model;

[0078] Specifically, in this embodiment, the point cloud digital surface model construction algorithm includes but is not limited to inverse distance weighted interpolation and Kriging interpolation, and it is necessary to select appropriate algorithms and parameters according to the point cloud density and terrain complexity;

[0079] Specifically, the crop canopy height information is a normalized surface height obtained by subtracting the ground point generated digital terrain model from the point cloud digital surface model in the first fused data; the ground point generated digital terrain model can be obtained by ground point interpolation, and the algorithm used in this embodiment includes but is not limited to TIN interpolation and ANUDEM interpolation method;

[0080] Step S302: extracting vegetation coverage information from the first fused data using a normalized difference vegetation index threshold method;

[0081] Step S303: spatially superimposing the crop canopy height information and the vegetation coverage information, and calculating the canopy coverage of each pixel in the first fused data based on the crop canopy height information and the vegetation coverage information;

[0082] Step S304: estimating the leaf area index using an empirical regression model based on the crop canopy height information and the normalized difference vegetation index information;

[0083] Specifically, the normalized difference vegetation index is calculated based on the reflectance of the red light and near-infrared bands, and can effectively distinguish between vegetation and non-vegetation areas; the selection of the normalized difference vegetation index threshold needs to be tested and adjusted according to the vegetation type and growth conditions of the study area, and in this embodiment, it is preferably between 0.2 and 0.5.

[0084] It should be noted that by extracting the structural parameters such as canopy height, canopy coverage and leaf area index, the three-dimensional growth state of the target crop can be quantitatively evaluated, which can better reflect the spatial structure, growth trend and health status of the target crop, and provide a scientific basis for crop growth monitoring and management. Specifically, the canopy height characterizes the vertical growth state of the crop and reflects the overall growth level and growth cycle of the crop; the canopy coverage describes the horizontal distribution characteristics of the crop and is used to evaluate the population density and light interception ability; the leaf area index is an important indicator for measuring the growth state of crop leaves, reflecting the photosynthesis potential and biomass accumulation; by comprehensively analyzing these parameters, the three-dimensional growth state of the target crop can be comprehensively characterized, including growth amount, coverage rate and health status, thereby achieving precise growth assessment and management.

[0085] Preferably, in step S4, estimating the physiological parameters of the target crop canopy specifically includes:

[0086] Step S401, based on the PROSAIL canopy reflectance model, simulate the radiation transmission process of the target crop canopy under different parameter conditions, set the input parameter range and step size of the PROSAIL canopy reflectance model according to prior knowledge, obtain a series of parameter-canopy reflectance pairs through forward simulation, and construct the parameter-canopy reflectance lookup table;

[0087] Step S402, implementing physiological parameter inversion through the lookup table, extracting the canopy reflectance in the first fused data, searching for a matching spectrum in the lookup table, obtaining a parameter inversion value corresponding to the canopy reflectance, completing single pixel parameter inversion, performing the physiological parameter inversion on the full image element in the first fused data, and obtaining a physiological parameter inversion map;

[0088] Specifically, the inversion map of physiological parameters including but not limited to chlorophyll content, leaf nitrogen content, etc., provides support for characterizing the biochemical condition of the target crop canopy; after obtaining the inversion map of physiological parameters, median filtering is used to eliminate outliers, and Kriging interpolation is used to improve the spatial continuity and smoothness of the inversion result.

[0089] Preferably, in step S5, the comprehensive analysis specifically includes: preprocessing and standardizing the structural parameters and the physiological parameters, wherein the standardization is implemented using the maximum and minimum values, and the conversion is specifically performed according to the following formula: ,in, is the original parameter value, and are the minimum and maximum values ​​of each of the structural parameter and the physiological parameter of the target crop, respectively, is the standardized parameter value, and its value range is .

[0090] Preferably, in step S5, constructing a quantitative crop growth assessment index includes:

[0091] Step S501: using a principal component analysis model, the standardized structural parameters and physiological parameters are subjected to dimensionality reduction and integration, and principal component variables reflecting crop growth are extracted. The principal component analysis model converts the original variables into linearly independent new variables, i.e., principal component variables, through orthogonal transformation, so that the first few principal components can retain as much variance information of the original data as possible;

[0092] Step S502, calculating the principal component score matrix F through the principal component analysis model, and using the first k columns of the principal component score matrix F as the crop growth assessment index;

[0093] Preferably, the principal component analysis model in this embodiment is as follows:

[0094] Features Crop parameter indicators, each indicator has observations, forming The original data matrix X = [ x 1 , x 2 , ..., x n ] T ;

[0095] The process of extracting the principal component variables by the principal component analysis model includes: Centering is done by subtracting the mean of each row to obtain a centralized matrix ; Calculate the centralization matrix The covariance matrix of : ; For the covariance matrix Perform eigenvalue decomposition: ,in, For the eigenvalues, is the corresponding eigenvector; sort the eigenvalues ​​from large to small: ; Take the front The eigenvector corresponding to the largest eigenvalue , forming the eigenvector matrix P =[ v 1 , v 2 , ..., v k ] ; The central matrix Projection onto the eigenvector matrix On the above, we get the principal component score matrix : ; The principal component score matrix Each row is the principal component score of an observed sample, and each column is a principal component. The principal component score reflects the position of the original variable in the new coordinate system; calculate the variance contribution rate of each principal component : ; The variance contribution rate Indicates The principal components can explain the covariance matrix ratio; take the previous The principal component with a cumulative variance contribution rate greater than 80% is used as the principal component variable;

[0096] Specifically, when calculating the principal component score, the principal component score coefficient matrix needs to be The calculation of the variance contribution rate of the different principal components Different, the weight of each principal component needs to be adjusted, and the principal component score coefficient matrix The calculation formula is: ; Wherein, P is the eigenvector matrix, is a diagonal matrix consisting of the square roots of the reciprocals of the eigenvalues;

[0097] The original data matrix With the principal component score coefficient matrix Multiply them together to get the principal component score matrix : ; The principal component score matrix Each column is a weighted value of the principal component score, and the number of rows is equal to the number of original observation samples, that is, the observation value , the principal component score matrix Before The columns are used as the crop growth assessment indicators.

[0098] Preferably, in step S5, based on the principal component analysis model, the crop growth model is introduced to dynamically evaluate the growth of the target crop. The crop growth model is a mathematical model that quantitatively describes crop growth and development and yield formation based on crop physiological and ecological processes, including but not limited to WOFOST. The crop growth model includes the following modules:

[0099] a. Phenological development module: simulates the changes in crop phenological stages, such as the emergence, tillering, heading, flowering, maturity and other key phenological periods and their duration; crop phenological development is mainly affected by environmental factors such as temperature and photoperiod;

[0100] Specifically, the phenological development module generates the following data of the target crop:

[0101] Developmental stage index : Used to characterize the progress of the crop’s current growth stage;

[0102] Status information of key phenological stages (such as germination, heading, flowering, maturity, etc.);

[0103] Cumulative heat: Effective accumulated temperature formula: ,in Indicates the cumulative heat from sowing to the current date, Indicates The average daily temperature, It represents the base temperature, below which crops will not grow;

[0104] Specifically, the development stage index formula of the target crop is: ,in Indicates the total accumulated temperature required to enter the current phenological stage.

[0105] b. Photosynthetic production module: simulates the photosynthesis process of crops and calculates the canopy's interception and conversion efficiency of light energy; photosynthetic production is mainly affected by factors such as light intensity, CO2 concentration, and leaf area index;

[0106] Specifically, the photosynthetic production module generates the following data of the target crop:

[0107] Total photosynthetic products : dry matter mass produced by photosynthesis;

[0108] Specifically, the photosynthetic product formula of the target crop is: ;in Indicates the light energy utilization efficiency (unit: ), It represents the photosynthetically active radiation intercepted by the canopy, and the formula is: , Represents photosynthetically active radiation (unit: ), represents the light interception coefficient; Represents leaf area index.

[0109] c. Material distribution module: simulates the distribution ratio of photosynthetic products among various organs of the plant, including but not limited to leaves, stems, roots, fruits, etc.; material distribution is regulated by crop genetic characteristics and growth stages;

[0110] Specifically, the data generated by the substance distribution module is as follows:

[0111] Organ dry matter distribution : Distribution of photosynthetic products among different organs (leaves, stems, roots, fruits);

[0112] Specifically, the dry matter distribution formula of the target crop is: ,in Indicates The dry matter distribution of each organ, Indicates The distribution coefficient of each organ varies with the developmental stage; for example, during the heading stage, the distribution coefficient of the fruit may be higher, while during the tillering stage, the distribution coefficient of the roots and leaves is higher.

[0113] d. Soil moisture module: simulates the movement and balance of soil moisture and calculates the amount of soil water available to crops; soil moisture dynamics are affected by factors such as precipitation, irrigation, evaporation, and plant transpiration;

[0114] Specifically, the data generated by the soil moisture module is as follows:

[0115] Water stress factor : Used to measure the degree to which soil moisture limits crop growth;

[0116] Specifically, the water balance formula of the target crop is: ,in represents the change in soil moisture, Indicates the amount of precipitation, Indicates the irrigation amount, It represents evapotranspiration, i.e. evaporation and transpiration. Indicates the amount of deep leakage;

[0117] Specifically, the water stress factor formula is: ,in Indicates the amount of water in the current soil available for the target crop (unit: mm). Indicates the theoretical maximum water requirement of the target crop under current meteorological conditions (unit: mm). Used to measure whether soil moisture meets the water requirements of target crops.

[0118] e. Soil nutrient module: simulates the transformation and circulation process of soil nutrients including but not limited to nitrogen, phosphorus, potassium, etc., and calculates the nutrient supply that can be absorbed by crops; soil nutrient dynamics are affected by factors such as fertilization, leaching, and absorption;

[0119] Specifically, the data generated by the soil nutrient module is as follows:

[0120] Nutrient stress factors : Used to measure the degree to which nutrient supply limits crop growth;

[0121] Specifically, the nutrient absorption formula of the target crop is: ,in Indicates the amount of nutrients absorbed by the target crop, It indicates the amount of nutrients available in the soil. Indicates the root absorption capacity and soil supply capacity;

[0122] Specifically, the nutrient stress factor formula is: ,in Indicates the nutrients in the soil that can be absorbed by the target crop (unit: kg / ha), Indicates the theoretical nutrient requirement of the target crop at the current growth stage (unit: kg / ha). Used to measure whether soil nutrients meet the growth requirements of target crops.

[0123] It should be noted that, through a series of mathematical equations in the crop growth model, the interaction and dynamic changes of environmental factors, physiological processes, agronomic management, etc. during the growth and development of crops are quantitatively described, and the growth state and yield level of crops under specific conditions are predicted; based on the first fusion data obtained in this embodiment, the initial parameters and process parameters of the crop growth model can be locally calibrated to improve the regional applicability of the crop growth model. The calibrated crop growth model can be used as an ideal reference for crop growth, that is, the potential growth state of the target crop under the current or specific variety, weather, soil, and management conditions;

[0124] Preferably, in step S5, the determination of whether the growth of the target crop is abnormal is to compare the monitored crop growth evaluation index with the ideal crop growth index simulated by the crop growth model, quantitatively determine whether the growth of the target crop is abnormal, and estimate the abnormal growth level of the target crop;

[0125] Preferably, the simulation result of the target crop is calculated by combining the data generated by the above modules of the principal component analysis model with the following formula: , including:

[0126] Integrating the stress factors of water and nutrients, the calculation formula is: ,in is a comprehensive stress factor, reflecting the minimum limitation of water and nutrients;

[0127] Dry matter accumulation: ,in represents the total dry matter mass after comprehensive stress effects, Indicates the dry matter allocated to each organ;

[0128] The final simulation results need to select different output indicators according to the specific research objectives. ,For example:

[0129] formula Used to predict the yield of the target crop;

[0130] formula Used to estimate the total dry matter mass of the target crop.

[0131] For example, the crop growth assessment index of the target crop in the i-th monitoring period is assumed to be , calculated by the principal component analysis model, represents the actual growth state of the target crop; the simulation result of the crop growth model of the target crop in the same period is , represents the ideal growth state of the target crop, and the relative difference between the actual growth state and the ideal growth state It is expressed as: ,in, It is used to indicate the abnormal degree of growth of the target crop in the i-th monitoring period. A positive value indicates that the growth is lower than the ideal state, and a negative value indicates that the growth is better than the ideal state.

[0132] Specifically, referring to the general law of crop growth and combining the knowledge of local agricultural experts, in this embodiment, the warning threshold of the abnormal growth level is set; preferably: When When the growth abnormality is detected, it is determined to be a severe abnormality; on the basis of determining the growth abnormality level, further analyzing the dominant factors causing the growth abnormality, considering factors including but not limited to water deficiency, nitrogen deficiency, pests and diseases, etc., to provide a reference for formulating crop management measures; in this embodiment, the dominant factor is preferably determined by establishing a correlation model between growth abnormality and environmental factors and management factors through machine learning methods such as random forests, and calculating the importance score of each factor; based on the growth abnormality level and the dominant factor, corresponding management decision suggestions are given.

[0133] In this embodiment, by extracting lidar point cloud data and multispectral satellite image data of the target crop and performing data preprocessing, the data are unified into the same coordinate system for data fusion, the structural parameters of the target crop canopy are extracted, the physiological parameters of the target crop canopy are estimated, the structural parameters and the physiological parameters are comprehensively analyzed, and quantitative crop growth assessment indicators are constructed. The ideal growth state of the target crop is simulated through a crop growth model, it is determined whether the growth of the target crop is abnormal, and an abnormal diagnosis result is generated. Finally, the growth warning information and management decision suggestions of the target crop are obtained, which provides a reference for formulating crop management measures and can significantly improve the accuracy, timeliness and quantification level of crop growth monitoring.

[0134] Those skilled in the art will readily appreciate other embodiments of the present invention after considering the specification and practicing the invention disclosed herein. This application is intended to cover any modification, use or adaptation of the present invention, which follows the general principles of the present invention and includes common knowledge or customary techniques in the art that are not disclosed by the present invention. The specification and examples are intended to be exemplary only, and the true scope and spirit of the present invention are indicated by the claims. It should be understood that the present invention is not limited to the precise structure described above and shown in the accompanying drawings, and various modifications and changes can be made without departing from the scope thereof. The scope of the present invention is limited only by the appended claims.

Claims

1. A crop monitoring method integrating laser radar and multispectral imaging, characterized in that: include: Extracting laser radar point cloud data and multispectral satellite image data of the target crop, and performing data preprocessing on the laser radar point cloud data and the multispectral satellite image data to obtain first point cloud data and first image data of the target crop respectively; Unifying the first point cloud data and the first image data into the same coordinate system through data registration, and performing data fusion to obtain first fused data of the target crop; Based on the first fusion data, extracting structural parameters of the canopy of the target crop to describe the three-dimensional growth state of the target crop, the structural parameters including canopy height, canopy coverage and leaf area index; Based on the first fused data, the physiological parameters of the canopy of the target crop are estimated by using a radiation transfer model inversion method to characterize the biochemical condition of the target crop, wherein the physiological parameters include chlorophyll content and leaf nitrogen content; Comprehensively analyzing the structural parameters and the physiological parameters, constructing a quantitative crop growth assessment index, simulating the ideal growth state of the target crop through a crop growth model, determining whether the growth of the target crop is abnormal, and generating an abnormal diagnosis result; Generating growth warning information and management decision suggestions for the target crop according to the abnormal diagnosis results; The step of judging whether the growth of the target crop is abnormal includes: assuming that the crop growth evaluation index of the target crop in the i-th monitoring period is , calculated by the principal component analysis model, represents the actual growth state of the target crop; the simulation result of the crop growth model of the target crop in the same period is , represents the ideal growth state of the target crop, and the relative difference between the actual growth state and the ideal growth state It is expressed as: ,in, It is used to indicate the abnormal degree of growth of the target crop in the i-th monitoring period, where a positive value indicates that the growth is lower than the ideal state, and a negative value indicates that the growth is better than the ideal state; The quantitative crop growth assessment index is constructed as follows: The standardized structural parameters and physiological parameters are subjected to dimensionality reduction and integration through a principal component analysis model to extract principal component variables reflecting crop growth, wherein the principal component analysis model converts the original variables into linearly independent new variables, namely, principal component variables, through orthogonal transformation; Calculate the principal component score matrix F by the principal component analysis model, and use the first k columns of the principal component score matrix F as the crop growth assessment index; The process of extracting the principal component variables by the principal component analysis model includes: Features Crop parameter indicators, each indicator has observations, forming The original data matrix ; For the original data matrix Centering is done by subtracting the mean of each row to obtain a centralized matrix ; Calculate the centralization matrix The covariance matrix of : ; For the covariance matrix Perform eigenvalue decomposition: ,in, For the eigenvalues, is the corresponding eigenvector; sort the eigenvalues ​​from large to small: ; Take the front The eigenvector corresponding to the largest eigenvalue , forming the eigenvector matrix ; The centralization matrix Projection onto the eigenvector matrix On the above, we get the principal component score matrix : ; The principal component score matrix Each row of is the principal component score of an observed sample, and each column is a principal component. The principal component score reflects the position of the original variable in the new coordinate system; Calculate the variance contribution of each principal component : ; The variance contribution rate Indicates The principal components can explain the covariance matrix proportion; Before The principal component with a cumulative variance contribution rate greater than 80% is used as the principal component variable.

2. The crop monitoring method integrating laser radar and multispectral imagery according to claim 1, characterized in that: The data fusion adopts the Gram-Schmidt functional transformation method to fuse the point cloud digital surface model after the data registration and the first image data in the feature space to obtain the first fused data, specifically including: Resampling the point cloud digital surface model to a coordinate system consistent with the first image data, and performing standardization processing on the point cloud digital surface model and the first image data respectively; Orthogonalizing the feature vectors of the point cloud digital surface model and the first image data by using a Gram-Schmidt functional transformation method, wherein the first principal component direction is selected as the first orthogonal direction during the orthogonal transformation; The orthogonally transformed feature vectors are weighted fused, and the feature vectors are converted into the original data space through the inverse Gram-Schmidt functional transformation method to obtain the fused point cloud digital surface model and the first image data, that is, the first fused data.

3. The crop monitoring method integrating laser radar and multispectral imagery according to claim 1, characterized in that: The step of extracting the structural parameters of the target crop canopy specifically includes: Using a point cloud digital surface model construction algorithm, a high-resolution point cloud digital surface model is generated from the first fused data, and a difference method is used to extract crop canopy height information from the high-resolution point cloud digital surface model; extracting vegetation coverage information from the first fused data using a normalized difference vegetation index threshold method; spatially superimposing the crop canopy height information and the vegetation coverage information, calculating the canopy coverage of each pixel using the height information and coverage information of the crop canopy height information, and calculating the canopy coverage of each pixel in the first fused data; Based on the crop canopy height information and the normalized difference vegetation index information, an empirical regression model is used to estimate the leaf area index.

4. The crop monitoring method integrating laser radar and multispectral imagery according to claim 1, characterized in that: The estimating the physiological parameters of the target crop canopy specifically includes: Based on the PROSAIL canopy reflectance model, the radiation transfer process of the target crop canopy under different parameter conditions is simulated, the input parameter range and step size of the PROSAIL canopy reflectance model are set according to prior knowledge, a series of parameter-canopy reflectance pairs are obtained through forward simulation, and a lookup table of the parameter-canopy reflectance is constructed; The physiological parameter inversion is realized through the lookup table, the canopy reflectance in the first fused data is extracted, and the matching spectrum is searched in the lookup table to obtain the parameter inversion value corresponding to the canopy reflectance, and the single pixel parameter inversion is completed. The physiological parameter inversion is performed on the full image element in the first fused data to obtain a physiological parameter inversion map.

5. The crop monitoring method integrating laser radar and multispectral imagery according to claim 1, characterized in that: The comprehensive analysis specifically includes: preprocessing and standardizing the structural parameters and the physiological parameters, wherein the standardization is implemented using the maximum and minimum values, and the conversion is specifically performed according to the following formula: ,in, is the original parameter value, and are the minimum and maximum values ​​of each of the structural parameter and the physiological parameter of the target crop, respectively, is the standardized parameter value, and its value range is .

6. The crop monitoring method integrating laser radar and multispectral imagery according to claim 1, characterized in that: Based on the principal component analysis model, the crop growth model is introduced to dynamically evaluate the growth of the target crop, and the crop growth model is a mathematical model based on crop physiological and ecological processes that quantitatively describes crop growth and development and yield formation, including WOFOST; The crop growth model includes the following modules: a phenological development module, a photosynthetic production module, a material distribution module, a soil moisture module and a soil nutrient module.

7. The crop monitoring method integrating laser radar and multispectral imagery according to claim 1, characterized in that: The data preprocessing specifically includes: For the laser radar point cloud data, a statistical filtering algorithm is used to remove noise points and outliers, and an octree-based spatial indexing method is used to perform spatial segmentation into point cloud units, the point cloud units are classified into ground points and non-ground points, ground points are extracted to generate a digital terrain model, and the non-ground point data is extracted as the canopy information of the target crop, wherein the first point cloud data includes the digital terrain model and the canopy information; For the multispectral satellite image data, the DN value in the multispectral satellite image data is converted into apparent reflectance through satellite radiation calibration parameters to eliminate the influence of sensor response differences, and based on the 6S atmospheric radiation transfer model and the atmospheric parameters when the satellite passes, the multispectral satellite image data is atmospherically corrected to eliminate the influence of atmospheric scattering and absorption, and the multispectral satellite image data after processing is used as the first image data.

8. The crop monitoring method integrating laser radar and multispectral imagery according to claim 1, characterized in that: The data registration adopts an image registration algorithm based on mutual information, specifically: converting the first point cloud data into a point cloud digital surface model in a raster format, and having the same resolution and data range as the first image data; extracting scale-invariant feature points from the point cloud digital surface model and the first image data respectively through a scale-invariant feature transformation algorithm, and calculating the initial transformation matrix between the point cloud digital surface model and the first image data through matching the scale-invariant feature points; based on the initial transformation matrix, using the Powell optimization algorithm to minimize the mutual information measure between the point cloud digital surface model and the first image data, to obtain a registration transformation matrix, and based on the registration transformation matrix, resampling the point cloud digital surface model to a coordinate system consistent with the first image data to complete the data registration.

Images