EOF 分析实例 使用ferret绘制结果图

EOF 分析实例 使用ferret绘制结果图 EOF 分析实例 使用ferret绘制结果 图 EOF 分析实例:使用ferret绘制结果图 2011年11月08日 下面使用ferret绘图分析前面eof分析的结果。使用的数据是cdo输出的结果eofs.nc、eigenvalues.nc，以及pc00000.nc-pc00010.nc。 我们使用循环，对每个变量分别绘制一张前4个eof的空间分布图，一张前4个pc的时间序列图，以及一张前10个eof的方差贡献图。 新建一个文件，命名为draweof.jnl。 1. 外层循环控制 let vps={"ul","ur","ll","lr"} let vars={"tsw","slp","ustrw","vstrw"} 这是定义字符串数组。ferret的字符串要放在双引号中，而字符串数组则用大括号定义。这里定义的两个字符串数组分别存放子图的位置和变量的名字，下面会用到。 repeat/range=1:4/name=nvar \ (......... ;\ ......... ;\ ) 这是最外层循环，遍历所有4个要分析的变量。ferret的循环用repeat命令完成，其range属性定义循环变量的范围，name属性定义循环变量的名称，后面接一对括号，里面是要循环执行的命令。 由于ferret脚本其实就是在ferret控制台输入的命令的记录，所以脚本文件中的语法与在控制台下一样，都是基于命令的，以行为单位。repeat命令后面的括号中如果要带有多行命令，要用分号隔开，这些命令如果要写在 。所以最终的格式如上面所示。下面的命令都写在多行上，需要在行尾加\ repeat的括号中。 let var=vars[i=`nvar`];\ 这条命令将变量名字符数组的当前元素赋值给一个字符串。在ferret中，任何变量都是定义在至少一维空间上的，第一维空间的下标变量固定为i(其他维分别为j、k、l)。我们刚才定义的字符串数组是一维的，它的第一个元素要用vars[i=1]来表示。类似的，一个4维变量的某个元素也可以用a[i=1,j=2,k=3,l=4]来表示。这种下标还支持切片:a[i=1:10,j=2,k=3,l=2:8]。用ijkl表示的空间与实际的地球坐标和时间没有关系，仅用来表示一个抽象的数学空间里的变量。而定义在4维物理空间(3维空间+1维时间)中的变量，除了抽象坐标职位，还可以用物理坐标xyzt来访问，如 a[x=100e:70w,y=20n,z=100,t=1-Jan-2011]访问变量a在北纬20度上东经100-西经70度的一条线上的值，层次为海面以下100m，时间是2011年1月1 日。 由于ferret的基本单位是终端中的命令，而命令不能接受变量做参数，只能接受键盘输入的明确常量(如a[i=1])。为了使命令能接受变量，引入了命令代换表达式。上面的`nvar`就是一个代换表达式。 代换表达式就是将用两个`符号(键盘上1左边)包围起来的一个表达式解析成数字。这个解析过程在命令执行之前发生，效果就相当于用户在键盘上手动输入了解析结果这个数字一样。 如`nvar`解析表达式"nvar"，这个表达式只是一个变量的名字，因此返回这个变量的值。比如当前是第一次循环，因此nvar=1，代换表达式的结果就是1。把这个数字放在它所在的命令里，这条命令就变成了 let var=vars[i=1];\ 然后ferret解析并执行这条命令。结果就好像用户明确输入了[i=1]一样。使用代换表达式，我们就可以在命令里使用变量了。这是ferret脚本中非常重要的一个功能。 现在就能理解，上面的命令是用代换表达式得到循环变量的当前值，用它取出变量名字符串数组的相应元素，并赋值给一个变量。因此var变量是一个保存着当前变量名的字符串。 2.2 use eofs.nc;\ 这条命令打开eofs.nc文件，其中包含cdo所计算的EOF。 2.3 set mode meta `var`_eofs.plt;\ 这里再次用了代换表达式，表达式内仍是一个变量的名字，因此返回这个变量的值，所以在第一次循环里代换之后这条命令变成了 set mode meta tsw_eofs.plt;\ meta是metafile的缩写。set mode metafile将metafile模式打开，这个模式将会把屏幕上绘制的图形以ferret自己的格式保存在一个图元文件(metafile)中，命令后面的文件名就是指定图元文件，如果不指定，则ferret使用默认文件名保存。每次绘图操作都会更新这个图元文件，直到下一条set mode meta命令指定了另外一个图元文件。 2.4 repeat/range=1:4/name=nvp \ (........ ;\ ........ ;\ ) 一个嵌套的循环，循环变量名为nvp，范围为1-4.这层循环用于遍历我们要绘制的4个eof。2.5-2.7的代码在这层循环中。 2.5 let vp=vps[i=`nvp`];\ 形式与2.1一样，用了代换表达式和下标索引。vp成为一个保存着子图位置的字符串，供下面使用。 2.6 set view `vp`;\ 代换表达式的结果为(第一次循环) set view ul;\ set view是set viewport的缩写，用于指定子绘图区的位置，类似于matlab或octave的subplot命令。ferret的viewport是一些提前定义好的区域，表示在整个大图上的一块子区域，并有名字。ferret内部定义好的viewport有upper、lower、left、right、ul、ur、ll、lr，分别表示上、下、左、右、左上、右上、左下、右下。如果还想使用其他区域可以自己定义。这里指定viewport为ul即左上。 注意这里有两个代换表达式位于字符串中，但这并不影响它的代换操作，它的优先级是非常高的。 fill命令用于绘制填色等值线图。l是抽象空间第四维(时间维)的下标变量，l=1限定第一个时间点。title属性设定图的标题，后面跟的变量名自然就是要绘制的变量。 常用的绘图命令还有plot、contour、shade等。 2.8 cancel view;\ 这一层嵌套循环结束了。使用cancel view命令来清空屏幕图像(不清空图元文件)。 2.9 sp Fprint -R -p portrait -o `var`_eofs.eps `var`_eofs.plt sp是spawn的缩写，用于调用系统命令。后面的Fprint就是一条系统命令，它是ferret提供的一个把图元文件转换成ps文件的小程序。Fprint命令的参数含义如下: -R 不对图元文件进行重命名，否则操作之后图元文件将被重命名为一 个不可预知的名字。 -p 纸张方向， portrait为肖像画，即长边竖直;landscape为风景画，即长边水平。 -o 输出文件的名字 后面的`var`_eofs.plt是输入文件的名字。 于是调用系统命令Fprint将tsw_eofs.plt转换成了tsw_eofs.eps，可以供latex使用，写文章比较方便。 3.1 set mode meta `var`_pcs.plt 打开metafile模式，输出图元文件为tsw_pcs.plt。 repeat/range=1:4/name=nvp \ (........... ;\ ........... ;\ ) 又一个嵌套循环，用于遍历4个pc。3.2-3.5在这个循环中。 3.3 use pc0000`nvp-1`;\ 代换表达式，生成结果为 use pc00000;\ 打开pc00000.nc文件，里面存放着第一个pc。 3.4 plot/l=1:252/color=blue/title="`var` PC`nvp`" `var`;\ 绘制时间序列图，color属性指定线的颜色。 3.5 plot/overlay/color=red/thick=3/title="1yr running mean" `var`[l=@sbx:12] ;\ 绘制滑动平均后的曲线。overlay属性指定当前图形叠加在原有图像之上，否则会冲掉原有图像。thick属性执行线宽。 [l=@sbx:12]是变量描述，指定在l维(时间维，L,不是1)上应用@sbx变换。ferret中以@开头的叫做变换，能够对指定维进行一些常用操作，如平均、求和等，很有用。@sbx变换是滑动平均(running mean)，后面的:12是它的参数，指定用12个点进行平均。即对于每个点，用其前面和后面共12个点的平均值作为这个点的新值。这里用12是为了得到年平均数据，因为原数据是每月一个的月平均数据。 嵌套循环结束。 清空图形，转换成ps文件。 use eigenvalues.nc;\ 指定输出文件，打开数据文件。 4.2 let/title="explained variance (%)" rel=`var`/`var`[l=@sum]*100;\ 虽然长，但其实只是一个赋值，用到了代换表达式和变换。为了清楚，写出代换后的结果: let/title="explained variance (%)" rel=tsw/tsw[l=@sum]*100;\ @sum变换返回变量在指定维上求和的结果，tsw[l=@sum]即将tsw变量在时间维上求和。不用担心求和后结果变成了一个数字，ferret执行变换时会保持其他维的形状。例如一个i=1:3,j=1:3的矩阵，有两个时间点l=1:2，则[l=@sum]之后的结果是一个3*3的矩阵，l=1。 上面的命令首先计算了所有eof对应的解释方差的和，然后分别除掉这个和，得到的就是相对的方差贡献了。 注意到赋值命令let中带了title属性。给变量指定的属性当绘制这个变量的图时，会自动成为图的标题。 4.3 plot/l=1:10/line/symbol/vs l,rel;\ 绘制前10个eof对应的方差贡献折线图。l=1:10指定为前10个;line和symbol属性分别指定绘制线条和符号，他们还可以接受赋值，用于设置线型和使用不同的符号。vs属性用于手动指定x轴变量，否则ferret会自 动选择合适的轴作为x轴。比如rel其实只有l维(或物理空间的t维)是变化的，其他维的长度都为1，所以ferret会自动选择t维的刻度作为x轴。这里使用了vs属性后，命令后面x轴和y轴变量都要给出，即x轴为l(1:10)，y轴为rel。 4.4 ppl xlab "NEOFs";\ ppl ylab "Explained Variance (percentage)";\ ppl plot;\ 设置坐标轴标记。ferret的绘图功能是使用了pplus软件包，所以图像的微调和定制都是通过ppl命令进行的。xlab和ylab分别设置x轴和y轴的标记。设置后需要用ppl plot刷新图像。 5. ferret脚本作为函数 这样程序就完成了。我们的程序循环4个变量，每个变量绘制前4个eof和pc，并绘制前10个eof的方差贡献。保存成draweof.jnl后可以在命令行中输入go draweof来运行。 程序中我们每次创建新图都要用set mode meta命令，转换图像时又需要用sp调用系统命令Fprint，属于重复行操作。为了更加方便，可以把这个操作写成脚本，每次需要绘图和转换时只要调用这个脚本就可以了，类似于其他语言(如matlab/octave、python等的函数)。正好作为练习。 ferret的脚本文件中没有类似于matlab/octave等语言的函数，在ferret中，“函数”指的是ferret提供的内部函数或者用户用fortran编写并编译为动态链接库的外部函数，外部函数在linux中一般为.so文件，在windows中为.dll文件。但ferret的脚本其实也可以接受参数，并被其他脚本(或命令行)调用，因此ferret的脚本本身就具有作为函数的能力。 向一个脚本传递参数，只要在go脚本文件名后依次用空格分隔各个参即可。 数 其中define sym定义了一个符号(symbol)名为pltfile，并在下面的命令中用($pltfile)调用它。ferret中的符号类似于c语言的“宏”(macro)，它仅仅是字符串的替换，会在命令运行之前把符号替换成它的内容。比如定义了一个名为a的符号: define sym a=hello 然后使用它: say ($a) ferret将会在解释并运行命令say之前，将($a)替换为hello，然后运行say hello.所以，ferret的符号完全等同于c语言的宏，它所作的是在命令运行之前，将符号替换为它的内容，好像用户手动输入了一样。符号的用法与变量类似，可以用于避免重复输入一样的字符，并且方便一次性更改多个地方。但它的本质只是字符替换，与代换表达式也不一样。 所以，其实这两条命令完全等于下面的一条命令: set mode meta .plt 表示传递给脚本的第一个参数。类似的、....都是使用参数，ferret最多支持99个参数。而$*则表示所有参数构成的字符串。 ferret在真正解释执行一条命令之前，其实会按顺序做三件事: 1> 将类似的参数表示替换为参数的值; 2> 运行代换表达式，并将结果放在表达式所在的位置; 3> 将符号(symbol)字符串替换为它的内容。 然后才真正形成了要运行的内容。 比如某文件name.jnl的内容为 let names={"john","tom","mike"} define sym sentence=hello, my name is `names[i=]`. say ($sentence) 则go name 2运行时，先将替换为其值2，成为 define sym sentence=hello,my name is `names[i=2]`. 然后运行代换表达式，成为 define sym sencence=hello,my name is tom. 然后替换符号，成为: say hello, my name is tom. 最后运行这条命令，终端输出hello, my name is tom. 理解了这些，我们上面写的epsinit.jnl也就很容易理解了，它的作用是设置图形输出到以其第一个参数命名的plt文件。 也就很清楚了。不过这里用到了参数的默认值，即%portrait%，如果没有给定第二个参数，那么第二个参数就为portrait。这给我们提供了很大的方便。 6. 运行方式 现在我们可以把draweof.jnl里面所有的set mode meta命令改为go epsinit脚本调用了，把所有的sp Fprint改为go epsprint。在ferret控制台输入go draweof.jnl运行。 默认情况下ferret控制台会打印出文件的内容，并展开循环的每一步，打印其完成了参数替换、代换表达式、符号替换之后的命令。这对于调试程序是很有用的，但如果不想每次都显示这么多信息，可以在ferret控制台输入 cancel mode verify 即关闭verify模式。 另外，除了在ferret控制台运行脚本，还可以在linux终端用命令 ferret -script draweof.jnl 运行脚本。这会打开ferret程序并执行draweof.jnl，执行后ferret仍处于运行中。-batch参数可以使ferret以批处理模式启动，不打开图形窗口，不显示verify信息，脚本运行完后即退出。比较适合用在shell脚本中实现自动化。 repeat/range=1:4/name=nvar \ (let var=vars[i=`nvar`];\ use eofs.nc;\ go epsinit `var`_eofs;\ repeat/range=1:4/name=nvp \ (let vp=vps[i=`nvp`];\ set view `vp`;\ fill/l=`nvp`/title="`var` EOF`nvp`" `var`;\ );\ cancel view;\ go epsprint `var`_eofs;\ ;\ go epsinit `var`_pcs;\ repeat/range=1:4/name=nvp \ (let vp=vps[i=`nvp`];\ set view `vp`;\ use pc0000`nvp-1`;\ plot/l=1:252/color=blue/title="`var` PC`nvp`" `var`;\ plot/overlay/color=red/thick=3/title="1yr running mean" `var`[l=@sbx:12] ;\ );\ cancel view;\ go epsprint `var`_pcs;\ ;\ go epsinit `var`_evs;\ use eigenvalues.nc;\ let/title="explained variance (%)" rel=`var`/`var`[l=@sum]*100;\ plot/l=1:10/line/symbol/vs l,rel;\ ppl xlab "NEOFs";\ ppl ylab "Explained Variance (percentage)";\ ppl plot;\ go epsprint `var`_evs;\ ) 前四个EOF的空间分布。显然第一EOF为正常模态，第二EOF为ENSO模态。 前4个PC的时间序列。显然第一个EOF为年变化，第二EOF即ENSO 的周期为4-7年。 前10个EOF的方差贡献分布图。显然赤道太平洋和北太平洋海域最强 的变化为年变化，而第二强的变化为ENSO。 特别声明: 1:资料来源于互联网，版权归属原作者 2:资料内容属于网络意见，与本账号立场无关 3 :如有侵权，请告知，立即删除。