Crop water and fertilizer stress early warning method based on energy residual error and unmanned aerial vehicle thermal infrared

Abstract

The crop water and fertilizer stress early warning method based on energy residual error and unmanned aerial vehicle thermal infrared belongs to the technical field of precision agriculture and crop stress remote sensing monitoring.The present method cannot identify early stress, cannot distinguish stress types and has unstable prediction performance.The present application solves the problems.The theoretical canopy temperature of each observation time is calculated based on the simplified SEBAL energy balance equation, and the energy residual deviation is obtained by comparing with the measured thermal infrared temperature, and the multi-dimensional diagnostic index is extracted to amplify the weak signal of early stress, the gradient amplitude analysis of the energy residual error index sequence in the time dimension is carried out by using the time sequence gradient stress response decoupling module, the deviation response is decomposed into the fast response water stress channel and the slow response nutrient stress channel by the adaptive gradient threshold, and the double-channel segmented features are formed by segment aggregation, and then the stress identification is carried out by inputting the water and fertilizer stress residual double-channel early warning network.The present application can be applied to crop water and fertilizer stress early warning.

Description

Technical Field

[0001] This invention belongs to the field of precision agriculture and remote sensing monitoring technology of crop stress, specifically involving a crop water and fertilizer stress early warning method based on energy surplus and UAV thermal infrared. Background Technology

[0002] Crops are frequently subjected to water and nutrient stress during their growth. Failure to identify stress types early and implement appropriate irrigation or fertilization measures can lead to significant yield reductions and quality deterioration. UAVs equipped with thermal infrared sensors can quickly acquire large-area crop canopy temperature distribution information, providing a convenient means for field-scale crop stress monitoring. However, thermal infrared signals during the early stages of crop stress are extremely weak: in early water shortages, plants preferentially regulate stomatal conductance while canopy temperature has not yet significantly increased; in early nutrient deficiency, photosynthetic rate decreases but the thermal effect is minimal. More importantly, water and nutrient stress exhibit strong coupling in thermal infrared representations—water shortage affects nutrient absorption, and nutrient deficiency affects stomatal behavior—making it almost impossible to directly distinguish the thermal infrared signals of the two types of stress in a single observation. Therefore, researching methods to separate the coupling effects of water and nutrient stress from thermal infrared remote sensing signals under complex environmental disturbances and achieve early warning has significant scientific value and engineering application implications.

[0003] Using drone-based thermal infrared remote sensing for early warning and identification of crop water and fertilizer stress, the commonly used technologies in the industry currently include the following categories:

[0004] I. Stress Identification Method Based on Canopy Temperature Absolute Threshold: This method determines whether crops are under stress by setting a fixed temperature threshold. While simple in principle and easy to implement, canopy temperature is strongly influenced by environmental factors such as solar radiation intensity, air temperature, wind speed, and humidity. The fixed threshold cannot adapt to normal temperature fluctuations under different meteorological backgrounds. In the early stages of stress, when canopy temperature deviates from the normal range by only a few tenths of a degree Celsius, environmental noise completely drowns out the stress signal, leading to a very high false negative rate. Furthermore, this method can only determine the presence of stress but cannot distinguish between water stress and nutrient stress.

[0005] II. A single-phase method based on the Crop Water Stress Index (CWSI) calculates a normalized stress index by constructing dry and wet reference surfaces. This method has some diagnostic capability under moderate to severe water stress conditions. However, the determination of the dry and wet reference surfaces depends on synchronously acquired high-precision meteorological data and the subjectivity of reference surface selection. In the early, mild stress stage, the CWSI change is within the range of instrument measurement noise and is difficult to detect reliably. Furthermore, CWSI essentially only reflects the degree of evapotranspiration deviation and lacks the ability to distinguish non-evapotranspiration thermal responses caused by nutrient stress, thus failing to identify the types of water stress and nutrient stress.

[0006] III. Indirect Inference Method Based on Temporal Changes in Vegetation Indices: This method uses the temporal decline trend of visible and near-infrared vegetation indices such as NDVI to infer the occurrence of stress. This method reflects changes in chlorophyll content and canopy structure, and exhibits a significant lag compared to the actual occurrence of stress. Detectable changes in vegetation indices typically appear several days to weeks after the stress event, failing to meet the timeliness requirements for early warning. Furthermore, changes in vegetation indices cannot directly distinguish between water stress and nutrient stress, both of which lead to chlorophyll degradation and canopy structure deterioration.

[0007] IV. A deep learning-based end-to-end classification method for thermal infrared images utilizes convolutional neural networks to directly learn a stress classification model from thermal infrared images. This method has a strong learning ability for precisely labeled and sufficient samples in the training data. However, obtaining precisely labeled early stress training samples under actual field conditions is extremely difficult because ground verification of early stress requires destructive sampling or long-term tracking confirmation. More importantly, existing deep models directly use the raw thermal infrared temperature as input and fail to use energy balance physical priors to separate the contribution of environmental factors to canopy temperature. This results in the features actually learned by the model being mixed with a large amount of meteorological noise rather than the stress response itself, leading to unstable prediction performance under different time periods and meteorological conditions.

[0008] In summary, existing UAV-based thermal infrared diagnostic technology for crop stress faces three core challenges: early stress thermal infrared signals are extremely weak and easily masked by environmental and meteorological noise, and existing methods cannot effectively amplify early deviation signals; water stress and nutrient stress exhibit a strong coupling effect in thermal infrared radiation, and existing feature extraction methods cannot separate their temporal response contributions; and deep learning-based methods show unstable predictive performance under different time periods and meteorological conditions. Therefore, it is urgent to propose a new early warning method for crop water and fertilizer stress to address these issues. Summary of the Invention

[0009] The purpose of this invention is to address the problems of existing methods being unable to identify early stress, unable to distinguish between water stress and nutrient stress, and having unstable predictive performance under different meteorological conditions. Therefore, this invention proposes a crop water and fertilizer stress early warning method based on energy surplus and UAV thermal infrared.

[0010] The technical solution adopted by this invention to solve the above-mentioned technical problems is: a crop water and fertilizer stress early warning method based on energy surplus and UAV thermal infrared, the method specifically includes the following steps:

[0011] Step 1: Obtain UAV thermal infrared time-series observation data YJNewData for the area to be monitored. The time-series observation data YJNewData includes observation records at TXNum time points, and each observation record at each time point includes the following fields:

[0012] HWGC field: Thermal infrared observation of canopy temperature;

[0013] KQWD field: Records the ambient temperature;

[0014] XDSD field: Records the relative humidity of the environment;

[0015] FSSL field: Records the ambient wind speed;

[0016] TRFS field: Records total solar radiation;

[0017] ZBZS field: The obtained Normalized Difference Vegetation Index (NDVI) value;

[0018] Step 2: Use the temporal gradient stress response decoupling module SXTJModel to extract the water stress channel segmentation features SFTrack and nutrient stress channel segmentation features YFTrack from the observed data;

[0019] Step 3: Establish the water and fertilizer stress residual difference dual-channel early warning network SFYCNet. The water and fertilizer stress residual difference dual-channel early warning network SFYCNet includes a water stress channel branch, a nutrient stress channel branch, and an MLP; wherein: the water stress channel branch and the nutrient stress channel branch both include an LSTM model;

