ncl画斜剖面

趁着今天顺利挂上了github，赶紧同步一下！

主要参考ncl官网中的例子：narr_5.ncl，但该例子多了对曲面网格数据的处理,首先利用了ESMF_regrid_with_weights进行网格化生成直线网格，如果用再分析资料个人觉得可以不用这一步。

主要分为两步，首先是大圆路径，主要使用ncl中gc_latlon这个函数，该函数中设置变量时有一个常数npts，为x轴的点数，根据需要设置，一般设置越大数据越密。接下来进行插值，使用双线性插值从直线网格插值到非结构化网格，使用linint2_points_Wrap这个函数。

在实际画图中，遇到的问题是最后画图设置的经纬度，因为地球是球面，所以经纬度不是均匀的，设置时需要画gc_latlon 函数对应的经纬度。在选择的经度范围跨越半球时，两点之间的方向不可控制，这是因为gc_latlon这个函数的问题，试了多次也没解决，只能选择不跨越半球了。
If the two specified points are exactly opposite one another on the globe, the code does not fail, but the direction of the great circle route will be somewhat unpredictable. Setting npts <= 2 has the effect of just calculating the great circle distance between the two input points. No lat/lon points will be interpolated in between.
最后出现的警告，猜测是由于网格的不均匀分布造成的，但对画图的结果好像没有影响。
warning:ScalarFieldSetValues: irregular coordinate array sfXArray non-monotonic: defaulting sfXArray
更新一下，最后出现的警告是因为我！写代码的时候忘记删掉cyclic="True"这一行设置了！所以坐标就会出现警告！（问题不大

;=================================================================

;----------------斜刨面

;  calculate great circle along transect

npts=600 ;The actual number of interpolated points is npts-2.

N1=npts-1

nLabels      = 5   ;x轴坐标想要的标签数目

fiCyclicX    =True ;注意数据是否循环

;选取的两个点

leftlat  =  45

rightlat =  80

leftlon  =  300

rightlon =  357.5

;attention:数据是从0-360 or -180-180

;lonFlip 可以使用该函数转换

dist = gc_latlon(leftlat,leftlon, rightlat,rightlon, npts,2) ;大圆路径

;插值 变量为Fx_pos_10_area

xlon=Fx_pos_10_area&lon

xlat=Fx_pos_10_area&lat

Fx_pos_10_cross = linint2_points_Wrap(xlon,xlat,Fx_pos_10_area,fiCyclicX,dist@gclon,dist@gclat,2)

;插值后赋予变量属性

Fx_pos_10_cross!0      = "time"

Fx_pos_10_cross!1      = "level"

Fx_pos_10_cross&level= Fx_pos_10_area&level

;-----对应的经纬度

latXsecUser  = dist@gclat      ; convenience

lonXsecUser  = dist@gclon

print (dist@gclat+"  "+dist@gclon )  ; print the lats/lons

;--------画图时x轴坐标设置

;数据从0-360

XBValues_pos    = toint( fspan(0,N1,nLabels) )

XBLabels_pos    = new(nLabels,"string")

do i=0,nLabels-1

x = lonXsecUser(XBValues_pos(i))

y = latXsecUser(XBValues_pos(i))

if(x.lt.180)then
XBLabels_pos(i) = sprintf("%3.1f", y)+"N"+"~C~"+sprintf("%3.1f", x)+"E"

else if(x.gt.180.and.x.le.360)then

x1=360-x

XBLabels_pos(i) = sprintf("%3.1f", y)+"N"+"~C~"+sprintf("%3.1f",x1)+"W"

else if(x.eq.180.or.x.eq.360)then

XBLabels_pos(i) = sprintf("%5.1f", y)+"N"+"~C~"+sprintf("%5.1f", x)
end if

end if

end if

end do

res@tmXBMode   = "Explicit"

res@tmXBValues = XBValues_pos

res@tmXBLabels = XBLabels_pos

图片类似这样吧！

ps：最近读文献才发现，看你研究方向热不热门，就看看高质量期刊新出的文章跟你的方向差不多的多不多；然后发现我的方向真！离大谱！（不是我爱科研

Ran Xue

ncl画斜剖面

RWS波源