A method and system for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model
By using the least squares-based GRACE atmospheric demixing model, and by separating parameters using inverse fast Fourier transform and trigonometric orthogonality, combined with Kronecker sign normalization, the residual error problem of the atmospheric and oceanic demixing model in GRACE satellite data was solved, achieving efficient and accurate gravity field inversion.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- SOUTHWEST JIAOTONG UNIV
- Filing Date
- 2026-03-13
- Publication Date
- 2026-05-26
Smart Images
Figure CN121835318B_ABST
Abstract
Description
Technical Field
[0001] This application relates to the field of atmospheric demixing technology, and more specifically, to a method and system for modeling the GRACE atmospheric demixing model based on least squares. Background Technology
[0002] The GRACE (Gravity Recovery and Climate Experiment) satellite system, as a crucial tool for monitoring changes in Earth's gravity field, demonstrates unique advantages in revealing dramatic geological activity and long-term material migration. Currently, GRACE and its successor, GRACE-FO (GRACE Follow-On), provide key data for studying global material migration processes such as earthquake deformation, ice sheet melting, and changes in terrestrial water storage by monitoring Earth's time-varying gravity field. However, satellite observation data contains superimposed high-frequency non-tidal mass variation signals from the atmosphere and ocean, with frequencies far exceeding sampling capabilities. Without correction, this undersampling will cause severe mixing effects, interfering with the accuracy of gravity field inversion.
[0003] Therefore, current data processing workflows generally employ Atmospheric and Oceanic Demixing (AOD) models to subtract such high-frequency signals. However, research indicates that residual errors inherent in AOD products have become one of the main factors limiting the accuracy of the GRACE / GRACE-FO time-varying gravity field, resulting in calculation results that have not yet fully met the expected accuracy of the mission design. These residual errors mainly stem from two aspects:
[0004] Atmospheric demixing error mainly stems from significant uncertainties in the input data of atmospheric models, especially in high-latitude regions where surface pressure observation errors are more pronounced.
[0005] High-frequency errors can alias into the gravity field solution through nonlinear mechanisms; while low-frequency errors in models (such as ECMWF Operational Analysis) such as monthly mean surface pressure field will gradually accumulate, affecting the estimation of long-term processes such as ice sheet mass changes.
[0006] To improve the accuracy and reliability of AOD products, researchers have developed various methods, such as building high-precision models based on datasets like ERA5 and CRA-40, or attempting to consider model uncertainties within the traditional least squares framework. However, existing methods still have the following significant drawbacks:
[0007] 1. Except for a few least squares-based methods, conventional numerical methods have difficulty in incorporating and handling the significant uncertainties of the input data itself during the solution process.
[0008] 2. Traditional least squares methods suffer from low computational efficiency due to their equation construction characteristics, and severely limit the maximum solvable spherical harmonic order.
[0009] 3. The non-tidal demixing error information currently provided is mostly an approximate simulation of the actual error, rather than a rigorous uncertainty estimate output synchronously during the solution process. Therefore, its application in subsequent gravity field inversion lacks a rigorous theoretical basis. Summary of the Invention
[0010] The purpose of this application is to provide a method and system for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model, which can improve computational efficiency while ensuring the computational accuracy of the atmospheric demixing model.
[0011] This application is implemented as follows:
[0012] In a first aspect, this application provides a method for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model, comprising the following steps:
[0013] S1, based on atmospheric vertical integral anomalies For the observed values, with the associated Legendre function Earth's radius Average density and load Love number Given the quantities, construct a system based on the spherical harmonic potential coefficients. and The observation equation has unknown parameters;
[0014] S2, Observations at each fixed latitude Performing an inverse fast Fourier transform along the longitude direction yields the Fourier coefficients of that latitude with respect to longitude. and ;
[0015] S3. Utilizing the orthogonality of trigonometric functions, integrate the observation equation from step S1 along the longitude direction to combine the coupled unknown parameters. and In a specific order The following separation is performed to obtain a result containing only a single unknown parameter. or The simplified equation;
[0016] S4. The simplified equation obtained in step S3 is discretized and approximated based on the discrete observation data of the same latitude and longitude grid. The orthogonality and properties of discrete trigonometric functions are used, combined with Kronecker symbols, to normalize the discretized equation and form a standard observation equation form suitable for indirect adjustment.
[0017] S5. Based on the standard observation equations obtained in step S4, the specific order is obtained using indirect adjustment. The spherical harmonic coefficients are then used to form the complete spherical harmonic coefficients. and .
[0018] Based on the first aspect, the observation equation for step S1 is:
[0019]
[0020] in, For the remaining latitude longitude place atmospheric vertical integral of order Excluding long-term average The result after that, For order, and These are the spherical harmonic potential coefficients to be determined. It is a fully normalized associated Legendre function. For the Earth's radius, The average density of the Earth, This is the load Love number.
[0021] Based on the first aspect, the Fourier coefficients in step S2 and The calculation formula is:
[0022]
[0023]
[0024] in, and It consists of the real and imaginary parts after the dimension-wise inverse Fast Fourier Transform.
[0025] Based on the first aspect, the specific steps of step S3 include:
[0026] Separating the unknown parameters using the orthogonality of trigonometric functions is achieved by multiplying both sides of equation (1) by a specific degree. Below Integrating, we have:
[0027]
[0028]
[0029] In the summation part on the right side of the equation, we have:
[0030]
[0031]
[0032]
[0033]
[0034] In a specific order The following equation simplifies to:
[0035]
[0036] .
[0037] Based on the first aspect, the specific steps of step S4 include:
[0038] For continuous functions Discrete processing uses equidistant points Approximating, we have:
[0039]
[0040] in, This represents the number of sampling points along the longitude direction.
[0041] Based on the cosine term, equation (10) can be written as:
[0042]
[0043] The discrete orthogonal sum on the right side of the equation is:
[0044]
[0045] By introducing Kronecker notation for normalization, we finally obtain:
[0046]
[0047]
[0048] in, If the left side of the equation is considered as the observed value, and the right side is the product of the coefficient and the unknown parameter, then equations (15) and (16) can be solved directly using the indirect adjustment principle; presented in matrix form:
[0049]
[0050]
[0051] in, and These are the real and imaginary parts of the original observations after performing an inverse Fast Fourier transform, and are related to a specific order. Relevant quantities:
[0052]
[0053] .
[0054] Based on the first aspect, the specific steps of step S5 include:
[0055] The process involves iterating over order l to construct the corresponding observation equations and form the normal equations for the current order l.
[0056] In the process of solving the normal equation, for a specific order Perform successive iterations, construct the equations shown in equations (17) and (18) respectively, and solve them using the indirect adjustment method to obtain the spherical harmonic coefficient fragments corresponding to the current order, and combine the spherical harmonic coefficient fragments of each order to form the complete spherical harmonic coefficients corresponding to the current order l;
[0057] Extract spherical harmonic coefficient fragments of order l from complete spherical harmonic coefficients;
[0058] Traverse the order l from 0 to the maximum truncation order, and combine the extracted spherical harmonic coefficient fragments of each order to obtain the complete set of spherical harmonic coefficients.
[0059] Secondly, this application provides a system for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model, comprising:
[0060] Module for constructing observation equations: using atmospheric vertical integral anomalies For the observed values, with the associated Legendre function Earth's radius Average density and load Love number Given the quantities, construct a system based on the spherical harmonic potential coefficients. and The observation equation has unknown parameters;
[0061] Longitude Fourier Transform Module: Observations at each fixed latitude. Performing an inverse fast Fourier transform along the longitude direction yields the Fourier coefficients of that latitude with respect to longitude. and ;
[0062] Unknown parameter separation module: Utilizing the orthogonality of trigonometric functions, the observation equation from step S1 is integrated along the longitude direction to separate the coupled unknown parameters. and In a specific order The following separation is performed to obtain a result containing only a single unknown parameter. or The simplified equation;
[0063] Equation Discretization and Normalization Module: The simplified equation obtained in step S3 is discretized and approximated based on discrete observation data of equal latitude and longitude grid. The orthogonality and properties of discrete trigonometric functions are used in combination with Kronecker symbols to normalize the discretized equation, forming a standard observation equation form suitable for indirect adjustment.
[0064] Spherical Harmonic Coefficient Output Module: Based on the standard observation equation obtained in step S4, a specific order is obtained using indirect adjustment. The spherical harmonic coefficients are then used to form the complete spherical harmonic coefficients. and .
[0065] Thirdly, this application provides an electronic device, comprising:
[0066] Memory, used to store one or more programs;
[0067] processor;
[0068] The above method is implemented when one or more programs are executed by the processor.
[0069] Fourthly, this application provides a computer-readable storage medium having a computer program stored thereon, which, when executed by a processor, implements the above-described method.
[0070] Compared with the prior art, this application has at least the following advantages or beneficial effects:
[0071] This application provides a method and system for modeling the GRACE atmospheric demixing model based on least squares. It surpasses the conventional numerical integration method in terms of computational efficiency, and is significantly superior to the conventional numerical integration method in terms of statistical optimality of output results and the error information that can be attached. It also completely overcomes the defect that the traditional least squares method cannot be used for high-order large-scale solutions. Attached Figure Description
[0072] To more clearly illustrate the technical solutions of the embodiments of this application, the accompanying drawings used in the embodiments will be briefly introduced below. It should be understood that the following drawings only show some embodiments of this application and should not be regarded as a limitation of the scope. For those skilled in the art, other related drawings can be obtained based on these drawings without creative effort.
[0073] Figure 1 The flowchart below illustrates a method for modeling an atmospheric demixing model based on the least squares GRACE atmospheric demixing model according to this application.
[0074] Figure 2A comparison of the order of geoid height differences using different demixing methods;
[0075] Figure 3 The graph shows a comparison of the residual power spectral density of KBRR after processing with different demixing methods.
[0076] Figure 4 This is a schematic diagram of the system structure for modeling a GRACE atmospheric demixing model based on least squares, as described in this application.
[0077] Figure 5 This is a schematic diagram of the structure of an electronic device according to this application.
[0078] icon:
[0079] 1. Module for constructing observation equations; 2. Module for Fourier transform of longitude direction; 3. Module for separating unknown parameters; 4. Module for equation discretization and normalization; 5. Module for outputting spherical harmonic coefficients; 6. Processor; 7. Memory; 8. Communication interface. Detailed Implementation
[0080] To make the objectives, technical solutions, and advantages of the embodiments of this application clearer, the technical solutions of the embodiments of this application will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of this application, and not all embodiments. The components of the embodiments of this application described and shown in the accompanying drawings can generally be arranged and designed in various different configurations.
[0081] The following detailed description of some embodiments of this application is provided in conjunction with the accompanying drawings. Unless otherwise specified, the various embodiments and features described below can be combined with each other.
[0082] Example
[0083] This application provides a method and system for modeling the GRACE atmospheric demixing model based on least squares, which can improve computational efficiency while ensuring the computational accuracy of the atmospheric demixing model.
[0084] This method is essentially a spherical harmonic analysis method used in the calculation of an atmospheric demixing model. The calculation process of the atmospheric demixing model can be generally divided into two main stages: the atmospheric vertical integration process and the spherical harmonic analysis process. The method in this application only involves the spherical harmonic analysis process, which will be analyzed in detail below.
[0085] Please refer to Figure 1 The method for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model includes the following steps:
[0086] S1, based on atmospheric vertical integral anomalies For the observed values, with the associated Legendre function Earth's radius Average density and load Love number Given the quantities, construct a system based on the spherical harmonic potential coefficients. and The observation equation has unknown parameters;
[0087] Specifically, the physical observation of atmospheric vertical integral anomaly is linked to the spherical harmonic coefficients to be determined through rigorous gravitational field theory (spherical harmonic synthesis formula), providing a physically correct and mathematically solvable model foundation for the entire inversion.
[0088] As one implementation method, the observation equation for step S1 is:
[0089]
[0090] in, For the remaining latitude longitude place atmospheric vertical integral of order Excluding long-term average The result after that, For order, and These are the spherical harmonic potential coefficients to be determined. It is a fully normalized associated Legendre function. For the Earth's radius, The average density of the Earth, This is the load Love number.
[0091] S2, Observations at each fixed latitude Performing an inverse fast Fourier transform along the longitude direction yields the Fourier coefficients of that latitude with respect to longitude. and ;
[0092] Specifically, traditional least squares methods require directly processing massive grid data to construct giant normal equations, resulting in extremely high computational complexity. This step, by performing FFT on the data for each latitude circle, transforms the convolution / correlation operations along the longitude direction into efficient spectral operations, reducing computational complexity and improving computational efficiency.
[0093] As one implementation method, the Fourier coefficients in step S2 and The calculation formula is:
[0094]
[0095]
[0096] in, and It consists of the real and imaginary parts after the dimension-wise inverse Fast Fourier Transform.
[0097] S3. Utilizing the orthogonality of trigonometric functions, integrate the observation equation from step S1 along the longitude direction to combine the coupled unknown parameters. and In a specific order The following separation is performed to obtain a result containing only a single unknown parameter. or The simplified equation;
[0098] Specifically, by utilizing the orthogonality of trigonometric functions, we can cleverly incorporate the information about... and The coupled two-dimensional global problem is decoupled into a series of independent one-dimensional problems of each order (dependent only on the latitude θ). This transforms the large, dense matrix into multiple independently solvable, much smaller strip or diagonal matrices, paving the way for efficient solutions.
[0099] As one implementation method, step S3 includes the following specific steps:
[0100] Separating the unknown parameters using the orthogonality of trigonometric functions is achieved by multiplying both sides of equation (1) by a specific degree. Below Integrating, we have:
[0101]
[0102]
[0103] In the summation part on the right side of the equation, we have:
[0104]
[0105]
[0106]
[0107]
[0108] In a specific order The following equation simplifies to:
[0109]
[0110] .
[0111] S4. The simplified equation obtained in step S3 is discretized and approximated based on the discrete observation data of the same latitude and longitude grid. The orthogonality and properties of discrete trigonometric functions are used, combined with Kronecker symbols, to normalize the discretized equation and form a standard observation equation form suitable for indirect adjustment.
[0112] Specifically, the continuous integral equations are adapted to discrete grid observation data. By introducing Kronecker notation for normalization, the orthogonality of the discretized trigonometric functions is strictly guaranteed, avoiding numerical ill-conditioning caused by discrete sampling and truncation errors, and ensuring the numerical stability of the final solution.
[0113] As one implementation method, step S4 includes the following specific steps:
[0114] For continuous functions Discrete processing uses equidistant points Approximating, we have:
[0115]
[0116] in, This represents the number of sampling points along the longitude direction.
[0117] Based on the cosine term, equation (10) can be written as:
[0118]
[0119] The discrete orthogonal sum on the right side of the equation is:
[0120]
[0121] By introducing Kronecker notation for normalization, we finally obtain:
[0122]
[0123]
[0124] in, If the left side of the equation is considered as the observed value, and the right side is the product of the coefficient and the unknown parameter, then equations (15) and (16) can be solved directly using the indirect adjustment principle; presented in matrix form:
[0125]
[0126]
[0127] in, and These are the real and imaginary parts of the original observations after performing an inverse Fast Fourier transform, and are related to a specific order. Relevant quantities:
[0128]
[0129] .
[0130] S5. Based on the standard observation equations obtained in step S4, the specific order is obtained using indirect adjustment. The spherical harmonic coefficients are then used to form the complete spherical harmonic coefficients. and .
[0131] Specifically, after efficiently organizing the data and constructing the equations in steps S2-S4, this step utilizes the least squares adjustment principle for solution. The purpose of this setup is: 1. To obtain the statistically optimal (minimum variance) estimate of the spherical harmonic coefficients under a given observation model and weights; 2. To simultaneously output the standard deviation or variance-covariance matrix of the coefficient solution, i.e., a quantitative estimate of the uncertainty. This is the core reason why this method outperforms conventional numerical integration methods; it not only calculates quickly and accurately, but also provides information on the accuracy of the calculation.
[0132] As one implementation method, step S5 includes the following specific steps:
[0133] The process involves iterating over order l to construct the corresponding observation equations and form the normal equations for the current order l.
[0134] In the process of solving the normal equation, for a specific order Perform successive iterations, construct the equations shown in equations (17) and (18) respectively, and solve them using the indirect adjustment method to obtain the spherical harmonic coefficient fragments corresponding to the current order, and combine the spherical harmonic coefficient fragments of each order to form the complete spherical harmonic coefficients corresponding to the current order l;
[0135] Extract spherical harmonic coefficient fragments of order l from complete spherical harmonic coefficients;
[0136] Traverse the order l from 0 to the maximum truncation order, and combine the extracted spherical harmonic coefficient fragments of each order to obtain the complete set of spherical harmonic coefficients.
[0137] In summary, considering that the traditional least squares method requires directly constructing and solving giant normal equations, its computational complexity is extremely high, severely limiting its practical application. This invention reduces computational complexity and achieves a breakthrough improvement in efficiency by introducing inverse fast Fourier transform to process longitude direction data. Simultaneously, it utilizes the orthogonality of trigonometric functions for parameter decoupling, decomposing the complex global coupling problem into a series of independent, low-dimensional sub-problems, further simplifying the solution scale. This enables the method of this invention to efficiently and stably solve high-order spherical harmonic models, solving the long-standing computational limitations of traditional methods. Furthermore, considering that conventional numerical integration methods cannot quantitatively assess and incorporate the uncertainty of input data during the solution process, this invention employs the indirect adjustment principle at the end of the efficient computational framework, simultaneously outputting the optimal spherical harmonic coefficient estimate and its corresponding standard error estimate. This provides a quantitative evaluation index for the accuracy of the demixed product, providing crucial prior information for error analysis and weighted processing in downstream gravity field inversion, significantly improving the rigor and reliability of the entire data chain processing. Furthermore, this invention incorporates normalization based on Kronecker notation during the discretization process, rigorously ensuring that the trigonometric functions on the discrete grid maintain their orthogonality. This measure effectively avoids numerical ill-conditioning that may result from discretization sampling and truncation errors.
[0138] To facilitate understanding of the method described in this application, different demixing methods are employed for gravity field inversion. Please refer to [link / reference needed]. Figure 2 , Figure 2 This is a comparison of the order of geoid height differences using different demixing methods. In the figure, the horizontal axis represents the spherical harmonic order; a higher order corresponds to a higher spatial resolution (especially for shorter wavelength high-frequency signals). The vertical axis represents the geoid height difference, which can be understood as the degree of inconsistency between the two products. The lower the value, the closer the results of the two methods are, and the higher the consensus. FFT-LS represents the least squares demixing method of this application, LS represents the traditional least squares demixing method, and NI represents the conventional numerical integration method. Higher-order geoid methods reflect the magnitude of the differences between pairwise products. Figure 2 As can be seen, the LS vs NI curves maintain a very low difference throughout, indicating that the results of the traditional least squares (LS) method and the conventional numerical integration (NI) method are basically consistent, representing the level of existing technology. However, the FFT-LS vs LS and FFT-LS vs NI curves, representing the comparison between the new and old methods, show significantly higher differences across the entire order range, especially in the high-frequency band corresponding to mid-to-high orders (>60 orders), and are also higher than the task baseline. Therefore, it can be concluded that the method in this application detects richer high-frequency true signals that were smoothed out or could not be separated by the old methods, achieving higher accuracy and lower noise in the high-frequency band when applied to gravity field inversion. Please refer to [reference needed]. Figure 3 , Figure 3 The results after processing using different demixing methods Figure 3 Furthermore, it is evident from the power spectral density of the inter-satellite distance variability (KBRR) residuals obtained from gravity field inversion after processing with different demixing methods that, in the high-frequency range (corresponding to the first 15 orbital frequency harmonics), the amplitude of the curve of the proposed method (FFT-LS) is significantly lower than that of the traditional methods (LS and NI). Therefore, the least squares demixing method of this application has a smaller amplitude of the inter-satellite distance variability residual in the high-frequency part when used for gravity field inversion, and the atmospheric demixing model calculated using this method can reduce the gravity field inversion error by suppressing high-frequency noise.
[0139] Please refer to Figure 4 This application also provides a system for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model, comprising:
[0140] Module 1 for constructing observation equations: using atmospheric vertical integral anomalies For the observed values, with the associated Legendre function Earth's radius Average density and load Love number Given the quantities, construct a system based on the spherical harmonic potential coefficients. and The observation equation has unknown parameters;
[0141] Longitude Fourier Transform Module 2: Observations for each fixed latitude Performing an inverse fast Fourier transform along the longitude direction yields the Fourier coefficients of that latitude with respect to longitude. and ;
[0142] Module 3 for separating unknown parameters: Utilizing the orthogonality of trigonometric functions, the observation equation from step S1 is integrated along the longitude direction to separate the coupled unknown parameters. and In a specific order The following separation is performed to obtain a result containing only a single unknown parameter. or The simplified equation;
[0143] Equation Discretization and Normalization Module 4: The simplified equation obtained in step S3 is discretized and approximated based on the discrete observation data of the same latitude and longitude grid. The orthogonality and property of discrete trigonometric functions are used in combination with Kronecker symbols to normalize the discretized equation, forming a standard observation equation form suitable for indirect adjustment.
[0144] Spherical Harmonic Coefficient Output Module 5: Based on the standard observation equation obtained in step S4, a specific order is obtained using indirect adjustment. The spherical harmonic coefficients are then used to form the complete spherical harmonic coefficients. and .
[0145] For a detailed implementation of the system based on the least squares GRACE atmospheric demixing model, please refer to the above-mentioned implementation of the method based on the least squares GRACE atmospheric demixing model. Further details will not be provided here.
[0146] Please refer to Figure 5 This application also provides an electronic device, including:
[0147] Memory 7 is used to store one or more programs;
[0148] Processor 6; Processor 6 is connected to memory 7 via communication interface 8;
[0149] When one or more programs are executed by processor 6, all or some of the above methods are implemented.
[0150] This application also provides a computer-readable storage medium having a computer program stored thereon, which, when executed by a processor 6, implements all or part of the methods described above.
[0151] It will be apparent to those skilled in the art that this application is not limited to the details of the exemplary embodiments described above, and that this application can be implemented in other specific forms without departing from the spirit or essential characteristics of this application. Therefore, the embodiments should be considered illustrative and non-limiting in all respects, and the scope of this application is defined by the appended claims rather than the foregoing description. Thus, all variations falling within the meaning and scope of equivalents of the claims are intended to be included within this application. No reference numerals in the claims should be construed as limiting the scope of the claims.
Claims
1. A method for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model, characterized in that, Includes the following steps: S1, based on atmospheric vertical integral anomalies For the observed values, with the associated Legendre function Earth's radius Average density and load Love number Given the quantities, construct a system based on the spherical harmonic potential coefficients. and The observation equation has unknown parameters; S2, the observed values for each fixed latitude Performing an inverse fast Fourier transform along the longitude direction yields the Fourier coefficients of that latitude with respect to longitude. and ; S3. Utilizing the orthogonality of trigonometric functions, integrate the observation equation from step S1 along the longitude direction to integrate the coupled unknown parameters. and In a specific order The following separation is performed to obtain a result containing only a single unknown parameter. or The simplified equation; S4. The simplified equation obtained in step S3 is discretized and approximated based on the discrete observation data of the same latitude and longitude grid. The orthogonality and properties of discrete trigonometric functions are used, combined with Kronecker symbols, to normalize the discretized equation and form a standard observation equation form suitable for indirect adjustment. S5. Based on the standard observation equations obtained in step S4, the specific order is obtained using indirect adjustment. The spherical harmonic coefficients are then used to form the complete spherical harmonic coefficients. and .
2. The method for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model according to claim 1, characterized in that, The observation equation for step S1 is: in, For the remaining latitude longitude place atmospheric vertical integral of order Excluding long-term average The result after that, For order, and These are the spherical harmonic potential coefficients to be determined. It is a fully normalized associated Legendre function. For the Earth's radius, The average density of the Earth, This is the load Love number.
3. The method for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model according to claim 2, characterized in that, Fourier coefficients in step S2 and The calculation formula is: in, and It consists of the real and imaginary parts after the dimension-wise inverse fast Fourier transform.
4. The method for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model according to claim 3, characterized in that, The specific steps of step S3 include: Separating the unknown parameters using the orthogonality of trigonometric functions is achieved by multiplying both sides of equation (1) by a specific degree. Below Integrating, we have: In the summation part on the right side of the equation, we have: In a specific order The following equation simplifies to: 。 5. The method for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model according to claim 4, characterized in that, The specific steps of step S4 include: For continuous functions Discrete processing uses equidistant points Approximating, we have: in, This represents the number of sampling points along the longitude direction. Based on the cosine term, equation (10) can be written as: The discrete orthogonal sum on the right side of the equation is: By introducing Kronecker notation for normalization, we finally obtain: in, If the left side of the equation is considered as the observed value, and the right side is the product of the coefficient and the unknown parameter, then equations (15) and (16) can be solved directly using the indirect adjustment principle; presented in matrix form: in, and These are the real and imaginary parts of the original observations after performing an inverse Fast Fourier transform, and are related to a specific order. Relevant quantities: 。 6. The method for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model according to claim 5, characterized in that, The specific steps of step S5 include: The process involves iterating over order l to construct the corresponding observation equations and form the normal equations for the current order l. In the process of solving the normal equations, for a specific order Perform successive iterations, construct the equations shown in equations (17) and (18) respectively, and solve them using the indirect adjustment method to obtain the spherical harmonic coefficient fragments corresponding to the current order, and combine the spherical harmonic coefficient fragments of each order to form the complete spherical harmonic coefficients corresponding to the current order l; Extract a spherical harmonic coefficient fragment of order l from the complete spherical harmonic coefficients; Traverse the order l from 0 to the maximum truncation order, and combine the extracted spherical harmonic coefficient fragments of each order to obtain the complete set of spherical harmonic coefficients.
7. A system for modeling atmospheric demixing based on the least squares GRACE atmospheric demixing model, characterized in that, include: Module for constructing observation equations: using atmospheric vertical integral anomalies For the observed values, with the associated Legendre function Earth's radius Average density and load Love number Given the quantities, construct a system based on the spherical harmonic potential coefficients. and The observation equation has unknown parameters; Longitude Fourier Transform Module: For the observed values at each fixed latitude Performing an inverse fast Fourier transform along the longitude direction yields the Fourier coefficients of that latitude with respect to longitude. and ; Unknown parameter separation module: Utilizing the orthogonality of trigonometric functions, the observation equation is integrated along the longitude direction to separate the coupled unknown parameters. and In a specific order The following separation is performed to obtain a result containing only a single unknown parameter. or The simplified equation; Equation Discretization and Normalization Module: The simplified equation is discretized and approximated based on discrete observation data of equal latitude and longitude grid, and the orthogonality and property of discrete trigonometric functions are used in combination with Kronecker notation to normalize the discretized equation to form a standard observation equation form suitable for indirect adjustment. Spherical harmonic coefficient output module: Based on the standard observation equation, a specific order is obtained using indirect adjustment. The spherical harmonic coefficients are then used to form the complete spherical harmonic coefficients. and .
8. An electronic device, characterized in that, include: Memory, used to store one or more programs; processor; When the one or more programs are executed by the processor, the method as described in any one of claims 1-6 is implemented.
9. A computer-readable storage medium having a computer program stored thereon, characterized in that, When the computer program is executed by a processor, it implements the method as described in any one of claims 1-6.