A mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT

Abstract

The present invention discloses a mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT. The method comprises the following steps: first, meteorological data is obtained as input to a WRF model, and the WRF model outputs the meteorological field of the study area; then, CO2 concentration is monitored by monitoring stations, and the CO2 concentration is sorted according to time to obtain an observation sequence; the observation sequence is subjected to outlier removal and window sliding average processing to obtain a CO2 elevated concentration sequence; then, the obtained meteorological field is used as input to the Flexwrf model and the STILT model to obtain a footprint field corresponding to the processed observation sequence; a priori total carbon flux is obtained, and a better posterior total carbon flux is obtained by Bayesian estimation in combination with the footprint field. The method of the present invention does not require manual intervention, that is, the CO2 flux of the surrounding area can be directly solved based on the observation values of the monitoring station, and the method has timeliness and higher inversion accuracy.

Description

Technical Field

[0001] The present invention belongs to the field of atmospheric technology, and in particular relates to a mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT. Background Art

[0002] Greenhouse gases, primarily carbon dioxide (CO2), are considered the most influential gas in global climate change. Most climatologists believe that human emissions of CO2 have led to a sharp rise in atmospheric CO2 concentrations over the past 200 years. With the accelerating pace of global industrialization, the continued increase in anthropogenic carbon emissions has had significant impacts on the global climate, economy, and ecology. This has garnered widespread attention from governments, research institutions, and the public worldwide, with the imposition of taxes on countries and businesses based on carbon emissions becoming a growing trend. Changes in atmospheric CO2 concentrations reflect the exchange of CO2 between the earth and the atmosphere, including CO2 exchange between vegetation and the atmosphere, as well as CO2 emissions from anthropogenic fossil fuel combustion. Cities, as major greenhouse gas emitters, contribute 70% of global anthropogenic CO2 emissions despite occupying only 2% of the land area. Therefore, quantitatively studying CO2 fluxes in urban areas and assessing their impact on atmospheric CO2 concentrations is crucial.

[0003] At present, there are two methods for calculating CO2 flux: one is the inventory method based on statistical data (such as energy consumption data), and the other is the inversion method based on monitoring data. The CO2 flux obtained by the inventory method is limited by the availability of statistical data, and can only obtain annual results or results for large areas (above the city level). The existing inversion rules are limited by the relatively extensive meteorological data (1°×1° or 0.5°×0.5°) and are difficult to carry out at the mesoscale. The existing inversion results often have a coarse spatial resolution. Therefore, there is an urgent need for an inversion method that can be carried out at the mesoscale to update the high temporal and spatial resolution CO2 flux in areas below the city level. Summary of the Invention

[0004] In view of the shortcomings of the existing technology, the present invention provides a mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT, which can quickly calculate CO2 flux results with high accuracy (3km or even 1km scale).

[0005] To achieve the above objectives, the technical solution proposed in the present invention is a mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT, including a WRF model, a FLEXPART model and a STILT model, comprising the following steps:

[0006] Step 1: Obtain meteorological data as input to the WRF model. The WRF model outputs a meteorological field with high temporal and spatial resolution near the study area. The spatial resolution of the meteorological field is related to the accuracy of CO2 flux inversion, and the temporal resolution of the meteorological field is related to the FLEXPART model and the STILT model.

[0007] Step 2: Use monitoring stations to monitor CO2 concentrations and sort the CO2 concentrations by time to obtain an observation sequence; remove outliers from the observation sequence and perform window sliding average processing to obtain a CO2 elevated concentration sequence; the site selection requirements for the monitoring station are: in principle, it should be as open as possible to avoid emissions from local anthropogenic and natural sources, and it is recommended to be set up in an open area downwind of the study area.

[0008] Step 3: Use the meteorological field obtained in step 1 as the input of the atmospheric transport models FLEXPART and STILT. According to the operating time and geographical location of the monitoring station, the footprint field corresponding to the observation sequence processed in step 2 is obtained based on the output of the atmospheric transport models FLEXPART and STILT.