[0020] The water stress channel segmentation feature SFTrack and the nutrient stress channel segmentation feature YFTrack are used as inputs to the water stress channel branch and the nutrient stress channel branch, respectively. The outputs of the water stress channel branch and the nutrient stress channel branch are then concatenated along the channel dimension. The concatenated result is then processed by an MLP, and the detection result of the area to be monitored is output through the MLP.

[0021] If water or nutrient stress exists, an early warning should be issued; otherwise, no action is necessary.

[0022] The beneficial effects of this invention are:

[0023] This invention constructs an energy residual deviation diagnostic index extraction module and a time-series gradient stress response decoupling module, and on this basis, establishes a dual-channel early warning network for water and fertilizer stress residual deviation. The energy residual deviation diagnostic index extraction module calculates the theoretical canopy temperature at each observation time based on the simplified SEBAL energy balance equation, and compares it with the measured thermal infrared temperature to obtain the energy residual deviation. From this, multi-dimensional diagnostic indicators including the original residual value, relative residual deviation, latent heat flux anomaly ratio, Bowen ratio deviation, and water vapor pressure difference sensitivity are extracted to amplify the weak early stress signals. The time-series gradient stress response decoupling module is then used to further analyze the energy residual deviation. The coupling module performs gradient magnitude analysis on the energy residual error index sequence at multiple time points in the time dimension. By using an adaptive gradient threshold, the deviation response is decomposed into a fast-response water stress channel and a slow-response nutrient stress channel, and the segments are aggregated to form dual-channel segmented features. Then, the dual-branch long short-term memory neural network model in the water and fertilizer stress residual error dual-channel early warning network is used to learn the temporal evolution patterns of the water channel and the nutrient channel respectively, and the coupling discrimination is completed in the fusion layer. This enables early warning and type identification of crop water stress and nutrient stress, and ensures the prediction performance under different meteorological conditions. Attached Figure Description

[0024] Figure 1 This is a flowchart of a crop water and fertilizer stress early warning method based on energy surplus and UAV thermal infrared radiation according to the present invention. Detailed Implementation

[0025] Specific implementation method one: Combining Figure 1 This embodiment describes a crop water and fertilizer stress early warning method based on energy surplus and UAV thermal infrared radiation. The method specifically includes the following steps:

[0026] Step 1: Obtain UAV thermal infrared time-series observation data YJNewData for the area to be monitored. The time-series observation data YJNewData includes observation records at TXNum time points, and each observation record at each time point includes the following fields:

[0027] HWGC field: Thermal infrared observation of canopy temperature, which is a temperature value in degrees Celsius;

[0028] The KQWD field records the ambient temperature, which is a temperature value in degrees Celsius.

[0029] The XDSD field records the ambient relative humidity, which is a humidity value in percentage form.

[0030] FSSL field: Records the ambient wind speed, which is a wind speed value in meters per second;

[0031] TRFS field: Records total solar radiation, a radiation intensity value in watts per square meter;

[0032] ZBZS field: The obtained Normalized Difference Vegetation Index (NDVI) value, which is a floating-point number between 0 and 1;

[0033] Step 2: Use the temporal gradient stress response decoupling module SXTJModel to extract the water stress channel segmentation features SFTrack and nutrient stress channel segmentation features YFTrack from the observed data;

[0034] Step 3: Establish the water and fertilizer stress residual difference dual-channel early warning network SFYCNet. The water and fertilizer stress residual difference dual-channel early warning network SFYCNet includes a water stress channel branch, a nutrient stress channel branch, and an MLP (Multilayer Perceptron). Among them, the water stress channel branch and the nutrient stress channel branch both include an LSTM (Long Short-Term Memory) neural network model.

[0035] The water stress channel segmentation feature SFTrack and the nutrient stress channel segmentation feature YFTrack are used as inputs to the water stress channel branch and the nutrient stress channel branch, respectively. The outputs of the water stress channel branch and the nutrient stress channel branch are then concatenated along the channel dimension. The concatenated result is then processed by an MLP, and the detection result of the area to be monitored is output through the MLP.

[0036] If water or nutrient stress exists, an early warning should be issued; otherwise, no action is necessary.

[0037] The energy deviation diagnostic index extraction module of this invention does not rely on a fixed temperature threshold or artificially set dry and wet reference surfaces. Instead, it uses the simplified SEBAL (Surface Energy Balance Algorithm for Land) energy balance physical equation to calculate the theoretical temperature benchmark of the crop canopy under the current meteorological conditions. The deviation between the measured temperature and the theoretical temperature is the pure stress response signal after removing environmental factors. It is not directly affected by fluctuations in solar radiation intensity, air temperature and wind speed. In the early stress stage, it amplifies the weak temperature anomaly of only a few tenths of a degree Celsius into an observable multidimensional diagnostic index deviation. It can accurately detect early stress signals without the need to simultaneously acquire dry and wet reference surfaces. The temporal gradient stress response decoupling module utilizes the physical fact that there is an essential difference in the gradient response rates of the two types of stress in the time dimension (i.e., water stress causes rapid stomatal adjustment, resulting in a rapid change in residual error in a short period of time, forming a fast response channel, while nutrient stress causes a gradual accumulation of residual error by slowly affecting the photosynthetic mechanism, forming a slow response channel). By using an adaptive gradient threshold, the two are separated into independent channels. Without the need for a large number of training samples with pre-labeled stress types, structured dual-channel feature separation can be achieved through physical priors, which significantly improves the accuracy of water and fertilizer stress type identification and the timeliness of early warning.

[0038] Specific Implementation Method Two: This implementation method is a further limitation of Specific Implementation Method One. The specific process of step two is as follows:

[0039] S101. The water stress channel segmentation feature SFTrack is a list of SegNum vectors, and each vector in the list is initialized as an NYCZDim-dimensional all-zero vector, where NYCZDim is the dimension of the energy residual error diagnostic index, and the default value of NYCZDim is 10.

[0040] The nutrient stress channel segmentation feature YFTrack is also a list of SegNum vectors, and each vector in the list is initialized as an NYCZDim dimensional all-zero vector;

[0041] Establish the residual index time series matrix SXTJFeaMat, and initialize the residual index time series matrix SXTJFeaMat as a matrix of all zeros with rows of TXNum and columns of NYCZDim;

[0042] Use the time-series observation data YJNewData as the input data SXTJInput of the time-series gradient stress response decoupling module SXTJModel, and initialize the decoupling time-series counter SXTJTCounter=1;

[0043] S102, Decoupling the timing, the first temporary variable SXTJTTemp1 = retrieves the SXTJTCounter observation record of SXTJInput, and uses the retrieved data as the input data NYCZInput of the energy residual deviation diagnostic index extraction module NYCZModel;

[0044] S103, the second temporary variable of the decoupling timing, SXTJTTemp2, is the result of processing the input data NYCZInput using the module NYCZModel, NYCZOutput.

[0045] Write SXTJTTemp2 into the SXTJTCounter row of the residual index time series matrix SXTJFeaMat;

