Why

相关分析我之前也写过，为什么还要有这篇文章？因为我希望将来用模块化的思维使用NCL。计算是一部分，画图是另一部分。这样每次一个routine需要做什么就调用之前写好的template就可以。思路上共3部分：计算；输出；leadlag以及合并。

load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"

begin
lead_lag_num=12
target_month=13

target_index_name="tp_sh_stat"
target_index_raw   = asciiread("~/project/iobm_plus/sst_eof/iobm_lsf/IOBM_monthly_index_"+target_month+".txt", -1, "float")
year               = ispan(1979,2014,1)
year_start         = 1981
year_end           = 2008
year_num           = year_end-year_start+1
target_index       = target_index_raw
target_index!0     = "year"
target_index&year  = year
;index              = dim_standardize_Wrap(target_index({year_start:year_end},{target_month}),0)
index              = dim_standardize_Wrap(dtrend(target_index({year_start:year_end}),False),0)
copy_VarCoords(target_index({year_start:year_end}),index)
printVarSummary(index)

file_list=systemfunc("ls /uwnd.mon.mean.1958-2014.nc")
a=addfile(file_list,"r")
sst_summer=short2flt(a->uwnd((year_start-1958)*12+target_month-lead_lag_num-1:(year_end-1958)*12+target_month-lead_lag_num-1:12,:,:,:))
printVarSummary(sst_summer)
sst_summer_dtrend=dtrend_n(sst_summer,True,0)
copy_VarCoords(sst_summer,sst_summer_dtrend)
coe=escorc(index,sst_summer_dtrend(lev|:,lat|:,lon|:,time|:))
copy_VarCoords(sst_summer_dtrend(0,:,:,:),coe)
printVarSummary(coe)
printMinMax(coe,False)

coe_ttest=coe
coe_ttest=coe@_FillValue
dim_coe_ttest=dimsizes(coe_ttest)
do i=0,dim_coe_ttest(0)-1
    do j=0,dim_coe_ttest(1)-1
      do k=0,dim_coe_ttest(2)-1
        n    = dimsizes(index) ; n=11
        df   = n-2         ; degree of freedom
        t    = coe(i,j,k)*sqrt((n-2)/(1-coe(i,j,k)^2))
        p    = student_t(t, df)
        coe_ttest(i,j,k)=p
      end do 
    end do 
end do 
copy_VarCoords(coe,coe_ttest)
printMinMax(coe_ttest,False)
printVarSummary(coe_ttest)

output_filename=get_script_prefix_name+"_"+target_month+"_"+lead_lag_num+".nc"
system("rm -f "+output_filename)
PDataOut=addfile("./"+get_script_prefix_name+"/"+output_filename,"c")
;***************************************************
; define file options [Version a033]
;***********************************************
setfileoption(PDataOut, "prefill", False)
setfileoption(PDataOut, "suppressclose", True)
setfileoption(PDataOut, "definemode", True)
;***********************************************
; assign file attributes
;***********************************************   
;*****************************************************                 
fAtt               = True                       
fAtt@title         = "correlation between GODAS ts and "+target_index_name+"monthly data"         
fAtt@source        = "GODAS station obs"
fAtt@Conventions   = "None"
fAtt@creation_date = systemfunc("date")                                         
fileattdef( PDataOut, fAtt )
;***********************************************
; predefine coordinate information               
;***********************************************
dimNames = (/"lev", "lat", "lon" /)                       
dimSizes = (/dimsizes(coe&lev),dimsizes(coe&lat),dimsizes(coe&lon)/)                         
dimUnlim = (/False,False,False/)                         
filedimdef(PDataOut, dimNames,dimSizes,dimUnlim)

filevardef(PDataOut, "lat", typeof(coe&lat), getvardims(coe&lat))
filevarattdef(PDataOut, "lat", coe&lat)

filevardef(PDataOut, "lon", typeof(coe&lon), getvardims(coe&lon) )
filevarattdef(PDataOut, "lon", coe&lon)     

