Soil surface radon exhalation rate correction method based on Kalman filtering algorithm

Abstract

The invention discloses a soil surface radon exhalation rate correction method based on a Kalman filtering algorithm. The method comprises a measurement process, a data processing process and a calculation process. In the measurement process, performing nonlinear fitting on radon concentration values measured by a radon measuring instrument at different time to obtain a measured soil surface radon exhalation rate value Jt; in the data processing process, environment noise in the measurement process and measurement noise carried by an instrument are set to be random signals in the Kalman filtering process, a corrected radon concentration value is obtained through a Kalman filtering algorithm, and the radon exhalation rate corrected through the Kalman filtering algorithm is obtained through nonlinear fitting. And comparing the corrected value, the measured value and the theoretical value of the soil surface radon exhalation rate Kalman filtering algorithm. According to the method, measurement errors caused by instruments and environments in the measurement process can be corrected, the accuracy of the radon exhalation rate on the soil surface is improved, and a basis is provided for subsequent application and treatment.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of soil surface radon exhalation rate measurement, in particular to a soil surface radon exhalation rate correction method based on Kalman filtering algorithm. BACKGROUND

[0002] Soil radon is usually the main source of indoor radon gas, accounting for about 70% to 90% of the total indoor radon. Radon is the second largest cause of lung cancer, second only to smoking. Therefore, monitoring the exhalation rate of soil radon is very important for evaluating indoor radon exposure risk.

[0003] Soil surface radon exhalation rate measurement mainly includes closed loop method and open loop method. The equipment for measuring soil surface radon exhalation rate by closed loop method mainly includes RAD7 radon detector, drying tube and accumulation chamber. The equipment for measuring soil surface radon exhalation rate by open loop method mainly includes RAD7 radon detector, drying tube and ventilation chamber. The traditional closed loop method and open loop method have system errors in measurement due to the influence of instrument noise and environmental noise, affecting the measurement accuracy.

[0004] Kalman filtering algorithm is an algorithm that uses linear system state equation to make optimal estimation of system state through system input and output observation data. It can effectively handle the uncertainty in linear system, and make optimal estimation of system state in the presence of noise by reasonably balancing prediction information and observation information. Kalman filtering algorithm has been widely used in signal processing field. Therefore, it is a very meaningful work to further process the measured data by Kalman filtering algorithm to improve the measurement accuracy. SUMMARY

[0005] The purpose of the present application is to overcome the above-mentioned shortcomings of the prior art and provide a soil surface radon exhalation rate correction method based on Kalman filtering algorithm.

[0006] The technical solution of the present application is: a soil surface radon exhalation rate correction method based on Kalman filtering algorithm, including measurement process, data processing process and calculation process.

[0007] The measurement process includes: first, using a radon detector to measure the radon concentration in the environment around the soil to be measured to obtain the initial radon concentration in the environment; second, placing the radon collection cover upside down on the surface of the soil to be measured, under the action of the built-in air pump of the RAD7 radon detector, the radon-containing air on the soil surface passes through the radon collection cover, the drying tube and the hose into the RAD7 radon detector, and the radon concentration in the radon collection cover is continuously measured in the RAD7 radon detector; the soil surface radon concentration is calculated by the RAD7 radon detector using the following formula: (1) Wherein, c(t) table C is the radon concentration in the radon collection cover, unit is Bq / m3 ;J t represents radon exhalation rate, unit is ; S represents the bottom area of the radon collection hood, unit is ; V represents the volume of the radon collection hood, hose and the internal pump of the RAD7 radon meter, unit is ; represents the effective decay constant, which is the sum of the decay coefficient and the leakage coefficient λ , the diffusion coefficient λ b and the cavity leakage coefficient λ leak , all units are .

[0008] Because the environmental radon concentration value at the initial moment is very low, approximately zero, it is ignored in calculation, at this time the solution of equation (1) is: (2) Substitute the radon concentration values measured by the RAD7 radon meter at multiple different times into equation (2) for nonlinear fitting to obtain the measured soil surface radon exhalation rate value J t .