[0046] S104. Set SXTJTCounter = SXTJTCounter + 1;

[0047] If SXTJTCounter is less than or equal to TXNum, proceed to S102; otherwise, proceed to S105.

[0048] S105. Calculate the time series gradient matrix SXTJGradMat based on the residual index time series matrix SXTJFeaMat.

[0049] S106. Calculate the gradient magnitude sequence SXTJGradMag based on the temporal gradient matrix SXTJGradMat, wherein the gradient magnitude sequence SXTJGradMag is a vector containing TXNum-1 elements;

[0050] S107. Calculate the adaptive gradient decoupling threshold SXTJThreshold:

[0051] Take the median of all elements in SXTJGradMag as SXTJThreshold. If SXTJThreshold is less than 0.001, set SXTJThreshold=0.001; otherwise, leave SXTJThreshold unchanged.

[0052] S108. Calculate the water channel feature matrix SXTJSFMat based on SXTJGradMag, SXTJThreshold, and SXTJFeaMat; calculate the nutrient channel feature matrix SXTJYFMat based on SXTJGradMag, SXTJThreshold, and SXTJFeaMat.

[0053] Then, based on the water channel feature matrix SXTJSFMat, the segmented features of the water stress channel SFTrack are obtained, and based on the nutrient channel feature matrix SXTJYFMat, the segmented features of the nutrient stress channel YFTrack are obtained.

[0054] The other steps and parameters are the same as in Specific Implementation Method 1.

[0055] Specific Implementation Method Three: This implementation method is a further limitation of Specific Implementation Method Two. The working process of the energy residual deviation diagnostic index extraction module NYCZModel is as follows:

[0056] S201. Initialize NYCZOutput as a zero-based vector with dimension NYCZDim;

[0057] S202. Obtain the first temporary variable NYCZTemp1 to the 18th temporary variable NYCZTemp18 based on the input data NYCZInput;

[0058] S203, set NYCZOutput[1]=NYCZTemp13 / 3.0, that is, use the original value of the residual difference to divide by the empirical temperature difference scale constant (3 degrees Celsius) for dimensional compression;

[0059] Among them, NYCZTemp13 is the 13th temporary variable, and NYCZOutput[1] is the first element in NYCZOutput;

[0060] S204. Set NYCZOutput[2]=NYCZTemp13 / (ABS(NYCZTemp12)+0.5), where the relative residual NYCZOutput[2] is the ratio of the residual of the original energy value to the absolute value of the theoretical canopy temperature;

[0061] Where ABS is the absolute value, NYCZTemp12 is the 12th temporary variable, and NYCZOutput[2] is the 2nd element in NYCZOutput;

[0062] S205, set NYCZOutput[3]=(1-NYCZTemp15 / (ABS(NYCZTemp10)+0.01))×0.8+0.1, NYCZOutput[3] is the latent heat flux anomaly ratio, the larger the value, the more the actual evaporation deviates from the reference value;

[0063] Among them, NYCZTemp15 is the 15th temporary variable, NYCZTemp10 is the 10th temporary variable, and NYCZOutput[3] is the 3rd element in NYCZOutput;

[0064] S206. Set NYCZOutput[4] to:

[0065]

[0066] Among them, NYCZTemp14 is the 14th temporary variable, NYCZTemp11 is the 11th temporary variable, and NYCZOutput[4] is the 4th element in NYCZOutput; NYCZOutput[4] is the Bowen ratio deviation term, which reflects the degree of deviation of the sensible heat and latent heat distribution relative to the unstressed reference.

[0067] S207, Set NYCZOutput[5]=NYCZTemp13×NYCZTemp18 / (QJYCBaseline×10+0.5); NYCZOutput[5] is the water vapor pressure difference sensitivity factor. The product of the residual error and the water vapor pressure difference (VPD) amplifies the deviation signal under arid conditions.

[0068] Among them, QJYCBaseline is the global residual deviation from the baseline, NYCZOutput[5] is the 5th element in NYCZOutput, and NYCZTemp18 is the 18th temporary variable;

[0069] S208, Set NYCZOutput[6]=NYCZTemp13×(1+0.2 / (NYCZTemp4+0.05));

[0070] NYCZOutput[6] is the wind speed correction residual. Under low wind speed conditions, the heat loss is reduced, which makes the residual amplified.

[0071] Among them, NYCZOutput[6] is the 6th element in NYCZOutput, and NYCZTemp4 is the 4th temporary variable;

[0072] S209, Set NYCZOutput[7]=NYCZTemp13 / (NYCZTemp7×0.001+0.3); NYCZOutput[7] is the radiation normalization residual, which eliminates the dimensional differences under different radiation intensities;

[0073] Among them, NYCZOutput[7] is the 7th element in NYCZOutput, and NYCZTemp7 is the 7th temporary variable. The radiation unit is converted to kilowatt by dividing by 1000.

[0074] S210, Set NYCZOutput[8]=NYCZTemp13×(0.4+0.6×NYCZTemp6)+0.05; NYCZOutput[8] is the vegetation-weighted residual difference, and the residual difference in high vegetation coverage areas is more diagnostically significant;

[0075] Among them, NYCZOutput[8] is the 8th element in NYCZOutput, and NYCZTemp6 is the 6th temporary variable;

[0076] S211. Set NYCZOutput[9]. NYCZOutput[9] is the global correction of the crown temperature difference gradient. Use the global mean difference to compensate for the systematic offset.

[0077]

[0078] Where QJHWGCMean is the global average canopy temperature, and QJKQWDMean is the global average air temperature;

[0079] S212. Set NYCZOutput

[10] , where NYCZOutput

[10] is the humidity correction residual. The Gaussian decay window suppresses observations far from the typical humidity range with 60% relative humidity as the center.

[0080]

[0081] Among them, NYCZOutput

[10] is the 10th element in NYCZOutput, NYCZTemp3 is the 3rd temporary variable, and QJFTXZCoeff is the wind temperature heat conduction correction coefficient;

[0082] If the value of NYCZOutput

[10] is greater than 2.0, then set NYCZOutput

[10] =2.0; otherwise, the value of NYCZOutput

[10] remains unchanged.

[0083] S213. Let the final output NYCZOutput = 1 / (1+exp(-NYCZOutput)), that is, perform Sigmoid compression on each element calculated from S203 to S212 to make the output within the range of 0 to 1, ensuring that the value is bounded.

[0084] S214. Use the final output NYCZOutput obtained in S213 as the result output of NYCZModel.

[0085] The other steps and parameters are the same as in Specific Implementation Method Two.

[0086] Specific Implementation Method Four: This implementation method is a further limitation of Specific Implementation Method Three. The specific process of S202 is as follows:

[0087] Step 1: Extract the HWGC field of NYCZInput as the first temporary variable for residual difference, NYCZTemp1;

[0088] Extract the KQWD field of NYCZInput as the second temporary variable for residual difference, NYCZTemp2;

[0089] Extract the XDSD field of NYCZInput and use it as the third temporary variable for the residual difference, NYCZTemp3;