[0009] Step 4: Obtain the prior total carbon flux and the CO2 elevation concentration series, and use Bayesian estimation in combination with the footprint field to obtain a better posterior total carbon flux.

[0010] In the above technical solution, further, step 1 is specifically as follows: first, the surface parameters are interpolated to select the study area, then the required meteorological parameters are extracted from the meteorological data, then the meteorological parameters are interpolated to the study area, and finally the WRF model is run to obtain the meteorological field of the study area.

[0011] Furthermore, in step 2, for the observation sequence within a window (recommended to be no less than one month), the non-outlier range is determined by calculating the IQR (interquartile range) within the window, thereby screening the non-outlier values in the observation sequence. The non-outlier range is as follows:

[0012] (Q1-1.5*IQR, Q3+1.5*IQR)

[0013] Where Q1 represents the lower quartile of the observation sequence, Q3 represents the upper quartile of the observation sequence, and IQR represents the difference between Q3 and Q1.

[0014] Observations that are not greater than the lower ten percentile (the values of the observation sequence are arranged from small to large, and the top 10%) are selected from the observation sequence that has been filtered for outliers. Fourteen days are used as a sliding window, and the lower ten percentile observations are taken as a sliding average to obtain the baseline sequence. The elevated CO2 concentration sequence is the difference between the observation sequence and the baseline sequence. Values less than or equal to 0 in the elevated concentration sequence are eliminated, and only values greater than 0 in the elevated concentration sequence are selected to enter the next step of inversion.

[0015] Furthermore, the implementation method for generating the footprint field in step 3 is to perform corresponding screening on the observation values of a certain monitoring station, and only use the footprints of the observed concentrations corresponding to non-outliers and positive uplift concentrations to enter the next inversion step.

[0016] Furthermore, the implementation of the posterior total carbon flux in step 4 is as follows:

[0017]

[0018] Where x is the a posteriori carbon flux to be solved; x a is the prior carbon flux; S a is the covariance matrix of the prior carbon flux; y obs is the rising concentration series of CO2; S o is the covariance matrix of the CO2 elevation concentration series; H is the footprint field obtained by the backward simulation of the atmospheric transport model in step 3. The cost function is minimized to obtain:

[0019] x=x a +S a H T (HS a H T +S o ) -1 (y obs -Hx a )

[0020] Furthermore, when using Bayesian estimation to invert CO2 flux in step 4, the covariance matrix S of the prior carbon flux should be considered. a The covariance matrix S of the CO2 concentration series o The off-diagonal elements of S correspond to spatial correlation and temporal correlation, respectively. a and S o The off-diagonal elements of are calculated as

[0021]

[0022] where σ ij represents the off-diagonal element in row i and column j of the covariance matrix, σ ii and σ jjrepresents the diagonal elements of the covariance matrix, L represents the temporal or spatial correlation coefficient, and d(i, j) represents the temporal or spatial distance.

[0023] On the other hand, the present invention provides a device for implementing the mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT, which includes the following modules:

[0024] The first module is used to obtain meteorological data as input to the WRF model, which outputs the meteorological field of the study area;

[0025] The second module is used to sort the CO2 concentration monitored by the monitoring station by time to obtain an observation sequence, and to remove outliers from the observation sequence and perform window sliding average processing to obtain a CO2 rising concentration sequence;

[0026] The third module is configured to use the meteorological field as input to the atmospheric transport models FLEXPART and STILT, and obtain the footprint field corresponding to the processed observation sequence based on the outputs of the atmospheric transport models FLEXPART and STILT according to the operating time and geographical location of the monitoring station;

[0027] The fourth module is used to obtain a better posterior total carbon flux using Bayesian estimation based on the prior total carbon flux and the CO2 elevation concentration series combined with the footprint field.