[0009] The data processing process includes: first, set the environmental noise in the measurement process and the measurement noise carried by the instrument itself as random signals and in the Kalman filtering process, assuming that they are white noise independent of each other and subject to normal distribution: (3) (4) Where, are the probabilities of noise, , the process excitation noise covariance matrix, and the observation noise covariance matrix.

[0010] Second, estimate the state variable of the discrete-time process using the Kalman filter, which is described by the following discrete random difference equation: (5) Where, is the estimated value at the last moment, and ​respectively, A is an n x n gain matrix, where n is the dimension of the state variable, B is the gain of optional control input

[0011] Define the observation variable , the measurement equation is obtained: (6) where, when the control function or the process excitation noise is zero, the n x n gain matrix A in the difference equation (5) linearly maps the previous time k-1 state to the current time k state; the n x l matrix B represents the gain of optional control input , where =1; the m x n matrix H in the measurement equation (6) represents the gain of the state variable to the measurement variable , where m is the number of observation variables.

[0012] The calculation process includes: first, define as the prior state estimate of the kth step under the condition that the state before the kth step is known, define as the posterior state estimate of the kth step when the measurement variable is known, and thus the calculation expressions of the prior estimation error and the posterior estimation error are: (7) (8) The calculation expression of the covariance of the prior estimation error is: (9) The calculation expression of the covariance of the posterior estimation error is: (10) The calculation expression of the posterior state estimate obtained from the prior state estimate and the observation variable is: (11) where K is the gain or mixing factor of the residual, which is an m x n matrix, so that the posterior estimation error covariance is minimized.

[0013] Simultaneous equations (8), (9), (10) and (11) are set The first derivative of K with respect to the derivative is zero, so that the value of K is solved, and the expression of K is: (12) ​Where R is the observation noise covariance matrix, and from formula (12), we know that the smaller the observation noise covariance R, the larger the residual gain K; in particular, when R approaches zero, we have: (13) On the other hand, prior estimation error covariance The smaller the value, the smaller the residual gain K; in particular, when As it approaches zero, we have: (14) As the measurement noise covariance R approaches zero, the measurement variable... The weight of is getting bigger and bigger, and Prediction The weight of the prior estimate covariance is decreasing; on the other hand, as the prior estimate covariance... Approaching zero, the measured variable The weight of is getting smaller and smaller, while predict Its weight is getting bigger and bigger.

[0014] The Kalman filter includes a time update equation and a measurement update equation. The time update equation is obtained from the difference equation (5): (15) Due to the excitation noise in equation (5) It is unobservable during the prediction process, and therefore reflected in the prior covariance matrix. middle.

[0015] From equations (11) and (12), the measurement update equation is obtained as follows: (16) First, calculate the gain of the Kalman filter. Secondly, measure the output to obtain Then, the posterior estimate of the state is generated according to equation (11), and finally, the state is estimated according to equation (16). Estimate the posterior covariance of the state, calculate the time update equation and the measurement update equation, set the posterior estimate obtained in the previous calculation as the prior estimate for the next calculation, and repeat the process.

[0016] In the formula (16) The estimated value is the radon concentration value after correction by the Kalman filter algorithm, and the radon release rate after correction by the Kalman filter algorithm is obtained by nonlinear fitting.

[0017] A further technical solution of the present invention is: the radon collection cover is placed upside down on the surface of the soil to be tested and inserted into the soil for 2-6 cm, and its surroundings are covered with soil to seal it.

[0018] A further technical solution of the present invention is that the dimensions of the state variables are time and radon concentration, and the number of observed variables is 12.

[0019] Compared with the prior art, the present invention has the following characteristics: This invention avoids measurement errors caused by instruments in existing technologies. The Kalman filter algorithm can correct measurement errors caused by instruments and the environment, making the final calculated radon exhalation rate on the soil surface more accurate.