[0090] Extract the FSSL field of NYCZInput as the fourth temporary variable of the residual, NYCZTemp4. If NYCZTemp4 is less than 0.1, set NYCZTemp4=0.1 to prevent the aerodynamic impedance from tending to infinity under zero wind speed conditions. If NYCZTemp4 is greater than or equal to 0.1, then NYCZTemp4 remains unchanged.

[0091] Extract the TRFS field of NYCZInput as the fifth temporary variable for residual difference, NYCZTemp5. If NYCZTemp5 is less than 1.0, set NYCZTemp5=1.0 to prevent division by zero under zero radiation conditions at night. If NYCZTemp5 is greater than or equal to 1.0, then NYCZTemp5 remains unchanged.

[0092] Extract the ZBZS field of NYCZInput as the 6th temporary variable of the residual NYCZTemp6;

[0093] Step 2: Calculate the 7th temporary variable NYCZTemp7 of the net radiation residual:

[0094] Where Albedo is the surface albedo, with a default value of 0.23; Emissivity is the surface emissivity, with a default value of 0.95; NYCZTemp7 represents the Stefan-Boltzmann constant, and NYCZTemp7 represents the simplified net radiation Rn.

[0095] Step 3: Calculate the 8th temporary variable NYCZTemp8 for the residual soil heat flux:

[0096]

[0097] Wherein, SHFRatio is the soil heat flux ratio coefficient, with a default value of 0.3; NYCZTemp8 is the soil heat flux G, and the higher the vegetation cover, the smaller the proportion of soil heat flux.

[0098] Step 4: Calculate the 9th temporary variable NYCZTemp9 for the aerodynamic impedance residual error:

[0099]

[0100] Where log is the natural logarithm function; NYCZTemp9 is the aerodynamic impedance ra; xtck is the drone hovering reference altitude (2.0 meters by default in this invention); zwccd is the crop surface roughness (0.05 micrometers by default in this invention); and 0.41 is the von Kármán constant.

[0101] Step 5: Calculate the 10th temporary variable NYCZTemp10 of the reference latent heat flux residual:

[0102]

[0103] Where PMKQ is the Penman reference evapotranspiration coefficient (calculated using the Penman method), and the default value of PMKQ is 0.65; NYCZTemp10 is the reference latent heat flux LE_ref under no-stress conditions;

[0104] Step 6: Calculate the 11th temporary variable NYCZTemp11 of the reference inductive heat flux residual:

[0105]

[0106] Where NYCZTemp11 is the reference induced heat flux H_ref under no-stress conditions;

[0107] Step 7: Calculate the 12th temporary variable, NYCZTemp12, of the theoretical canopy temperature residual:

[0108]

[0109] Wherein, NYCZTemp12 is the theoretical canopy temperature based on energy balance without stress, 1.225 is the air density (in kilograms per cubic meter), and 1004 is the specific heat capacity of air at constant pressure (in joules per kilogram per Kelvin).

[0110] Step 8: Calculate the residual of the original energy value, the 13th temporary variable NYCZTemp13:

[0111]

[0112] Among them, NYCZTemp13 is the core energy residual deviation. A positive value indicates that the measured temperature is higher than the theoretical temperature, which implies the presence of stress.

[0113] Step 9: Calculate the 14th temporary variable NYCZTemp14, representing the actual sensible heat flux residual:

[0114]

[0115] Step 10: Calculate the 15th temporary variable NYCZTemp15 for the actual latent heat flux residual:

[0116]

[0117] Step 11: Calculate the 16th temporary variable NYCZTemp16 of the saturated water vapor pressure residual difference:

[0118]

[0119] Where exp represents an exponential function with the natural constant e as its base;

[0120] Step 12: Calculate the 17th temporary variable NYCZTemp17 of the actual water vapor pressure difference:

[0121]

[0122] Step 13: Calculate the residual error of the vapor pressure difference (VPD) using the 18th temporary variable NYCZTemp18.

[0123] .

[0124] The other steps and parameters are the same as in Specific Implementation Method 3.

[0125] Specific Implementation Method Five: This implementation method further defines Specific Implementation Method Four. The global canopy temperature mean QJHWGCMean, global air temperature mean QJKQWDMean, global residual deviation from baseline QJYCBaseline, and wind temperature heat conduction correction coefficient QJFTXZCoeff are calculated based on the UAV thermal infrared time series training dataset WRJCDataSet and saved for training and subsequent actual monitoring. The calculation method for the global canopy temperature mean QJHWGCMean, global air temperature mean QJKQWDMean, global residual deviation from baseline QJYCBaseline, and wind temperature heat conduction correction coefficient QJFTXZCoeff is as follows:

[0126] S301. Obtain the UAV thermal infrared time series training dataset WRJCDataSet. The total number of training samples in the dataset WRJCDataSet is WRJCTZNum, and the number of observation times for each training sample is TXNum.

[0127] S302. Establish a global statistical accumulation variable QJLJTemp1=0 and a global statistical counting variable QJLJTemp2=0;

[0128] S303. Establish a global statistical sample counter QJYBCounter=1;

[0129] S304. Establish a global statistical time counter QJSKCounter=1;

[0130] S305. Set the global statistical temporary variable QJTemp1 to retrieve the HWGC field of the QJSKCounter observation record of the QJYBCounter sample in WRJCDataSet;

[0131] Let QJLJTemp1 = QJLJTemp1 + QJTemp1;

[0132] QJLJTemp2=QJLJTemp2+1;

[0133] S306. Set QJSKCounter = QJSKCounter + 1;

[0134] If QJSKCounter is less than or equal to TXNum, proceed to S305; otherwise, proceed to S307.

[0135] S307. Let QJYBCounter = QJYBCounter + 1;

[0136] If QJYBCounter is less than or equal to WRJCTZNum, then proceed to S304; otherwise, proceed to S308.

[0137] S308. Calculate the global canopy temperature mean value QJHWGCMean=QJLJTemp1 / QJLJTemp2;

[0138] Calculate the global temperature mean QJKQWDMean: Calculate the arithmetic mean of the KQWD field for all observation times of all samples in WRJCDataSet. The resulting arithmetic mean is QJKQWDMean.

[0139] S309. Calculate the global residual deviation from the baseline QJYCBaseline: For all samples in WRJCDataSet at all observation times, calculate the absolute value of the HWGC field minus the KQWD field, and then take the median of all absolute values ​​as QJYCBaseline.

[0140] If QJYCBaseline is less than 0.01, then set QJYCBaseline=0.01; otherwise, QJYCBaseline remains unchanged.

[0141] S310, Calculate the wind temperature heat conduction correction coefficient QJFTXZCoeff:

[0142] QJFTXZCoeff=tanh(QJHWGCMean-QJKQWDMean) / (QJYCBaseline+0.1)

[0143] Where tanh is the hyperbolic tangent function; QJFTXZCoeff characterizes the degree of correction of the canopy and air temperature differences relative to the residual baseline at the global scale.

[0144] The other steps and parameters are the same as in Specific Implementation Method Four.

[0145] In addition, each sample in the training set contains the following fields:

[0146] WRID field: Sample ID;

[0147] The TXSeqData field specifically includes HWGC, KQWD, XDSD, FSSL, TRFS, and ZBZS;