[0028] On the other hand, the present invention provides an electronic device comprising one or more processors; a memory for storing one or more programs; and when the one or more programs are executed by the one or more processors, the one or more processors implement the mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT.

[0029] On the other hand, the present invention provides a computer-readable storage medium having computer instructions stored thereon, characterized in that the computer instructions are used to enable a computer to execute the steps of the mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT.

[0030] In the present invention, meteorological data are obtained as the input of the WRF model, and the WRF model obtains the meteorological field with high temporal and spatial resolution near the study area through interpolation; the obtained meteorological field is used as the input of the Flexwrf model and the STILT model to obtain the footprint field corresponding to each observation value; the prior total carbon flux (EDGAR, ODIAC) is obtained, and the Bayesian estimation is used to obtain a better posterior total carbon flux and the spatial distribution of the posterior flux.

[0031] The present invention has the following advantages:

[0032] (1) No human intervention is required, and the CO2 flux in the surrounding area can be directly solved based on the observation values of the monitoring station.

[0033] (2) Compared with the statistics of urban CO2 emission inventory, this method is timely.

[0034] (3) Compared with the CO2 inversion results of coarse grids, this method provides higher inversion accuracy. BRIEF DESCRIPTION OF THE DRAWINGS

[0035] Figure 1 is a flow chart of an embodiment of the present invention;

[0036] Figure 2 are the real observation sequence, the priori simulated concentration, and the posterior simulated concentration of the embodiment of the present invention;

[0037] Figure 3 These are the a priori flux space distribution and the a posteriori flux space distribution of the embodiment of the present invention. DETAILED DESCRIPTION

[0038] The technical solution of the present invention is described in detail below with reference to the accompanying drawings and embodiments.

[0039] In an embodiment of the present invention, anthropogenic carbon emission fluxes can be obtained from the EDGAR or ODIAC emission inventories. Based on the geographical location and operating time of the monitoring stations, the WRF-FLEXPART model and the WRF-STILT model are used to obtain the footprint field, and the posterior flux can be obtained using Bayesian inversion.

[0040] in,

[0041] EDGAR: emissions database for global atmospheric research, global atmospheric research emissions database;

[0042] ODIAC: Open-Data Inventory for Anthropogenic Carbon dioxide, Open Data Inventory for Anthropogenic Carbon dioxide;

[0043] See also Figure 1 The embodiment of the present invention provides a mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT, wherein the monitoring station data is provided by a ground monitoring network. WRF-FLEXPART and WRF-STILT are both Lagrangian transmission models. The implementation process of the CO2 flux inversion method includes the following steps:

[0044] (1) Generate meteorological fields, including obtaining CFSR (Climate Forecast System Reanalysis) data or FNL (Final Operation Global Analysis) data within a time period as the initial meteorological field; and use the WRF (The Weather Research and Forecasting Model) model to interpolate the initial meteorological field in space and time.

[0045] In this embodiment, CFSR meteorological data (or FNL data) are obtained as the input of the WRF model with a time resolution of 6 hours and a spatial resolution of 0.5°. The WRF model can be used to obtain a meteorological field with high temporal and spatial resolution near the study area, with a time resolution of up to 1 hour and a spatial resolution of up to 3 km. The meteorological field obtained by the WRF model will be provided as input to the FLEXPART model and the STILT model.

[0046] In specific implementation, the WRF model can be preferentially set to double-layer nesting of 12km and 3km (or triple-layer nesting of 12km, 3km and 1km), and output NC format meteorological fields (for driving FLEXPART) and ARL format meteorological fields (for driving the STILT model).

[0047] The FLEXPART model and the STILT model are driven by the high-precision meteorological field output by the WRF model. In the embodiment, the WRF model adopts a two-layer nested bidirectional feedback setting, and the spatial resolution of the two-layer nested simulation area is set to 12 and 3 km, and the number of grid points in the east-west and north-south directions is set to 120*120 and 141*141.