filevardef(PDataOut,"lev",typeof(coe&lev), getvardims(coe&lev) )
filevarattdef(PDataOut,"lev", coe&lev)

filevardef(PDataOut,"coe",typeof(coe), getvardims(coe) )
filevarattdef(PDataOut,"coe", coe)

filevardef(PDataOut,"coe_ttest",typeof(coe_ttest), getvardims(coe_ttest) )
filevarattdef(PDataOut,"coe_ttest", coe_ttest)
;***********************************************
; terminate define mode: not necessary but for clarity
;***********************************************
setfileoption(PDataOut, "definemode", False)
;***********************************************
; write data values to predefined locations 
; (/ .../) operator transfer values only                     
;***********************************************   
PDataOut->lat    = (/ coe&lat /)                                             
PDataOut->lon    = (/ coe&lon /)
PDataOut->lev   = (/ coe&lev /)
PDataOut->coe   = (/ coe /)
PDataOut->coe_ttest   = (/ coe_ttest /)
end

接下来用nco的命令ncecat将之前算的文件合并即可。之所以不直接都写在ncl里，一是考虑到高内聚低冗杂的模块化思想，二是考虑到可以降低内存方面的负载。至于为什么不在bash里优雅地使用循环。。好吧，长时间不用已经不想写bash了，觉得bash这种语言可读性不高，简单地句子分分钟明白（写几层循环过几天看整个人都不好了），而且本人Vim玩的很熟，上面地leadlag也就是录个宏搞定地，写下面这个东西也是几个正则表达式就搞定了。