[0148] JBLabel field: The stress category label of the sample, specifically including four categories: normal, water stress only, nutrient stress only, and water-fertilizer coupled stress. 1 indicates normal, 2 indicates water stress only, 3 indicates nutrient stress only, and 4 indicates water-fertilizer coupled stress.

[0149] Similarly, the SFTrackData and YFTrackData fields for each sample can be obtained;

[0150] SFTrackData field: Segmentation characteristics of moisture stress channels;

[0151] YFTrackData field: Nutrient stress channel segmentation characteristics.

[0152] The water stress channel segmentation features and nutrient stress channel segmentation features of the training set samples are used as inputs to the water and fertilizer stress residual difference dual-channel early warning network SFYCNet, and the stress category label of the samples is used as the output target. When training the water and fertilizer stress residual difference dual-channel early warning network SFYCNet using the training set samples, the cross-entropy loss function and Adam optimizer are used. The default training epochs are set to 200 epochs and the default learning rate is set to 0.001, so that the trained SFYCNet has the decision-making ability to predict the crop water and fertilizer stress category based on the input dual-channel segmentation feature sequence.

[0153] Specific Implementation Method Six: This implementation method is a further limitation of Specific Implementation Method Five. The specific process of S105 is as follows:

[0154] Step 1) Initialize the decoupling gradient counter SXTJGCounter=1;

[0155] Step 2) Let the first temporary decoupled gradient vector SXTJGTemp1 = the NYCZDim dimension difference vector obtained by subtracting the SXTJGCounter row from the SXTJGCounter+1 row of SXTJFeaMat;

[0156] Write SXTJGTemp1 into the SXTJGCounter row of SXTJGradMat;

[0157] Step 3) Set SXTJGCounter = SXTJGCounter + 1;

[0158] If SXTJGCounter is less than or equal to TXNum-1, proceed to step 2; otherwise, proceed to step S106.

[0159] The other steps and parameters are the same as in Specific Implementation Method 5.

[0160] Specific Implementation Method Seven: This implementation method is a further limitation of Specific Implementation Method Six. The specific process of S106 is as follows:

[0161] Step (1): Set the decoupling amplitude counter SXTJMCounter = 1;

[0162] Step (2): Let the first temporary variable for decoupling amplitude, SXTJMTemp1, be:

[0163] Take out each element in the SXTJMCounter row of SXTJGradMat, calculate the square of each element, sum all the squares, and then take the square root of the sum. The result of the square root operation is SXTJMTemp1.

[0164] Step (3): Set SXTJMTemp1 as the SXTJMCounter-th element of SXTJGradMag;

[0165] Let SXTJMCounter = SXTJMCounter + 1;

[0166] If SXTJMCounter is less than or equal to TXNum-1, proceed to step (2); otherwise, proceed to step S107.

[0167] The other steps and parameters are the same as in Specific Implementation Method Six.

[0168] Specific Implementation Method Eight: This implementation method is a further limitation of Specific Implementation Method Seven. The specific process of S108 is as follows:

[0169] Step ①: Initialize the water channel feature matrix SXTJSFMat=TXNum rows and NYCZDim columns, and set all elements in the matrix to 0.

[0170] Initialize the nutrient channel feature matrix SXTJYFMat=TXNum rows and NYCZDim columns, with all elements in the matrix having a value of 0;

[0171] Step 2: Set the first row of SXTJSFMat to the first row of SXTJFeaMat multiplied by 0.5; set the first row of SXTJYFMat to the first row of SXTJFeaMat multiplied by 0.5;

[0172] Step ③: Set the decoupling counter SXTJFCounter = 1;

[0173] Step 4: Decouple and separate the first temporary variable SXTJFTemp1, and retrieve the SXTJFCounter-th element of SXTJGradMag;

[0174] The second temporary variable, SXTJFTemp2, is used to calculate the water channel weight decoupling separation.

[0175] SXTJFTemp2=1 / (1+exp(-8.0×(SXTJFTemp1 / SXTJThreshold-1.0)))

[0176] When SXTJFTemp2 is close to 1, it indicates that the current gradient is a fast response and belongs to the water channel.

[0177] The third temporary variable, SXTJFTemp3, is used to calculate the nutrient channel weight decoupling separation.

[0178] SXTJFTemp3=1-SXTJFTemp2

[0179] Step 5: Multiply each element of the SXTJFCounter+1 row of SXTJSFMat by SXTJFTemp2.

[0180] Let each element of the SXTJFCounter+1 row of SXTJYFMat be multiplied by SXTJFTemp3;

[0181] Step 6: Set SXTJFCounter = SXTJFCounter + 1;

[0182] If SXTJFCounter is less than or equal to TXNum-1, proceed to step ④; otherwise, proceed to step ⑦.

[0183] Step 7: Let the number of time points contained in each observation data segment be SegLen (an integer greater than or equal to 1), satisfying SegLen=Round(TXNum / SegNum), where Round is the rounding process, and SegNum is the segment size of the observation data. When the result of TXNum / SegNum is not an integer, the number of time points contained in the last segment is less than SegLen.

[0184] Set the decoupling segment counter SXTJSegCounter to 1;

[0185] Step 8: Let the first temporary variable of the decoupling segment, SXTJSegTemp1, be equal to the row from (SXTJSegCounter-1)×SegLen+1 to the Min(SXTJSegCounter×SegLen, TXNum) row of SXTJSFMat. Calculate the mean of the elements in the same column of each row to obtain an NYCZDim dimensional vector.

[0186] And set the SXTJSegCounter row of the water stress channel segmentation feature SFTrack to SXTJSegTemp1;

[0187] Where Min represents the smaller of the two numbers;

[0188] Let the second temporary variable of the decoupling segment, SXTJSegTemp2, be equal to the row from (SXTJSegCounter-1)×SegLen+1 to the Min(SXTJSegCounter×SegLen, TXNum) of SXTJYFMat. Then, calculate the mean of the elements in the same column of each row to obtain an NYCZDim dimensional vector.

[0189] And set the SXTJSegCounter row of the nutrient stress channel segmentation feature YFTrack to SXTJSegTemp2;

[0190] Step 9: Set SXTJSegCounter = SXTJSegCounter + 1;

[0191] If SXTJSegCounter is less than or equal to SegNum, proceed to step ⑧; otherwise, proceed to step ⑩.

[0192] Step 10: Perform arctan compression on each element in SFTrack to obtain the final SFTrack;

[0193] Specifically, the Arctan compression method is as follows:

[0194] x'=arctan(x) / (π / 2)

[0195] Where arctan is the arctangent function, x is any element in SFTrack, x' represents the arctan compression result of x, and π is pi;

[0196] Each row in YFTrack is subjected to exponential decay weighting to obtain the final YFTrack;

[0197] Specifically, for the v-th element in any row, the exponential decay weighted method is as follows:

[0198] y'=y×exp(-0.05v)

[0199] Where y is the v-th element in a row, and y' represents the exponential decay weighted result corresponding to y; the value of v ranges from 1 to NYCZDim. The larger the indicator dimension number, the stronger the decay, so that the lower dimension number (such as the original value of the residual, the relative residual, and other core indicators) can obtain higher weight in the nutrient channel.

[0200] The other steps and parameters are the same as in Specific Implementation Method Seven.