[0048] Specifically, before running the WRF model, the following three pre-processing operations are required. These three functions can be completed by the three subroutines of WPS (WRF Pre-Processing System) during implementation.

[0049] geogrid.exe: Its main function is to define the model projection, regional range, and nesting relationship. The implementation example is used to interpolate surface parameters to select the study area.

[0050] ungrib.exe: Its main function is to extract the required meteorological parameters from the meteorological data in grib format (such as CFSR meteorological data, FNL meteorological data, etc.);

[0051] metgrid.exe: The main function is to interpolate meteorological parameters to the study area.

[0052] (2) Using monitoring stations to monitor CO2 concentrations, and sorting the CO2 concentrations by time to obtain an observation sequence. The observation sequence is subjected to outlier screening and window averaging to obtain a CO2 elevation concentration sequence. The specific method of the embodiment is:

[0053] For a time series within a month, the non-outlier range is determined by calculating the IQR (interquartile range) within the window, thereby screening non-outlier values in the time series. The non-outlier range is as follows:

[0054] (Q1-1.5*IQR, Q3+1.5*IQR)

[0055] Where Q1 represents the lower quartile of the observation sequence, Q3 represents the upper quartile of the observation sequence, and IQR represents the difference between Q3 and Q1.

[0056] Observations that are not greater than the lower ten percentile (the values of the observation sequence are arranged from small to large, and the top 10%) are selected from the observation sequence that has been filtered for outliers. Fourteen days are used as a sliding window, and the lower ten percentile observations are taken as a sliding average to obtain the baseline sequence. The elevated CO2 concentration sequence is the difference between the observation sequence and the baseline sequence. Values less than or equal to 0 in the elevated concentration sequence are eliminated, and only values greater than 0 in the elevated concentration sequence are selected to enter the next step of inversion.

[0057] (3) The generated meteorological field is used as the input of the FLEXPART model and the STILT model. According to the operating time and geographical location of the monitoring station, the footprint field corresponding to the observation sequence processed in step (2) is obtained based on the output of the atmospheric transport model FLEXPART and STILT:

[0058] The monitoring station is a ground monitoring network, and the inversion time of the atmospheric transport model can be set to 10 days.

[0059] The FLEXPART model is a Lagrangian transport and dispersion model suitable for simulating large-scale atmospheric transport processes. In addition to transport and turbulent diffusion, it can also simulate dry and wet deposition, decay, and linear chemistry; it can be used in either forward or backward mode. The specific principle is to release particles to simulate the backward motion of gases under various meteorological conditions. The footprint weight is quantified by calculating the total number of particles within a certain altitude in a certain upstream region and the residence time of each particle.

[0060] The STILT model is a transport model based on Lagrangian random walk theory. It links upstream and source (sink) fluxes with concentration changes at the observation point using footprint weights. The specific principle is to simulate the backward motion of gas driven by turbulence and average wind direction by releasing a large number of air particles backward. The footprint weight is quantified by calculating the total number of particles within a certain height in the boundary layer of a certain upstream region and the residence time of each particle.

[0061] The purpose of combining the results of the two models for calculation is to reduce the calculation deviation caused by the results of a single model.

[0062] (4) Obtain the prior total carbon flux (i.e., the prior CO2 flux, mainly anthropogenic carbon flux) and the CO2 elevation concentration sequence, and use Bayesian estimation in combination with the footprint field to obtain a better posterior total carbon flux, that is, to invert the posterior CO2 flux, providing more accurate data for further research on small-scale carbon emissions and other issues.

[0063] In an embodiment, anthropogenic carbon emission fluxes can be obtained from the EDGAR or ODIAC emission inventory. Based on the prior CO2 flux, the prior CO2 flux covariance matrix, the CO2 elevated concentration series, the covariance matrix of the CO2 elevated concentration series, and the footprint field, an optimized and more refined CO2 posterior flux result and spatial distribution can be obtained.

