CN117669270B - 一种对微站组网数据进行时空一致性订正的方法 - Google Patents
一种对微站组网数据进行时空一致性订正的方法 Download PDFInfo
- Publication number
- CN117669270B CN117669270B CN202410130460.XA CN202410130460A CN117669270B CN 117669270 B CN117669270 B CN 117669270B CN 202410130460 A CN202410130460 A CN 202410130460A CN 117669270 B CN117669270 B CN 117669270B
- Authority
- CN
- China
- Prior art keywords
- concentration
- data
- time
- station
- space
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 72
- 230000006855 networking Effects 0.000 title claims abstract description 34
- 239000003344 environmental pollutant Substances 0.000 claims abstract description 89
- 231100000719 pollutant Toxicity 0.000 claims abstract description 88
- 238000012544 monitoring process Methods 0.000 claims abstract description 75
- 238000009792 diffusion process Methods 0.000 claims abstract description 44
- 238000004364 calculation method Methods 0.000 claims abstract description 38
- 238000012545 processing Methods 0.000 claims abstract description 25
- 230000008569 process Effects 0.000 claims abstract description 22
- 238000005315 distribution function Methods 0.000 claims abstract description 12
- 238000000528 statistical test Methods 0.000 claims abstract description 6
- 238000012937 correction Methods 0.000 claims description 67
- 230000006870 function Effects 0.000 claims description 29
- 238000004088 simulation Methods 0.000 claims description 26
- 238000003908 quality control method Methods 0.000 claims description 17
- 238000005070 sampling Methods 0.000 claims description 17
- 239000000356 contaminant Substances 0.000 claims description 16
- 239000013598 vector Substances 0.000 claims description 16
- 238000004422 calculation algorithm Methods 0.000 claims description 9
- 238000001914 filtration Methods 0.000 claims description 8
- 238000010586 diagram Methods 0.000 claims description 7
- 239000000779 smoke Substances 0.000 claims description 7
- 230000008859 change Effects 0.000 claims description 6
- 230000002159 abnormal effect Effects 0.000 claims description 5
- 238000005259 measurement Methods 0.000 claims description 5
- 238000007476 Maximum Likelihood Methods 0.000 claims description 4
- 239000011159 matrix material Substances 0.000 claims description 4
- 238000010606 normalization Methods 0.000 claims description 4
- 238000002790 cross-validation Methods 0.000 claims description 3
- 230000007613 environmental effect Effects 0.000 claims description 3
- 238000011835 investigation Methods 0.000 claims description 3
- 230000001419 dependent effect Effects 0.000 claims description 2
- 230000008571 general function Effects 0.000 claims description 2
- 230000033764 rhythmic process Effects 0.000 claims description 2
- 239000000126 substance Substances 0.000 claims description 2
- 230000009897 systematic effect Effects 0.000 abstract description 5
- 239000007789 gas Substances 0.000 description 5
- 239000012855 volatile organic compound Substances 0.000 description 4
- 238000004458 analytical method Methods 0.000 description 3
- 230000010354 integration Effects 0.000 description 3
- 241000282414 Homo sapiens Species 0.000 description 2
- 239000000809 air pollutant Substances 0.000 description 2
- 231100001243 air pollutant Toxicity 0.000 description 2
- 230000007547 defect Effects 0.000 description 2
- 238000001514 detection method Methods 0.000 description 2
- 230000007774 longterm Effects 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 238000007405 data analysis Methods 0.000 description 1
- 238000012217 deletion Methods 0.000 description 1
- 230000037430 deletion Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000007689 inspection Methods 0.000 description 1
- 230000002452 interceptive effect Effects 0.000 description 1
- 238000010801 machine learning Methods 0.000 description 1
- 238000012423 maintenance Methods 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000002093 peripheral effect Effects 0.000 description 1
- 239000002243 precursor Substances 0.000 description 1
- 238000003672 processing method Methods 0.000 description 1
- 238000012216 screening Methods 0.000 description 1
- 230000011218 segmentation Effects 0.000 description 1
- 230000001360 synchronised effect Effects 0.000 description 1
- 238000012876 topography Methods 0.000 description 1
- 238000004800 variational method Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/15—Correlation function computation including computation of convolution operations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q50/00—Information and communication technology [ICT] specially adapted for implementation of business processes of specific business sectors, e.g. utilities or tourism
- G06Q50/10—Services
- G06Q50/26—Government or public services
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Business, Economics & Management (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Data Mining & Analysis (AREA)
- Mathematical Analysis (AREA)
- Computational Mathematics (AREA)
- Tourism & Hospitality (AREA)
- Mathematical Physics (AREA)
- General Engineering & Computer Science (AREA)
- Human Resources & Organizations (AREA)
- Health & Medical Sciences (AREA)
- Strategic Management (AREA)
- Primary Health Care (AREA)
- Computing Systems (AREA)
- Marketing (AREA)
- General Health & Medical Sciences (AREA)
- Economics (AREA)
- Algebra (AREA)
- General Business, Economics & Management (AREA)
- Educational Administration (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- Development Economics (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明提供一种对微站组网数据进行时空一致性订正的方法，使用监测数据处理子模块处理微站监测数据，进而对监测数据做四维订正，结合模拟区域的经纬度，进行大气污染物排放扩散模型及其建立求解过程，对模拟区域内风场空间进行同化处理，对污染物日浓度变化幅度参量进行假设，并将各个站点的污染物日浓度参量对于时间作时间分布函数，最终在大气污染物排放总量反推与统计检验计算模块运行中拟合出大气污染物浓度与路径积分长度函数的关系。本发明用于修正微站单站时间序列的***误差，消除微站组网监测数据的站别差异，提高时空代表性。
Description
技术领域
本发明属于测量数据处理技术领域，涉及一种对微站组网数据进行时空一致性订正的方法，用于对气体和大气颗粒物的浓度时空分布进行参数化订正及数值建模。
背景技术
目前存在的针对工业园区大气排放的微站监测项目，主要是通过使用固定式和移动式等多种不同形式的微站组网，实现对一定区域范围内空气质量网格化的在线监测和实时数据集成的常规综合分析，此项技术不仅可以用来判断一次污染物、二次污染物及其前体物的浓度和浓度变化趋势，对城市内的工业园区及工厂复合污染的成因和机理也有合理的推导和认证，从而更好的服务于大气科学研究。
微站监测项目采取了多种形式的微站组网，微站仪器虽然具有价格低廉维护简易等优势，但是与标准大气在线监测站相比，存在由于不同点位仪器之间由于校准调试不一致或长期运行中非定向漂移导致***误差较大、数据一致性较差、监测采样时间较短和分析数据存在较大不确定性等一系列问题。此外，在复杂的建筑分布和地形地物等条件下进行大范围布点，存在不同点位受到周边环境非均一性影响而造成数据时空代表性差异的现象。以上问题对微站数据组网的分析的应用可靠性和广泛性均造成了影响。
因此开发一种快速订正方案，对微站单站时间序列进行***误差修正，消除微站组网监测数据站别差异，进而提高其时空代表性，就显得愈发重要。
发明内容
发明目的：为了克服现有技术中存在的不足，本发明提供一种具有高可靠性和高实时性的针对微站组网数据进行快速订正及初步处理的方法，用于解决大气颗粒物浓度检测中误差订正的问题。
技术方案：为解决上述技术问题，本发明提供一种对微站组网数据进行时空一致性订正的方法，使用监测数据处理子模块处理微站监测数据，进而对监测数据做四维订正；结合模拟区域的经纬度，进行大气污染物排放扩散模型的建立与求解；同化处理模拟区域内风场空间，假设污染物日浓度变化幅度参量近似守恒，并将各个站点的污染物日浓度参量作时间分布函数；最终在大气污染物排放总量反推与统计检验计算模块运行中拟合出大气污染物浓度与路径积分长度函数的关系；所述监测数据处理子模块包括浓度监测采集粗质控子模块、站点微环境影响因子的变分插值子模块、浓度时空分布的四维同化子模块和基于大气扩散模拟的松弛逼近子模块。
具体地，所述浓度监测采集数据的粗质控模块用于消除直接采样得到的数据集中的随机数据野点。监测数据通过园区周边的监测监控点采集获取，检测和监控点可提供同步多点、时间连续的大气污染物成分的浓度采样数据。但是由于复杂工况、仪器鲁棒性和近距离偶发排放干扰和其他环境偶然因素的影响，直接采样得到的数据集中至少存在随机数据野点，故使用本模块消除上述影响因素。
具体地，浓度信息块T i,j,k 是独立空间格点的污染物浓度观测值，表示为：
T i,j,k =T true（i,j,k） +ε i,j,k (1)
式中T true(i,j,k) 表示对应空间格点上的真值，i,j,k是表示空间三个维度上的脚标，其中ε i,j,k 表示同类测量值中的未知的误差；
基于观测采样的独立性假设，忽略不同微站的采样之间彼此的干扰，令不同空间点位上的同一时间的采样数据之间的相关性为零：
E(ε i,j,k )=0 （2）
E(ε i,j,k )为同类测量值中的未知的误差的数学期望值，令其等于0；
同时对于每一台特定采样仪器，本文件中定义污染物浓度的时序观测方差为：
E(ε 2)=δ 2 (3)
ε 2为污染物浓度的时序观测值，δ 2 为污染物浓度的时序观测方差；
对有关对特定点位上的真实浓度观测值T true 的统计信息作出假设：
T true =δ(i,j,k)*T i,j,k (4)
由上述公式（4）推导得：δ(i,j,k)=1 (5)
由上述公式（5）最终得到：
(6)
式中为浓度真值的时序观测方差，/>同类测量值中未知误差的平方；
通过上述公式（6），根据已知的污染物浓度观测值，确定该时刻特定区域范围内的污染物浓度观测值的理论统计方差限，由此对所有采样数据进行误差判定，当某值落在误差限之外时，将该数值当作野点数据删除处理。
具体地，所述微环境影响因子的变分插值子模块用于消除微环境影响因子，或者对不同站点的监测数据进行订正后使得该影响因子保持一致；
在以往的理论研究中，对于监测布点的选取需满足以下条件——边地形足够平坦开阔，地理环境和周边物体对该点的大气流动，以及周边小微排放源对于该点浓度的影响要足够小，同时在整个监测网中的所有站点，以上各种因素均需保持对污染物监测影响的一致性。
站点微环境影响因子的变分插值子模块的核心思路是通过求取观测值的时空分布函数的变分极小值来得到对该变量的最佳估计。这个关于时空分布的变分函数我们定义为污染物浓度估计值T和其周边若干采样值之间的引力权重函数，如下式所示：
(7)
式中T i,j,k 为独立空间格点的污染物浓度观测值，ε i,j,k 为同类测量值中的未知的误差，T为本***中定义的污染物浓度估计值，J(T)为区域周边采样值之间的引力权重，N为采样值的个数；
应用黎曼变分理论证明上述公式（7）中：
当T= T true 时， (8)
式中T true(i,j,k) 表示对应空间格点上的真值；
以大气污染物扩散的高斯分布为例，某一站点位置受微环境影响后得到的浓度为：
（9）
式中C i,j,k (T||T true )为某站点位置受微环境影响后的污染物浓度；
而其污染物浓度的真值似然值计算公式为：
（10）
式中为污染物浓度的真值似然值；
推得污染物浓度的极大似然值为：
（11）
式中const为根据具体污染物种类和排放节律所设定的常数， H为空间的观测算子，表示数组规模为n的观测变量T的极大似然值。
通过对空间中所有观测点位当前时刻的观测值进行公式（11）的似然函数变分极小值计算，得出对于后验概率极大的污染物浓度观测真值。
具体地，所述浓度时空分布四维同化子模块用于减少计算误差，所述计算误差包括逐站产生随时间流逝向特定方向增加的***误差，以及不同站点之间表示为非特定方向且分别逐渐增加的***误差。
由于存在不同点位的仪器之间的人工校调、机器型号和软硬件性能非定向飘移等多方面因素，会导致逐站产生随时间流逝向特定方向增加的***误差，以及不同站点之间表示为非特定方向且分别逐渐增加的***误差。本***中的浓度时空分布四维同化子模块用于消除上述影响，减少计算误差。
选取a,b,c三个格点位置，定义此时的污染物浓度观测场为：
浓度观测场： （12）
式中X obs 为浓度观测场，为a点的浓度观测场，/>为b点的浓度观测场,为c点浓度观测场；
观测污染物浓度时周边环境的背景场为：
（13）
X bg 为观测污染物浓度时周边环境的背景场，为在a点观测污染物浓度时周边环境的背景场，/>为在b点观测污染物浓度时周边环境的背景场，/>为在c点观测污染物浓度时周边环境的背景场；
浓度观测向量为：
（14）
Y obs 为浓度观测向量，为某一个格点的浓度观测向量，/>为其他格点的浓度观测向量；
则环境背景场、浓度观测向量、浓度观测场之间应有下式关系：
（15）
上式中H为该空间的观测算子；其中h 1a、 h 1b 和h 1c 分别为在a,b和c点的基于本站点历史浓度数据获得的背景场影响因子；h 2a、 h 2b 和h 2c 分别为在a,b和c点的基于本站以外的周边其他站点历史浓度数据获得的背景场影响因子；
同时背景误差协方差矩阵为：（16）
式中B为背景误差协方差矩阵、b ab 表示a与b两个站点处背景误差的协方差因子，其它类同。b aa 表示a站点处背景误差的自协方差因子，其它类同。
最终得到空间浓度的观测算子、浓度观测向量、背景误差之间存在泛函数关系：
（17）
式中J obs 为区域周边采样值之间的引力权重系数，t 0 为初始时刻、t为任一观测时刻，L(t||t 0 ) T 为针对观测变量T（T为污染物浓度估计值，是多种大气污染物的监测浓度）根据初始时刻标定数据得到的对当前时刻观测值的时间订正函数，表示在任一空间位置（i,j,k）处关于观测量变T的空间订正系数、R为对于该H系数做有关距离的二阶订正的修正函数，其取值与位置有关；/>表示在任一位置(i,j,k)处浓度空间订正值与背景浓度观测向量之间的差。
令公式（17）= 0作为控制函数，求该函数的极小值，求得任一站点对初始时刻旧的污染物浓度观测值与终止时刻污染物浓度新的观测值之间的时空订正值，进而利用其在每个观测值上分别进行求解和订正，最终实现本子模块的订正目标函数。
具体地，所述基于大气扩散模拟的松弛逼近计算子模块，其松弛方程为：
（18）
上式中c为污染物浓度，τc为松弛尺度值，为势函数。其中松弛尺度τ的选项需要建立在对大气扩散方程的反复迭代的基础上，并依赖于求解变量本身的时空分布特征。若τ值太小，则动力过程的时间调整时间过短从而高估排放点源的数量；反之若τ值过大，则订正误差会在松弛逼近的过程中逐渐积累，从而造成计算结果可靠性的下降。通过上述的公式（18）反复迭代计算对松弛尺度值τ进行调整。/>为风速向量，/>为浓度的空间梯度项，cobs为扩散浓度的观测值，ctrue为扩散浓度的理想真值。
具体地，采用监测数据的四维订正来消除上述误差，是对原始观测浓度数据进行空间滤波处理、时序特征的多点一致性订正，以及根据布点空间相对位置和复杂地物分布特征对扩散路径上的阻抗差异进行的归一化订正。
具体地，通过园区地图定位确定原站点的经纬度坐标(lon,lat)，是将原站点的经纬度坐标数据(lon,lat)经过卫星比对、实地考察、远程走访的方法，对地图上园区位置进行标点，得到较精准的站点分布数据集。
具体地，所述监测数据处理子模块对原始监测风场数据的空间分辨率进行判别和质控，对有明显区别的不同情况进行归一化交叉验证处理，监测数据的采集与质控完成后的进行大气污染物排放和扩散模型的建立与求解；
进行大气污染物排放和扩散模型的建立与求解，是对监测数据处理子模块中的风场数据进行同化，对原始监测风场数据——包括站点实时测得的风速大小，风场方向，湍流强度的测算质量、测站站点的数量及数据的空间分辨率进行判别和质控。对有明显区别的不同情况进行归一化交叉验证处理。
若测站数量较多即测站数量在20个及以上，且园区提供原始风场的数据质量较高，使用二维空间内的散点插值算法在某一实时案例中使用matlab软件中的Scatteredinterolant函数，或是基于其他数学语言的高阶拉格朗日插值算法对根据原始风场数据构建的原始风场进行多点插值处理，该过程即为模拟区域内的风场空间同化。
若风场测站数量有限即测站数量在10个及以内，则对各个站点的数据质控处理过程进行动态权重分配，根据单一化风场的模拟结果比较选择处最优的权重系数。得出最优的单一化风场模拟结果，更加贴合实际园区或工厂周围的风场情况。
若风场测站仅为一个，则直接使用该站点测值生成模拟空间内均匀单一的风场分布。
监测数据采集与质控模块运行过程中对原始浓度数据进行基本的质量控制处理，剔除错误数据、空白数据、及明显不符合物理规律的数据。
具体地，大气扩散相关气象和化学系数计算过程，是根据监测数据采集和质控模块运行结束后得出的风场同化结果，同时依据烟羽扩散的物理模型，在受同源影响的各个站点区域范围内，依据模拟尺度区域较小的假设，设定一近似守恒量——污染物日浓度（变化幅度）参量计算方法如下
（19）
式中Λ为污染物日浓度参量，用于判断某站点当日的污染源排放情况；ΔX表示站点污染物日浓度的浮动幅度值，表示站点污染物日浓度均值；计算过程结束后将各个站点的日浓度参量Λ相对于时间t做分布函数图；
上述计算过程使用的烟羽扩散模型为连续点源排放的烟羽扩散,该模型假定某种大气污染物浓度分布服从正态分布的经验模式。烟羽以平均风速向下风向输送,忽略沿水平风输送方向的湍流扩散,垂直方向和水平横向扩散，由扩散参数描述。
上述烟羽扩散模型中正态型分布的烟羽浓度仅在平稳和均匀的条件下才能成立,并不符合实际大气的情形。但由于公式的形式和参数比较简单,仍广泛应用于扩散计算中,可近似地适用于开阔平坦地形上风速不太小且定常情况下的小尺度扩散计算,也可以根据地理条件、气象条件和排放源几何形状的特殊性,转换成各种变型模式。
其中将各个站点的日浓度参量Λ的时间t的分布函数进行空间卡曼滤波。根据滤波修正函数或者沿时间维滑动平均的逐站点修正因子，实时监测原始数据叠加修正因子后生成四维数据订正结果；通过与一个基于卡曼滤波方法的阈值进行判断修正因子在时间维度上的稳定性，在某一实施案例中，经判断后修正因子在时间维度上不稳定，则判定站点在该时刻的浓度数据为异常，并对其剔除。
对上述园区内经过初步筛查处理后的污染物浓度数据使用贝塞尔插值函数进行插值处理，结合精准的站点分布数据集得到模拟区域的浓度格点数据分布
对上述插值后的污染物浓度数据生成订正函数，用于分析浓度贝塞尔插值后的区域浓度分布状况，从而确定观测仪器出现校准误差的站点。
对于观测仪器出现校准误差的站点，根据误差站点周边经判定后划为正常站点的大气污染物浓度数据，根据风场数据加权后计算出的浓度年度平均值生成订正参数
具体地，根据浓度订正参数对观测仪器出现校准误差的站点进行订正，计算方法如下：
（20）
式中X 0* 表示对特定站点根据年均值生成的浓度订正参数，X 0 表示订正前浓度值，X 1 、X 2 、X 3 表示各时刻根据订正参数得到的浓度订正值；D 1 、D 2 、D 3 均为计算系数。根据使用上述函数订正后的浓度空间分布进行扩散模态的分析，确定计算园区或工厂内大气污染物排放源强所使用的数值模拟方法。
在大气污染物排放总量反推与统计检验计算模块运行过程中，根据监测数据采集与质控模块运行结束后得出的风场同化结果，对不同区域浓度分布平行于主要的风力扩散方向进行遍历得到污染物的主扩散路径，即采用变分法取极值，拟合大气污染物浓度与路径积分长度函数关系。
发明原理：根据排放源强和大气扩散方程进行浓度时空分布的计算，是偏微分方程领域的基本算法之一，已有了很多成熟稳定而快速的方法。然而从浓度时空分布反推排放源强的反过程推导计算，则一直是数值计算中常见的非适定问题，无法通过上述大气扩散方程的后向求解直接获得，是现有数值计算中存在的技术难题。使用基于伴随算子方法的反推方法时常消耗巨量的算力，同时仍旧需要对排放源强分布的进行强假设，对计算精度和实时性能均造成严重影响。
本算法提供了另一种关于扩散步长的强假设，对排放源的位置和下风向烟羽的时空分布特征进行计算，并据此再运行大气扩散模式，同时对基于观测和扩散模拟的浓度格点化时空分布函数进行逐步长的松驰逼近（Nudging），最终可以在算力较少的影响下得到基本符合实际情况的基于采样监测数据的浓度时空格点化分布数据。解决以往计算过程中存在的问题。
有益效果：本发明通过使用四个算法子模块处理微站监测数据，进而对监测数据做四维订正，结合模拟区域的经纬度，进行大气污染物排放扩散模型的建立与求解，对模拟区域内风场空间进行同化处理，对污染物浓度变化幅度参量进行假设，并将各个站点的日浓度参量对于时间作时间分布函数，最终在大气污染物排放总量反推与统计检验计算模块运行中拟合出大气污染物浓度与路径积分长度函数的关系；
该方法过程对微站单站时间序列进行***误差修正，消除微站组网监测数据站别差异，进而提高其时空代表性。这种具有高可靠性和高实时性的针对微站组网数据进行快速订正及初步处理的方法，可以进一步提高现有和未来城市和重点局地区域微站组网数据的应用可靠性，进而实现对特定区域空气质量监测数据网格化的和高时空精度化，可以单独用于多站位监测数据分析，也可以作为基于气体浓度监测数据快速计算排放源强方法的前置关键步骤。
附图说明
图1是本发明实施例的执行流程图；
图2是本发明使用的自然邻域插值方法使用的流程图；
图3为本发明实施例中微站组网数据订正后进行NO 2 浓度监测数据快速计算排放源强结果示意图；
图4为本发明实施例中微站组网数据订正后进行PM 2.5 浓度监测数据快速计算排放源强结果示意图；
图5为本发明实施例中微站组网数据订正后进行PM 10 浓度监测数据快速计算排放源强结果示意图；
图6为本发明实施例中微站组网数据订正后进行SO 2 浓度监测数据快速计算排放源强结果示意图；
图7为本发明实施例中微站组网数据订正后进行VOC浓度监测数据快速计算排放源强结果示意图。
具体实施方式
本实施例选取地形地物分布复杂的城市和工业园区内，区域内建筑和其他地物的遮挡导致站点局地时空代表性的不足，同时周边排放对于站点的影响往往不可忽略；不同站点周边的微气象和排放环境很难保持一致，这些问题将造成微站数据组网和监测信息时空分布数据的整合过程往往存在显著的不确定性，故采用本发明的方法进行订正。
具体选取约15km×8km的区域，原始数据为2020.1.1至2020.12.31的气象要素和若干污染气体与颗粒物的多站点连续浓度监测数据，时间分辨率为１小时。站点包括两个空气站、29个微站及敏感站。主要记录了TVOC和NO x 的污染物浓度数据；
在本示例中根据实地考察、卫星对比、远程走访的方式对站点在奥维互动地图软件中进行地图标点定位。在同一地图软件中可以确保各个站点的相对经纬度坐标精准，最终得到较精准得订正后站点经纬度数据。
在本示例中区域较小，风力数据有限，风场站点唯一，认为区域内风场为均匀平流，以空气中心站为标准，处理得到单一化风场。
运行本发明的监控数据采集和质控模块，收集园区观测站点的原始风场数据。通过园区观测站点获得原始风场数据，以及NO 2 、PM 2.5 、PM 10 、SO 2 、VOC污染物浓度数据，并对所述原始数据进行初步的统计和整理。
如图1所示，本方法首先进行对监测数据进行四维订正的步骤，包括浓度监测采集数据的粗质控模块，站点微环境影响因子的变分插值模块、浓度时空分布四维同化模块和基于大气扩散模拟的松弛逼近模块，监测数据订正后，进行大气污染物排放扩散模型及其建立求解过程，对模拟区域内风场空间进行同化处理，对污染物浓度变化幅度参量进行假设，并将各个站点的日浓度参量对于时间作时间分布函数，最终在大气污染物排放总量反推与统计检验计算模块运行中拟合出大气污染物浓度与路径积分长度函数的关系。
在本实施例中将获得的站点测得的原始大气污染物浓度数据中的仪器异常值进行查找替换，并进行去风频化长时间尺度数据整合，通过机器学习生成订正权重函数，并对浓度进行集中订正。最终生成最后的大气污染物浓度订正数据。
利用烟羽扩散得物理模型进行该区域大气污染物浓度场的进行判定和分配。
其中对逐站做守恒量Λ值的多源分配进行定标。当站点的Λ值符合多源分配结果时，则保留该站，否则舍弃。
如图2所示，通过订正后站点的经纬度数据，风场同化结果，以及经过多源分配定标后的站点浓度数据对单个排放源进行自然邻点插值方法，生成烟羽浓度空间分布的格点数据。再对多排放源进行叠加，得到多源浓度分布。其中自然邻点插值过程中首先要应用Delaunay三角网格算法对水平空间进行非组织化网格分割，然后对逐个网格应用Voronoi面元算法计算空间权重系数，用于下一步生成烟羽浓度空间分布。
最终得到模拟区域内污染物浓度分布的格点数据，如图3-图7所示，包括实验区域内站点01-站点31，NO 2 、PM 2.5 、PM 10 、SO 2 、VOC的浓度场分布结果以及排放强度月变化趋势结果。
本发明中阐述的方法***适用于在大气污染物集中排放区域内，对通过微站组网监测***得到的大气成分浓度监测数据进行实时订正。
本发明构建了基于在线气象场和大气成分浓度场监控数据，对微站组网的大气成分实时监测数据进行采集和时空一致性订正的数值工程处理方法，来消除现有不同点位仪器之间由于校调不一致或者长期运行中存在的非定向漂移的***误差，减小复杂地理环境对结果的监测影响，并解决由于扩散路径的阻抗差异的存在而导致的问题。
该算法适用于在工业园区或者大气污染物集中排放区域内，通过微站组网监测***得到的大气成分浓度监测数据的实时订正。并基于在线气象场和大气成分浓度场监控数据，对微站组网的大气成分实时监测数据进行采集和时空一致性订正的数值工程处理。
经本实施例的实地检验及实验比对，证明本发明的方法明显提升了工业园区排放精准源头治理的水平。协调人类生产生活和自然生态，为人与自然的和谐相处贡献科技力量。
Claims (9)
1.一种对微站组网数据进行时空一致性订正的方法，使用监测数据处理子模块处理微站监测数据，进而对监测数据做四维订正；结合模拟区域的经纬度，进行大气污染物排放扩散模型的建立与求解；同化处理模拟区域内风场空间，假设污染物日浓度变化幅度参量近似守恒，并将各个站点的污染物日浓度参量作时间分布函数；最终在大气污染物排放总量反推与统计检验计算模块运行中拟合出大气污染物浓度与路径积分长度函数的关系；所述监测数据处理子模块包括浓度监测采集粗质控子模块、站点微环境影响因子的变分插值子模块、浓度时空分布的四维同化子模块和基于大气扩散模拟的松弛逼近子模块；其特征在于：所述浓度监测采集粗质控子模块用于消除直接采样得到的数据集中的随机数据野点；
浓度信息块T i,j,k 是独立空间格点的污染物浓度观测值，表示为：
T i,j,k =T true（i,j,k） +ε i,j,k (1)
式中T true(i,j,k) 表示对应空间格点上的真值，i,j,k是表示空间三个维度上的脚标，ε i,j,k 表示同类测量值中的待检验误差；
令不同空间点位上同一时间的采样数据之间的相关性为零；
E(ε i,j,k )=0 （2）
E(ε i,j,k )是同类测量值中的待检验误差的数学期望值；
同时对于每一台采样仪器，定义污染物浓度的时序观测方差为：
E(ε 2)=δ 2 (3)
ε 2为污染物浓度的时序观测值，δ 2 为污染物浓度的时序观测方差；
对有关点位上的真实浓度观测值T true 的统计信息作出假设：
T true =δ(i,j,k)*T i,j,k (4)
由公式（4）推导得：δ(i,j,k)=1 (5)
由上述公式（5）最终得到：
(6)
式中为浓度真值的时序观测方差，/>为同类测量值中待检验误差的平方；
通过上述公式（6），根据已知的污染物浓度观测值，确定区域范围内的污染物浓度观测值的理论统计方差限，由此对所有采样数据进行误差判定，将落在误差限之外的数值当作野点数据删除处理。
2.根据权利要求1所述的一种对微站组网数据进行时空一致性订正的方法，其特征在于：所述微环境影响因子的变分插值子模块用于消除微环境影响因子，或者对不同站点的监测数据进行订正后使得微环境影响因子保持一致；
求取观测值的时空分布函数的变分极小值，所述时空分布的变分函数是污染物浓度估计值T和其周边若干采样值之间的引力权重函数：
(7)
式中T i,j,k 为独立空间格点的污染物浓度观测值，T为污染物浓度估计值，J(T)为区域周边采样值之间的引力权重， N为采样值的个数；
应用黎曼变分理论证明公式（7）中：
当T= T true 时， (8)
式中J为区域内某一站点周边采样值之间的引力权重；
以大气污染物扩散的高斯分布，某一站点位置受微环境影响后得到的浓度为：
（9）
式中C i,j,k (T||T true )为某站点位置受微环境影响后的污染物浓度；
污染物浓度的真值似然值计算公式为：
（10）
式中为污染物浓度的真值似然值；
推得污染物浓度的极大似然值为：
（11）
式中const为根据具体污染物种类和排放节律所设定的常数， H为空间的观测算子，表示数组规模为n的T的极大似然值；
通过对空间中所有观测点位当前时刻的观测值进行公式（11）的似然函数变分极小值计算，得出对于后验概率极大的污染物浓度观测真值。
3.根据权利要求2所述的一种对微站组网数据进行时空一致性订正的方法，其特征在于：
所述浓度时空分布四维同化子模块用于减少计算误差，所述计算误差包括逐站产生随时间流逝而逐渐积累的***误差，以及不同站点之间发生方式和频率各自不同但是分别逐渐积累增加的***误差；
选取a,b,c三个格点位置，定义此时的污染物浓度观测场为：
（12）
式中X obs 为浓度观测场，为a点的浓度观测场，/>为b点的浓度观测场, />为c点浓度观测场；
观测污染物浓度时周边环境的背景场为：
（13）
X bg 为观测污染物浓度时周边环境的背景场，为在a点观测污染物浓度时周边环境的背景场，/>为在b点观测污染物浓度时周边环境的背景场，/>为在c点观测污染物浓度时周边环境的背景场；
浓度观测向量为：
（14）
Y obs 为浓度观测向量，为一个格点的浓度观测向量，/>为其他格点的浓度观测向量；
则环境背景场、浓度观测向量、浓度观测场之间应有下式关系：
（15）
式中h 1a、 h 1b 和h 1c 分别为在a,b和c点的基于本站点历史浓度数据获得的背景场影响因子；h 2a、 h 2b 和h 2c 分别为在a,b和c点的基于本站以外的周边其他站点历史浓度数据获得的背景场影响因子；
同时背景误差协方差矩阵为：（16）
式中B为背景误差协方差矩阵、b ab 表示a与b两个站点处背景误差的协方差因子；b aa 表示a站点处背景误差的自协方差因子；
最终得到空间浓度的观测算子、浓度观测向量、背景误差之间存在泛函数关系：
（17）
式中J obs 为区域周边采样值之间的引力权重系数，t 0 为初始时刻、t为任一观测时刻，L(t||t 0 ) T 为针对T根据初始时刻标定数据得到的对当前时刻观测值的时间订正函数，表示在任一空间位置i,j,k处关于T的空间订正系数、R为对于H系数做有关距离的二阶订正的修正函数，其取值与位置有关；/>表示在任一位置i,j,k处浓度空间订正值与背景浓度观测向量之间的差；
令公式（17）= 0作为控制函数，求该函数的极小值，求得任一站点对初始时刻旧的污染物浓度观测值与终止时刻污染物浓度新的观测值之间的时空订正值，进而利用其在每个观测值上分别进行求解和订正，最终实现本子模块的订正目标函数。
4.根据权利要求3所述的一种对微站组网数据进行时空一致性订正的方法，其特征在于：
所述基于大气扩散模拟的松弛逼近计算子模块，其松弛方程为：
（18）
上式中c为污染物浓度，τc为松弛尺度值，为势函数，/>为浓度的空间梯度项，为风速向量、c obs 为扩散浓度的观测值、c true 为扩散浓度的理想真值，Δ为考虑松弛之后的时间变量、f为地转偏向力的水平分量。
5.根据权利要求1所述的一种对微站组网数据进行时空一致性订正的方法，其特征在于：对监测数据做四维订正，是对原始观测浓度数据进行空间滤波处理、时序特征的多点一致性订正，以及根据布点空间相对位置和复杂地物分布特征对扩散路径上的阻抗差异进行的归一化订正。
6. 根据权利要求1所述的一种对微站组网数据进行时空一致性订正的方法，其特征在于：模拟区域的经纬度，是对园区地图定位确定原站点的经纬度坐标 (lon,lat)，经过卫星比对、实地考察、远程走访的方法，对地图上园区位置进行标点，得到精准的站点分布数据集。
7.根据权利要求1所述的一种对微站组网数据进行时空一致性订正的方法，其特征在于：所述监测数据处理子模块对原始监测风场数据的空间分辨率进行判别和质控，对有明显区别的不同情况进行归一化交叉验证处理，监测数据的采集与质控完成后的进行大气污染物排放和扩散模型的建立与求解；
所述原始监测风场数据包括站点实时测得的风速大小、风场方向、湍流强度的测算质量、测站站点的数量及数据；
当测站数量在20个及以上时，使用二维空间内的散点插值算法对风场数据的同化；
当风场测站数量在10个及以内时，则对各个站点的数据质控处理过程进行动态权重分配，根据单一化风场的模拟结果比较选择处最优的权重系数，得出最优的单一化风场模拟结果；
若风场测站仅为一个，则直接使用该站点测值生成模拟空间内均匀单一的风场分布；
监测数据采集与质控模块运行过程中对原始浓度数据进行质量控制处理，剔除错误数据、空白数据及明显不符合物理规律的数据。
8. 根据权利要求1所述的一种对微站组网数据进行时空一致性订正的方法，其特征在于：还包括大气扩散相关气象和化学系数计算过程，是根据监测数据采集和质控模块运行结束后得出的风场同化结果，同时依据烟羽扩散的物理模型，在受同源影响的各个站点区域范围内，依据模拟尺度区域较小的假设，设定污染物日浓度变化幅度参量Λ为近似守恒量：
（19）
式中ΔX表示站点污染物日浓度的浮动幅度值，表示站点污染物日浓度均值；Λ用于判断某站点当日的污染源排放情况；计算过程结束后将各个站点的污染物日浓度变化幅度参量Λ相对于时间t做分布函数图， 将分布函数进行空间卡曼滤波，根据滤波修正函数或者沿时间维滑动平均的逐站点修正因子，实时监测原始数据叠加修正因子后生成四维数据订正结果；通过与一个基于卡曼滤波方法的阈值，判断修正因子在时间维度上的稳定性，当判断修正因子在时间维度上不稳定时，则判定站点的浓度数据为异常，并对其剔除。
9. 根据权利要求1所述的一种对微站组网数据进行时空一致性订正的方法，其特征在于：根据浓度订正参数对观测仪器出现校准误差的站点进行订正：
（20）
式中X 0* 表示站点根据年均值生成的浓度订正参数，X 0 表示订正前浓度值，X 1 、X 2 、X 3 表示各时刻根据订正参数得到的浓度订正值；D 1 、D 2 、D 3 均为计算系数；
根据订正后的浓度空间分布进行扩散模态的分析，确定计算园区或工厂内大气污染物排放源强。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202410130460.XA CN117669270B (zh) | 2024-01-31 | 2024-01-31 | 一种对微站组网数据进行时空一致性订正的方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202410130460.XA CN117669270B (zh) | 2024-01-31 | 2024-01-31 | 一种对微站组网数据进行时空一致性订正的方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN117669270A CN117669270A (zh) | 2024-03-08 |
| CN117669270B true CN117669270B (zh) | 2024-04-19 |
Family
ID=90068371
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202410130460.XA Active CN117669270B (zh) | 2024-01-31 | 2024-01-31 | 一种对微站组网数据进行时空一致性订正的方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN117669270B (zh) |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111537023A (zh) * | 2020-05-13 | 2020-08-14 | 浙江大学 | 一种工业园区大气污染物扩散模拟与溯源方法 |
| CN116258101A (zh) * | 2023-05-09 | 2023-06-13 | 南京信息工程大学 | 一种基于气体浓度监测数据快速计算排放源强的方法 |
| CN116930423A (zh) * | 2023-07-25 | 2023-10-24 | 北京工业大学 | 一种空气质量模型模拟效果的自动验证评估方法及*** |
-
2024
- 2024-01-31 CN CN202410130460.XA patent/CN117669270B/zh active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111537023A (zh) * | 2020-05-13 | 2020-08-14 | 浙江大学 | 一种工业园区大气污染物扩散模拟与溯源方法 |
| CN116258101A (zh) * | 2023-05-09 | 2023-06-13 | 南京信息工程大学 | 一种基于气体浓度监测数据快速计算排放源强的方法 |
| CN116930423A (zh) * | 2023-07-25 | 2023-10-24 | 北京工业大学 | 一种空气质量模型模拟效果的自动验证评估方法及*** |
Also Published As
| Publication number | Publication date |
|---|---|
| CN117669270A (zh) | 2024-03-08 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN112905560B (zh) | 一种多源时空大数据深度融合的空气污染预测方法 | |
| CN110046771B (zh) | 一种pm2.5浓度预测方法与装置 | |
| CN109784552B (zh) | 一种基于Re-ESF算法的空间变系数PM2.5浓度估算模型的构建方法 | |
| CN109916788B (zh) | 一种区分不同区域排放变化和气象条件变化对pm2.5浓度影响的方法 | |
| CN112699959B (zh) | 基于能量泛函模型的多源多尺度降水数据融合方法和装置 | |
| CN110232471B (zh) | 一种降水传感网节点布局优化方法及装置 | |
| CN114280695A (zh) | 一种空气污染物监测预警方法及云平台 | |
| CN112800603B (zh) | 一种基于集合最优插值算法的大气环境数据同化方法 | |
| CN111353634B (zh) | 一种基于模式初始场同化的二氧化硫排放源优化方法 | |
| CN109858686A (zh) | 一种基于EnKF的地面排放清单反演优化方法 | |
| Ma et al. | Statistical spatial-temporal modeling of ambient ozone exposure for environmental epidemiology studies: A review | |
| Wang et al. | Inverse modeling of black carbon emissions over China using ensemble data assimilation | |
| CN116223395A (zh) | 近地表痕量气体浓度反演模型及反演方法 | |
| CN114118613A (zh) | 一种焦化土壤污染空间分布预测优化方法及*** | |
| CN114077798B (zh) | 估算小尺度区域大气污染物网格浓度的方法 | |
| CN116258101B (zh) | 一种基于气体浓度监测数据快速计算排放源强的方法 | |
| CN114724647A (zh) | 一种高时空分辨率co2通量反演***及方法 | |
| CN113156395B (zh) | 气溶胶激光雷达数据融合方法及*** | |
| CN113834902B (zh) | 一种基于四维变分同化的二氧化硫排放源反演方法 | |
| CN117669270B (zh) | 一种对微站组网数据进行时空一致性订正的方法 | |
| CN116466368B (zh) | 基于激光雷达和卫星资料的沙尘消光系数廓线估算方法 | |
| CN115840908B (zh) | 基于lstm模型的微波链路构建pm2.5三维动态监测场的方法 | |
| CN114371519B (zh) | 一种基于无偏差亮温的地基微波辐射计大气温湿廓线反演方法 | |
| CN114219176A (zh) | 一种多源算法集成的太阳辐射预报***及方法 | |
| CN114757103A (zh) | 基于时空光梯度增强机的地表o3浓度估算方法 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |