wstats_mod.F90

路径

LMDZ.MARS\libf\phymars\wstats_mod.F90

所属目录 / 模块

libf/phymars

文件定位

wstats_mod.F90 实现传统统计诊断文件 stats.nc 的生成流程。它不是逐时间步保存瞬时场，而是在每个 sol 固定保存 istime=12 个本地时刻槽位，对调用方传入的 2D/3D 字段累加和、平方和，最后由 mkstats 原地改写为平均值和总标准差。physiq_mod.F 中的 WSTATS 段按固定顺序调用 wstats 输出地表、辐射、动力、尘埃、水循环、CO2 和 tracer 诊断，并在 lastcall.and.callstats 时调用 mkstats 完成归一化。

模块级开关 callstats 由 conf_phys.F 读取 callstats 配置键初始化。若 callstats=.false.，wstats 立即返回；若运行目录存在 stats.def，则只有列在文件中的变量会进入 stats.nc。源码只处理结构化经纬网格，grid_type==unstructured 时直接返回。

定义的符号

符号	类型	行号	作用
`wstats_mod`	module	1	保存 `stats.nc` 统计输出开关、采样状态和相关例程。
`callstats`	saved logical	9	控制是否生成 `stats.nc`，由 `conf_phys` 读取配置。
`istats`	saved integer	13	每隔多少个物理时间步累加一次统计，由 `inistats` 计算。
`istime`	integer parameter	17	每个 sol 的统计时刻槽位数量，固定为 12。
`count(istime)`	saved integer array	19	每个时刻槽位累计了多少个样本。
`wstats`	subroutine	24	按变量逐次累加 `stats.nc` 中的和与平方和字段。
`inistats`	subroutine	479	创建 `stats.nc`，写经纬度、高度、hybrid sigma 坐标和采样间隔。
`inivar`	subroutine	628	新变量首次创建时写入初始样本。
`def_var_stats`	subroutine	749	在 `stats.nc` 中定义变量并写 `title/units` 属性。
`mkstats`	subroutine	813	结束时把累加量原地转换为平均值和标准差。

依赖的模块

use 模块	only 列表	用途	待确认
`mod_phys_lmdz_para`	`is_mpi_root`, `is_master`, `klon_mpi_begin`	限制 NetCDF 写入到 master/root；`klon_mpi_begin` 在当前源码中未使用。	`klon_mpi_begin` 是否历史遗留待确认。
`mod_phys_lmdz_transfert_para`	`gather`, `bcast`	并行收集字段，广播 `stats.def` 筛选信息。	-
`mod_grid_phy_lmdz`	`klon_glo`, `Grid1Dto2D_glo`, `nbp_lon`, `nbp_lat`, `nbp_lev`, `grid_type`, `unstructured`	网格尺寸、结构化网格检查和物理网格到经纬网格重排。	-
`vertical_layers_mod`	`ap`, `bp`, `aps`, `bps`, `preff`, `pseudoalt`, `presnivs`	`inistats` 写垂直坐标和 pseudo-altitude。	-
`nrtype`	`pi`	经纬度弧度转度。	-
`time_phylmdz_mod`	`daysec`, `dtphys`	计算每 sol 物理步数和 `istats`。	-
`regular_lonlat_mod`	`lon_reg`, `lat_reg`	写 `stats.nc` 经纬度坐标。	-

调用的关键例程

被调用例程	所在模块 / 文件	调用位置	作用
`inistats`	本文件	`wstats` 首次有效调用	创建 `stats.nc` 并写固定坐标。
`inivar`	本文件	新变量首次创建	为平均量和平方和变量写第一帧样本。
`def_var_stats`	本文件	`inistats` 与 `wstats`	定义坐标或诊断变量并写属性。
`abort_physic`	外部错误处理	多处错误分支	`stats.def` 过长、时间步不整除、NetCDF 读写失败时中止。
`gather`, `Grid1Dto2D_glo`	并行/网格模块	并行 `dim=2/3` 路径	把局地物理网格字段收集并重排为全局经纬网格。
NetCDF `NF_*`	NetCDF F77 接口	多处	创建、打开、关闭、定义、读写 `stats.nc`。

输入

输入	来源	类型 / 维度	单位	含义
`callstats`	`conf_phys.F` / `callphys.def`	logical	-	是否启用 `stats.nc`。
`stats.def`	运行目录可选文件	text file	-	若存在，仅输出列出的变量。
`ngrid`	`wstats` 调用方	integer	-	当前物理网格点数。
`nom`, `titre`, `unite`	`wstats` 调用方	character	-	变量名、标题和单位。
`dim`	`wstats` 调用方	integer 2/3	-	诊断变量维度码；源码实际只处理 2D 和 3D。
`px(ngrid,nbp_lev)`	`wstats` 调用方	real array	variable-dependent	待累加诊断场。
`daysec`, `dtphys`	时间模块	real	seconds	用于计算一个 sol 有多少物理步，以及每个统计槽位间隔。

输出

输出	去向	类型 / 维度	单位	含义
`stats.nc`	当前运行目录	NetCDF file	-	统计日诊断文件。
`Time`	`stats.nc`	`istime` 个写入点	hours	`mkstats` 写入 `i*24/istime` 的日内时刻。
`<nom>`	`stats.nc`	2D 或 3D + `Time`	`unite`	运行结束后为样本平均值；结束前暂为样本和。
`<nom>_sd`	`stats.nc`	2D 或 3D + `Time`	`unite`	运行结束后为标准差；结束前暂为平方和。

共享状态与副作用

callstats、istats、count、firstvar、firstcall、indx、step 和 stats.def 筛选状态均为 SAVE / THREADPRIVATE。
第一次有效 wstats 调用以 NF_CLOBBER,NF_64BIT_OFFSET 创建 stats.nc，会覆盖同名文件。
wstats 对每个变量反复打开、读、写、关闭 stats.nc；新变量首次出现时先创建 <nom> 和 <nom>_sd 两个变量。
mkstats 会原地覆盖 stats.nc 中的累加量，把 <nom> 改成均值，把 <nom>_sd 改成标准差。
grid_type==unstructured 直接返回，不生成统计输出。

核心逻辑

wstats 首先检查 callstats 和 grid_type，关闭统计或非结构网格时直接返回。
首次检查 stats.def：master 读取变量名单，再广播 stats_def/n_name_def/name_def。
firstcall 阶段确定 firstvar。若存在 stats.def，第一个被列入的变量才会成为 firstvar；否则第一个调用变量就是 firstvar。
首次有效调用执行 inistats：检查 daysec/dtphys 是整数物理步数，且能被 istime=12 整除，计算 istats=nsteppd/istime，创建 stats.nc 和固定坐标。
每次回到 firstvar 时递增 step；只有 mod(step,istats)==0 时才继续累加。
按 stats.def 二次筛选变量；未列入变量返回。
对 2D/3D 字段进行并行 gather 或串行经纬重排，得到带冗余经度点的 dx2/dx3。
若变量尚不存在，定义 <nom> 与 <nom>_sd，用 inivar 写第一帧并返回。
若变量已存在，读取当前槽位累加量，更新 count(indx)，累加 sum += dx 和 square += dx**2。
最后 mkstats 遍历所有变量和 istime 槽位，写 Time，并把 sum/count 与 sqrt(max(0,square/count-mean**2)) 写回文件。

伪代码

wstats(ngrid, nom, titre, unite, dim, px):
  if not callstats or unstructured grid:
    return
  on first stats.def check:
    master reads optional stats.def
    broadcast selected names
  on first valid variable:
    choose firstvar
    call inistats()
    allocate mean/sd/dx buffers
  if nom == firstvar:
    step += 1
  if step is not on an istats boundary:
    return
  if stats.def exists and nom is not listed:
    return
  gather/recast px to lon-lat dx field
  open stats.nc
  if variable does not exist:
    define nom and nom_sd
    initialize first sample with inivar()
    close and return
  update count and current time-slot index
  read previous sum and square sum
  write sum + dx and square + dx**2
  close stats.nc

mkstats():
  increment last count slot
  for each variable pair in stats.nc:
    for each of 12 time slots:
      mean = sum / count
      sd = sqrt(max(0, square / count - mean**2))
      overwrite stats.nc values

参与的主题流程

主题	参与方式
传统统计诊断	生成 `stats.nc`，供离线分析和后处理工具读取统计日字段。
CO2 / 水循环 / 尘埃诊断	`physiq` 在相关分支调用 `wstats` 输出柱积分、混合比、光学厚度和微物理诊断。
配置系统	`conf_phys` 读取 `callstats`，`stats.def` 在运行目录中进一步筛选输出变量。

写法特点

stats.nc 中变量在运行结束前不是最终均值，而是累加和/平方和；必须执行 mkstats 后才是页面语义上的统计量。
Time 不是模拟绝对时间，而是统计日内 12 个槽位的小时值。
2D 变量写入使用 3 个维度，源码把 start/sizes 的第 4 个元素设为 0；这是为了复用固定长度数组，需依赖 NetCDF F77 调用只读取实际维数。
firstvar 绑定调用顺序；若 stats.def 排除了 physiq 中靠前变量，firstvar 会变成第一个被允许的变量。

复现要点

daysec/dtphys 必须接近整数，且 nsteppd 必须能被 istime=12 整除，否则 inistats 直接中止。
mkstats 必须在最后一次物理调用后执行；源码在 physiq_mod.F 的 lastcall.and.callstats 分支调用。
若使用 stats.def，必须保持每个物理步的首个允许变量一致，否则 step/indx/count 槽位推进会错位。
当前源码实际只处理 dim=2 和 dim=3；不要把 0D/1D 字段直接交给 wstats。

待确认

klon_mpi_begin 在 wstats 的 use mod_phys_lmdz_para 中导入但未见使用，可能是历史遗留。
mkstats 中 count(istime)=count(istime)+1 用来补最后一步计数；该补偿对所有调用序列是否都严格正确，需结合 physiq 最后一次调用顺序运行核验。