[0064] The specific method for calculating the CO2 posterior flux is:

[0065] x=x a +S a H T (HS a H T +S o ) -1 (y obs -Hx a )

[0066] Where x is the CO2 posterior flux to be solved; x a is the CO2 prior flux; S a is the covariance matrix of CO2 prior flux; y obs is the rising concentration series of CO2; S o is the covariance matrix of the CO2 elevation concentration series; H is the footprint field obtained by the backward simulation of the atmospheric transport model in step (3).

[0067] When using Bayesian estimation to invert CO2 flux, the covariance matrix S of the CO2 prior flux should be considered. a The covariance matrix S of the CO2 concentration series oThe off-diagonal elements correspond to spatial correlation and temporal correlation respectively. The calculation formula of off-diagonal elements is

[0068]

[0069] where σ ij represents the off-diagonal element in row i and column j of the error matrix, σ ii and σ jj represents the diagonal elements of the error matrix, L represents the temporal or spatial correlation coefficient, and d(i, j) represents the temporal or spatial distance.

[0070] Example

[0071] See also Figure 2 This example uses data from the Anji County Monitoring Station in Zhejiang Province (30.64°N, 119.67°E) in June 2023, and uses the EDGAR inventory as the prior anthropogenic carbon flux. Since EDGAR provides the annual average CO2 flux, it is necessary to multiply it by the time scale factor corresponding to June as the prior anthropogenic carbon flux for June. The prior simulated concentration and the posterior simulated concentration of CO2 at this location are calculated through step (4). The concentration simulation of the posterior flux (green line) is closer to the observed time series (blue line) than the concentration simulation of the prior flux (red line). See Figure 3 Figure 2 shows the CO2 flux optimization results before and after June 2023 near Anji County, Zhejiang Province. The left side shows the EDGAR results before optimization, and the right side shows the results using the example of the present invention. The results on the right side of the present invention mainly optimize the flux on roads in the central and southeastern parts of Anji County, which may be due to the high industrial output and heavy traffic volume in that month.

[0072] In specific implementation, the method proposed in the technical solution of the present invention can be automatically run by those skilled in the art using computer software technology. System devices that implement the method, such as computer-readable storage media that store the corresponding computer program of the technical solution of the present invention and computer equipment that runs the corresponding computer program, should also be within the scope of protection of the present invention.

[0073] In some possible embodiments, a device based on the WRF-FLEXPART model and WRF-FLEXPART inversion of mesoscale CO2 flux is provided, comprising the following modules:

[0074] The first module is used to obtain meteorological data as input to the WRF model, which outputs the meteorological field of the study area;

[0075] The second module is used to sort the CO2 concentration monitored by the monitoring station by time to obtain an observation sequence, and to remove outliers from the observation sequence and perform window sliding average processing to obtain a CO2 rising concentration sequence;

[0076] The third module is configured to use the meteorological field as input to the atmospheric transport models FLEXPART and STILT, and obtain the footprint field corresponding to the processed observation sequence based on the outputs of the atmospheric transport models FLEXPART and STILT according to the operating time and geographical location of the monitoring station;

[0077] The fourth module is used to obtain a better posterior total carbon flux using Bayesian estimation based on the prior total carbon flux and the CO2 elevation concentration series combined with the footprint field.

[0078] In some possible embodiments, an electronic device is provided, comprising one or more processors; a memory for storing one or more programs; and when the one or more programs are executed by the one or more processors, the one or more processors implement the mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT.

[0079] In some possible embodiments, a computer-readable storage medium is provided, on which computer instructions are stored, wherein the computer instructions are used to cause a computer to execute the steps of the mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT.

[0080] The specific embodiments described herein are merely illustrative of the spirit of the present invention. Persons skilled in the art may make various modifications, additions, or substitutions to the described embodiments without departing from the spirit of the present invention or exceeding the scope of the appended claims.

Claims