ncecat -u lead_lag "${0%.*}"_1_-12.nc "${0%.*}"_1_-11.nc "${0%.*}"_1_-10.nc "${0%.*}"_1_-9.nc "${0%.*}"_1_-8.nc "${0%.*}"_1_-7.nc "${0%.*}"_1_-6.nc "${0%.*}"_1_-5.nc "${0%.*}"_1_-4.nc "${0%.*}"_1_-3.nc "${0%.*}"_1_-2.nc "${0%.*}"_1_-1.nc "${0%.*}"_1_0.nc "${0%.*}"_1_1.nc "${0%.*}"_1_2.nc "${0%.*}"_1_3.nc "${0%.*}"_1_4.nc "${0%.*}"_1_5.nc "${0%.*}"_1_6.nc "${0%.*}"_1_7.nc "${0%.*}"_1_8.nc "${0%.*}"_1_9.nc "${0%.*}"_1_10.nc "${0%.*}"_1_11.nc "${0%.*}"_1_12.nc                           "${0%.*}"_1.nc
ncecat -u lead_lag "${0%.*}"_2_-12.nc "${0%.*}"_2_-11.nc "${0%.*}"_2_-10.nc "${0%.*}"_2_-9.nc "${0%.*}"_2_-8.nc "${0%.*}"_2_-7.nc "${0%.*}"_2_-6.nc "${0%.*}"_2_-5.nc "${0%.*}"_2_-4.nc "${0%.*}"_2_-3.nc "${0%.*}"_2_-2.nc "${0%.*}"_2_-1.nc "${0%.*}"_2_0.nc "${0%.*}"_2_1.nc "${0%.*}"_2_2.nc "${0%.*}"_2_3.nc "${0%.*}"_2_4.nc "${0%.*}"_2_5.nc "${0%.*}"_2_6.nc "${0%.*}"_2_7.nc "${0%.*}"_2_8.nc "${0%.*}"_2_9.nc "${0%.*}"_2_10.nc "${0%.*}"_2_11.nc "${0%.*}"_2_12.nc                           "${0%.*}"_2.nc
ncecat -u lead_lag "${0%.*}"_3_-12.nc "${0%.*}"_3_-11.nc "${0%.*}"_3_-10.nc "${0%.*}"_3_-9.nc "${0%.*}"_3_-8.nc "${0%.*}"_3_-7.nc "${0%.*}"_3_-6.nc "${0%.*}"_3_-5.nc "${0%.*}"_3_-4.nc "${0%.*}"_3_-3.nc "${0%.*}"_3_-2.nc "${0%.*}"_3_-1.nc "${0%.*}"_3_0.nc "${0%.*}"_3_1.nc "${0%.*}"_3_2.nc "${0%.*}"_3_3.nc "${0%.*}"_3_4.nc "${0%.*}"_3_5.nc "${0%.*}"_3_6.nc "${0%.*}"_3_7.nc "${0%.*}"_3_8.nc "${0%.*}"_3_9.nc "${0%.*}"_3_10.nc "${0%.*}"_3_11.nc "${0%.*}"_3_12.nc                           "${0%.*}"_3.nc
ncecat -u lead_lag "${0%.*}"_4_-12.nc "${0%.*}"_4_-11.nc "${0%.*}"_4_-10.nc "${0%.*}"_4_-9.nc "${0%.*}"_4_-8.nc "${0%.*}"_4_-7.nc "${0%.*}"_4_-6.nc "${0%.*}"_4_-5.nc "${0%.*}"_4_-4.nc "${0%.*}"_4_-3.nc "${0%.*}"_4_-2.nc "${0%.*}"_4_-1.nc "${0%.*}"_4_0.nc "${0%.*}"_4_1.nc "${0%.*}"_4_2.nc "${0%.*}"_4_3.nc "${0%.*}"_4_4.nc "${0%.*}"_4_5.nc "${0%.*}"_4_6.nc "${0%.*}"_4_7.nc "${0%.*}"_4_8.nc "${0%.*}"_4_9.nc "${0%.*}"_4_10.nc "${0%.*}"_4_11.nc "${0%.*}"_4_12.nc                           "${0%.*}"_4.nc
ncecat -u lead_lag "${0%.*}"_5_-12.nc "${0%.*}"_5_-11.nc "${0%.*}"_5_-10.nc "${0%.*}"_5_-9.nc "${0%.*}"_5_-8.nc "${0%.*}"_5_-7.nc "${0%.*}"_5_-6.nc "${0%.*}"_5_-5.nc "${0%.*}"_5_-4.nc "${0%.*}"_5_-3.nc "${0%.*}"_5_-2.nc "${0%.*}"_5_-1.nc "${0%.*}"_5_0.nc "${0%.*}"_5_1.nc "${0%.*}"_5_2.nc "${0%.*}"_5_3.nc "${0%.*}"_5_4.nc "${0%.*}"_5_5.nc "${0%.*}"_5_6.nc "${0%.*}"_5_7.nc "${0%.*}"_5_8.nc "${0%.*}"_5_9.nc "${0%.*}"_5_10.nc "${0%.*}"_5_11.nc "${0%.*}"_5_12.nc                           "${0%.*}"_5.nc
ncecat -u lead_lag "${0%.*}"_6_-12.nc "${0%.*}"_6_-11.nc "${0%.*}"_6_-10.nc "${0%.*}"_6_-9.nc "${0%.*}"_6_-8.nc "${0%.*}"_6_-7.nc "${0%.*}"_6_-6.nc "${0%.*}"_6_-5.nc "${0%.*}"_6_-4.nc "${0%.*}"_6_-3.nc "${0%.*}"_6_-2.nc "${0%.*}"_6_-1.nc "${0%.*}"_6_0.nc "${0%.*}"_6_1.nc "${0%.*}"_6_2.nc "${0%.*}"_6_3.nc "${0%.*}"_6_4.nc "${0%.*}"_6_5.nc "${0%.*}"_6_6.nc "${0%.*}"_6_7.nc "${0%.*}"_6_8.nc "${0%.*}"_6_9.nc "${0%.*}"_6_10.nc "${0%.*}"_6_11.nc "${0%.*}"_6_12.nc                           "${0%.*}"_6.nc
ncecat -u lead_lag "${0%.*}"_7_-12.nc "${0%.*}"_7_-11.nc "${0%.*}"_7_-10.nc "${0%.*}"_7_-9.nc "${0%.*}"_7_-8.nc "${0%.*}"_7_-7.nc "${0%.*}"_7_-6.nc "${0%.*}"_7_-5.nc "${0%.*}"_7_-4.nc "${0%.*}"_7_-3.nc "${0%.*}"_7_-2.nc "${0%.*}"_7_-1.nc "${0%.*}"_7_0.nc "${0%.*}"_7_1.nc "${0%.*}"_7_2.nc "${0%.*}"_7_3.nc "${0%.*}"_7_4.nc "${0%.*}"_7_5.nc "${0%.*}"_7_6.nc "${0%.*}"_7_7.nc "${0%.*}"_7_8.nc "${0%.*}"_7_9.nc "${0%.*}"_7_10.nc "${0%.*}"_7_11.nc "${0%.*}"_7_12.nc                           "${0%.*}"_7.nc
ncecat -u lead_lag "${0%.*}"_8_-12.nc "${0%.*}"_8_-11.nc "${0%.*}"_8_-10.nc "${0%.*}"_8_-9.nc "${0%.*}"_8_-8.nc "${0%.*}"_8_-7.nc "${0%.*}"_8_-6.nc "${0%.*}"_8_-5.nc "${0%.*}"_8_-4.nc "${0%.*}"_8_-3.nc "${0%.*}"_8_-2.nc "${0%.*}"_8_-1.nc "${0%.*}"_8_0.nc "${0%.*}"_8_1.nc "${0%.*}"_8_2.nc "${0%.*}"_8_3.nc "${0%.*}"_8_4.nc "${0%.*}"_8_5.nc "${0%.*}"_8_6.nc "${0%.*}"_8_7.nc "${0%.*}"_8_8.nc "${0%.*}"_8_9.nc "${0%.*}"_8_10.nc "${0%.*}"_8_11.nc "${0%.*}"_8_12.nc                           "${0%.*}"_8.nc
ncecat -u lead_lag "${0%.*}"_9_-12.nc "${0%.*}"_9_-11.nc "${0%.*}"_9_-10.nc "${0%.*}"_9_-9.nc "${0%.*}"_9_-8.nc "${0%.*}"_9_-7.nc "${0%.*}"_9_-6.nc "${0%.*}"_9_-5.nc "${0%.*}"_9_-4.nc "${0%.*}"_9_-3.nc "${0%.*}"_9_-2.nc "${0%.*}"_9_-1.nc "${0%.*}"_9_0.nc "${0%.*}"_9_1.nc "${0%.*}"_9_2.nc "${0%.*}"_9_3.nc "${0%.*}"_9_4.nc "${0%.*}"_9_5.nc "${0%.*}"_9_6.nc "${0%.*}"_9_7.nc "${0%.*}"_9_8.nc "${0%.*}"_9_9.nc "${0%.*}"_9_10.nc "${0%.*}"_9_11.nc "${0%.*}"_9_12.nc                           "${0%.*}"_9.nc
ncecat -u lead_lag "${0%.*}"_10_-12.nc "${0%.*}"_10_-11.nc "${0%.*}"_10_-10.nc "${0%.*}"_10_-9.nc "${0%.*}"_10_-8.nc "${0%.*}"_10_-7.nc "${0%.*}"_10_-6.nc "${0%.*}"_10_-5.nc "${0%.*}"_10_-4.nc "${0%.*}"_10_-3.nc "${0%.*}"_10_-2.nc "${0%.*}"_10_-1.nc "${0%.*}"_10_0.nc "${0%.*}"_10_1.nc "${0%.*}"_10_2.nc "${0%.*}"_10_3.nc "${0%.*}"_10_4.nc "${0%.*}"_10_5.nc "${0%.*}"_10_6.nc "${0%.*}"_10_7.nc "${0%.*}"_10_8.nc "${0%.*}"_10_9.nc "${0%.*}"_10_10.nc "${0%.*}"_10_11.nc "${0%.*}"_10_12.nc  "${0%.*}"_10.nc
ncecat -u lead_lag "${0%.*}"_11_-12.nc "${0%.*}"_11_-11.nc "${0%.*}"_11_-10.nc "${0%.*}"_11_-9.nc "${0%.*}"_11_-8.nc "${0%.*}"_11_-7.nc "${0%.*}"_11_-6.nc "${0%.*}"_11_-5.nc "${0%.*}"_11_-4.nc "${0%.*}"_11_-3.nc "${0%.*}"_11_-2.nc "${0%.*}"_11_-1.nc "${0%.*}"_11_0.nc "${0%.*}"_11_1.nc "${0%.*}"_11_2.nc "${0%.*}"_11_3.nc "${0%.*}"_11_4.nc "${0%.*}"_11_5.nc "${0%.*}"_11_6.nc "${0%.*}"_11_7.nc "${0%.*}"_11_8.nc "${0%.*}"_11_9.nc "${0%.*}"_11_10.nc "${0%.*}"_11_11.nc "${0%.*}"_11_12.nc  "${0%.*}"_11.nc
ncecat -u lead_lag "${0%.*}"_12_-12.nc "${0%.*}"_12_-11.nc "${0%.*}"_12_-10.nc "${0%.*}"_12_-9.nc "${0%.*}"_12_-8.nc "${0%.*}"_12_-7.nc "${0%.*}"_12_-6.nc "${0%.*}"_12_-5.nc "${0%.*}"_12_-4.nc "${0%.*}"_12_-3.nc "${0%.*}"_12_-2.nc "${0%.*}"_12_-1.nc "${0%.*}"_12_0.nc "${0%.*}"_12_1.nc "${0%.*}"_12_2.nc "${0%.*}"_12_3.nc "${0%.*}"_12_4.nc "${0%.*}"_12_5.nc "${0%.*}"_12_6.nc "${0%.*}"_12_7.nc "${0%.*}"_12_8.nc "${0%.*}"_12_9.nc "${0%.*}"_12_10.nc "${0%.*}"_12_11.nc "${0%.*}"_12_12.nc  "${0%.*}"_12.nc
ncecat -u month "${0%.*}"_1.nc "${0%.*}"_2.nc "${0%.*}"_3.nc "${0%.*}"_4.nc "${0%.*}"_5.nc "${0%.*}"_6.nc "${0%.*}"_7.nc "${0%.*}"_8.nc "${0%.*}"_9.nc "${0%.*}"_10.nc "${0%.*}"_11.nc "${0%.*}"_12.nc "${0%.*}".nc 