[0020] The detailed structure of the present invention will be further described below with reference to the accompanying drawings and specific embodiments. Attached Figure Description

[0021] Appendix fig. 1 This is a schematic diagram of a soil surface radon exhalation rate measurement device based on the Kalman filter algorithm. Appendix fig. 2 A fitted curve of actual radon measurements at the soil surface diffusion site; Appendix fig. 3 A fitted curve of radon concentration measurements after correction using the Kalman filter algorithm; Appendix fig. 4 This is a fitted curve of the theoretical value of diffusion on the soil surface. Detailed Implementation

[0022] Example 1, Appendix figs. 1-4 A method for correcting the radon exhalation rate on soil surface based on the Kalman filter algorithm, including the measurement process, data processing process and calculation process.

[0023] The measurement process includes: First, using a RAD7 radon meter 1, the radon concentration in the environment surrounding the soil to be tested is measured to obtain the initial radon concentration in the environment. Next, the radon collection hood 3 is inverted onto the surface of the soil to be tested and inserted 2-6 cm into the soil, and its surroundings are covered with soil to seal it. Under the action of the built-in air pump of the RAD7 radon meter 1, radon-containing air from the soil surface enters the RAD7 radon meter 1 through the radon collection hood 3, the drying tube 2, and the flexible tube. The radon concentration inside the radon collection hood 3 is continuously measured within the RAD7 radon meter 1. fig. 1 As shown. The radon concentration at the soil surface is calculated using the following formula from the RAD7 radon meter 1: (1) in, c(t) table The radon concentration inside the radon enclosure is shown in units of... ; J t This indicates the radon release rate, in units of... ; S This represents the base area of ​​the radon collection enclosure, in units of... ;V V represents the volume of the radon collection chamber, the hose and the internal pump of the RAD7 radon meter, in ; Veff represents the effective decay constant, which is the decay coefficient and the leakage coefficient λ , the diffusion coefficient λ b and the leakage coefficient of the chamber λ leak , all in .

[0024] Since the ambient radon concentration value at the initial moment is very low, approximately zero, it is neglected in the calculation, and the solution of equation (1) at this time is: (2) The radon concentration values measured by the RAD7 radon meter at multiple different times are substituted into equation (2) for nonlinear fitting to obtain the measured soil surface radon exhalation rate value J t . The soil surface radon exhalation rate values at different measurement times J t are fitted to obtain the fitting curve diagram of the radon concentration values measured at the soil surface diffusion site as shown in the accompanying fig. 2 .

[0025] After obtaining the radon concentration value inside the radon collection chamber, according to previous measurement experience, we understand that the environmental noise fluctuation interval during the measurement operation is 0.05 to 1; at the same time, the measurement noise carried by the measurement instrument is also in the range of 0.05 to 1. Therefore, the measured soil surface radon exhalation rate value J t exists a certain error with the true soil surface radon exhalation rate value, so it is necessary to further process the measured soil surface radon exhalation rate value.

[0026] The data processing process includes: first, set the environmental noise during the measurement process and the measurement noise carried by the instrument itself as the random signals and in the Kalman filtering process, assuming that they are white noise independent of each other and subject to normal distribution: (3) (4) wherein, are the probability density functions of the noise, is the process excitation noise covariance matrix, is the observation noise covariance matrix.

[0027] ​Second, the Kalman filter is used to estimate the state variables of the discrete-time process , which is described by the following discrete stochastic difference equation: (5) where is the estimate at the previous time step k-1, and are the excitation noise and observation noise at the previous time step k-1, respectively, A is an n x n gain matrix, where n is the dimension of the state variable, and B is the gain of the optional control input In this example, n = 2, representing time and radon concentration, respectively.

[0028] The observation variable is defined, and the measurement equation is obtained: (6) where, when the control function or the process excitation noise is zero, the n x n gain matrix A in the difference equation (5) linearly maps the state at the previous time step k-1 to the state at the current time step k; the n x l matrix B represents the gain of the optional control input , where = 1; the m x n matrix H in the measurement equation (6) represents the gain of the state variable to the measurement variable , where m is the number of observation variables. In this example, m = 12.

[0029] The calculation process includes: first, define as the prior state estimate at the kth step under the condition that the state before the kth step is known, and define as the posterior state estimate at the kth step when the measurement variable is known, so the calculation expressions of the prior estimation error and the posterior estimation error are: (7) (8) The calculation expression of the covariance of the prior estimation error is: (9) The calculation expression of the covariance of the posterior estimation error is: (10) The calculation expression of the posterior state estimate obtained from the prior state estimate and the observation variable is: (11) where K is the residual gain or mixing factor, which is an mxn matrix, such that the posterior estimation error covariance is minimized.

[0030] Solving equations (8), (9), (10) and (11) simultaneously, and setting The first derivative of K with respect to K is zero, thus the value of K is solved as (12) where R is the observation noise covariance matrix. From equation (12), it can be seen that the smaller the observation noise covariance R, the larger the residual gain K; in particular, when R tends to zero, we have (13) On the other hand, the prior estimation error covariance is smaller, the residual gain K is smaller; in particular, when tends to zero, we have (14) Another interpretation of the gain K is that as the measurement noise covariance R tends to zero, the weight of the measurement variable is getting larger and larger, while the weight of the prediction is getting smaller and smaller; on the other hand, as the prior estimation covariance tends to zero, the weight of the measurement variable is getting smaller and smaller, while the weight of the prediction is getting larger and larger. A Kalman filter estimates the process state using a feedback control approach, including a time update equation and a measurement update equation. The time update equation is used to advance the numerical values of the current state variable and the error covariance estimate accurately and timely, thus constructing a priori prediction for the subsequent time state. The measurement update equation, on the other hand, assumes the function of feedback, i.e., combining the a priori prediction with the latest measurement variable to generate a more accurate posteriori estimate. In this sense, the time update equation can also be regarded as a prediction equation, while the measurement update equation can be regarded as a correction equation.

[0031] The time update equation is obtained from the difference equation (5) as

[0032] (15) Since the excitation noise in equation (5) is unobservable in the prediction process, it is reflected in the prior covariance matrix , which is the covariance matrix of the excitation noise.

[0033] ​​​From equations (11) and (12), the measurement update equation is obtained as follows: (16) The first step in the measurement update equation is to calculate the gain of the Kalman filter. Secondly, measure the output to obtain Then, the posterior estimate of the state is generated according to equation (11), and finally, the state is estimated according to equation (16). Estimate the posterior covariance of the state. After calculating the time update equation and the measurement update equation, the entire process is repeated. The posterior estimate obtained in the previous calculation is used as the prior estimate for the next calculation.

[0034] In the formula (16) The estimated value is the radon concentration value corrected by the Kalman filter algorithm, and a nonlinear fitting is performed on it to obtain the result shown in the attached figure. fig. 3 The radon extraction rate curve shown is the result of correction using the Kalman filter algorithm.

[0035] To evaluate the effectiveness of the corrected radon emanation rate value, the radon concentration values ​​at the soil surface were measured and fitted to obtain the values ​​shown in the attached figure. fig. 2 The radon release rate values ​​shown J t Furthermore, by setting appropriate environmental noise fluctuation values ​​and the measurement noise fluctuation values ​​built into the measuring instrument, the theoretical radon concentration value of the soil surface is obtained, and nonlinear fitting is performed to obtain the results as shown in the attached figure. fig. 4 The theoretical radon release rate curve is shown. This is compared with the measured radon release rate values ​​at the soil surface. J t The effectiveness of the Kalman filter algorithm is evaluated by assessing the closeness between the corrected values ​​obtained after setting environmental noise and measurement noise using the Kalman filter algorithm and the theoretical values.