1. A mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT, characterized by: The following steps are involved: Step 1: Obtain meteorological data as input to the WRF model, which outputs the meteorological field of the study area. Step 2: Monitor the CO2 concentration using a monitoring station and sort the CO2 concentration by time to obtain an observation sequence; remove outliers from the observation sequence and perform window sliding average processing to obtain a CO2 elevation concentration sequence; Step 3: Using the meteorological field obtained in step 1 as the input of the atmospheric transport models FLEXPART and STILT, and obtaining the footprint field corresponding to the observation sequence processed in step 2 based on the output of the atmospheric transport models FLEXPART and STILT according to the operating time and geographical location of the monitoring station; Step 4: Based on the prior total carbon flux and the CO2 elevation concentration series, the footprint field is combined with Bayesian estimation to obtain a better posterior total carbon flux; In step 4, the calculation formula of the posterior total carbon flux is as follows: x=x a +S a H T (HS a H T +S o ) -1 (y obs -Hx a ) Where x is the posterior carbon flux to be solved; x a is the prior carbon flux; S a is the covariance matrix of the prior carbon flux; y obs is the rising concentration series of CO2; S o is the covariance matrix of the CO2 elevation concentration series; H is the footprint field; The covariance matrix S of the prior carbon flux a The covariance matrix S of the CO2 concentration series o The off-diagonal elements of are calculated as: Among them, σ ij represents the off-diagonal element in the i-th row and j-th column of the matrix; σ ii and σ jj Represents the diagonal elements of the matrix; L represents the time or space correlation coefficient, corresponding to the matrix S a and matrix S o ; d(i,j) represents the time or space distance, corresponding to the matrix S a and matrix S o .

2. The mesoscale CO2 flux inversion method based on WRF-FLEXPART and WRF-STILT according to claim 1, characterized in that: In step 2, the observation sequence is subjected to outlier removal and window sliding average processing to obtain a CO2 elevated concentration sequence. The specific method is: For an observation sequence within a window, the non-outlier range is determined by calculating the IQR within the window, thereby screening non-outlier values in the observation sequence. The non-outlier range is as follows: (Q1-1.5*IQR, Q3+1.5*IQR) Among them, Q1 represents the lower quartile value of the observation sequence, Q3 represents the upper quartile value of the observation sequence, and IQR represents the difference between Q3 and Q1; Observations that are no greater than the lower ten percentile in the observation sequence after outliers are removed are screened out. Fourteen days are used as a sliding window, and the lower ten percentile observations are taken as a sliding average to obtain the baseline sequence. The difference between the observation sequence and the baseline sequence is the CO2 elevated concentration sequence. Values less than or equal to 0 in the CO2 elevated concentration sequence are eliminated.

3. A mesoscale CO2 flux inversion device based on WRF-FLEXPART and WRF-STILT, used to perform the method of claim 1 or 2, characterized in that: include The first module is used to obtain meteorological data as input to the WRF model, which outputs the meteorological field of the study area; The second module is used to sort the CO2 concentration monitored by the monitoring station by time to obtain an observation sequence, and to remove outliers from the observation sequence and perform window sliding average processing to obtain a CO2 rising concentration sequence; The third module is configured to use the meteorological field as input to the atmospheric transport models FLEXPART and STILT, and obtain the footprint field corresponding to the processed observation sequence based on the outputs of the atmospheric transport models FLEXPART and STILT according to the operating time and geographical location of the monitoring station; The fourth module is used to obtain a better posterior total carbon flux using Bayesian estimation based on the prior total carbon flux and the CO2 elevation concentration series combined with the footprint field.

4. An electronic device, characterized in that: include: one or more processors; a memory for storing one or more programs; When the one or more programs are executed by the one or more processors, the one or more processors implement the method according to claim 1 or 2.

5. A computer-readable storage medium having computer instructions stored thereon, characterized in that: The computer instructions are used to enable a computer to execute the steps of the method according to claim 1 or 2.

Images