[0201] Specific Implementation Method Nine: This implementation method further defines Specific Implementation Method Eight. The processing procedure of the water and fertilizer stress residual difference dual-channel early warning network SFYCNet is as follows:

[0202] Step 3: 1. Use the segmented features of the water channel as the input of the water stress channel branch. Within the water stress channel branch, process the segmented features of the water channel through the LSTM model and extract the hidden state of the last time step of the LSTM model as the output of the water stress channel branch. The output of the water stress channel branch is a 32-dimensional hidden state vector.

[0203] The nutrient channel segmentation features are used as the input to the nutrient stress channel branch. Within the nutrient stress channel branch, the nutrient channel segmentation features are processed through an LSTM model, and the hidden state of the last time step of the LSTM model is extracted as the output of the nutrient stress channel branch. The output of the nutrient stress channel branch is a 32-dimensional hidden state vector.

[0204] Step 32: Concatenate the outputs of the water stress channel branch and the nutrient stress channel branch along the dimensional direction to form a 64-dimensional fusion vector;

[0205] Step 33: Use MLP to complete the water and fertilizer stress detection of the splicing results.

[0206] The other steps and parameters are the same as in Specific Implementation Method 8.

[0207] Specific Implementation Method Ten: This implementation method further defines Specific Implementation Method Nine. The working process of the MLP is as follows:

[0208] Step 331: Pass the splicing result through the first hidden layer. The first hidden layer maps the 64-dimensional splicing result to 32 dimensions. Then, pass the mapping result of the first hidden layer through the first ReLU activation function layer (the ReLU function sets negative values ​​to 0 and keeps positive values ​​unchanged).

[0209] Step 332: Use the output of the first ReLU activation function layer as the input of the second hidden layer. The second hidden layer maps the 32-dimensional features to 16-dimensional features, and the mapping result of the second hidden layer is passed through the second ReLU activation function layer.

[0210] Step 3: The output of the second ReLU activation function layer is passed through a fully connected layer. The 16-dimensional features are mapped to the SFJBNum dimension (SFJBNum represents the number of categories, which is 4 in this invention). The output of the fully connected layer is passed through a softmax activation function layer. The probability distribution of each stress category is obtained through the softmax activation function layer. The category with the highest probability is output as the detection result.

[0211] The other steps and parameters are the same as in Specific Implementation Method Nine.

[0212] The present invention takes corresponding processing measures based on the detection results:

[0213] If there are no obvious signs of water or nutrient stress at present, it is recommended to maintain the existing water and fertilizer management plan;

[0214] If the output result indicates water stress, and the canopy energy surplus shows a rapid upward trend, it is recommended to irrigate in time to replenish water to prevent the stress from worsening.

[0215] If the output result indicates nutrient stress and the residual energy deficit in the canopy shows a slow response and gradual accumulation trend, it is recommended to apply an appropriate amount of fertilizer as soon as possible to restore photosynthetic activity.

[0216] If the output result indicates water and fertilizer coupled stress, and the canopy energy surplus shows dual-channel anomalies of fast and slow response, it is recommended to simultaneously irrigate and top-dress to synergistically alleviate the compound stress.

[0217] The above examples of the present invention are merely illustrative of the computational model and process of the present invention, and are not intended to limit the implementation of the present invention. Those skilled in the art will recognize that other variations or modifications can be made based on the above description. It is impossible to exhaustively list all possible implementations here. Any obvious variations or modifications derived from the technical solutions of the present invention are still within the scope of protection of the present invention.

Claims

