天气诊断分析报告.docx

1、天气诊断分析报告第二章 资料处理和客观分析 2.1 资料处理随着气象观测手段的发展和现代化，得到的资料数量增多，门类拓广，有常规站点 观测，有非定点海洋观测，有飞机观测，有气象卫星观测，有定时观测和非定常观测等 等，怎样使用这些资料呢？一般来说，对气象资料的要求有两方面：一是可靠性，二是 便于使人所周知。气象测站的分布是不规则的，因此我们只能得到这些不规则点上的气 象资料，但是数值预报中的网格点是规则的，因而资料无法直接使用。另外，无论是用 穿孔纸带或用电信号的形式将气象电报直接输入电子计算机， 都要首先按照专门的程序 进行译码、检查、整理。因为气象电报的内容，是按照气象电码格式编发的，而它的

2、形 式又是按照邮电电码格式编发的，从观测、编码、发报，到传递、转换、接收等，在每 个工序和环节上，都存在着出错的可能性。因此，我们所接收到的气象电码，不可避免 地存在着一些错误或不妥之处。 所以，要正确使用这些气象资料， 必须经过必要的处理。 比较简单的资料处理可分为以下几个方面：2.1.1记录错情判断在利用接收到的气象资料之前，首先对资料要作错情判断。一般的做法是根据不同 等压面上各种要素值的大小，给出相应略大于其最大值的一个数作为其上限值；也给出 相应的略小于其最小值的一个数作为其下限值。例如在我国范围内的各测站，冬季 500 百帕层上的位势高度最小值不超过 500位势什米，我们则取 50

3、0 作为其下限值；最大值 不超过 600 位势什米，我们就取 600作为其上限值。然后利用比较大小的子程序，由计 算机对每一组数据进行判断，凡是大于上限值或小于下限值的记录，我们就认为它是错 误的，予以舍掉，作为缺测记录。另外还可利用气象要素在时间变化上的连续性和空间 分布上的连续性，来判断一个气象要素记录是否错误。对于错误记录，可用下面介绍的 补缺测记录的方法，另外补一个值。2.1.2 补缺测或漏传记录一般可把缺测或漏传记录的测点看作是一个网格点， 然后由下节介绍的客观分析方法，利用周围已有的测站记录，插一个值补上2.1.3实测风矢量的分解气象台站观测到的风场资料，是一个既有大小又有方向的风

4、矢量，为便于该资料的 利用经常将实测风分解为东西和南北两个分量。 分别用u、v表示，并规定：u向东为正, V向北为正。其数量值分别由下式计算：u Msi n( /180 ) (21)v Vcos( /180 )图2.1实测风矢量的分解这里同为实测风速值，为测风报告发布的风向度数。女口，实测风为西南风，为240,风速10m/s，算得：u 10 sin(240 /180 ) 8.7(m/s)v 10 cos(240 /180 ) 5(m/s)2.1.4风场订正诊断分析一般都是在有限区域内进行的，多数都采用正方形网格。客观分析后所得 到的网格点上的u、v分量值，并不处处与网格区的x,y轴平行，因此还

5、必须进行风 向订正。因为只有在基线上的网格点，其东西、南北方向与x,y方向一致，其他网格 点上东西、南北方向与x,y方向总有一个偏差角，这显然会给计算带来误差，特别是 当计算范围取得较大时，边缘的网格上，这种风向的误差显得更加突出。不进行适当订 正是不行的。如图2.2,先考虑在基线以西的某一网格点 A, N是北极，NO为基线，NA和CA 分别为经过A点的经线和纬线，MA和LA分别和这个正方形网格系统的 X轴及丫轴平 行。假设A点的风速在经纬的分量分别为 u和V。而在网格的X,Y方向的分量分别为u ，由于A点所在的经线不与基线相平行，故u、v分量和u 、v 分量彼此也不平行,而是有一夹角a ,

6、a =Z ANO。由图可看出，它们之间有如下的换算关系：图2.2风场订正示意图在基线以西的网格 AO0，所以 =sin-1（AO/AN）0 ;在基线以东的网格，AO0，所以0;在基线上的网格，a = 0，则u =u、v =v。即不必订。设在某网格点上，u = v= 10米/秒，a = 45由（2.2）式订正后u 10 cos45 10 sin45 14（米 /秒）v 10 cos45 10 sin45 0（米 / 秒）2.1.5平滑和滤波气象观测资料，总存在着各种各样的误差。比如由气象仪器安装不标准等带来的，非偶然性误差（器差）和由工作人员在观测、编码、发收报等造成的偶然性误差，以及 将要素值

7、内插到网格点上时，产生的舍入和插值误差等等。无疑，这些误差都将会影响 计算的结果，为了减少误差的影响，通常在计算之前先对原始资料等进行平滑和过滤 （滤 波），滤掉那些次要的小的天气意义的东西，而保留和突出主要的量，或者，为了研究 的需要滤去资料中某些波长的量，而保留与问题有关的量。现分别介绍如下：(1) 一维平滑算子这是最简单的平滑算子。利用同一直线上三点的资料，又称三点平滑算子j 1 2fj (1 S)fj 2(2.3)这里f j表示第j点平滑后的值，fj表示第j点平滑前的值。S为平滑系数(可正可负)该平滑算子对j点对称，其权重除j+1，j和j-1点外均为0。对函数f (x)可展成富氏级数，

8、在fj点，可写成:(2.4)fj(x) C Aeik(Xj )这里C为常数，A为波动的振幅。k 2 / L为波数，L为波长，LJ为位相。同样在Xj 1和Xj 1点，函数f (x)可写成比较(2.4) (2.6)式，平滑后的波相未变，改变的只是波的振幅，平滑后的振幅为：A A1 S(1 cosk x) (2.7)令r=LAZA ，称响应函数，即平滑波幅同原波幅之比，表示平滑后的变化。显然, R= 1,表示平滑后波幅一样，R1，表示平滑后使原波幅被放大。由(2.7式知：R(k,s) A/ A 1 S(1 cosk x)或者 |R(L,S) 1 2S?si n2( x/训 (2.8)可见对于固定的网

9、格距，响应函数 R只与波数k (波长L)以及平滑系数S有关。由于0W sin2( x/L) 1若希望平滑后使原波动衰减，以致滤掉(但不希望出现反 位相情况)，由只须要有0W R 1,于是由(2.6式知：0 S1，则必有平滑系数S2 x的波，可以使其波幅有所恢复。比如取 L=6 x时：一次平滑得 R1 (6A x,1/2)=0.75二次平滑得 R1-2=R(6A x,1/2) R (6 x,-1/2)=0.94二者相比，二次平滑使该波幅恢复了 19%。表明这种平滑对保留长波是有益的。(2)二维平滑算子对于平面的问题，须进行二维空间的平滑，把一维推广到二维有两种处理方法：将计算的场先分别在X方向和

10、丫方向进行平滑，然后取平均，即这里用到的是i,j点及其前、后、左、右共五个点的资料，故称为五点平滑格式将场先在一个方向平滑，然后再在另一方向上平滑，即-nr S 2 S2 2仆 fi,j fi,j T(1 S)辂怙 (2.13)2 4其中： 2fij 同(2.12),2*fj=f i+1,j+1 +f i+1,j-1 + fi-1,j+1 + f i-1,j-1-4fi,j o 这里用的是1,J点及其前后,左右及前后点的左右点(式左右点的前后点)共九个点的资料，故称为九点平滑格式。用Kx，Ky，Lx, Ly分别表示X，丫方向上的波数和波长，其平面波的表示形式可写成C Aei(kxXkY丫)或C

11、 Ae(2X/Lx 2 丫山)其中A为振幅。以之代入(2.12),(2.13易得其相应的响应函数1 2 2R五占 一Rx Ry 1 Ssin (kx X/2) sin 你丫 Y/2)21 Ssin2( X/ Lx) sin2( Y/ Ly) (2.14)R九点 RX ?Ry 1 2Ssin2(kX X /2) ? 1 2Ssin2(kY Y/2)1 2Ssiri( X/LX)?1 2Ssirf( 丫/L丫) (2.15)2.1.6尺度分离实际的大气运动，包含了各种尺度的天气系统，为着研究的需要，经常要将实际的 扰动，分离成不同尺度的波。应用富里叶级数将扰动展成不同尺度的波是常用的分离状 的谱分

12、析方法，但是它只能用于沿整个纬圈的半球性分析，不便考察某一特定地区和某 些特定尺度系统等之间的关系。而利用上节所讲的平滑滤波技术同时可以满足这些方面 的要求，同时又有着明显的天气意义。为了分离不同尺度的波，须事先设计一个合适的滤波器。对于任意天气变量A ,可以写成：A= ? +(A-?)o ?表示平滑后的A场，显然：(A- ?) 表示原始场与平滑场的差。如果在作平滑的过程中，滤去高频波(短波) ，保留低频波(长波)，这种平滑就称为一个低通滤波器，它的作用是让低频波通过。表示场的低频部分，即大尺度运动。而(A- ?)部分，为该场的高颇(短波)部分，即小尺度系统。对 (A- ?)进行运算，也就相当

13、于一个高通滤波器。这里问题的关键，是如何选择滤波器，如何选择平滑系统和平滑次数，要看对具体问题的分析。比如要想分离出短波(波长L为500km左右)可以选用三点平滑算子,=1/2，连续平滑三次，就得到单波响应函数： R1-3=(1-sin2( x/L)3。如果计算网格距d取为100km,在？场中波长为5d(500km以下的波，衰减了70%以上。保留了 30%以下，而对波长L=20d(2000km)以上的波，竟保留了 95%以上。这表明，经过如此过滤之后，在 ？场中主要是波长为2000公里以上的天气波。而在(A- ?)场中，主要是波长在500km以下的次天气尺度的短波系统。如果研究的是面(二维)上

14、的问题，即可采用二维平滑，如用九点平滑算子:S2 S 2 S 2Ai,j A,j 2(1 S) Ai,j 4 A,j其中2Ai,j Ai 1,j Ai 1,j Ai,j 1 Ai,j 1 4A,jAi,j Ai 1,j 1 Ai 1,j 1 Ai 1,j 1 Ai 1, j 1 4 Ai,j(2.16)将实测风场，分解为u、v分量，取d= 100km, S=1/2，u、v分别用(2.15式连续平滑三次，然后再合成，即得到相应的流场。此外，还可以根据这种方法，设计只允许某一范围的波通过，而其它波不能通过的所谓带通滤波器，以分析研究某种波长的波的活动情况，这里变不再介绍了。 2.2客观分析诊断分析

15、一般所需要的资料是网格点上的，而常规的气象观测资料是在固定地点(地 面和高空观测站)和固定时间观测到的。为了由这些离散的分布不规则的资料计算出某些 物理量，从原则上必须得到每一观测变量在时、空上呈连续分布的场。但实际上并不需 要观测值(如位势高度)在(X,Y)平面上连续变化的分布。因为一般可用有限差分方法来 计算所需要的导数和梯度，这时只需要把空间上分布不均匀的台站资料内插到规则分布 的网格点上就行了。为了得到网格上的资料，可采用两种方法进行内插：一种是主观内 插法，即手工分析各种气象要素场的等值线，然后按网格点读取格点数，这种方法叫主观分析；另一种方法是根据直接联系格点值与台站值的方程，从数

16、值上(用计算机)进行内插。这种方法叫客观分析。所用的方程或函数可以是不同的。常 用的有限元、多项式、样条等，数值天气预报中还常使用逐步订正法、最优插值法、谱 方法、变分法等。另外，曲面拟合方法用于台站稀少的地区的物理量计算，也是比较好 的方法。本节我们主要介绍有限元法、多项式法和逐次订正法。客观分析方案可进行两种类型分析：一种是向量场分析，例如风场，其中所处理的 资料不公有量值，还有方向；另一种处理量只有量值的资料，如温度场、湿度场等，这 是标量场分析。2.2.1有限元素法图2.5三角形区域单元将我们要进行物理量计算的区域内的所有测站划分成有限个三个角形单元， 如图2.5就是其中一个单元，A、

17、B、C表示三个邻近的测站，0是I ABC|内的一个网格点，若三 个测站的要素值在一个平面上，根据有限元素法，它们可按其位置坐标展开以下的线性(2.17)关系式:SALiL2XAL3YASBLiL2XBL3YBScLiL2XCL3YC其中(XA,YA), (XB,YB)，(XC,YC)，分别为A、B、C三测站相对于任意坐标系的位置坐标。SA，SB, SC分别为测站所测得的要素值，这些都是已知的。因此，若有三组观 测资料，根据线性代数中的克莱姆法则或高斯消去法等方法，就可求解出方程组的三个 系数 L1、L2、L3。对于I ABC|内的任意一个风格点0,已知它在同一坐标系里的位置坐标（Xo,Y。），

18、则表示0点的要素值So可由下式求得：S0=Li+L2Xo+L3丫0 （2.18）以上介绍的是有限元素法的一般原理，要用它来内插出所有网格点上的要素值。还 要借助于电子计算机。具体做法是：先任取一网格点（i,j）并以此为圆心，以一个网格距为 扫瞄半径画圆，然后依次计算各测站到该网格点（i,j）的距离。凡是距离小于或等于扫瞄半 径的测站被入选，大于的被舍去，这样一个过程。称为第一次“扫瞄” 。这次扫瞄后，若圆内测站数小于3,则可加大扫瞄半径（取1.5个网格距或两个网格距），按上述的过程 再扫瞄一次，若圆内测站仍小于 3个，可继续加大扫瞄半径，并重复上述过程，一直到 圆内测站数大于或等于3个为止。如

19、果是等于3个，就可由这3个测站的资料，构成三 元一次线性方程组。用高斯一亚当消去法或用其它的方法求出三个系数 Li、L2、L3。，将这三个系数网格点（i,j）的位置坐标代入（2.17式，便可内插出该网格点上的要素值。如果 圆内的测站数大于3个，则再调用一个比较距离的子程序，选取离网格点（i,j）最近的一个 测站。禾用其要素值，按上述方法，内插出该网格点上相应的要素值。网格点（i,j）的要素值被内插出后，再进行下一个循环，同样用上述方法内插出下一个 网格点的要素值，这样依次循环下去，直至内插出所有网格点上的要素值为止。实际上，由于观测误差等原因，通常希望有比 3个测站还要多的资料，这样求得的系数

20、，可以减少观测随机误差的影响。为了解决这个问题，常采用最小二乘法来确定系数，由此得到（2.18）式:式中n是以扫瞄半径所画圆内的测站数，n4,解方程组（2.19求得系数 也、亘和口， 然后再利用公式（2.17）求得圆内格点上的要素值2.2.2多项式法多项式法的原理是寻找一个由多项式表示的的曲面来逼近网格点周围区域各测站 实测的气象要素值。如果这个曲面被找到，它即可代表该气象要素在这一网格点附近的空间分布状态，从而可求得网格点上的要素值，具体做法如下：假如某一等压面上的位 势高度分布Z(x,y)可以用一个m次多项式表示，其表达式为次幂项数p=0P=1P=2P=3?P=m? 缶式中X，丫是该等压面

21、上的两个位置坐标，p，q分别是它们的次数，ak是xp-qyq的系 数。如果想用一个二次曲面去逼近高度场，那么可取 m = 2;如果想用一个三次曲面去逼近高度场，那么可取 3;若想用一个更高次曲面去逼近高度场，那么 m取更大的数值。现在我们所面临的问题是如何运用网格点周围测站已知的一组实测高度值来确定 系数ak，使得(2.20)式能最佳地逼近网格附近的高度场。对于给定的一组观测资料来说，坐标(Xi,yi)(i=1,2,.n是已知的，实测高度值Zb也是 已知的，m值是按我们的要求确定的，只有 Q是未知数，即待确定的。当(2.19)式中m = 2时，二次多项式的一般形式为2 2Z(x, y) a1

22、a?x a3y aqX a5xy aey | (2.21)当3时，三次多项式的一般形式为Z(x, y) a1 a?x a3y aqX2 axy aey23 2 2 3 (.)a?x agx y agxy ay(2.21)式中含有六个待定系数，至少要有六组观测值才能求解。 (2.2)式中含有10个待定系数，至少需要110组观测资料。m与待定系数个数 M的关系为，M 2（m 1）?（m 2）。实际上，如有限元素一节中所述，为了减少观测随机误差的影响，得到最佳的逼近，常采用最小乘 二法来确定系数 刨。如果所考虑的网格点附近有 n个测站，n个测站实测高度值为Z0b （i 1,2,.呵，根据各测站的坐标

23、|（3）|，由多项式可以再得相应的高度值，某测 站i的实测值与计算值的差为|Zi0b乙|，称为残差，最小二乘法就是取残差的平方和为 最小。残差平方和为：n(2.23)E(aj) (Zb Zi)2i 1上式为m元二次函数，自变量为（a）,残差平方和为最小的条件是它的一阶偏导数为零， 即:整理得:代如（2.24得:上式为求解M个未知数 包的线性方程组（若 m = 2,则M = 6;若m= 3, M = 10）, 形式与（2.19式相同，解这个方程组，即可确定系数 凶，从而最后确定我们所要寻求的 多项式，并用此内插出网格点上的要素值。在电子计算机上实现的方法和步骤，与有限 元素法基本相同。2.2.3

24、逐步订正法逐步订正法就是将格点周围站记录与终点值进行比较， 用格点周围不同半径范围内各测站的观测数值情形与估计值之差的加权平均作修正量逐步对其订正， 最终使格点分析值与周围测站记录相比达到完全合理为止。下面以位势高度的格点分析值的计算为例，简要介绍一下逐步订正法。以格点为中心，以R为半径划一个圆，设落在圆内的观测值有 N个（图2.6）,我们将根据这N个观测值来订正估值。将初估场（预备场）的格点值| d|内插到某一个观测点上，例如i点，设内插出的高 度值为，若1点的高度观测值为 ，则观测值与内插值的差可以算出：i i i对圆内所有N个观测值都算出它们与各自内插值的差 1，I2，. N I，并计算

25、这些差值的加权平均N / NC w 订 w （2.28）i 1 / i 1就是来自高度场的订正值。式中 回为权重系数。画有很多种形式，但最主 要的特征是和格点到1点的距离的平方成反比，其物理意义是，观测值对格点值的影响 程度随着它们之间距离的平方而减小。显然，这种公重系数是各向同性的，它在平面上 的分布是一个个同心圆。例如可以取|w （R2）/（R2）|，式中且是格点到1点的 距离，R是扫描半径。也可以取为 W 1（1 和，是经验系数，随不同层次，不同 要素而异。如果在求取格点的高度分析值时还要使用风场的观测资料， 则还需要求出来自风场的订正值Cv0。Z/ r/ 0、01 o 4图2.6圆形区

26、域d|在1点展开:设所要求取的格点分析值为口，则有| d d d假定在图2.5所示的圆内位势高度是线性变化的，利用台劳展式将即:f v x u y id i xyddx iy i对圆内所有有测风报告的站点加权平均，便得到来自风场的订正值MW i f(v X u y)i(2.30)i 1M di 1C a CM为圆内有测风报告的站点个数，显然 MN于是，总订正值为:(2.31)匾可以取匕和aj为加权系数，可根据需要加以调整。比如，对于风场好用的地区, 得大一些，反之亦然，但U和刨之间必须满足以下关系:a av 1由此得到格点的分析值为d d C| (2.32)这是对于影响半径R进行的第一次许订正

27、，然后逐渐缩小R进行重复订正。每次订 正后的分析场用作下一次订正的估计值。在连续订正过程中，缩小影响半径可以去掉估 计场的大尺度误差，使分析场越来越逼近观测结果。在逐次订正中还可以采用不同形式的权重系数。2.2.4拉格朗日垂直插值法气象资料在垂直方向上往往是不等距的，为计算方便起见，常常需要将资料(观测 结果)进行垂直插值。图2.7线性插值多项式最简单的插值是线性插值。设对应变量X有函数值丫，对点I XoXjX?,Xn|所对 应的函数值为丫0,丫，丫2,Yn |，取Y=P(X)线性插值，就是在点|Xo,Xl上有函数值丫， 两点作直线于是可得线性插值。Y 丫0 (丫！ Yo)(X XoWX! X。) R(x) (2.33)也就是一次插值多项式。显然，由此所示得的丫值其误差较大(如图2.7)。常用的是二次插值多项式，为：Y Pn(X) Lo(X)Y0 Li(X)Yi Ln(X)Yn (2.35)式中的系数，又是一个n阶多项式由X Xo时，丫 丫0 ； X X1时，丫 丫；X Xn时，丫 Yn ;可知|Li(X)应满足条件:Li(Xj) 0 j (2.36)0 J I考虑Lo(Xj)，除J=0，Lo(Xj) =1夕卜，其它点Lo(Xj)=O，故可将Lo(Xj)表示为多项式:Lo(X) C(X X1)(X X2)(X Xn)| (2