A subsurface soil moisture inversion method and system based on P-band radiometer
Through the combination of multi-angle observation based on P-band radiometer and multiple model algorithms, dynamic adjustment of weights and iterative optimization, high-precision inversion of deeper soil moisture is achieved, solving the limitations of traditional methods for monitoring shallow soil moisture.
Patent Information
- Application Number
- CN202411757462.8
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2024-12-03
- Publication Date
- 2025-02-28
- Estimated Expiration
- 2044-12-03
AI Technical Summary
The existing passive microwave remote sensing method is sensitive to surface roughness and vegetation moisture content in soil moisture monitoring, and can only observe shallow soil depths (0-5cm) and cannot effectively invert deeper soil moisture.
The subsurface soil moisture inversion method based on P-band radiometer is adopted, and the soil moisture inversion is inverted through various models and algorithms such as dielectric hybrid model, Choudhury module and radiation transmission equation, combined with multi-angle observation data and optimization algorithms, the weights and iterative optimization are dynamically adjusted to invert soil moisture.
It breaks the limitations of traditional methods and can invert deeper soil layers (about 10-15cm), reduces sensitivity to vegetation and surface roughness, and improves the accuracy of soil moisture inversion and model adaptability.
Smart Images

Figure CN119246561B_ABST
Abstract
Description
Technical Field
[0001] The invention belongs to the field of hydrological remote sensing science and discloses a subsurface soil moisture inversion method and system based on a P-band radiometer. Background Art
[0002] Soil moisture is a key factor in maintaining plant growth and ecological balance. The most commonly used method is passive microwave remote sensing, which monitors soil moisture based on the correlation between passive microwave sensors and specific wavelengths. However, passive microwave remote sensing methods have limitations because the inversion process of the L-band radiometer used in most passive microwave remote sensing methods is extremely sensitive to surface roughness and the moisture content of vegetation covering the ground, and it can currently only observe shallow soil depths (0-5cm).
[0003] Therefore, further research and innovation are needed to solve the above problems existing in the prior art. Summary of the invention
[0004] The purpose of the invention is to provide a sub-surface soil moisture inversion method and system based on a P-band radiometer to solve the above-mentioned problems existing in the prior art.
[0005] The technical solution is a subsurface soil moisture inversion method based on a P-band radiometer, which includes the following steps:
[0006] S1. Collect observation data and preprocess them to build initial parameters for the preconfigured dielectric hybrid model; the observation data includes raw brightness temperature data, soil temperature data and initial roughness parameter set from a P-band radiometer at different incident angles;
[0007] S2. Based on the initial parameters of the dielectric hybrid model and the preset initial value of soil moisture, the relative dielectric constant of the soil is calculated, and based on the relative dielectric constant and the initial roughness parameter set, the soil emissivity at each incident angle is calculated;
[0008] S3, based on soil temperature data and relative dielectric constant, call the preconfigured Choudhury module to calculate the effective soil temperature;
[0009] S4. Based on the soil emissivity and effective temperature, the simulated brightness temperature is calculated using the zero-order approximate model of the radiation transfer equation; the brightness temperature difference between the original brightness temperature data and the simulated brightness temperature is calculated, and based on the brightness temperature difference, a cost function based on multi-angle data is constructed;
[0010] S5. Based on the cost function, the optimization algorithm is called to update the soil moisture and roughness parameters to obtain the optimized parameters, and the weight of each angle data in the cost function is dynamically adjusted according to the performance of different angles;
[0011] S6. Based on the optimized parameters, repeat steps S2 to S5, iteratively perform soil moisture inversion calculation until a preset convergence condition is met, and obtain an inverted soil moisture result.
[0012] A subsurface soil moisture inversion system based on P-band radiometer, comprising:
[0013] at least one processor; and,
[0014] a memory communicatively connected to at least one of the processors; wherein,
[0015] The memory stores instructions that can be executed by the processor, and the instructions are used to be executed by the processor to implement the sub-surface soil moisture inversion method based on the P-band radiometer.
[0016] Beneficial effects: The present invention proposes a sub-surface soil moisture inversion method and system based on a P-band radiometer, which not only breaks the limitations of traditional methods, detects deeper soil layers through the P-band, reduces the sensitivity to vegetation and surface roughness, and improves the accuracy of soil moisture inversion; it also ensures the adaptability and flexibility of the model through dynamic adjustment of weights and iterative optimization. BRIEF DESCRIPTION OF THE DRAWINGS
[0017] Figure 1 It is a flow chart of the present invention.
[0018] Figure 2 This is a flow chart of step S1 of the present invention.
[0019] Figure 3 This is a flow chart of step S2 of the present invention.
[0020] Figure 4 This is a flow chart of step S3 of the present invention.
[0021] Figure 5 This is a flow chart of step S4 of the present invention.
[0022] Figure 6 This is a flow chart of step S5 of the present invention.
[0023] Figure 7 This is a flow chart of step S6 of the present invention. DETAILED DESCRIPTION
[0024] like Figure 1 As shown, this application proposes a subsurface soil moisture inversion method based on a P-band radiometer, comprising the following steps:
[0025] S1. Collect observation data and preprocess them to build initial parameters for the preconfigured dielectric hybrid model; the observation data includes raw brightness temperature data, soil temperature data and initial roughness parameter set from a P-band radiometer at different incident angles;
[0026] S2. Based on the initial parameters of the dielectric hybrid model and the preset initial value of soil moisture, the relative dielectric constant of the soil is calculated, and based on the relative dielectric constant and the initial roughness parameter set, the soil emissivity at each incident angle is calculated;
[0027] S3, based on soil temperature data and relative dielectric constant, call the preconfigured Choudhury module to calculate the effective soil temperature;
[0028] S4. Based on the soil emissivity and effective temperature, the simulated brightness temperature is calculated using the zero-order approximate model of the radiation transfer equation; the brightness temperature difference between the original brightness temperature data and the simulated brightness temperature is calculated, and based on the brightness temperature difference, a cost function based on multi-angle data is constructed;
[0029] S5. Based on the cost function, the optimization algorithm is called to update the soil moisture and roughness parameters to obtain the optimized parameters, and the weight of each angle data in the cost function is dynamically adjusted according to the performance of different angles;
[0030] S6. Based on the optimized parameters, repeat steps S2 to S5, iteratively perform soil moisture inversion calculation until a preset convergence condition is met, and obtain an inverted soil moisture result.
[0031] This embodiment combines the original brightness temperature data, soil temperature data and initial roughness parameter set of the P-band radiometer, and comprehensively calculates the relative dielectric constant, soil emissivity and effective temperature of the soil through multiple models and algorithms such as dielectric mixing model, Kirchhoff's law and radiation transfer equation, and finally achieves accurate inversion of soil moisture. Among them, the P-band (about 40cm) has a longer wavelength, which can invert deeper (about 10-15cm) soil moisture compared to the traditional L-band (about 21cm). At the same time, the longer wavelength of the P-band also reduces the impact of surface roughness and vegetation on the inversion accuracy, thereby improving the accuracy of soil moisture inversion.
[0032] like Figure 2 As shown, according to one aspect of the present application, the basic soil data is collected and preprocessed in step S1 further as follows:
[0033] S11, collecting observation data, including original brightness temperature data, soil temperature data and initial roughness parameter set at different incident angles;
[0034] S12, preprocessing the raw brightness temperature data at each angle, including removing outliers, unit conversion, and normalization or standardization, to obtain preprocessed raw brightness temperature data;
[0035] S13, performing multi-angle adaptive fusion based on the pre-processed original brightness temperature data to obtain adaptively fused brightness temperature data;
[0036] The brightness temperature data after adaptive fusion is:
[0037] Adapted TB (t)=∑ i=1 n α i (t)·Preprocessed TBθi (t)+γ∑ i=1 n (1-α i (t))·Temporal Filtered _ TBθi (t);
[0038] Among them, γ is a smoothing factor, which is used to balance the impact of the current time and the previous and next time data. TB Represents the brightness temperature data after adaptive fusion, Preprocessed TBθi (t) represents the incident angle θ i and the preprocessed brightness temperature data at time t, Temporal Filtered _ TBθi (t) represents the incident angle θ i and the brightness temperature data after time filtering at time t, n represents the number of angles, α i (t) represents the weight of brightness temperature data.
[0039] In this example, brightness temperature data (TB) is obtained from a P-band radiometer, including brightness temperature values at different polarization directions (horizontal polarization and vertical polarization) and different incident angles (θ), and TB data is read from the raw data file using a Python script. Soil temperature data (ST) is prepared, including the surface temperature (T surf ) and deep soil temperature (T deep ), data can be obtained from a weather station or temperature probe. Initialize the roughness parameters, including H RP (height correlation length ratio), N RP (roughness index) and Q R (polarization mixing parameter), the initial value can be set based on experience, such as H RP∈[-1,1], N RP ∈[-1,1], Q R ∈[0,1]. Set the parameters of the Dielectric Mixing Model. These parameters are used to relate soil moisture to the dielectric constant. The parameters of the Mironov model, the Dobson model, or the Peplinski model can be used.
[0040] In one embodiment of the present application, a P-band radiometer is used to collect raw brightness temperature data (TBθ1, TBθ2, ..., TBθn) at different incident angles (θ1, θ2, ..., θn); and corresponding soil temperature data and an initial roughness parameter set are collected simultaneously.
[0041] The raw brightness temperature data at each angle are preprocessed, including removal of outliers, unit conversion and normalization, to obtain the preprocessed brightness temperature data (TBθ1', TBθ2', ..., TBθn').
[0042] Multi-angle data matching is performed based on the following formula:
[0043] ;
[0044] Among them, α i (t) is the dynamic weight of the i-th angle at time t, satisfying ∑ i=1n α i (t) = 1.
[0045] Dynamic weight α i The calculation of (t) is based on the following optimization model:
[0046]
[0047] Where ΔTBθi(t) = |TBθi'(t) - TBar(t)|, TBar(t) is the average brightness temperature at all angles, and β is the adjustment parameter.
[0048] In one embodiment of the present application, the brightness temperature data after adaptive fusion of multiple angles is:
[0049] Adapted TB (t)=∑ i=1 n α i (t)·Preprocessed TBθi (t)+γ∑ i=1 n (1-α i (t))·Temporal Filtered _TBθi (t);
[0050] The weight update formula is:
[0051] α i (t) = 1 / Z · exp(-∑ τ=1 T (|Preprocessed TBθi (t) - Simulated TB (t)| / σ i (t))²);
[0052] where Z is the normalization constant, σ i (t) is the standard deviation of the brightness temperature error at the i-th angle, T is the time window size, Simulated TB (t) represents the simulated brightness temperature data at time t.
[0053] This application can effectively reduce the noise and error in single-angle observation data and improve the accuracy of brightness temperature data by fusing multi-angle observation data. Observation data from different angles can complement each other to make up for the shortcomings of single-angle observation, thereby obtaining more accurate brightness temperature data. Time filtering processing is introduced to balance the data influence of the current moment and the previous and next moments through the smoothing factor γ. This processing method can effectively reduce random fluctuations and outliers in the data, making the changes in brightness temperature data in the time series smoother and continuous, and improving the stability and reliability of the data.
[0054] The adaptive fusion method is based on the weight coefficient α of the observation data at different angles. i (t) Weighted fusion can dynamically adjust the contribution of data from each angle. This adaptive method can automatically optimize the fusion strategy according to the actual observation situation and improve the accuracy and consistency of the fused brightness temperature data. In the preprocessing stage, the original brightness temperature data can be effectively cleaned and standardized by removing outliers, unit conversion, and normalization or standardization. These preprocessing steps can eliminate outliers and inconsistencies in the data and provide a high-quality data foundation for subsequent fusion processing.
[0055] Multi-angle observation can provide more comprehensive surface information, and observation data from different angles can reflect different characteristics of the surface. By fusing multi-angle data, the changes in surface brightness temperature can be captured more comprehensively, and the spatial resolution and information richness of the data can be improved. This method of adaptively fusing brightness temperature data from multiple angles can be widely used in remote sensing, meteorology, environmental monitoring and other fields. By improving the accuracy and stability of brightness temperature data, it can provide more reliable data support for surface temperature monitoring, climate change research, agricultural management, etc.
[0056] This embodiment can improve the accuracy and reliability of brightness temperature data through technical means such as data precision improvement, time smoothing, adaptive fusion, outlier processing and multi-angle observation, and has broad application prospects and important technical value.
[0057] like Figure 3 As shown, according to one aspect of the present application, step S2 is further:
[0058] S21, reading the initial parameters of the dielectric hybrid model, combining the initial value of soil moisture, and calculating the relative dielectric constant of the soil;
[0059] S22. Combining the relative dielectric constant and the incident angle, the ideal surface reflectivity is calculated using the Fresnel equation;
[0060] S23, reading the initial roughness parameter set, combining with the ideal surface reflectivity, calling the pre-configured HQN model to calculate the actual reflectivity taking into account the roughness effect;
[0061] S24. Based on the actual reflectivity, Kirchhoff's law is used to calculate the soil emissivity at each angle.
[0062] In one embodiment of the present application, the Fresnel equation is used to calculate the reflectivity Γ* of an ideal smooth surface. This calculation is based on the relative dielectric constant ε of the soil and the observed incident angle θ.
[0063]
[0064] where Γ H *For horizontal polarization, that is, to calculate the reflectivity when the electric field vector of the electromagnetic wave is parallel to the incident plane; Γ V *For vertical polarization, that is, to calculate the reflectivity when the electric field vector of the electromagnetic wave is perpendicular to the incident plane.
[0065] Applying the HQN model to calculate the reflectivity Γ considering the roughness effect P :
[0066] ;
[0067] Use Kirchhoff's law to calculate the soil emissivity eP: eP = 1-Γ P .
[0068] This embodiment can accurately calculate the relative dielectric constant of the soil by reading the initial parameters of the dielectric hybrid model and the initial value of soil moisture. It provides accurate basic data for subsequent reflectivity and emissivity calculations, and improves the accuracy of the overall method. The ideal surface reflectivity is calculated using the Fresnel equation, which can accurately describe the physical process of electromagnetic wave reflection on the soil surface. The Fresnel equation takes into account the influence of the incident angle and the dielectric constant, so that the calculation results are more in line with the actual situation. By calling the preconfigured HQN model, combining the initial roughness parameter set and the ideal surface reflectivity, the actual reflectivity considering the roughness effect is calculated. The HQN model can effectively simulate the influence of soil surface roughness on reflectivity and improve the accuracy of reflectivity calculation. Soil emissivity calculation based on Kirchhoff's law: Based on the actual reflectivity, the soil emissivity at each angle is calculated using Kirchhoff's law. Kirchhoff's law points out that under thermal equilibrium conditions, the emissivity of an object is equal to its absorptivity, which can accurately reflect the radiation characteristics of the soil and improve the reliability of soil emissivity calculation. By calculating the soil emissivity at each angle, the fusion of multi-angle observation data can be achieved, which can fully utilize the advantages of multi-angle observation data and improve the accuracy and stability of soil moisture inversion. By comprehensively considering the dielectric properties, surface roughness and radiation properties of the soil, and combining multiple models and equations, it can more comprehensively reflect the physical properties of the soil and improve the accuracy and reliability of soil moisture inversion.
[0069] This embodiment can improve the accuracy and reliability of soil moisture inversion by accurately calculating the soil dielectric constant, ideal surface reflectivity, actual reflectivity and soil emissivity.
[0070] like Figure 4 As shown, according to one aspect of the present application, step S3 is further:
[0071] S31. Calculating a temperature gradient based on soil temperature data; wherein the soil temperature data includes surface temperature and deep soil temperature;
[0072] S32, based on the temperature gradient and relative dielectric constant, calling the pre-configured Lv model to calculate the CT parameters;
[0073] S33, based on the soil temperature data and CT parameters, calling the pre-configured Choudhury module to calculate the effective soil temperature;
[0074] In one embodiment of the present application, the Choudhury scheme is used to calculate the effective soil temperature T eff s :T eff s =T depth +CT (T surf -Tdepth ), where T surf and T depth It is the physical temperature of the surface and deep soil; CT parameter determines the surface and deep temperatures to the effective soil temperature T eff s In the field of soil science and remote sensing, CT parameters usually refer to parameters related to temperature gradients and are calculated using the Lv scheme: , where exp( ) represents the exponential function, Δx represents the thickness of the soil layer, and λ0 represents the wavelength of the electromagnetic wave.
[0075] This embodiment can accurately reflect the vertical distribution of soil temperature by calculating the temperature gradient based on the surface temperature and deep soil temperature data. The temperature gradient is a key parameter in the soil heat conduction process, which can provide information on the heat transfer inside the soil and provide basic data for the subsequent CT parameter calculation. The CT parameter is an important parameter calculated based on the temperature gradient and relative dielectric constant. The relative dielectric constant reflects the electromagnetic properties of the soil, while the temperature gradient reflects the thermal properties of the soil. By combining these two parameters, the physical state of the soil can be more accurately described and the accuracy of soil moisture inversion can be improved. By calling the pre-configured Choudhury module, the effective soil temperature is calculated based on the soil temperature data and the CT parameter. The Choudhury module is an advanced calculation model that can comprehensively consider the thermal and electromagnetic properties of the soil and provide a more accurate soil temperature estimate. The effective soil temperature is a key parameter in the soil moisture inversion and can directly affect the accuracy of the inversion result. By integrating multi-level data such as surface temperature, deep soil temperature, temperature gradient, relative dielectric constant and CT parameters, the physical and thermal properties of the soil can be fully reflected. This multi-level data fusion method can improve the accuracy and stability of soil moisture inversion and reduce the errors caused by a single data source.
[0076] This embodiment can improve the accuracy and reliability of soil moisture inversion through the application of multi-level data fusion and advanced computing models.
[0077] like Figure 5 As shown, according to one aspect of the present application, step S4 is further:
[0078] S41, based on soil emissivity and effective temperature, the simulated brightness temperature is calculated using the radiation transfer equation;
[0079] S42, calculating a brightness temperature difference value based on the original brightness temperature data and the simulated brightness temperature;
[0080] S43, obtaining current optimization parameters, and calculating parameter variation based on the current optimization parameters and the initial roughness parameter set;
[0081] S44, constructing a cost function based on multi-angle data based on brightness temperature difference values and parameter changes;
[0082] The cost function is: J(t) =∑ i=1 n α i (t) (TB θi , obs (t) - TB θi ,sim (t))² + λ∑ j=1 m (P j (t)- P j,prior (t))²;
[0083] TB θi , obs (t) represents the time t and the incident angle θ i The original brightness temperature data under TB θi ,sim (t) represents the time t and the incident angle θ i The simulated brightness temperature under j (t) represents the estimated value of the jth parameter at time t, P j,prior (t) represents the prior value of the jth parameter at time t.
[0084] This embodiment can accurately simulate the brightness temperature values of soil under different conditions by using the radiation transfer equation to calculate the simulated brightness temperature based on the soil emissivity and effective temperature. The radiation transfer equation can take into account the physical properties and environmental conditions of the soil, so that the simulated brightness temperature is closer to the actual observed value. The brightness temperature difference value is calculated based on the original brightness temperature data and the simulated brightness temperature, and the deviation between the simulated brightness temperature and the actual observed brightness temperature can be effectively evaluated. It can help identify and correct errors in the simulation process and improve the accuracy of the simulation results. By obtaining the current optimization parameters and calculating the parameter change based on these parameters and the initial roughness parameter set, the model parameters can be dynamically adjusted. In this way, the model can be continuously optimized according to the actual observation data to improve the adaptability and accuracy of the model. By constructing a cost function based on multi-angle data, the advantages of multi-angle observation data can be fully utilized. The cost function includes the difference between the original brightness temperature data and the simulated brightness temperature data at different angles, as well as the difference between the parameter estimation value and the prior value. By minimizing the cost function, the effective fusion of multi-angle data can be achieved, and the accuracy of soil moisture inversion can be improved. The second term in the cost function takes into account the difference between the parameter estimate and the prior value. By introducing the regularization term, the overfitting problem in the parameter estimation process can be prevented, and the reliability and stability of the parameter estimation can be improved. This soil moisture inversion method based on multi-angle data fusion and dynamic optimization has broad application prospects. This embodiment can improve the accuracy and reliability of soil moisture inversion by accurately simulating brightness temperature, calculating brightness temperature difference values, dynamically adjusting optimization parameters and multi-angle data fusion.
[0085] like Figure 6 As shown, according to one aspect of the present application, step S5 is further:
[0086] S51, based on the cost function, using an optimization algorithm to calculate the parameter update amount;
[0087] S52, based on the parameter update amount, updating the soil moisture and roughness parameters to obtain the optimized parameters, and in the optimization process, using an adaptive weight adjustment algorithm to dynamically adjust the weight of each angle in the cost function according to the performance of different angles;
[0088] Dynamically adjust the weights as follows:
[0089] α i ,new (t) = exp(-J θi (t)) / ∑ k=1 n exp(-J θk (t));
[0090] Where exp( ) represents the exponential function, J represents the cost function, and θ i and θ k Represents different incident angles.
[0091] This embodiment can dynamically adjust the model parameters to make it more consistent with the actual observation data by using the optimization algorithm to calculate the parameter update amount based on the cost function, so as to effectively reduce the model error and improve the accuracy of soil moisture inversion. In the optimization process, the adaptive weight adjustment algorithm is used to dynamically adjust the weight of each angle in the cost function according to the performance of different angles, so that the angle data with better performance occupies a larger proportion in the fusion process, improve the accuracy and stability of multi-angle data fusion, and make full use of the advantages of multi-angle observation data. Data from different angles can complement each other and make up for the shortcomings of single angle observation, so as to obtain more accurate and comprehensive soil moisture information. By optimizing parameters and adaptive weight adjustment, the accuracy of soil moisture inversion can be improved. The optimized parameters are more consistent with the actual situation, and the weights after adaptive adjustment make data fusion more reasonable, thereby improving the reliability of the inversion results. The adaptive weight adjustment algorithm can dynamically adjust the weight according to the performance of the actual observation data, and has strong robustness. Even if there is a certain amount of noise or outliers in the observation data, the algorithm can still effectively adjust the weight to reduce the impact on the inversion result. This embodiment can significantly improve the accuracy and reliability of soil moisture inversion through parameter optimization and adaptive weight adjustment.
[0092] like Figure 7 As shown, according to one aspect of the present application, step S6 is further:
[0093] S61, based on the optimized parameters, repeating steps S2 to S5, performing soil moisture inversion calculation, and obtaining an inverted soil moisture result; until the cost function converges or reaches the maximum number of iterations, outputting the final inverted soil moisture result;
[0094] S62, reading the soil moisture observation data, comparing the final inverted soil moisture result with the soil moisture observation data, calculating the error statistics, and evaluating the inversion accuracy based on the error statistics.
[0095] In one embodiment of the present application, based on the optimization parameters, steps S2 to S5 are repeated until the cost function converges or the maximum number of iterations is reached to obtain the final inverted soil moisture result; the actual soil moisture observation data is read from the field sensor, and the error between the inverted soil moisture result and the observation data is calculated, such as the correlation coefficient (R), the root mean square error (RMSE) and the mean absolute error (MAE). Based on the error statistics, the accuracy of the inversion result is evaluated: if the error is within an acceptable range, the inversion result is considered valid; otherwise, the model parameters or optimization algorithm may need to be adjusted.
[0096] This embodiment can effectively improve the accuracy of soil moisture inversion by repeatedly performing soil moisture inversion calculation based on optimized parameters and adjusting parameters in each iteration until the cost function converges or reaches the maximum number of iterations, and can make full use of the difference between the observed data and the model to gradually approach the true soil moisture value. After the inversion calculation is completed, the accuracy of the inversion result can be intuitively evaluated by comparing the final inverted soil moisture result with the actual observed data and calculating the error statistics. The calculation of the error statistics helps to identify the deviations and deficiencies in the inversion process, thereby providing a basis for further optimizing the inversion algorithm. Evaluating the inversion accuracy based on the error statistics can provide feedback for the optimization of the inversion algorithm. By analyzing the error statistics, the performance of the inversion algorithm under different conditions can be found, so as to improve and optimize it in a targeted manner, and improve the robustness and applicability of the algorithm. By setting the maximum number of iterations, infinite loops in the inversion calculation process can be avoided and the calculation efficiency can be improved. At the same time, the iterative calculation method based on the optimized parameters can reduce unnecessary calculations and improve the overall calculation efficiency under the premise of ensuring accuracy. By continuously optimizing the inversion algorithm and parameters, the applicability of the model under different soil types and environmental conditions can be improved. It is not only applicable to soil moisture inversion in a specific area, but can also be extended to a wider range of areas and different soil types, improving the versatility and practicality of the model. This embodiment not only improves the accuracy and computational efficiency of soil moisture inversion by optimizing parameters, iterative calculation, error evaluation and other methods, but also provides strong support for the optimization of the inversion algorithm and the applicability of the model.
[0097] According to another aspect of the present application, it also includes:
[0098] S71, obtaining a parameter set of the current iteration, including soil moisture, roughness parameters, and dielectric hybrid model parameters, and based on the parameter set of the current iteration, using a numerical differentiation method, calculating a partial derivative of the cost function with respect to each parameter to obtain a gradient vector;
[0099] S72, reading the velocity vector of the previous iteration, the current learning rate and the momentum factor, and updating the velocity vector using the momentum gradient descent method; calculating the parameter increment based on the updated velocity vector; applying the parameter increment to the parameter set of the current iteration to obtain an updated parameter set;
[0100] S73, obtaining the current number of iterations and the preset random perturbation strength, calculating the current perturbation strength based on the current number of iterations and the random perturbation strength, and generating random perturbations using Gaussian distribution; adding the random perturbations to the updated parameter set to obtain a final new parameter set;
[0101] S74, obtaining the current cost function value and the cost function value of the previous iteration, and comparing them, if the current cost function value is less than the cost function value of the previous iteration, increasing the learning rate; otherwise, reducing the learning rate; outputting the updated learning rate;
[0102] S75, obtaining the original brightness temperature data and simulated brightness temperature of each current incident angle, and calculating the current brightness temperature difference of each angle; based on the current brightness temperature difference, using an exponential function to update the angle weight, normalizing the angle weight, and obtaining an updated angle weight;
[0103] S76, based on the current cost function value, the cost function value of the previous iteration and the preset convergence threshold, calculate the change of the cost function; if the change of the cost function is less than the convergence threshold, the marking algorithm has converged and is ready to exit the loop; otherwise, continue iterating;
[0104] S77. Store the new parameter set, the current cost function value, the updated velocity vector, the updated learning rate and the updated angle weight as input data for the next iteration until a preset number of iterations is reached.
[0105] In the embodiment of the present application, the parameter set P of the current iteration is obtained. current , including soil moisture SM and roughness parameters (H RP , Q R , N RP ). Using numerical differentiation methods, calculate the partial derivatives of the cost function J with respect to each parameter and obtain the gradient vector grad. Read the velocity vector v of the previous iteration prev , the current learning rate η and momentum factor β. Apply the momentum gradient descent method to update the velocity vector: v = β * v prev + (1 -β) * grad. According to the updated velocity vector v, calculate the parameter increment: delta P = -η* v. Apply parameter increment to current parameter set: P new =P current + delta P , get the updated parameter set P new .
[0106] Get the current number of iterations iter and the preset random perturbation intensity σ. According to the iteration progress, calculate the current perturbation intensity: σ current =σ* (1 - iter / max iter ). Generate random disturbances using Gaussian distribution: noise ~ N(0, σ current ²). Add random perturbations to the updated parameter set: P new =P new+ noise, and get the final new parameter set P final .
[0107] Get the current cost function value J current and the cost function value J of the previous iteration prev Compare J current and J prev , if J current < J prev , then increase the learning rate: η new =η* 1.05; otherwise, reduce the learning rate: η new =η*0.95. Update the learning rate η to η new .
[0108] Get the observed brightness temperature (original brightness temperature data) TB for each current incident angle θ obs (θ) and simulated brightness temperature TB sim (θ). Calculate the brightness temperature difference at each angle: diff(θ) = |TB obs (θ) - TB sim (θ)|. Update the weights using an exponential function: w new (θ) = w(θ) * exp(-α * diff(θ)), where α is the preset adjustment factor. Normalize the new weight to get the updated angle weight w updated (θ).
[0109] Get the current cost function value J current , the cost function value J of the previous iteration prev And the preset convergence threshold ε. Calculate the change in the cost function: delta J = |J current -J prev |. Compare delta J With ε, if delta J < ε, the algorithm is marked as converged and ready to exit the loop; otherwise, it continues to iterate.
[0110] The parameter set P of this iteration final , cost function value J current , velocity vector v, learning rate η new and the angle weight w updated (θ) is stored as the input data for the next iteration. Increase the iteration counter: iter = iter + 1.
[0111] Compare the current number of iterations iter with the maximum number of iterations max iter If iter > max iterOr the algorithm has been marked as converged, the main loop ends; otherwise, return to the beginning to continue the next round of iterations.
[0112] Through this optimization main loop, the algorithm continuously updates soil moisture and surface roughness parameters, while dynamically adjusting the learning strategy to adapt to the complex parameter space.
[0113] This embodiment utilizes gradient information, random search and adaptive strategies to effectively improve the accuracy and robustness of sub-surface soil moisture inversion using a P-band radiometer.
[0114] In another embodiment of the present application, for the cost function, the parameters are initialized: the original brightness temperature data TB of the P-band radiometer at different incident angles (20°-60°) are obtained. obs (θ); Initialize soil moisture SM and roughness parameter set (H RP ,Q R , N RP ); Set the initial learning rate η=0.01, momentum factor β=0.9, random perturbation intensity σ=0.1; define the maximum number of iterations max iter and the convergence threshold ε.
[0115] For the main loop, the current iteration number iter is from 1 to the maximum iteration number max iter In each cycle, the gradient is calculated: for each parameter p (including soil moisture SM and roughness parameter (H RP , Q R , N RP )) Using numerical differentiation method, we calculate the partial derivatives of cost function J with respect to each parameter and get the gradient vector grad p :grad p = ΨJ / Ψp (using numerical differentiation or automatic differentiation), where Ψ represents the partial derivative; for each parameter p, update: v p = β * v p + (1 - β) * grad p ; p = p - η * v p + N(0, σ²) , where v p represents the momentum term of parameter p, β represents the attenuation coefficient of the momentum term, and grad p represents the gradient of the cost function J with respect to the parameter p, η represents the learning rate, N(0, σ²) represents Gaussian noise with a mean of 0 and a variance of σ²; Gaussian noise is added to the parameters. Dynamically adjust the learning rate: If the current cost function value J current <The cost function value J of the previous iteration prev, then increase the learning rate, η=η* 1.05, if the cost decreases, slightly increase the learning rate; otherwise η=η* 0.95, if the cost increases, slightly reduce the learning rate. Adaptive adjustment of angle weights: For each angle θ:w(θ) = w(θ) * exp(-α * |TB obs (θ) - TB sim (θ)|), where exp( ) represents the exponential function, α represents the parameter for adaptively adjusting the angle weight, and w(θ) is normalized so that its sum is 1. Update the random perturbation intensity σ = σ * (1 - iter / max iter ), and gradually reduce the disturbance with iteration. Check convergence. If the difference between the current cost function value and the cost function value of the previous iteration is less than the preset convergence threshold ε, stop the iteration. Output the optimal soil moisture SM and roughness parameters.
[0116] In another embodiment of the present application, the soil of 0-20 cm is divided into 4 layers, each layer is 5 cm. We assign an initial soil moisture value to each layer, for example: Layer 1 (0-5 cm): 0.20 m 3 / m 3 ; Layer 2 (5-10cm): 0.22 m 3 / m 3 ; Layer 3 (10-15cm): 0.24 m 3 / m 3 ; Layer 4 (15-20cm): 0.26 m 3 / m 3 . Using these initial values and the forward model (including the radiation transfer model and the dielectric model) to calculate the theoretical brightness temperature, the theoretical brightness temperature obtained is 290K. Compare the observed value and the model value: the actual observed brightness temperature is 285K, and the calculated difference is: 285K - 290K = -5K. Construct a cost function: Define a cost function J, for example: J = (TB obs -TB model )² + λ * R(SM), where R(SM) is a regularization term used to maintain the smoothness of the soil moisture distribution. Calculate the Jacobian matrix and calculate the sensitivity of brightness temperature to soil moisture in each layer. For example: ΨTB / ΨSM1 = -50 K / (m 3 / m 3 );ΨTB / ΨSM2 = -40 K / (m 3 / m 3 );ΨTB / ΨSM3 = -30 K / (m 3 / m 3 );ΨTB / ΨSM4= -20 K / (m 3 / m 3). Use an optimization algorithm (such as the Levenberg-Marquardt algorithm) to calculate the updated amount of soil moisture. For example: ΔSM1 = 0.02 m 3 / m 3 ; ΔSM2 = 0.015 m 3 / m 3 ; ΔSM3 = 0.01 m 3 / m 3 ; ΔSM4 = 0.005 m 3 / m 3 Applying the update, the new soil moisture estimate is: Layer 1: 0.20 + 0.02 = 0.22 m 3 / m 3 ; Layer 2: 0.22 +0.015 = 0.235 m 3 / m 3 ; Layer 3: 0.24 + 0.01 = 0.25 m 3 / m 3 ; Layer 4: 0.26 + 0.005 = 0.265 m 3 / m 3 The theoretical brightness temperature is recalculated using the new soil moisture estimate and compared to the observed value, and this process is repeated until the termination condition is met (such as the number of iterations reaches the upper limit, or the cost function changes very little). When the iteration ends, the final soil moisture estimate is the inversion result.
[0117] The forward model of this embodiment is used to calculate the theoretical brightness temperature under a given soil moisture distribution. The optimized parameterization scheme is reflected in the soil stratification method, the choice of regularization terms, and possible constraints (such as the reasonable range of soil moisture). The core of this embodiment is to make the brightness temperature predicted by the model as close as possible to the observed brightness temperature by continuously adjusting the soil moisture estimation, while considering the physical rationality of the soil moisture distribution.
[0118] In another embodiment of the present application, the soil effective temperature is calculated:
[0119] T eff = T deep + (T surf -T deep )·[1 - exp(-∑ i=1 n W i ·B i )] + ∑ i=1 n T i ·W i ·(1-exp(-B i ));
[0120] Where T deep is the deep soil temperature, T surf is the surface soil temperature, T i is the soil temperature of the i-th layer, B i =Δx i ·(4π / λ)·(ε i '' / 2sqrt(ε i ')), Δx i is the thickness of the i-th layer, λ is the wavelength of the microwave radiation, ε i '' and ε i ' are the imaginary and real parts of the dielectric constant of the i-th layer of soil, W i is a newly introduced weight factor, which indicates the effect of the i-th layer on the overall T eff The contribution of i=1 n W i = 1.
[0121] The cost function is: J = Σ(TB obs -TB sim )²+λ1Σ(SM - SM prior )²+λ2Σ(H RP -H RPprior )²+λ3Σ(Q R -Q Rprior )²+λ4Σ(N RP -N RPprior )²;
[0122] TB sim is the simulated brightness temperature, SM prior , H RPprior , Q Rprior and N RPprior is the prior information, λ1, λ2, λ3 and λ4 are weight coefficients.
[0123] In another embodiment of the present application, step S5 may also be:
[0124] S5a, based on the initial value of soil moisture, construct an initial population, assign an initial strategy usage probability vector to each individual in the initial population, and the initial strategy usage probability vector includes at least the SH0 algorithm, the WOA algorithm and the DE algorithm;
[0125] S5b, substitute the simulated brightness temperature into the cost function and calculate the fitness;
[0126] S5c, based on fitness, the population is divided into elite layer, exploration layer and diversity layer, and the initial strategy is called to use the probability vector, specifically: the improved SHO algorithm is applied to the elite layer, the soil moisture and roughness parameters are adjusted, and the local search solution set is generated; the improved WOA algorithm is applied to the exploration layer, the soil moisture and roughness parameters are adjusted in a larger range, and the exploration solution set is generated; the improved DE algorithm is applied to the diversity layer, and the soil moisture and roughness parameters are adjusted at the same time through mutation, crossover and selection operations to generate a diversity solution set;
[0127] S5d, merging the local search solution set, the exploration solution set and the diversity solution set into a new population, repeating step S5b for each individual in the new population, and calculating a new fitness value; comparing the fitness of the new population with that of the original population, and selecting better individuals to form an updated population;
[0128] S5e, calculate the improvement ΔI of each initial strategy in reducing the cost function value; use the softmax function to update the initial strategy usage probability vector SP new = softmax(SP + α*ΔI); update the probability vector SP and learning rate α; according to the changing trend of the cost function during the update process, adaptively adjust the key parameters of the SHO, WOA and DE algorithms, as well as the weight coefficients in the cost function;
[0129] S5f, storing the roughness parameter corresponding to the optimal soil moisture solution in the updated population into the memory bank; if the memory bank is full, replacing the record with the largest cost function value;
[0130] Calculate the standard deviation of the soil moisture values in the updated population. If the standard deviation is less than a preset threshold or reaches the maximum number of iterations, the algorithm terminates; otherwise, return to step S5c to continue iteration;
[0131] The solution with the smallest cost function value is selected from the memory library as the final soil moisture inversion result; at the same time, the corresponding roughness parameter, strategy usage probability and key parameter settings are output.
[0132] According to another aspect of the present application, a subsurface soil moisture inversion system based on a P-band radiometer comprises:
[0133] at least one processor; and,
[0134] a memory communicatively connected to at least one of the processors; wherein,
[0135] The memory stores instructions that can be executed by the processor, and the instructions are used to be executed by the processor to implement the sub-surface soil moisture inversion method based on P-band radiometer described in any of the above embodiments.
[0136] The preferred embodiments of the present invention are described in detail above; however, the present invention is not limited to the specific details in the above embodiments. Within the technical concept of the present invention, various equivalent transformations can be made to the technical solutions of the present invention, and these equivalent transformations all belong to the protection scope of the present invention.
Claims
1. The subsurface soil moisture inversion method based on P-band radiometer is characterized by: The steps include: S1, collect observation data and preprocess them to build initial parameters for the preconfigured dielectric hybrid model; The observation data include original brightness temperature data, soil temperature data and initial roughness parameter set from P-band radiometer at different incident angles; S2. Based on the initial parameters of the dielectric hybrid model and the preset initial value of soil moisture, the relative dielectric constant of the soil is calculated, and based on the relative dielectric constant and the initial roughness parameter set, the soil emissivity at each incident angle is calculated; S3, based on soil temperature data and relative dielectric constant, call the preconfigured Choudhury module to calculate the effective soil temperature; S4, based on soil emissivity and effective temperature, the simulated brightness temperature is calculated using the zero-order approximation model of the radiation transfer equation; Calculate the brightness temperature difference between the original brightness temperature data and the simulated brightness temperature, and construct a cost function based on multi-angle data based on the brightness temperature difference; S5. Based on the cost function, the optimization algorithm is called to update the soil moisture and roughness parameters to obtain the optimized parameters, and the weight of each angle data in the cost function is dynamically adjusted according to the performance of different angles; S6, based on the optimized parameters, repeating steps S2 to S5, iteratively performing soil moisture inversion calculation until a preset convergence condition is met, and obtaining an inverted soil moisture result; In step S1, basic soil data is collected and preprocessed as follows: S11, collecting observation data, including original brightness temperature data, soil temperature data and initial roughness parameter set at different incident angles; S12, preprocessing the raw brightness temperature data at each angle, including removing outliers, unit conversion, and normalization or standardization, to obtain preprocessed raw brightness temperature data; S13, performing multi-angle adaptive fusion based on the pre-processed original brightness temperature data to obtain adaptively fused brightness temperature data; The brightness temperature data after adaptive fusion is: Adapted TB (t)=∑ i=1 n α i (t)·Preprocessed TBθi (t)+γ∑ i=1 n (1-α i (t))·Temporal Filtered _ TBθi (t); Among them, γ is a smoothing factor, which is used to balance the impact of the current time and the previous and next time data. TB Represents the brightness temperature data after adaptive fusion, Preprocessed TBθi (t) represents the incident angle θ i and the preprocessed brightness temperature data at time t, Temporal Filtered _ TBθi (t) represents the incident angle θ i and the brightness temperature data after time filtering at time t, n represents the number of angles, α i (t) represents the weight of brightness temperature data; Step S3 is further as follows: S31. Calculating a temperature gradient based on soil temperature data; wherein the soil temperature data includes surface temperature and deep soil temperature; S32, based on the temperature gradient and relative dielectric constant, calling the pre-configured Lv model to calculate the CT parameters; S33, based on the soil temperature data and CT parameters, calling the pre-configured Choudhury module to calculate the effective soil temperature; The preconfigured Choudhury modules are: T eff s =T depth +CT (T surf -T depth ); Among them, T surf and T depth It is the physical temperature of the surface and deep soil; CT parameter determines the surface and deep temperatures to the effective soil temperature T eff s Contribution ratio; , where exp( ) represents the exponential function, Δx represents the thickness of the soil layer, and λ0 represents the wavelength of the electromagnetic wave; Step S4 is further as follows: S41, based on soil emissivity and effective temperature, the simulated brightness temperature is calculated using the radiation transfer equation; S42, calculating a brightness temperature difference value based on the original brightness temperature data and the simulated brightness temperature; S43, obtaining current optimization parameters, and calculating parameter variation based on the current optimization parameters and the initial roughness parameter set; S44, constructing a cost function based on multi-angle data based on brightness temperature difference values and parameter changes; The cost function is: J(t) =∑ i=1 n α i (t) (TB θi , obs (t) - TB θi,sim (t))² + λ∑ j=1 m (P j (t) -P j,prior (t))²; TB θi , obs (t) represents the time t and the incident angle θ i The original brightness temperature data under TB θi,sim (t) represents the time t and the incident angle θ i The simulated brightness temperature under j (t) represents the estimated value of the jth parameter at time t, P j,prior (t) represents the prior value of the jth parameter at time t.
2. The subsurface soil moisture inversion method based on P-band radiometer according to claim 1 is characterized in that: Step S2 is further as follows: S21, reading the initial parameters of the dielectric hybrid model, combining the initial value of soil moisture, and calculating the relative dielectric constant of the soil; S22, combining relative dielectric constant and incident angle, using Fresnel equation to calculate ideal surface reflectivity; S23, reading the initial roughness parameter set, combining with the ideal surface reflectivity, calling the pre-configured HQN model to calculate the actual reflectivity taking into account the roughness effect; S24. Based on the actual reflectivity, Kirchhoff's law is used to calculate the soil emissivity at each angle.
3. The subsurface soil moisture inversion method based on P-band radiometer according to claim 2 is characterized in that: Step S5 is further as follows: S51, based on the cost function, using an optimization algorithm to calculate the parameter update amount; S52, based on the parameter update amount, updating the soil moisture and roughness parameters to obtain the optimized parameters, and in the optimization process, using an adaptive weight adjustment algorithm to dynamically adjust the weight of each angle in the cost function according to the performance of different angles; Dynamically adjust the weights as follows: α i,new (t) = exp(-J θi (t)) / ∑ k=1 n exp(-J θk (t)); Where exp( ) represents the exponential function, J represents the cost function, and θ i and θ k Represents different incident angles.
4. The subsurface soil moisture inversion method based on P-band radiometer according to claim 3 is characterized in that: Step S6 is further as follows: S61, based on the optimized parameters, repeating steps S2 to S5, performing soil moisture inversion calculation, and obtaining an inverted soil moisture result; until the cost function converges or reaches the maximum number of iterations, outputting the final inverted soil moisture result; S62, reading the soil moisture observation data, comparing the final inverted soil moisture result with the soil moisture observation data, calculating the error statistics, and evaluating the inversion accuracy based on the error statistics.
5. The subsurface soil moisture inversion method based on P-band radiometer according to claim 4 is characterized in that: Also includes: S71, obtaining a parameter set of the current iteration, including soil moisture, roughness parameters, and dielectric hybrid model parameters, and based on the parameter set of the current iteration, using a numerical differentiation method, calculating a partial derivative of the cost function with respect to each parameter to obtain a gradient vector; S72, reading the velocity vector of the previous iteration, the current learning rate and the momentum factor, and updating the velocity vector using the momentum gradient descent method; and calculating the parameter increment based on the updated velocity vector; Apply the parameter increment to the parameter set of the current iteration to obtain an updated parameter set; S73, obtaining the current number of iterations and the preset random perturbation strength, calculating the current perturbation strength based on the current number of iterations and the random perturbation strength, and generating random perturbations using Gaussian distribution; adding the random perturbations to the updated parameter set to obtain a final new parameter set; S74, obtaining the current cost function value and the cost function value of the previous iteration, and comparing them, if the current cost function value is less than the cost function value of the previous iteration, increasing the learning rate; otherwise, reducing the learning rate; outputting the updated learning rate; S75, obtaining original brightness temperature data and simulated brightness temperature at each current incident angle, and calculating the current brightness temperature difference at each angle; Based on the current brightness temperature difference, the angle weight is updated using an exponential function, and the angle weight is normalized to obtain the updated angle weight; S76, calculating a change in the cost function based on the current cost function value, the cost function value of the previous iteration, and a preset convergence threshold; If the change in the cost function is less than the convergence threshold, the marking algorithm has converged and is ready to exit the loop; otherwise, continue iterating; S77. Store the new parameter set, the current cost function value, the updated velocity vector, the updated learning rate and the updated angle weight as input data for the next iteration until a preset number of iterations is reached.
6. A subsurface soil moisture inversion system based on a P-band radiometer, characterized in that: include: at least one processor; as well as, a memory communicatively connected to at least one of the processors; wherein, The memory stores instructions that can be executed by the processor, and the instructions are used to be executed by the processor to implement the sub-surface soil moisture inversion method based on a P-band radiometer as described in any one of claims 1 to 5.
Citation Information
Patent Citations
Method for assimilating remote sensing data of soil humidity in watershed scale
CN102034027A
Method for detecting radio frequency interference of multi-source passive microwave data
CN108736990A