之前主要用的是再分析资料,它们基本上是规整的经纬度网格。
最近在看PDO对亚洲季风的影响,器测数据的长度不够,就挑了MPI-ESM的piControl来对比着看。
然后,在计算速度势的时候,遇到了点麻烦。
GrADS有内置函数可以去算速度势,其中fish这一系列函数应该是迭代法解的,而shfilt这一系列函数是用球谐函数解的。
但是这俩族函数都只能用在规整的经纬度坐标上,这对于再分析资料来说没啥问题,可是模式输出经常是在谱网格上。
所以直接用的话,会出现这之类的报错:
ga-> d fish_chi(ua.1,va.2)
Error from UADV (VADV): Nonlinear U/V grids
Both xdef and ydef should be linear
Operation Error: Error from madvu function
Operation Error: Error from fish function
Operation Error: Error from fish_chi function
以及
ga-> d sh_chi(ua.1,va.2)
Error from UADV (VADV): Nonlinear U/V grids
Both xdef and ydef should be linear
Operation Error: Error from madvu function
Operation Error: Error from sh_fish function
Operation Error: Error from sh_chi function
GrADS倒是可以去做插值,比如lterp这个函数,其具体说明可以在这个网址上找到http://cola.gmu.edu/grads/gadoc/gadocindex.html (这里有GrADS很多命令的解释)
但是lterp需要弄俩ctl,也就是模板的ctl也得有,稍微有点麻烦。
所以我一开始是用MATLAB的interp2函数插值到经纬度坐标,然后再写成dat和ctl,让GrADS算速度势。
(不直接用MATLAB解泊松方程,是因为这是个球面上的泊松方程,跟直角坐标不一样。但我记得气象家园有帖子说怎么用MATLAB计算速度势辐散风,这是另一种解决问题的方法啦)
fish的介绍里,提到了它要求是全球且均匀的经纬度网格,也提到了如果资料不满足条件,可以用re这个函数去插值到规整的经纬度网格。感觉这样做会比较简单,因为不需要在不同软件(GrADS和MATLAB)之间倒腾来倒腾去。
re的用法挺简单的,例如说我想从T62(192x94)的高斯网格插值到2.5x2.5的经纬度网格,这一句就够了
define var25=re(varT62 , 2.5 )
这个函数还有一些参数可以设定,截个图吧,具体的不说了,看它的解释文件就好了
插值之后还有一个小坑,计算速度势的时候,直接用fish_chi(uu,vv)是不行的,其中uu和vv是插值之后的风场。会报错如下:
ga-> d fish_chi(uu,vv)
Error in gagchk: axis sizes are not the same
Error from UADV (VADV, WADV): Incompatable grids
Dimension ranges unequal
Operation Error: Error from madvu function
Operation Error: Error from fish function
Operation Error: Error from fish_chi function