【X-STILT模型第一期】X-STILT 模型概述

目录X-STILT 模型概述支持的观测平台与气体物种GitHub 仓库中的内置核心脚本/功能X-STILT 模型的下载安装一、 下载与安装模型 (Download and install model)二、 前置条件与数据准备 (Prerequisites)1. 依赖卫星观测的柱浓度模拟 (For SATELLITE-dependent column simulation)2. 无卫星依赖的理想模拟 (For IDEAL column simulation)3. 气象场数据 (Meteorological fields) - 绝对必选4. (可选) 计算柱浓度增量 (For calculating column enhancements)5. (可选) 执行传输误差分析 (For performing transport error analyses)6. (可选) 执行排放误差分析 (For performing emission error analyses)获取“柱浓度足迹Obtain column footprint一、 核心概念什么是“柱浓度足迹 (Column Footprint)”二、 运行流程run_xstilt.r 中的 6 个关键步骤确定背景柱浓度Determine background XCO2一、 核心概念M3 背景计算法与前向模拟二、 运行流程compute_bg.r 中的逻辑开关 (Logical Flags)估算水平和垂直的传输误差Transport Errors参考X-STILTColumn-Stochastic Time-Inverted Lagrangian Transport model是一个专门用于大气微量气体“柱浓度”的大气传输模型。它在地学、大气科学和碳循环研究中非常重要。GitHub -X-STILT: Column-Stochastic Time-Inverted Lagrangian Transport modelX-STILT 模型概述支持的观测平台与气体物种模型目前可以与多种先进的卫星和地基观测数据结合使用OCO-2 3 (XCO2)美国宇航局的轨道碳观测者卫星专门高精度测量大气中的二氧化碳柱浓度 (XCO2)。TROPOMI (CO, CH4, NO2)欧洲空间局的对流层观测仪X-STILT 可以结合它来模拟一氧化碳 (CO)、甲烷 (CH4) 以及二氧化氮 (NO2)。其中处理 NO2 的特殊模块被称为 STILT-NOx因为 NO2 的化学性质活泼需要特殊的化学传输计算。TCCON总碳柱观测网。这是一个全球性的地基傅里叶变换光谱仪网络用于验证卫星数据。X-STILT 也支持模拟 TCCON 站点观测到的 CO2, CO, N2O, CH4 柱浓度对应脚本 run_xstilt_tccon.r。EM27/SUN一种便携式的地基光谱仪目前相关代码还在开发/建设中。Ideal column simulations (理想柱浓度模拟)即使没有特定的真实卫星数据也可以输入经纬度和时间进行理想化的模拟计算不过这种情况下不会应用卫星的平均核函数 AK即 Averaging Kernel它代表卫星对不同高度气体的敏感度。GitHub 仓库中的内置核心脚本/功能这段代码库包含了几个用于不同科学任务的 R 语言脚本1、计算后向轨迹与“足迹” (run_xstilt.r)功能从一个大气柱向后推算轨迹backward trajectories并计算“足迹column footprint”。什么是足迹 在大气科学中“足迹”代表了某一次观测对上风向地表排放的敏感度。有了足迹地图就知道某地的高浓度气体是由哪些区域的工厂或城市排放吹过来的。2、计算前向轨迹与背景值 (compute_bg.r)功能从站点周围的区域向前模拟轨迹forward-time trajectories并结合卫星观测数据计算出大气的“背景值background”。意义要知道一个城市的排放量必须先知道干净空气的浓度背景值然后用观测到的总浓度减去背景值得到的就是该城市的“增强量enhancement”。3、估算不确定性 / 误差 (run_xstilt.r)功能估算风场wind和行星边界层高度PBL的误差并把这些气象误差转化为最终计算出的 XCO2 的误差。这对于科学研究中评估结果的可靠性至关重要。4、模拟不同卫星在特定探空点的浓度 (run_sim_multi.r)功能针对 OCO-2/3 和 TROPOMI 卫星的具体观测点soundings正向模拟出二氧化碳、一氧化碳、甲烷和二氧化氮的柱浓度数值。X-STILT 模型的下载安装一、 下载与安装模型 (Download and install model)X-STILT 并不是一个完全从零开始编写的代码库它的底层依赖于 STILT-R version 2 和 HYSPLITv5美国 NOAA 开发的权威轨迹模型。在 GitHub 上这些底层模型被作为“子模块 (submodule)”嵌套在 X-STILT 中。STILT模型的安装可参考另一博客-【STILT模型实操第1期】STILT模型安装环境配置及安装cd/data/wanzhou/STILT/gitclone--recursivehttps://github.com/uataq/X-STILT.git–recursive 的作用是不仅下载 X-STILT 自己的代码还会自动一并下载所有关联的子模块代码。只需要打开核心脚本 run_xstilt.r 运行它的前几行代码。它会自动检测并安装所需的 R 包前提是你的 R 语言版本必须大于 3.5.0。同时它会调用系统的 Fortran 编译器自动编译底层的 HYSPLITv5 fortran 核心程序sourceactivate stilt_envcd/data/wanzhou/STILT/X-STILT/ Rscript run_xstilt.r错误总结1、库包安装package ‘geosphere’# 方式1condainstall-cconda-forge r-geosphere# 方式2install.packages(geosphere)二、 前置条件与数据准备 (Prerequisites)运行 X-STILT 模拟大气传输和气体浓度需要输入大量的外部数据。这部分根据你的研究目的分为了“必选”和“可选”两类1. 依赖卫星观测的柱浓度模拟 (For SATELLITE-dependent column simulation)如果想结合真实的卫星数据来算比如研究某天某卫星飞过某城市时拍到的画面所需数据下载 OCO-2 Lite 级别数据或 TROPOMI Level 2 级别数据。参数设置需要在代码中修改对应的传感器类型 (obs_sensor)、数据版本 (obs_ver)、观测气体种类 (obs_species) 以及数据存放路径 (oco_path或trp_path)。核心机制平均核函数与垂直权重模型会读取卫星文件里的Averaging Kernels (AK, 平均核函数)。由于卫星对大气层不同高度的气体敏感度不同AK 代表了这种敏感度。X-STILT 会自己计算气压权重函数把不同高度释放的空气块三维的足迹结合卫星的 AK按垂直方向压缩最终生成一个等效的“二维柱浓度足迹column footprints”。2. 无卫星依赖的理想模拟 (For IDEAL column simulation)如果不想用真实的卫星数据只想做理论研究比如假设某年某月某日在特定地点有个观测点所需数据准备一个包含接收器观测点信息的.csv文件至少要有经度 (lon)、纬度 (lat) 和时间 (time)。示例代码库里提供了一个名为receptor_demo.csv的模板供参考。3. 气象场数据 (Meteorological fields) - 绝对必选无论哪种模拟都需要知道“风是怎么吹的”这依赖于气象场数据。格式要求必须是ARL 格式NOAA 空气资源实验室使用的特定格式。数据来源可以直接从 NOAA 网站下载现成的 ARL 气象数据如 NARR, GDAS。如果想使用自己跑的WRF (天气研究与预报模型)的高精度气象结果必须先用工具把 WRF 的输出格式转换为 ARL 格式。详细使用教程可参考另一博客-【STILT模型实操第2期】运行 STILT 模型提供 WRF 示例。4. (可选) 计算柱浓度增量 (For calculating column enhancements)如果不仅想知道“风往哪吹”还想算出“污染源到底让这里的二氧化碳增加了多少 ppm”所需数据需要输入排放清单 (Emission inventories)。默认推荐模型默认使用ODIAC一个全球高分辨率的化石燃料 CO2 排放清单。ODIAC的详细介绍可参考另一博客-【CO2数据集】全球1 km分辨率化石燃料 CO₂ 排放数据集 ODIAC。需要下载 1km 分辨率的 ODIAC.tif格式文件并在代码中指定odiac_path。5. (可选) 执行传输误差分析 (For performing transport error analyses)模型算出来的风和实际的风总有误差。如果想评估这种误差如何影响你的结果风场误差数据需要下载NOAA 探空仪气象气球数据 (Radiosonde data)作为实际风速风向的参考。要求下载 FSL 格式且风速单位选择 0.1 m/s (tenths of m s-1)。背景浓度数据需要下载Carbon-Tracker (CT-NRT)数据。这是一种全球碳同化系统它能提供除了人类排放之外的全球二氧化碳大背景浓度即自然环境的总 CO2。6. (可选) 执行排放误差分析 (For performing emission error analyses)不同的科学家做的排放清单本身也是有差异的这会导致结果的不确定性操作方式你可以同时引入多种“自下而上Bottom-up”的排放清单比如FFDAS和EDGAR。通过集合对比不同的清单来评估先验排放数据本身的误差对最终模拟结果的影响。获取“柱浓度足迹Obtain column footprint一、 核心概念什么是“柱浓度足迹 (Column Footprint)”在模型中足迹的单位是[ppm / (umol m-2 s-1)]。物理意义它代表了“源-受体敏感度source-receptor sensitivities”。简单来说它是一张敏感度地图告诉你地面上某个特定区域每排放一定量1 μ mol m − 2 s − 1 1 \, \mu\text{mol m}^{-2} \text{s}^{-1}1μmol m−2s−1的气体会导致卫星观测点受体的大气柱浓度增加多少 ppm。数学意义它本质上是一个“雅可比矩阵Jacobian Matrices”即浓度增量结果与排放通量原因之间的偏导数矩阵。有了这张“足迹地图”只要把它和真实的排放清单比如工厂、城市的排放数据相乘就能计算出观测点实际增加的温室气体浓度。二、 运行流程run_xstilt.r 中的 6 个关键步骤用户主要通过修改主脚本run_xstilt.r来初始化模型和参数。流程分为以下 6 个步骤STEPSTEP 1: 选择卫星过境数据 (Select a satellite overpass)需要告诉模型去模拟哪一次卫星观测。方式可以手动输入具体信息也可以让模型自动搜索。自动搜索机制模型会在设定的一个大空间范围例如2 ∘ × 2 ∘ 2^\circ \times 2^\circ2∘×2∘和一个核心近场范围例如0.6 ∘ × 0.6 ∘ 0.6^\circ \times 0.6^\circ0.6∘×0.6∘内寻找所有包含有效探测点soundings的卫星过境事件。STEP 2: 选择模拟类型 (Modify logical flags)这里有多个“开关”逻辑变量设为T代表 True/开启F代表 False/关闭需要根据研究目的打开相应的开关run_trajec T或run_foot T【核心功能】计算后向轨迹并生成经过垂直加权的柱浓度足迹。run_hor_err T进行水平传输误差分析前提是已经计算了风场误差。run_wind_err T计算模型风场相对于探空气球Radiosonde数据的误差。run_ver_err T进行垂直传输误差分析。做法是在 STEP 4 中设定比例因子人为地把“混合层高度边界层高度”调高或调低看这对结果有多大影响。run_sim T【模拟最终浓度】基于足迹来模拟化石燃料产生的 CO2 柱浓度增量FFCO2 XCO2。注意运行这一步时需要已经算好了足迹并且要把前面算轨迹和足迹的开关run_trajec和run_foot关掉以节省计算资源。run_emiss_err T进行先验排放清单的误差分析同样需要先算好足迹。STEP 3: 配置柱状接收器参数 (Modify parameters for placing column receptors)因为卫星看的是一整根“大气柱”模型需要把这根柱子切分成很多层来模拟。需要设置这根柱子的最大高度、相邻两个垂直层之间的间隔米、每一层释放的粒子数量以及接收器的位置。近场优化X-STILT 默认会在距离污染源或观测点较近的“近场区域NEAR-FIELD”放置更密集的接收器因为近距离的扩散过程对浓度影响最直接、也最剧烈需要更高的精度。STEP 4: 修改气象场设置 (Modify meteorological fields)指定模型使用哪里的气象数据如 ARL 格式的风场。这部分配置与底层的 STILT 模型通用可以参考 STILT 的官方文档。STEP 5: 设置足迹的水平分辨率和空间范围 (Modify horizontal resolution and spatial extent)定义最终生成的“敏感度地图”有多大、网格有多细。多重足迹生成X-STILT 非常灵活基于同一套粒子轨迹它可以同时生成多套不同格式的足迹比如将小时数据积分汇总的足迹或者是逐小时展开的动态足迹。如果需要生成第二套或第三套足迹只需在此处修改参数即可无需重新跑耗时的轨迹计算。STEP 6: 设置 SLURM 并行计算 (Set up SLURM parallel computing)由于大气拉格朗日轨迹计算尤其是释放大量粒子时非常消耗计算资源模型内置了对 SLURM 调度系统的支持。允许在高性能计算集群超算中心上进行并行计算大幅缩短运行时间。确定背景柱浓度Determine background XCO2这段内容详细讲解了 X-STILT 模型中一个极其关键的步骤如何确定背景柱浓度Determine background XCO2。在温室气体研究中卫星如 OCO-2 或 TROPOMI观测到的是大气中的总浓度。为了算出某个城市或工厂到底排放了多少温室气体科学家必须用“总浓度”减去“干净空气的浓度”这个干净空气的浓度就是背景值Background。模型使用核心脚本compute_bg.r来完成这项工作。以下是对该部分内容的详细拆解一、 核心概念M3 背景计算法与前向模拟M3 (Method 3) 特定过境背景值这是 Wu et al. (2018) 论文中提出的一种动态背景计算方法。因为风向每天都在变卫星每次过境时的背景浓度也是不一样的所以需要为每一次观测单独计算背景值。计算原理寻找城市羽流 Urban Plume模型从城市或污染源中心位置随着时间推移向前forward-time释放虚拟的空气块代表污染物。这些空气块会随风飘散形成一条“污染带”这在大气科学中被称为城市羽流Urban Plume。模型使用二维核密度估计2D kernel density这种统计方法在地图上画出这条羽流的清晰边界。既然知道了“脏空气羽流”在哪里那么卫星在羽流之外周边拍到的数据自然就是干净的“背景值”了。二、 运行流程compute_bg.r中的逻辑开关 (Logical Flags)与前面算足迹类似需要在compute_bg.r脚本中修改以下几个关键的T / F开关来控制模型跑哪些步骤run_trajec T计算前向轨迹划定污染范围功能强制模型从站点周围的一个小方块区域由参数box.len控制大小开始计算前向轨迹。机制在这个模式下模型会按照你设定的时间间隔参数dtime.sep源源不断地释放空气块以此来模拟污染物的持续排放和扩散过程。run_trajec T run_for_err T带有误差的前向轨迹计算功能在计算前向轨迹的同时加入风场误差分析。意义风预报往往是不准确的这会导致算出来的“羽流”位置有偏差。打开这个开关模型会引入一个风场误差参数siguverr通过模拟风向的随机波动来评估羽流位置的不确定性。run_wind_err T计算模型风场误差功能将模型使用的气象风场数据与实际的气象气球探空仪 radiosonde观测数据进行对比算出模型风场的具体误差到底有多大。需要事先下载探空仪数据。run_bg T【核心步骤】计算最终的背景值功能结合之前算出的前向轨迹羽流位置和实际的卫星数据OCO-2/3 或 TROPOMI正式计算背景浓度。自动绘图该脚本内置的函数会自动生成可视化图表。图上会标出哪些是羽流被污染的增强区域并假设羽流外围上/下/左/右即北/南/东/西方向的数据为干净的背景区域。数据存储如果设置了writeTF T开启写入最终估算出的背景数值将作为结果保存在你设定的路径file.path(store_path, fn)中供后续计算排放量使用。估算水平和垂直的传输误差Transport Errors在大气模拟中最大的不确定性来源往往是“风”。如果气象模型提供的风速或风向哪怕只有一点点偏差随着时间的推移模型算出的“污染源”位置就会与实际位置产生巨大的偏离。这段内容就是教你如何量化这种“风的偏差”对最终二氧化碳浓度XCO2结果的影响。一、 核心概念X-STILT 的误差估算方法通常科学家评估误差的方法是“模式集合model ensembles”也就是用好几套不同的气象数据比如欧洲的、美国的、日本的分别跑一遍模型看看结果差多少。但 X-STILT 采用了不同的方法误差传播Error Propagation。它的原理是不依赖多套气象数据而是直接将真实世界中风场u-v风速的随机误差基于特定的空间相关长度和时间尺度注入到模型中。模型会将这些风场误差一步步“传播propagates”从粒子释放的高度推算到观测点的位置最终转化为卫星过境时 XCO2 浓度的误差。这种方法最早由 Lin 和 Gerbig (2005)提出现在被扩展应用到了“柱浓度”的研究中。这张图直观地展示了“风场误差”在后向轨迹模拟中的威力。图中展示了从阿拉伯半岛的一个观测点白色圆点 Receptor向后倒推 48 小时的粒子轨迹。橙色点WRF-STILT NO errors这是没有加入风场误差的“标准轨迹”。可以看到橙色粒子形成了几条相对集中的、清晰的“带状”路径。模型“认为”气流是沿着这些清晰的路径从非洲大陆吹过红海到达观测点的。蓝色点WRF-STILT Transport errors这是加入了风场误差后的轨迹。模型在每一步计算时都加上了随机的风速/风向扰动。可以观察到两个重要现象扩散与发散Dispersion相比于集中的橙色点蓝色点分布得非常散形成了一大片“云”。这代表了现实中风的不确定性——真实的污染源可能并不在橙色那条细线上而是可能位于广阔的蓝色区域内的任何地方。误差随时间累积请看时间标记-12hrs, -24hrs, -36hrs, -48hrs。在 -12 小时的时候蓝点和橙色点还靠得比较近但退到 -48 小时的时候蓝点的分布范围已经极其广阔了。这说明追溯的时间越长风场误差导致的源区定位不确定性就越大。图片注记Caption补充图注提到这里只画了 3 公里以下的粒子。对于每一个时间步长处于更高垂直高度的粒子比靠近地面的粒子更偏西因为高空风速通常更快。二、 代码执行逻辑根据文本描述如果你想在自己的研究中复现图片中的效果并计算误差需要通过修改 run_xstilt.r 中的开关分两步走1、生成带有误差的轨迹蓝色点首先必须在 STEP 2 中将 run_hor_err 设置为 T (True)。这会指示模型在计算轨迹时引入风场误差统计数据从而生成如上图所示的大范围发散的“蓝色轨迹集”。2、计算最终的浓度误差XCO2 errors生成完这些蓝色轨迹后你需要关闭 run_hor_err然后开启 run_sim T。此时模型会对比两套轨迹无误差的橙色轨迹 vs 有误差的蓝色轨迹并结合排放清单最终计算出这两种轨迹会导致观测点Receptor的 XCO2 浓度产生多少 ppm 的偏差。参考1、论文-J2005-Accounting for the effect of transport errors on tracer inversions