[0036] In this embodiment, a RAD7 radon meter 1, a drying tube 2, and a radon collection hood 3 were used to measure the radon release rate of soil in a certain location. Based on measurement experience, the environmental noise fluctuation value and the measurement noise fluctuation value of the instrument itself were both set to 0.1 during the measurement process, thereby obtaining the radon concentration on the soil surface. The measured values, the corrected values ​​after Kalman filtering, and the theoretical values ​​of radon concentration on the soil surface are shown in Table 1.

[0037] Table 1. Measured values ​​of radon concentration at soil surface, corrected values ​​after Kalman filter algorithm correction, and theoretical values ​​of radon concentration at soil surface. (Unit is) ) soil surface radon concentration value measurement value kalman filter algorithm correction after correction value soil surface radon concentration value theoretical value 0 0 9.3 143.5 94 53.0 93.8 93 116.6 102.6 99 188.3 272.9 205 263.1 529.5 405 339.1 396.8 399 415.4 396.8 397 491.7 790.2 640 567.8 876.2 785 643.7 577.1 656 719.4 540.8 584 794.8 The relevant parameters for the nonlinear fitting curve of the measured value of radon concentration on soil surface, the corrected value after Kalman filter algorithm correction, and the theoretical value of radon concentration on soil surface are shown in Table 2.

[0038] Table 2. Measured values ​​of radon exhalation rate from soil surface, corrected values ​​after Kalman filter algorithm correction, and soil surface radon exhalation rate. Relevant parameters fitted from the theoretical value of radon exudation rate parameter soil surface radon exhalation rate value measurement value kalman filter algorithm correction after correction value soil surface radon exhalation rate value theoretical value J (Bqm -2 s -1 )]]> ​ 0.012 ± 0.0047 0.009 ± 0.002 0.008 ± 1.3 x 10 -5 ]] (s -1 )]]> (5.08 ± 4.77) x 10 -4 ]] (1.98 ± 2.75) x 10 -4 ]] 1.5 x 10 -5 ± 1.6 x 10 -5 ]]> [R 2 ]]> 0.76 0.89 0.99 From Table 1, Table 2 and Appendix figs. 2-4 As can be seen, the corrected value after Kalman filtering is closer to the theoretical value of radon exhalation rate from the soil surface. Furthermore, the experimental results show that after correction using the Kalman filtering algorithm, both the radon exhalation rate from the soil surface and the corrected value are more consistent. Or the goodness of fit R 2 Both values ​​are closer to the theoretical value, indicating that Kalman filtering is of great significance for correcting the radon exhalation rate on the soil surface and obtaining more accurate measurement results.

Claims

1. A method for correcting soil surface radon exhalation rate based on Kalman filtering algorithm, comprising a measurement process, a data processing process and a calculation process; The measurement procedure comprises: First, the radon concentration in the environment around the soil to be measured is measured by using a radon measuring instrument to obtain the initial environmental radon concentration; second, the radon collection cover is inverted and placed on the surface of the soil to be measured, under the action of the built-in air pump of the RAD7 radon measuring instrument, the radon-containing air on the soil surface passes through the radon collection cover, the drying tube and the hose into the RAD7 radon measuring instrument, and the radon concentration in the radon collection cover is continuously measured in the RAD7 radon measuring instrument; the soil surface radon concentration calculated by the RAD7 radon measuring instrument is as follows: (1) wherein, c(t) represents the soil surface radon concentration at time t represents the radon concentration in the radon collection hood, in ; J t represents the radon exhalation rate, in ; S represents the bottom area of the radon collection hood, in ; V represents the volume of the radon collection hood, the hose, and the internal pump of the RAD7 radon detector, in ; represents the effective decay constant, which is the sum of the decay coefficient and the leakage coefficient λ , the diffusion coefficient λ b and the chamber leakage coefficient λ leak , all in ; Since the initial environmental radon concentration value is very low, approximately zero, it is ignored during calculation, and the solution of equation (1) at this time is: (2) The radon concentration values at different times measured by the RAD7 radon meter are substituted into equation (2) for nonlinear fitting to obtain the measured soil surface radon exhalation rate value J t ; The data processing process includes: firstly setting the environmental noise in the measurement process and the measurement noise carried by the instrument itself as random signals in the Kalman filtering process and , assuming that they are white noises independent of each other and subject to normal distribution: (3) (4) wherein are noise, is the probability density function, is the process excitation noise covariance matrix, is the observation noise covariance matrix; Second, the state variables of the discrete-time process are estimated using a Kalman filter The discrete-time process is described by the following discrete stochastic difference equation: (5) wherein is the estimate of the previous time step, and are the excitation noise and observation noise of the previous time step, respectively, A is a gain matrix of n x n order, where n is the dimension of the state variable, B is a gain of the optional control input . Defining the observation variable , resulting in the measurement equation: (6) where the control function or the process excitation noise is zero, the n x n gain matrix A in the difference equation (5) linearly maps the state at the previous time instant k-1 to the state at the current time instant k; the n x l matrix B represents the gain of the optional control input where = 1; the m x n matrix H in the measurement equation (6) represents the gain of the state variables on the measurement variables where m is the number of observation variables; The computation process includes first defining the prior state estimate for the kth step given the k-1th step state and the measurement variable , the posterior state estimate for the kth step, the prior estimate error and the posterior estimate error are computed as follows: (7) (8) The covariance calculation expression of the prior estimation error is: (9) The covariance calculation expression of the posterior estimation error is: (10) The calculation expression of the posterior state estimation obtained from the prior state estimation and the observation variable is: (11) Wherein, K is the residual gain or mixing factor, which is an m×n order matrix, so that the posterior estimation error covariance is minimized; Simultaneous equations (8), (9), (10) and (11) are set up The first derivative of K is zero, so that the value of K is solved, and the expression of the value of K is: (12) Wherein, R is the observation noise covariance matrix, and it can be known from formula (12) that the smaller the observation noise covariance R, the greater the residual gain K; in particular, when R tends to zero, there is: (13) On the other hand, the a priori estimate error covariance The smaller the residual gain K is; in particular, when Tends to zero, one has: (14) As the measurement noise covariance R approaches zero, the measurement variable... The weight of is getting bigger and bigger, and Prediction The weight of the prior estimate covariance is decreasing; on the other hand, as the prior estimate covariance... Approaching zero, the measured variable The weight of is getting smaller and smaller, while predict Its weight is getting bigger and bigger; The Kalman filter includes a time update equation and a measurement update equation, and the time update equation obtained from the difference equation (5) is: (15) Due to the excitation noise in equation (5) which is not observable in the prediction process, it is therefore reflected in the prior covariance matrix ​ The measurement update equation obtained from equation (11) and equation (12) is: (16) First, the gain of the Kalman filter is computed Second, the measurement output is obtained Then, the posterior estimate of the state is generated according to equation (11), and finally, the posterior covariance of the state is estimated according to equation (16) The time update equation and the measurement update equation are computed, and the posterior estimate obtained in the last computation is set as the prior estimate for the next computation, and the process is repeated again. In the formula (16) The estimated value is the radon concentration value after the Kalman filtering algorithm correction, and the radon exhalation rate after the Kalman filtering algorithm correction is obtained by nonlinear fitting.

2. The method for correcting soil surface radon exhalation rate based on Kalman filtering algorithm according to claim 1, characterized in that: The radon collection cover is inverted and placed on the surface of the soil to be measured and inserted into the soil 2-6 cm, and the surrounding soil is covered to seal it.

3. The method for correcting soil surface radon exhalation rate based on Kalman filtering algorithm according to claim 1, characterized in that: The dimensions of the state variable are time and radon concentration, Wherein, the number of observation variables is 12.

Images