1. A crop water and fertilizer stress early warning method based on energy surplus and UAV thermal infrared radiation, characterized in that, The method specifically includes the following steps: Step 1: Obtain UAV thermal infrared time-series observation data YJNewData for the area to be monitored. The time-series observation data YJNewData includes observation records at TXNum time points, and each observation record at each time point includes the following fields: HWGC field: Thermal infrared observation of canopy temperature; KQWD field: Records the ambient temperature; XDSD field: Records the relative humidity of the environment; FSSL field: Records the ambient wind speed; TRFS field: Records total solar radiation; ZBZS field: The obtained Normalized Difference Vegetation Index (NDVI) value; Step 2: Use the temporal gradient stress response decoupling module SXTJModel to extract the water stress channel segmentation features SFTrack and nutrient stress channel segmentation features YFTrack from the observed data; The specific process of step two is as follows: S101, the water stress channel segmentation feature SFTrack is a list of SegNum vectors, and each vector in the list is initialized as an NYCZDim dimensional all-zero vector; The nutrient stress channel segmentation feature YFTrack is also a list of SegNum vectors, and each vector in the list is initialized as an NYCZDim dimensional all-zero vector; Establish the residual index time series matrix SXTJFeaMat, and initialize the residual index time series matrix SXTJFeaMat as a matrix of all zeros with rows of TXNum and columns of NYCZDim; Use the time-series observation data YJNewData as the input data SXTJInput of the time-series gradient stress response decoupling module SXTJModel, and initialize the decoupling time-series counter SXTJTCounter=1; S102, Decoupling the timing, the first temporary variable SXTJTTemp1 = retrieves the SXTJTCounter observation record of SXTJInput, and uses the retrieved data as the input data NYCZInput of the energy residual deviation diagnostic index extraction module NYCZModel; S103, the second temporary variable of the decoupling timing, SXTJTTemp2, is the result of processing the input data NYCZInput using the module NYCZModel, NYCZOutput. Write SXTJTTemp2 into the SXTJTCounter row of the residual index time series matrix SXTJFeaMat; S104. Set SXTJTCounter = SXTJTCounter + 1; If SXTJTCounter is less than or equal to TXNum, proceed to S102; otherwise, proceed to S105. S105. Calculate the time series gradient matrix SXTJGradMat based on the residual index time series matrix SXTJFeaMat. The specific process of S105 is as follows: Step 1) Initialize the decoupling gradient counter SXTJGCounter=1; Step 2) Let the first temporary decoupled gradient vector SXTJGTemp1 = the NYCZDim dimension difference vector obtained by subtracting the SXTJGCounter row from the SXTJGCounter+1 row of SXTJFeaMat; Write SXTJGTemp1 into the SXTJGCounter row of SXTJGradMat; Step 3) Set SXTJGCounter = SXTJGCounter + 1; If SXTJGCounter is less than or equal to TXNum-1, then proceed to step 2; otherwise, proceed to step S106. S106. Calculate the gradient magnitude sequence SXTJGradMag based on the temporal gradient matrix SXTJGradMat. The specific process of S106 is as follows: Step (1): Set the decoupling amplitude counter SXTJMCounter = 1; Step (2): Let the first temporary variable for decoupling amplitude, SXTJMTemp1, be: Take out each element in the SXTJMCounter row of SXTJGradMat, calculate the square of each element, sum all the squares, and then take the square root of the sum. The result of the square root operation is SXTJMTemp1. Step (3): Set SXTJMTemp1 as the SXTJMCounter-th element of SXTJGradMag; Let SXTJMCounter = SXTJMCounter + 1; If SXTJMCounter is less than or equal to TXNum-1, proceed to step (2); otherwise proceed to S107. S107. Calculate the adaptive gradient decoupling threshold SXTJThreshold: Take the median of all elements in SXTJGradMag as SXTJThreshold. If SXTJThreshold is less than 0.001, set SXTJThreshold=0.001; otherwise, leave SXTJThreshold unchanged. S108. Calculate the water channel feature matrix SXTJSFMat based on SXTJGradMag, SXTJThreshold, and SXTJFeaMat; calculate the nutrient channel feature matrix SXTJYFMat based on SXTJGradMag, SXTJThreshold, and SXTJFeaMat. Then, based on the water channel feature matrix SXTJSFMat, the water stress channel segmentation feature SFTrack is obtained, and based on the nutrient channel feature matrix SXTJYFMat, the nutrient stress channel segmentation feature YFTrack is obtained. The specific process of S108 is as follows: Step ①: Initialize the water channel feature matrix SXTJSFMat=TXNum rows and NYCZDim columns, and set all elements in the matrix to 0. Initialize the nutrient channel feature matrix SXTJYFMat=TXNum rows and NYCZDim columns, with all elements in the matrix having a value of 0; Step 2: Set the first row of SXTJSFMat to the first row of SXTJFeaMat multiplied by 0.5; set the first row of SXTJYFMat to the first row of SXTJFeaMat multiplied by 0.5; Step ③: Set the decoupling counter SXTJFCounter = 1; Step 4: Decouple and separate the first temporary variable SXTJFTemp1, and retrieve the SXTJFCounter-th element of SXTJGradMag; The second temporary variable, SXTJFTemp2, is used to calculate the water channel weight decoupling separation. SXTJFTemp2=1 / (1+exp(-8.0×(SXTJFTemp1 / SXTJThreshold-1.0))) The third temporary variable, SXTJFTemp3, is used to calculate the nutrient channel weight decoupling separation. SXTJFTemp3=1-SXTJFTemp2 Step 5: Multiply each element of the SXTJFCounter+1 row of SXTJSFMat by SXTJFTemp2. Let each element of the SXTJFCounter+1 row of SXTJYFMat be multiplied by SXTJFTemp3; Step 6: Set SXTJFCounter = SXTJFCounter + 1; If SXTJFCounter is less than or equal to TXNum-1, proceed to step ④; otherwise, proceed to step ⑦. Step 7: Let the number of time segments contained in each observation data segment be SegLen, satisfying SegLen=Round(TXNum / SegNum), where Round is the rounding process and SegNum is the segment size of the observation data. Set the decoupling segment counter SXTJSegCounter to 1; Step 8: Let the first temporary variable of the decoupling segment, SXTJSegTemp1, be equal to the row from (SXTJSegCounter-1)×SegLen+1 to the Min(SXTJSegCounter×SegLen, TXNum) row of SXTJSFMat. Calculate the mean of the elements in the same column of each row to obtain an NYCZDim dimensional vector. And set the SXTJSegCounter row of the water stress channel segmentation feature SFTrack to SXTJSegTemp1; Where Min represents the smaller of the two numbers; Let the second temporary variable of the decoupling segment, SXTJSegTemp2, be equal to the row from (SXTJSegCounter-1)×SegLen+1 to the Min(SXTJSegCounter×SegLen, TXNum) of SXTJYFMat. Then, calculate the mean of the elements in the same column of each row to obtain an NYCZDim dimensional vector. And set the SXTJSegCounter row of the nutrient stress channel segmentation feature YFTrack to SXTJSegTemp2; Step 9: Set SXTJSegCounter = SXTJSegCounter + 1; If SXTJSegCounter is less than or equal to SegNum, proceed to step ⑧; otherwise, proceed to step ⑩. Step 10: Perform arctan compression on each element in SFTrack to obtain the final SFTrack; perform exponential decay weighting on each row in YFTrack to obtain the final YFTrack. The working process of the energy residual deviation diagnostic index extraction module NYCZModel is as follows: S201. Initialize NYCZOutput as a zero-based vector with dimension NYCZDim; S202. Obtain the first temporary variable NYCZTemp1 to the 18th temporary variable NYCZTemp18 based on the input data NYCZInput; The specific process of S202 is as follows: Step 1: Extract the HWGC field of NYCZInput as the first temporary variable for residual difference, NYCZTemp1; Extract the KQWD field of NYCZInput as the second temporary variable for residual difference, NYCZTemp2; Extract the XDSD field of NYCZInput and use it as the third temporary variable for the residual difference, NYCZTemp3; Extract the FSSL field of NYCZInput as the fourth temporary variable NYCZTemp4 for the residual. If NYCZTemp4 is less than 0.1, set NYCZTemp4=0.

1. If NYCZTemp4 is greater than or equal to 0.1, then NYCZTemp4 remains unchanged. Extract the TRFS field of NYCZInput as the fifth temporary variable of the residual, NYCZTemp5. If NYCZTemp5 is less than 1.0, set NYCZTemp5=1.