最后生成地文件大概就是这样：

Variable: f
Type: file
filename:   coe_uwnd_iobm
path:   coe_uwnd_iobm.nc
   file global attributes:
      creation_date : Wed Feb 22 20:00:53 CST 2017
      Conventions : None
      source : GODAS station obs
      title : correlation between GODAS ts and tp_sh_statmonthly data
      history : Fri Feb 24 07:28:23 2017: ncecat -u month coe_uwnd_iobm_1.nc coe_uwnd_iobm_2.nc coe_uwnd_iobm_3.c
Fri Feb 24 07:27:39 2017: ncecat -u lead_lag coe_uwnd_iobm_1_-12.nc coe_uwnd_iobm_1_-11.nc coe_uwnd_iobm_1_-10.nc
      nco_openmp_thread_number : 1
   dimensions:
      month = 12  // unlimited
      lat = 145
      lon = 288
      lev = 37
      lead_lag = 25
   variables:
      float lat ( lat )
         long_name :    latitude
         GridType : Cylindrical Equidistant Projection Grid
         units :    degrees_north
         Dj :   1.25
         Di :   1.25
         Lo2 :  -1.25
         La2 :  -90
         Lo1 :   0
         La1 :  90

      float lon ( lon )
         long_name :    longitude
         GridType : Cylindrical Equidistant Projection Grid
         units :    degrees_east
         Dj :   1.25
         Di :   1.25
         Lo2 :  -1.25
         La2 :  -90
         Lo1 :   0
         La1 :  90

      integer lev ( lev )
         long_name :    isobaric level
         units :    hPa

      float coe ( month, lead_lag, lev, lat, lon )
         _FillValue :   9.96921e+36

      float coe_ttest ( month, lead_lag, lev, lat, lon )
         _FillValue :   9.96921e+36