雷达降水估测中变分方法的改进研 究 慕熙昱 1 徐琪 2 刘 国庆 31.江苏省气象科学研究所 210008 南京2.民航华东空中交通管理局江苏分局 210000 南京3.南京工业大学 210000 南京摘 要: 雷达 作为 定量 测量 降水 的工 具, 在探 测降 水分 布和 结构 上比 较准 确, 但是 在精 度上 还是 低于 雨量计 本文 利 用变分 法, 综 合多普 勒天 气 雷达数 据及 雨 量计数 据进 行 定量降 水估 测 ,其本 质是 通 过调整模 式解以 达到 模 式解与 观测 数 据之间 最小 二 乘条件 最 开 始在利 用变 分 方法校 正雷 达 降雨场 的过 程 中,对变 分方程 中的 系 数是经 验选 取 的,但 是这 种 经验关 系有 很 大的局 限性 在对变 分方 程 求解过 程中 还 有一个问 题就是边界的处理 在对变分方程求解的过程中, 采用交叉验证法 交叉验证( cross-validation) 理论 是:对 某一 点 来说, 用其 余若 干个观 察值 对 该点进 行估 计, 每一点 都重 复 这个过 程, 用 估计值 与真 实 值进行 对比得 误差 参 量,最 后进 行处 理,使 这些 误差 参量达 到要 求 。
该方 法基 于 定理: 椭圆 型 偏微分 方程 中 ,区域 内的解 连续 依 赖于边 界 首 先利用 变分 方 程可以 求出 雨 量与边 界值 关 系的表 达式 ; 然后根 据已 有 的雨量计 数据, 利用 最 优化方 法可 以 求得边 界值 将边界 值代 入 变分方 程, 则 方程可 解, 就 求得每 个离 散 点上的雨 量这 种方 法 的优点 在于 不 需要对 边界 条 件进行 假定 , 改进了 以往 运 用变分 法进 行 雷达定 量降 水 估测时存 在的问题 通过对交叉验证方法对变分方程求解, 以及将求解后的雨量进行拟合, 得到新的 Z R 关系 利 用 实 测 雨 量 计数 据 对 新 的 Z R 关 系 进 行 检 验 结 果 表 明 交 叉验 证 方 法 能 够在 不 假 定 边 界条 件 的 情 况 下 , 对 变 分 方 程 求解 , 解 的 状 况比 较 符 合 实 际情 况 对 层 状云 降 水 来 说 ,经 典 的 Z R 关 系 能够 很 好 的 反 映 降 水 的 情 况 , 拟合 后 的 结 果 与直 接 应 用 经 典 Z R 关 系 得 到 的 雨量 没 有 较 大 差异 , 两 者 与 实际 情 况 的 符 合 程度都很高;对于积层混合云来说,经过拟合后的降水结果明显优于直接用 Z R 关系计算的结果。
关键词:变分方法,降水估测,交叉验证一、引言雷达作为定量测量降水的工具, 在探测降水分布和结构上比较准确, 但是在精度上还是 低于雨量计 尽管双偏振雷达的出现使雷达估测降水精度大大增加, 但是由于国内目前双偏 振雷达研究还处于刚起步阶段, 并且目前国内布网雷达基本上不具有双偏振功能, 因此, 如 何有效地利用雨量计校正雷达估测降水来提高雷达测雨精度是非常必要的 为了用雨量计校 正 雷 达 降 雨 , 国 内 外 科 研 工 作 者 提 出 了 很 多 种 方 法 , 如 平 均 法 ( Wilson, 1979) , 变 分 方 法(Ninom iya and Akiyama,1 978;S un Shouxiang et al.,1993;赵坤 等,2001 ) 、卡 尔曼滤波方法(A hnert, 1986;D akSyangLin and Krajewski,1 989) Daley 曾明确指出伴随变分和 Kalman 滤波是大气数据同化的未来发展方向如今他的 预见已被越来越多的工作所验证 变分方法通过非线性模式解与不同时次观测资料集的全局 调整达到同化的目的二 、 研 究方 法 的 提 出 在 雷 达 定 量 估 测 降 水 领 域 ,近 年 来 的 主 要 发 展 在 变 分 分 析 方 案 , 其 本 质 是 通 过 调 整 模 式解以达到模式解与观测数据之间最小二乘条件。
在对变分方程求解过程中的一个问题就是边 界的处理 已往的研究中对于边界的处理方法是假定边界 这种方法从理论上来说存在误差 为解决这个问题,本文尝试用交叉验证法对变分方程求解,以改善雷达定量测量降水结果 三 、 详 细算法首先考虑任一矩形区域 Ω , 有边界 , 在 Ω 区域内有 Ng 个雨量站, 即为集合 在 Ω 资助项目:江苏省自然科学基金 BK2010599本文 部分内容投稿至其它学术期刊1区域内建立如下变分:⎪⎧⎪⎫⎡⎛ ⎞ 22R ⎛ ⎞ ⎤R 2 R R dxdy2⎥ R RJ ⎢ ⎨Ω ⎪⎩ ⎬⎜ x ⎟⎜ y ⎟r g⎣⎢⎝⎠ ⎝⎠⎦⎥⎪⎭其中 Rr 表示由经典 Z-R 关系估测出来的降雨量, Rg 是地面面雨量站的实测降雨量第 一 项 约 束 表 示 R 在 Ω 区 域 内 的 光 滑 性 , 第 二 项 约 束 表 示 R 在 整 个 区 域 内 与 雷 达 估 测 降 水 的误差,第三项约束表示 R 在整个 Ω 区域内与雨量计测的降水的误差; 、 、 0, 分 别为各项的权重系数。
为使求得的 R 在整个区域内光滑并且整体误差达到最小, 也 就是使得泛函 J 在边界 围 成的区域 Ω 中达到最小,即⎧⎪⎫⎪⎡⎛ ⎞ 22R ⎛ ⎞ ⎤R 2 R R dxdy 02⎥ R R J ⎢ ⎨Ω ⎪⎩ ⎬⎜ x ⎟⎜ y ⎟r g⎢⎣⎝⎠ ⎝⎠⎥⎦ ⎭上式中 为变分算子相应于此式的 Euler 方程为 :⎛ 2 R 2 R ⎞ ⎜ ⎟ (R Rr ) (R Rg ) 0 ⎝ x2 y2 ⎠在 整 个 区 域 Ω 内 对 上 式 采 用 五 点 差 分 方 案 离 散 其 中 Ri, j R( xi , y j ) , 令 ( i, j )ΩA 2 , B R R , Δ x 为 i 方向步长, Δ y 为 j 方向的步长将上2Δ x2Δ y2i, j ri , j gi, j述离散方程中未知变量 Ri, j 提取出来表示为: R F其中:( 1)⎡ D1⎢D2 ⎢⎢⎢⎢⎣D2 D1 O⎤⎥⎥⎥⎥⎢ D2OD2,在 中OD1 D2D2 ⎥D1 ⎥⎦M M ⎡ ⎤⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥ Δ y2A⎡⎢⎤⎥⎥⎥⎥⎥⎥⎥⎥ ⎥⎢ Δx2A ⎢⎢⎢⎢⎢⎢ Δx2 Δ y2⎢ Δ y2⎢⎢D1 ⎢⎢⎢⎢⎢⎢⎢⎣, D2 ⎢O⎢⎢⎢⎢⎢⎢⎣O O Δ y2OΔxA 2Δ y2 ⎥。
⎥⎥⎥⎦ N N Δ y2⎥A Δx2 ⎦N N TR [ R1,1 ,R1,2 ,L,R1,N ,L,Ri,1 ,Ri,2 ,L,Ri,N ,RM ,1 ,LRM ,N ]M N ,右端项 F B C E2B [ B1 ,B2 ,LBM 1 ,BM ]M , Bi [ Bi,1 ,Bi,2 ,LBi,N 1 ,Bi,N ]N , Bi, j Rri, j Rgi, j T T T TC [ C1 ,C2 ,LCM 1 ,CM ]M C1 [ R0,1 ,R0,2 ,LR0,N 1 ,R0,N ]N Δ x2, C2 C3 L CM 1 0 , TCM [ RM 1,1 ,RM 1,2 ,LRM 1,N 1 ,RM 1,N ]N Δ x2 T TE [ E1 ,E2 ,LEM 1 ,EM ]M , Ei [ Ri,0 ,0,L0,Ri,N 1 ]N Δ y2可 写 为 : F Fn Fb , 其 中 Fn B 表 示 与 边 界 无 关 的 部 分 , Fb C E 表 示 与 边 界 相 关 的部分。
那么(1)式写为 : R Fn Fb ,推出:R 1 F F (2)n b(2) 式中 H 1 、 F 中所有元素已知, 只有 F 含有未知的与边界相关的项, 所以问题的n b关键在于求解 Fb 在椭圆型偏微分方程中有定理:椭圆型偏微分方程中,区域内的解连续依赖于边界设向量 R 中有 Ng 个已知量( 即集合 ),即集合 中, R Rg ,即( k ) R k H 1 k ( F ( 3) F ) k R g n b其中 R k 和 H 1 k 分别表示向量 R 中的第 k 个元 素和矩阵 H 1 中的第 k 行k H 1 k F H 1 k F(3)式变化 为: R ( 4)g n b在(4)式中 ,左边的两项都为已知,右边包含未知的边界项考 虑 0 ~ M 10 ~ N 1 的 矩 形 区 域 , Fb 中 共 有 2M 2N 4 个 变 量 , 即 边 界 点 共2M 2N 4 个, 但是 4 个 顶点不参与运算 , 不在边界点向量 中体现, 所以 (L ) 表示边界点向量, 其中 L 2M 2N 。
以 R(1,0)作为 的第一个点 (1) , 沿着边界顺时针方向增加为 ( 2 ) …… (L ) 可建立 Fb 与 的对应关系建立泛函: ( R( i,j ) R ( i,j ))2 min ,利用上述计算结果带入此公式,得到:g(i, j ) (H 1(F ( i N j ) F ( i N j )) R ( i, j ))2 min 此式中 H 1 、 F 、 R 均为已n b g n g(i, j )知, 仅 Fb 中含有未知量 , 即边界上的值 对边界上的每点, 即 (i ) 求变分, 得到结果 如下所示: l ( i,j )) R( i, j ) l 1~ L 2( R( i,j ) R 0g l(i, j ) l ( H 1( F (i N j ) F (i N j )) R ( i N j ))H 1 Fb ( i N j ) 0 (5)l 1~ L n b g l(i, j )因 F C E 中仅含变量 ,所以 Fb 的矩阵可写成一个仅含 0 和 1 的稀疏 矩阵。
b3F b定义: 为计算方便,将所有二维变量矩阵都写成一维列向量,并只NM L NM L取与 对应的 Ng 个元素i, j ) ⎧⎪R( i, j ) R( tr ) , Rg ( i, j ) Rg ( tr ) ,其中 , r 1 ~ Ng t ⎨r⎪⎩ tr i* N j 1(5)式为: l ( H 1(t )( F F ) R (t )) ( HN。