0. If NYCZTemp5 is greater than or equal to 1.0, NYCZTemp5 remains unchanged. Extract the ZBZS field of NYCZInput as the 6th temporary variable of the residual NYCZTemp6; Step 2: Calculate the 7th temporary variable NYCZTemp7 of the net radiation residual: Where Albedo is the surface albedo; Emissivity is the surface emissivity; Step 3: Calculate the 8th temporary variable NYCZTemp8 for the residual soil heat flux: Wherein, SHFRatio is the soil heat flux ratio coefficient; Step 4: Calculate the 9th temporary variable NYCZTemp9 for the aerodynamic impedance residual error: Where log is the natural logarithm function, xtck is the drone hovering reference altitude, and zwccd is the crop surface roughness; Step 5: Calculate the 10th temporary variable NYCZTemp10 of the reference latent heat flux residual: Where PMKQ is the Penman reference evapotranspiration coefficient; Step 6: Calculate the 11th temporary variable NYCZTemp11 of the reference inductive heat flux residual: Step 7: Calculate the 12th temporary variable, NYCZTemp12, of the theoretical canopy temperature residual: Step 8: Calculate the residual of the original energy value, the 13th temporary variable NYCZTemp13: Step 9: Calculate the 14th temporary variable NYCZTemp14, representing the actual sensible heat flux residual: Step 10: Calculate the 15th temporary variable NYCZTemp15 for the actual latent heat flux residual: Step 11: Calculate the 16th temporary variable NYCZTemp16 of the saturated water vapor pressure residual difference: Where exp represents an exponential function with the natural constant e as its base; Step 12: Calculate the 17th temporary variable NYCZTemp17 of the actual water vapor pressure difference: Step 13: Calculate the 18th temporary variable NYCZTemp18, representing the residual water vapor pressure difference: S203, Set NYCZOutput[1]=NYCZTemp13 / 3.0; Among them, NYCZTemp13 is the 13th temporary variable, and NYCZOutput[1] is the first element in NYCZOutput; S204. Set NYCZOutput[2]=NYCZTemp13 / (ABS(NYCZTemp12)+0.5); Where ABS is the absolute value, NYCZTemp12 is the 12th temporary variable, and NYCZOutput[2] is the 2nd element in NYCZOutput; S205, Set NYCZOutput[3]=(1-NYCZTemp15 / (ABS(NYCZTemp10)+0.01))×0.8+0.1; Among them, NYCZTemp15 is the 15th temporary variable, NYCZTemp10 is the 10th temporary variable, and NYCZOutput[3] is the 3rd element in NYCZOutput; S206. Set NYCZOutput[4] to: Among them, NYCZTemp14 is the 14th temporary variable, NYCZTemp11 is the 11th temporary variable, and NYCZOutput[4] is the 4th element in NYCZOutput; S207, Set NYCZOutput[5]=NYCZTemp13×NYCZTemp18 / (QJYCBaseline×10+0.5); Among them, QJYCBaseline is the global residual deviation from the baseline, NYCZOutput[5] is the 5th element in NYCZOutput, and NYCZTemp18 is the 18th temporary variable; S208, Set NYCZOutput[6]=NYCZTemp13×(1+0.2 / (NYCZTemp4+0.05)); Among them, NYCZOutput[6] is the 6th element in NYCZOutput, and NYCZTemp4 is the 4th temporary variable; S209, Set NYCZOutput[7]=NYCZTemp13 / (NYCZTemp7×0.001+0.3); Among them, NYCZOutput[7] is the 7th element in NYCZOutput, and NYCZTemp7 is the 7th temporary variable; S210, Set NYCZOutput[8]=NYCZTemp13×(0.4+0.6×NYCZTemp6)+0.05; Among them, NYCZOutput[8] is the 8th element in NYCZOutput, and NYCZTemp6 is the 6th temporary variable; S211, Set NYCZOutput[9]; Where QJHWGCMean is the global average canopy temperature, and QJKQWDMean is the global average air temperature; S212, Set NYCZOutput[10]; Among them, NYCZOutput[10] is the 10th element in NYCZOutput, NYCZTemp3 is the 3rd temporary variable, and QJFTXZCoeff is the wind temperature heat conduction correction coefficient; If the value of NYCZOutput[10] is greater than 2.0, then set NYCZOutput[10] = 2.0; otherwise, the value of NYCZOutput[10] remains unchanged. S213. Let the final output NYCZOutput = 1 / (1 + exp(-NYCZOutput)); S214. Use the final output NYCZOutput obtained in S213 as the result output of NYCZModel; Step 3: Establish the water and fertilizer stress residual difference dual-channel early warning network SFYCNet. The water and fertilizer stress residual difference dual-channel early warning network SFYCNet includes a water stress channel branch, a nutrient stress channel branch, and an MLP; wherein: the water stress channel branch and the nutrient stress channel branch both include an LSTM model; The water stress channel segmentation feature SFTrack and the nutrient stress channel segmentation feature YFTrack are used as inputs to the water stress channel branch and the nutrient stress channel branch, respectively. The outputs of the water stress channel branch and the nutrient stress channel branch are then concatenated along the channel dimension. The concatenation result is then processed by an MLP, and the detection result of the monitored area is output through the MLP. The detection result includes: no signs of water stress or nutrient stress, water stress, nutrient stress, and water-fertilizer coupling stress. If water or nutrient stress exists, an early warning should be issued; otherwise, no action is necessary.

2. The crop water and fertilizer stress early warning method based on energy surplus and UAV thermal infrared as described in claim 1, characterized in that, The calculation methods for the global canopy temperature mean QJHWGCMean, the global air temperature mean QJKQWDMean, the global residual deviation from the baseline QJYCBaseline, and the wind-temperature heat conduction correction coefficient QJFTXZCoeff are as follows: S301. Obtain the UAV thermal infrared time series training dataset WRJCDataSet. The total number of training samples in the dataset WRJCDataSet is WRJCTZNum, and the number of observation times for each training sample is TXNum. S302. Establish a global statistical accumulation variable QJLJTemp1=0 and a global statistical counting variable QJLJTemp2=0; S303. Establish a global statistical sample counter QJYBCounter=1; S304. Establish a global statistical time counter QJSKCounter=1; S305. Set the global statistical temporary variable QJTemp1 to retrieve the HWGC field of the QJSKCounter observation record of the QJYBCounter sample in WRJCDataSet; Let QJLJTemp1 = QJLJTemp1 + QJTemp1; QJLJTemp2=QJLJTemp2+1; S306. Set QJSKCounter = QJSKCounter + 1; If QJSKCounter is less than or equal to TXNum, proceed to S305; otherwise, proceed to S307. S307. Let QJYBCounter = QJYBCounter + 1; If QJYBCounter is less than or equal to WRJCTZNum, then proceed to S304; otherwise, proceed to S308. S308. Calculate the global canopy temperature mean value QJHWGCMean=QJLJTemp1 / QJLJTemp2; Calculate the global temperature mean QJKQWDMean: Calculate the arithmetic mean of the KQWD field for all observation times of all samples in WRJCDataSet. The resulting arithmetic mean is QJKQWDMean. S309. Calculate the global residual deviation from the baseline QJYCBaseline: For all samples in WRJCDataSet at all observation times, calculate the absolute value of the HWGC field minus the KQWD field, and then take the median of all absolute values ​​as QJYCBaseline. If QJYCBaseline is less than 0.01, then set QJYCBaseline=0.01; otherwise, QJYCBaseline remains unchanged. S310, Calculate the wind temperature heat conduction correction coefficient QJFTXZCoeff: QJFTXZCoeff=tanh(QJHWGCMean-QJKQWDMean) / (QJYCBaseline+0.1) Where tanh is the hyperbolic tangent function.

3. The crop water and fertilizer stress early warning method based on energy surplus and UAV thermal infrared as described in claim 2, characterized in that, The processing procedure of the water and fertilizer stress residual difference dual-channel early warning network SFYCNet is as follows: Step 3:

1. Use the segmented features of the water channel as the input of the water stress channel branch. Within the water stress channel branch, process the segmented features of the water channel through the LSTM model and extract the hidden state of the last time step of the LSTM model as the output of the water stress channel branch. The nutrient channel segmentation features are used as the input to the nutrient stress channel branch. Within the nutrient stress channel branch, the nutrient channel segmentation features are passed through the LSTM model, and the hidden state of the last time step of the LSTM model is extracted as the output of the nutrient stress channel branch. Step 32: Segment the outputs of the water stress channel branch and the nutrient stress channel branch along the dimensional direction; Step 33: Use MLP to complete the water and fertilizer stress detection of the splicing results.

4. The crop water and fertilizer stress early warning method based on energy surplus and UAV thermal infrared as described in claim 3, characterized in that, The working process of the MLP is as follows: Step 331: Pass the splicing result through the first hidden layer, and then pass the mapping result of the first hidden layer through the first ReLU activation function layer; Step 332: Use the output of the first ReLU activation function layer as the input of the second hidden layer. The second hidden layer maps the 32-dimensional features to 16-dimensional features, and the mapping result of the second hidden layer is passed through the second ReLU activation function layer. Step 3: Pass the output of the second ReLU activation function layer through a fully connected layer, then pass the output of the fully connected layer through a softmax activation function layer. Obtain the probability distribution of each stress category through the softmax activation function layer, and output the category with the highest probability as the detection result.

Images