<div dir="ltr"><div>I did it. Thanks. <br></div><div><br></div><div>I have been trying many times with this data to get olr and zonal wind anomaly and finally I got both min and max data values. <br></div><div><br></div><div><b>This is the result when I use printMinMax statement :</b> <br></div><div><br></div><div>Variable: olr<br>Type: double<br>Total Size: 5504200 bytes<br> 688025 values<br>Number of Dimensions: 3<br>Dimensions and sizes: [time | 365] x [lat | 13] x [lon | 145]<br>Coordinates: <br> time: [ 0..524160]<br> lat: [ -15.. 15]<br> lon: [ 0.. 360]<br>Number Of Attributes: 1<br> _FillValue : -999000000<br>ncl 19> printMinMax(olr,0)<br><b>(0) min=-185.5106410148043 max=77.65257890949087</b><br>ncl 20> print("------------------------------")<br>(0) ------------------------------<br>ncl 21> printVarSummary(u200)<br><br>Variable: u200<br>Type: double<br>Total Size: 5504200 bytes<br> 688025 values<br>Number of Dimensions: 3<br>Dimensions and sizes: [time | 365] x [lat | 13] x [lon | 145]<br>Coordinates: <br> time: [ 0..524160]<br> lat: [ -15.. 15]<br> lon: [ 0.. 360]<br>Number Of Attributes: 1<br> _FillValue : -999000000<br>ncl 22> printMinMax(u200,0)<br><b>(0) min=-48.11269687686956 max=57.21568789939967</b><br>ncl 23> print("------------------------------")<br>(0) ------------------------------<br>ncl 24> printVarSummary(u850)<br><br>Variable: u850<br>Type: double<br>Total Size: 5504200 bytes<br> 688025 values<br>Number of Dimensions: 3<br>Dimensions and sizes: [time | 365] x [lat | 13] x [lon | 145]<br>Coordinates: <br> time: [ 0..524160]<br> lat: [ -15.. 15]<br> lon: [ 0.. 360]<br>Number Of Attributes: 1<br> _FillValue : -999000000<br>ncl 25> printMinMax(u850,0)<br><b>(0) min=-31.22578327383305 max=41.24650545117836</b></div><div><br></div><div>I want to run mjoclivar_14.ncl (<a href="https://www.ncl.ucar.edu/Applications/mjoclivar.shtml">https://www.ncl.ucar.edu/Applications/mjoclivar.shtml</a>), but I got is always both min and max data value -999000000. Could somebody help me making multivariate EOF with this file :</div><div><br></div><div>OLR anomaly data:<div class="gmail_chip gmail_drive_chip" style="width:396px;height:18px;max-height:18px;background-color:rgb(245,245,245);padding:5px;color:rgb(34,34,34);font-family:arial;font-style:normal;font-weight:bold;font-size:13px;border:1px solid rgb(221,221,221);line-height:1"><a href="https://drive.google.com/file/d/1bH67Flh4XgfYfSfsrBRPd0F7Gwq9qce5/view?usp=drive_web" target="_blank" style="display:inline-block;max-width:366px;overflow:hidden;text-overflow:ellipsis;white-space:nowrap;text-decoration:none;padding:1px 0;border:none"><img style="vertical-align: bottom; border: none;" src="https://ssl.gstatic.com/docs/doclist/images/icon_10_generic_list.png"> <span dir="ltr" style="color:rgb(17,85,204);text-decoration:none;vertical-align:bottom">olranom.nc</span></a><img src="//ssl.gstatic.com/ui/v1/icons/common/x_8px.png" style="opacity: 0.55; cursor: pointer; float: right; top: -1px; display: none;"></div></div><div>200 hPa zonal wind anomaly data :<div class="gmail_chip gmail_drive_chip" style="width:396px;height:18px;max-height:18px;background-color:#f5f5f5;padding:5px;color:#222;font-family:arial;font-style:normal;font-weight:bold;font-size:13px;border:1px solid #ddd;line-height:1"><a href="https://drive.google.com/file/d/1Y_JmVgDatlxeQ2VsPSD4g-dh6Cc0XNna/view?usp=drive_web" target="_blank" style="display:inline-block;max-width:366px;overflow:hidden;text-overflow:ellipsis;white-space:nowrap;text-decoration:none;padding:1px 0;border:none"><img style="vertical-align: bottom; border: none;" src="https://ssl.gstatic.com/docs/doclist/images/icon_10_generic_list.png"> <span dir="ltr" style="color:#15c;text-decoration:none;vertical-align:bottom">uanom.200.nc</span></a><img src="//ssl.gstatic.com/ui/v1/icons/common/x_8px.png" style="display: none; opacity: 0.55; cursor: pointer; float: right; top: -1px;"></div></div><div>850 hPa zonal wind anomaly data:<div class="gmail_chip gmail_drive_chip" style="width:396px;height:18px;max-height:18px;background-color:#f5f5f5;padding:5px;color:#222;font-family:arial;font-style:normal;font-weight:bold;font-size:13px;border:1px solid #ddd;line-height:1"><a href="https://drive.google.com/file/d/1TwaQL7zF8J3wIyCs8mmnQC1lqimPFsCV/view?usp=drive_web" target="_blank" style="display:inline-block;max-width:366px;overflow:hidden;text-overflow:ellipsis;white-space:nowrap;text-decoration:none;padding:1px 0;border:none"><img style="vertical-align: bottom; border: none;" src="https://ssl.gstatic.com/docs/doclist/images/icon_10_generic_list.png"> <span dir="ltr" style="color:#15c;text-decoration:none;vertical-align:bottom">uanom.850.nc</span></a><img src="//ssl.gstatic.com/ui/v1/icons/common/x_8px.png" style="display: none; opacity: 0.55; cursor: pointer; float: right; top: -1px;"></div></div><div><br></div><div>This is my script for mjoclivar_14 :</div><div>ncl 0> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/diagnostics_cam.ncl"<br>dim_rmvmean_n( olr, 0) ; (time,lon)<br> u850 = dimncl 1> ;******************************************************<br>ncl 2> ;<br>ncl 3> ; mjoclivar_14.ncl<br>ncl 4> ;<br>ncl 5> ;***********************************************************<br>ncl 6> ; Combined EOFs<br>ncl 7> ; Latest Update: July, 2016: Eun-Pa Lim; Bureau of Meteorology, Australia<br>ncl 8> ;***********************************************************<br>ncl 9> ;;<br>ncl 10> ;; The following are automatically loaded from 6.2.0 onward<br>ncl 11> ;;load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl" <br>ncl 12> ;;load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl" <br>ncl 13> ;;load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl" <br>ncl 14> <br>ncl 15> undef("read_rename")<br>ncl 16> function read_rename(f[1]:file, varName[1]:string \<br>var_olr = dim_avg_n_Wrap( var_olr , 0) <br> zavg_var_u850 = dim_avg_n_Wrap( var_u850, 0)<br> zavg_var_u200 = dim_avg_n_Wrap(ncl 16> ,iStrt[1]:integer, iLast[1]:integer \<br>ncl 16> ,latS[1]:numeric , latN[1]:numeric )<br>ncl 17> ; Utility to force specific named dimensions<br>ncl 18> ; This is done for historical reasons (convenience) <br>ncl 19> begin<br>ncl 20> work = f->$varName$(iStrt:iLast,{latS:latN},:) ; (time,lat,lon)<br>ncl 21> work!0 = "time" ; CAM model names<br>ncl 22> work!1 = "lat"<br>ncl 23> work!2 = "lon"<br>ncl 24> return(work)<br>ncl 25> end<br>ncl 26> ; =========================> MAIN <==============================<br>ncl 27> begin<br>ncl 28> neof = 2<br>ncl 29> <br>ncl 30> latS = -15<br>ncl 31> latN = 15<br>ncl 32> <br>mncl 33> ymdStrt = 20180101 ; start yyyymmdd<br>ncl 34> ymdLast = 20181231 ; last <br>ncl 35> <br>yncl 36> yrStrt = ymdStrt/10000<br>ncl 37> yrLast = ymdLast/10000<br>ncl 38> <br>ancl 39> pltDir = "./" ; plot directory<br>ncl 40> pltType = "png" <br>ncl 41> pltName = "mjoclivar" <br>ncl 42> <br>ncl 43> diri = "./" ; input directory<br>ncl 44> <br>oncl 45> filolr = "<a href="http://olranom.nc">olranom.nc</a>"<br>ncl 46> filu200 = "<a href="http://uanom.200.nc">uanom.200.nc</a>"<br>ncl 47> filu850 = "<a href="http://uanom.850.nc">uanom.850.nc</a>"<br>ncl 48> <br>ncl 49> ;************************************************<br>ncl 50> ; create BandPass Filter<br>ncl 51> ;************************************************<br>oncl 52> ihp = 2 ; bpf=>band pass filter<br>ncl 53> nWgt = 201<br>ncl 54> sigma = 1.0 ; Lanczos sigma<br>ncl 55> fca = 1./100.<br>ncl 56> fcb = 1./20.<br>ncl 57> wgt = filwgts_lanczos (nWgt, ihp, fca, fcb, sigma )<br>ncl 58> <br>*ncl 59> ;***********************************************************<br>ncl 60> ; Find the indices corresponding to the start/end times<br>ncl 61> ;***********************************************************<br>ncl 62> f = addfile (diri+filolr , "r") <br>ncl 63> TIME = f->time ; days since ...<br>ncl 64> <br>ncl 65> YMD = cd_calendar(TIME, -2) ; entire (time,6)<br>ncl 66> <br>ncl 67> iStrt = ind(YMD.eq.ymdStrt) ; index start<br>ncl 68> iLast = ind(YMD.eq.ymdLast) ; index last <br>ncl 69> delete([/ TIME, YMD /])<br>ncl 70> <br> ncl 71> ;***********************************************************<br>ncl 72> ; Read anomalies<br>ncl 73> ;***********************************************************<br>ncl 74> <br>encl 75> work = read_rename(f,"olranom",iStrt,iLast,latS,latN) ; (time,lat,lon) <br>ncl 76> OLR = dim_avg_n_Wrap(work, 1) ; (time,lon)<br>ncl 77> delete(work)<br>ncl 78> <br>ncl 79> f = addfile (diri+filu850 , "r") <br>ncl 80> work = read_rename(f,"uanom",iStrt,iLast,latS,latN) ; (time,lat,lon) <br>ncl 81> U850 = dim_avg_n_Wrap(work, 1) ; (time,lon)<br>ncl 82> delete(work)<br>ncl 83> <br>ncl 84> f = addfile (diri+filu200 , "r") <br>ncl 85> work = read_rename(f,"uanom",iStrt,iLast,latS,latN) ; (time,lat,lon) <br>ncl 86> U200 = dim_avg_n_Wrap(work, 1) ; (time,lon)<br>ncl 87> <br>0ncl 88> dimw = dimsizes( work )<br>ncl 89> ntim = dimw(0)<br>ncl 90> nlat = dimw(1)<br>)ncl 91> mlon = dimw(2)<br>uncl 92> delete(work)<br>ncl 93> <br>oncl 94> lon = OLR&lon <br>ncl 95> time = OLR&time <br>ncl 96> date = cd_calendar(time, -2) ; yyyymmdd<br>ncl 97> <br>encl 98> ;************************************************<br>ncl 99> ; Apply the band pass filter to the original anomalies<br>ncl 100> ;************************************************<br>ncl 101> olr = wgt_runave_n_Wrap ( OLR, wgt, 0, 0) ; (time,lon)<br>ncl 102> u850 = wgt_runave_n_Wrap (U850, wgt, 0, 0)<br>ncl 103> u200 = wgt_runave_n_Wrap (U200, wgt, 0, 0)<br>ncl 104> <br>*ncl 105> ;************************************************<br>ncl 106> ; remove temporal means of band pass series: *not* necessary <br>ncl 107> ;************************************************<br>ncl 108> olr = dim_rmvmean_n( olr, 0) ; (time,lon)<br>ncl 109> u850 = dim_rmvmean_n(u850, 0)<br>ncl 110> u200 = dim_rmvmean_n(u200, 0)<br>ncl 111> <br>nncl 112> ;************************************************<br>ncl 113> ; Compute the temporal variance at each lon<br>ncl 114> ;************************************************<br>ncl 115> var_olr = dim_variance_n_Wrap( olr, 0) ; (lon)<br>ncl 116> var_u850 = dim_variance_n_Wrap(u850, 0)<br>ncl 117> var_u200 = dim_variance_n_Wrap(u200, 0)<br>ncl 118> <br>ncl 119> ;************************************************<br>ncl 120> ; Compute the zonal mean of the temporal variance<br>ncl 121> ;************************************************<br>ncl 122> zavg_var_olr = dim_avg_n_Wrap( var_olr , 0) <br>ncl 123> zavg_var_u850 = dim_avg_n_Wrap( var_u850, 0)<br>ncl 124> zavg_var_u200 = dim_avg_n_Wrap( var_u200, 0)<br>ncl 125> <br>rncl 126> ;************************************************<br>ancl 127> ; Normalize by sqrt(avg_var*)<br>ncl 128> ;************************************************<br>ncl 129> zsqrt = where(zavg_var_olr.gt.0,sqrt(zavg_var_olr),zavg_var_olr@_FillValue)<br>ncl 130> olr = olr/zsqrt ; (time,lon)<br>;ncl 131> z850sqrt =where(zavg_var_u850.gt.0,sqrt(zavg_var_u850),zavg_var_u850@_FillValue)<br>ncl 132> u850 = u850/z850sqrt<br>ncl 133> z200sqrt= where(zavg_var_u200.gt.0,sqrt(zavg_var_u200),zavg_var_u200@_FillValue)<br>)ncl 134> u200 = u200/z200sqrt<br>ncl 135> <br>cncl 136> ;************************************************<br>ncl 137> ; Combine the normalized data into one variable<br>ncl 138> ;************************************************<br>ncl 139> cdata = new ( (/3*mlon,ntim/), typeof(olr), getFillValue(olr))<br>ncl 140> do ml=0,mlon-1<br>ncl 141> cdata(ml ,:) = (/ olr(:,ml) /)<br>ncl 142> cdata(ml+ mlon,:) = (/ u850(:,ml) /)<br>ncl 143> cdata(ml+2*mlon,:) = (/ u200(:,ml) /)<br>ncl 144> end do<br>ncl 145> <br>=ncl 146> ;************************************************<br>ncl 147> ; Compute **combined** EOF; Sign of EOF is arbitrary<br>ncl 148> ;************************************************<br>ncl 149> eof_cdata = eofunc(cdata , neof, False) ; (neof,3*mlon)<br>ncl 150> print("==============")<br>ncl 151> printVarSummary(eof_cdata)<br>ncl 152> printMinMax(eof_cdata, True)<br>ncl 153> <br>)ncl 154> eof_ts_cdata = eofunc_ts(cdata,eof_cdata,False) ; (neof,3*ntim)<br>ncl 155> print("==============") <br>ncl 156> printVarSummary(eof_ts_cdata)<br>ncl 157> printMinMax(eof_ts_cdata, True)<br>ncl 158> <br>ncl 159> ;************************************************<br>ncl 160> ; For clarity, explicitly extract each variable. Create time series <br>ncl 161> ;************************************************<br>ncl 162> <br>tncl 163> nvar = 3 ; "olr", "u850", "u200"<br>ncl 164> ceof = new( (/nvar,neof,mlon/), typeof(cdata), getFillValue(cdata))<br>ncl 165> <br>dncl 166> do n=0,neof-1<br>ncl 167> ceof(0,n,:) = eof_cdata(n,0:mlon-1) ; olr<br>ncl 168> ceof(1,n,:) = eof_cdata(n,mlon:2*mlon-1) ; u850<br>ncl 169> ceof(2,n,:) = eof_cdata(n,2*mlon:) ; u200<br>ncl 170> end do<br>ncl 171> <br>(ncl 172> ceof!0 = "var"<br>ncl 173> ceof!1 = "eof"<br>ncl 174> ceof!2 = "lon" <br>ncl 175> ceof&lon = olr&lon<br>ncl 176> <br>incl 177> ceof_ts = new( (/nvar,neof,ntim/), typeof(cdata), getFillValue(cdata))<br>ncl 178> ceof_ts(0,:,:) = eofunc_ts_Wrap( olr(lon|:,time|:),ceof(0,:,:),False) ; (0,neof,ntim)<br>ncl 179> ceof_ts(1,:,:) = eofunc_ts_Wrap(u850(lon|:,time|:),ceof(1,:,:),False) ; (1,neof,ntim)<br>ncl 180> ceof_ts(2,:,:) = eofunc_ts_Wrap(u200(lon|:,time|:),ceof(2,:,:),False) ; (2,neof,ntim)<br>ncl 181> <br>tncl 182> ;**********************************************t*<br>ncl 183> ; Add code contributed by Marcus N. Morgan, Florida Institute of Technology; Feb 2015<br>ncl 184> ; Calculate % variance (pcv_ )accounted for by OLR, U850 and U200<br>ncl 185> ;************************************************<br>ncl 186> <br>ncl 187> pcv_eof_olr = new(neof,typeof(ceof))<br>ncl 188> pcv_eof_u850 = new(neof,typeof(ceof))<br>ncl 189> pcv_eof_u200 = new(neof,typeof(ceof))<br>ncl 190> <br>ncl 191> do n=0,neof-1<br>ncl 192> pcv_eof_olr(n) = avg((ceof(0,n,:)*sqrt(ceof@eval(n)))^2)*100<br>ncl 193> pcv_eof_u850(n) = avg((ceof(1,n,:)*sqrt(ceof@eval(n)))^2)*100<br>ncl 194> pcv_eof_u200(n) = avg((ceof(2,n,:)*sqrt(ceof@eval(n)))^2)*100<br>ncl 195> ;;print("pcv: neof="+(n+1)+": "+pcv_eof_olr(n)+" "+pcv_eof_u850(n)+" "+pcv_eof_u200(n))<br>ncl 196> end do<br>ncl 197> <br>ncl 198> ;************************************************<br>ncl 199> ; Change sign of EOFs for spatial structures of PC1 and PC2 <br>ncl 200> ; to represent convection over the tropical Indian Ocean and the tropical western Pacific Ocean, respectively <br>ncl 201> ; (Ad hoc approach) <br>ncl 202> ;************************************************<br>ncl 203> <br>ncl 204> imax_olr_eof1 = maxind(ceof(0,0,:)) <br>ncl 205> imax_olr_eof2 = maxind(ceof(0,1,:)) <br>ncl 206> <br>)ncl 207> lonmax_eof1 = ceof&lon(imax_olr_eof1) ; longitude of max value (i.e. suppressed convection)<br>ncl 208> lonmax_eof2 = ceof&lon(imax_olr_eof2)<br>ncl 209> <br>sncl 210> if (lonmax_eof1.ge.60 .and. lonmax_eof1.lt.180) then ; Change the sign of EOF1 <br>ncl 211> ceof(:,0,:) = -ceof(:,0,:) ; if OLR is positive<br>ncl 212> ceof_ts(:,0,:) = -ceof_ts(:,0,:) ; over the tropical Indian Ocean<br>ncl 213> eof_cdata(0,:) = -eof_cdata(0,:)<br></div><div>ncl 214> eof_ts_cdata(0,:) = -eof_ts_cdata(0,:)<br>ncl 215> end if<br>ncl 216> <br>encl 217> if (lonmax_eof2.ge.120 .and. lonmax_eof2.lt.180) then ; Change the sign of EOF2<br>ncl 218> ceof(:,1,:) = -ceof(:,1,:) ; if OLR is positive <br>ncl 219> ceof_ts(:,1,:) = -ceof_ts(:,1,:) ; over the tropical western Pacific Ocean<br>ncl 220> eof_cdata(1,:) = -eof_cdata(1,:)<br>ncl 221> eof_ts_cdata(1,:) = -eof_ts_cdata(1,:)<br>ncl 222> end if<br>ncl 223> <br>ncl 224> print("==============")<br>ncl 225> printVarSummary(eof_cdata)<br>ncl 226> printMinMax(eof_cdata, True)<br>nncl 227> <br>ncl 228> ;************************************************<br>ncl 229> ; Compute cross correlation of each variable's EOF time series at zero-lag<br>ncl 230> ;************************************************<br>ncl 231> r_olr_u850 = escorc(ceof_ts(0,:,:) , ceof_ts(1,:,:) ) ; (neof)<br>ncl 232> r_olr_u200 = escorc(ceof_ts(0,:,:) , ceof_ts(2,:,:) )<br>ncl 233> r_u850_u200 = escorc(ceof_ts(1,:,:) , ceof_ts(2,:,:) )<br>ncl 234> <br>ncl 235> print("==============")<br>ncl 236> do n=0,neof-1<br>ncl 237> print("neof="+n \<br>ncl 237> +" r_olr_u850=" +sprintf("%4.3f",r_olr_u850(n)) \<br>ncl 237> +" r_olr_u200=" +sprintf("%4.3f",r_olr_u200(n)) \<br>ncl 237> +" r_u850_u200="+sprintf("%4.3f",r_u850_u200(n)) )<br>ncl 238> end do<br>ncl 239> print("==============")<br>ncl 240> <br>ncl 241> ;************************************************<br>ncl 242> ; Compute cross correlation of the multivariate EOF; EOF 1 vs EOF 2<br>ncl 243> ;************************************************<br>ncl 244> <br>ncl 245> mxlag = 25<br>ncl 246> rlag_01 = esccr(eof_ts_cdata(0,:),eof_ts_cdata(1,:), mxlag) ; (N,mxlag+1)<br>ncl 247> rlag_10 = esccr(eof_ts_cdata(1,:),eof_ts_cdata(0,:), mxlag) ; (N,mxlag+1)<br>ncl 248> ccr_12 = new ( (/2*mxlag+1/), float) <br>ncl 249> <br>ncl 250> ccr_12(mxlag:) = rlag_10(0:mxlag) <br>ncl 251> ccr_12(0:mxlag) = rlag_01(::-1) ; reverse order<br>ncl 252> ;;print(ccr_12)<br>ncl 253> <br>ncl 254> <br>ncl 255> ;************************************************<br>ncl 256> ; Normalize the multivariate EOF 1&2 component time series<br>ncl 257> ; Compute (PC1^2+PC2^2): values > 1 indicate "strong" periods<br>encl 258> ;************************************************<br>ncl 259> eof_ts_cdata(0,:) = eof_ts_cdata(0,:)/stddev(eof_ts_cdata(0,:))<br>ncl 260> eof_ts_cdata(1,:) = eof_ts_cdata(1,:)/stddev(eof_ts_cdata(1,:))<br>ncl 261> <br>ncl 262> mjo_ts_index = eof_ts_cdata(0,:)^2 + eof_ts_cdata(1,:)^2 <br>ncl 263> mjo_ts_index_smt = runave(mjo_ts_index, 91, 0) ; 91-day running mean<br> ncl 264> <br>ncl 265> nGood = num(.not.ismissing(mjo_ts_index)) ; # non-missing<br>ncl 266> nStrong = num(mjo_ts_index .ge. 1.0)<br>ncl 267> print("nGood="+nGood+" nStrong="+nStrong+" nOther="+(nGood-nStrong))<br>ncl 268> <br>rncl 269> ;************************************************<br>uncl 270> ; Write PC results to netCDF for use in another example.<br>ncl 271> ;************************************************<br>ncl 272> mjo_ts_index!0 = "time"<br>t<br>ncl 273> mjo_ts_index&time = time <br>ncl 274> mjo_ts_index@long_name = "MJO PC INDEX" <br><br><br>;****ncl 275> mjo_ts_index@info = "(PC1^2 + PC2^2)" <br>ncl 276> <br>ncl 277> PC1 = eof_ts_cdata(0,:)<br>ncl 278> PC1!0 = "time"<br>********************************************<br>; MJO "strong" index <br>;********************ncl 279> PC1&time = time<br>ncl 280> PC1@long_name = "PC1"<br>ncl 281> PC1@info = "PC1/stddev(PC1)"<br>ncl 282> <br>oncl 283> PC2 = eof_ts_cdata(1,:)<br>ncl 284> PC2!0 = "time"<br>ncl 285> PC2&time = time<br>ncl 286> PC2@long_name = "PC2"<br>ncl 287> PC2@info = "PC2/stddev(PC2)"<br>ncl 288> <br> ncl 289> diro = "./"<br>ncl 290> filo = "MJO_PC_INDEX.nc"<br>ncl 291> system("/bin/rm -f "+diro+filo) ; remove any pre-existing file<br>ncl 292> ncdf = addfile(diro+filo,"c") ; open output netCDF file<br>ncl 293> ; make time an UNLIMITED dimension <br>ncl 294> filedimdef(ncdf,"time",-1,True) ; recommended for most applications<br>ncl 295> ; output variables directly<br>ncl 296> ncdf->MJO_INDEX = mjo_ts_index <br>ncl 297> ncdf->PC1 = PC1 <br>ncl 298> ncdf->PC2 = PC2 <br>ncl 299> <br>ncl 300> ;------------------------------------------------------------<br>ncl 301> ; PLOTS<br>ncl 302> ;------------------------------------------------------------<br>Pncl 303> <br>ncl 304> yyyymmdd = cd_calendar(time, -2)<br>ncl 305> yrfrac = yyyymmdd_to_yyyyfrac(yyyymmdd, 0.0)<br>ncl 306> delete([/ yrfrac@long_name, lon@long_name /])<br>ncl 307> <br>ncl 308> day = ispan(-mxlag, mxlag, 1)<br>ncl 309> ;day@long_name = "lag (day)"<br>ncl 310> <br>ncl 311> pltPath = pltDir+pltName<br>ncl 312> <br>ncl 313> wks = gsn_open_wks(pltType,pltPath)<br>ncl 314> gsn_define_colormap(wks,"default")<br>ncl 315> plot = new(3,graphic) <br>ncl 316> <br>ncl 317> ;************************************************<br>ncl 318> ; Multivariate EOF plots<br>ncl 319> ;************************************************<br>ncl 320> rts = True<br>ncl 321> rts@gsnDraw = False ; don't draw yet<br>ncl 322> rts@gsnFrame = False ; don't advance frame yet<br>ncl 323> rts@gsnScale = True ; force text scaling <br>ncl 324> <br>ncl 325> rts@vpHeightF = 0.40 ; Changes the aspect ratio<br>ncl 326> rts@vpWidthF = 0.85<br>ncl 327> rts@vpXF = 0.10 ; change start locations<br>ncl 328> rts@vpYF = 0.75 ; the plot<br>ncl 329> rts@xyLineThicknesses = (/2, 2, 2/)<br>ncl 330> rts@xyLineColors = (/"black","red","blue"/)<br>ncl 331> rts@gsnYRefLine = 0. ; reference line <br>ncl 332> rts@trXMaxF = max(lon)<br>ncl 333> rts@trXMinF = min(lon)<br>ncl 334> <br>ncl 335> rts@pmLegendDisplayMode = "Always" ; turn on legend<br>ncl 336> rts@pmLegendSide = "Top" ; Change location of <br>ncl 337> rts@pmLegendParallelPosF = 1.16 ; move units right<br>ncl 338> rts@pmLegendOrthogonalPosF = -0.50 ; move units down<br>ncl 339> rts@pmLegendWidthF = 0.15 ; Change width and<br>ncl 340> rts@pmLegendHeightF = 0.15 ; height of legend.<br>ncl 341> rts@lgLabelFontHeightF = 0.0175<br>ncl 342> <br>ncl 343> <br>ncl 344> rtsP = True ; modify the panel plot<br>ncl 345> ; rtsP@gsnMaximize = True ; large format<br>ncl 346> rtsP@gsnPanelMainString = "Multivariate EOF: 15S-15N: "+yrStrt+"-"+yrLast <br>ncl 347> <br>ncl 348> do n=0,neof-1<br>ncl 349> rts@xyExplicitLegendLabels = (/"OLR: "+sprintf("%4.1f", pcv_eof_u200(n)) +"%" \<br>ncl 349> ,"U850: "+sprintf("%4.1f", pcv_eof_u850(n))+"%" \<br>ncl 349> ,"U200: "+sprintf("%4.1f", pcv_eof_olr(n))+"%" /)<br>ncl 350> rts@gsnLeftString = "EOF "+(n+1)<br>ncl 351> rts@gsnRightString = sprintf("%3.1f",ceof@pcvar(n)) +"%"<br>ncl 352> plot(n) = gsn_csm_xy (wks,lon,ceof(:,n,:),rts)<br>ncl 353> end do<br>ncl 354> gsn_panel(wks,plot(0:1),(/2,1/),rtsP) ; now draw as one plot<br>ncl 355> <br>ncl 356> ;-----------------------------------------<br>ncl 357> ; The following doesn't work with some older versions of NCL<br>ncl 358> ; With old versions, the user must delete each individually.<br>ncl 359> ;----------------------------------------- <br>ncl 360> delete([/ rts@xyExplicitLegendLabels, rts@pmLegendDisplayMode \<br>ncl 360> , rts@xyLineThicknesses , rts@gsnLeftString \<br>ncl 360> , rts@gsnRightString , rts@xyLineColors \<br>ncl 360> , rts@trXMaxF , rts@trXMinF /] )<br>ncl 361> <br>ncl 362> lag = ispan(-mxlag,mxlag,1)<br>ncl 363> lag@long_name = "lag (days)"<br>ncl 364> <br>ncl 365> plot(0) = gsn_csm_xy (wks, lag ,ccr_12,rts)<br>ncl 366> rtsP@gsnPanelMainString = "Cross Correlation: Multivariate EOF: 15S-15N: " \<br>ncl 366> + yrStrt+"-"+yrLast <br>ncl 367> rtsP@gsnPaperOrientation = "portrait" ; force portrait<br>ncl 368> gsn_panel(wks,plot(0),(/1,1/),rtsP) ; now draw as one plot<br>ncl 369> <br>ncl 370> ;************************************************<br>ncl 371> ; MJO "strong" index <br>ncl 372> ;************************************************<br>ncl 373> rts@gsnYRefLine = 1.0<br>ncl 374> rts@gsnYRefLineColor = "black"<br>ncl 375> rts@xyMonoDashPattern = True<br>ncl 376> rts@xyLineColors = (/"black", "blue"/)<br>ncl 377> rts@xyLineThicknesses = (/1, 2/)<br>ncl 378> rts@pmLegendDisplayMode = "Always" ; turn on legend<br>ncl 379> rts@pmLegendWidthF = 0.12 ; Change width and<br>ncl 380> rts@pmLegendHeightF = 0.10 ; height of legend.<br>ncl 381> rts@pmLegendParallelPosF = 0.86 ; move units right<br>ncl 382> rts@pmLegendOrthogonalPosF = -0.40 ; move units down<br>ncl 383> rts@xyExplicitLegendLabels = (/"daily", "91-day runavg" /)<br>ncl 384> <br>ncl 385> mjo_ind_plt = new ( (/2,ntim/), typeof(mjo_ts_index))<br>ncl 386> mjo_ind_plt(0,:) = mjo_ts_index<br>ncl 387> mjo_ind_plt(1,:) = (/ mjo_ts_index_smt /)<br>ncl 388> plot(0) = gsn_csm_xy(wks, yrfrac,mjo_ind_plt,rts)<br>ncl 389> <br>ncl 390> rtsP@gsnPanelMainString = "MJO Index: (PC1^2+ PC2^2) : 15S-15N: "+yrStrt+"-"+yrLast <br>ncl 391> gsn_panel(wks,plot(0),(/1,1/),rtsP) ; now draw as one plot<br>ncl 392> <br>ncl 393> end<br><br><br> sub drveof: ier,jopt= 10 0 returned from vcmssm/crmssm<br>warning:eofunc: 10 eigenvectors failed to converge<br>(0) ==============<br><br>Variable: eof_cdata<br>Type: double<br>Total Size: 6960 bytes<br> 870 values<br>Number of Dimensions: 2<br>Dimensions and sizes: [2] x [435]<br>Coordinates: <br>Number Of Attributes: 5<br> _FillValue : -999000000<br> method : no transpose<br> matrix : covariance<br> pcvar : ( -9.99e+08, -9.99e+08 )<br> eval : ( -999000000, -999000000 )<br>(0) <br>(0) min=-999000000 max=-999000000<br>warning:eofunc_ts: 2 eigenvectors failed to converge<br>(0) ==============<br><br>Variable: eof_ts_cdata<br>Type: double<br>Total Size: 5840 bytes<br> 730 values<br>Number of Dimensions: 2<br>Dimensions and sizes: [2] x [365]<br>Coordinates: <br>Number Of Attributes: 3<br> _FillValue : -999000000<br> matrix : covariance<br> ts_mean : ( -999000000, -999000000 )<br>(0) <br>(0) min=-999000000 max=-999000000<br><b>warning:eofunc_ts: 2 eigenvectors failed to converge<br>warning:eofunc_ts: 2 eigenvectors failed to converge<br>warning:eofunc_ts: 2 eigenvectors failed to converge<br>warning:avg: Entire input array contains missing values, can't compute average<br>warning:avg: Entire input array contains missing values, can't compute average<br>warning:avg: Entire input array contains missing values, can't compute average<br>warning:avg: Entire input array contains missing values, can't compute average<br>warning:avg: Entire input array contains missing values, can't compute average<br>warning:avg: Entire input array contains missing values, can't compute average<br>fatal:Subscript out of range, error in subscript #0<br>fatal:An error occurred reading unknown<br>fatal:["Execute.c":8635]:Execute: Error occurred at or near line 207</b></div><div><br></div><div>Thank you. <br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Fri, Dec 20, 2019 at 1:15 PM Barry Lynn <<a href="mailto:barry.h.lynn@gmail.com">barry.h.lynn@gmail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr">Hi:<div><br></div><div>Could you please define what you mean by "missing." Did you PrntMinMax(XXXX,False)?</div><div><br></div><div><br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Fri, Dec 20, 2019 at 1:20 AM Rahpeni Fajarianti via ncl-talk <<a href="mailto:ncl-talk@ucar.edu" target="_blank">ncl-talk@ucar.edu</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div>Hi Dennis, Thank you for your time. <br></div><div><br></div><div>I am going to get OLR and (850 hPa and 200 hPa) zonal wind anomaly for year 2018 using mjoclivar_2.ncl (<a href="https://www.ncl.ucar.edu/Applications/Scripts/mjoclivar_2.ncl" target="_blank">https://www.ncl.ucar.edu/Applications/Scripts/mjoclivar_2.ncl</a>) to create netcdf file as input data for mjoclivar_14. <br></div><div><br></div><div>The dataset I used were OLR daily mean (<a href="https://www.esrl.noaa.gov/psd/data/gridded/data.interp_OLR.html" target="_blank">https://www.esrl.noaa.gov/psd/data/gridded/data.interp_OLR.html</a>) and u wind daily mean (<a href="https://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis.pressure.html" target="_blank">https://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis.pressure.html</a>). <br></div><div><br></div><div>When I run the script and I got the result, the wind anomaly data is missing. I used print statement to make sure the values.</div><div>can anybody help to fix this? Thank you<br></div><div><br></div><div><b>This is what I got for OLR anomaly:</b></div><div>Variable: f2<br>Type: file<br>filename: olr.day.anomalies.2018<br>path: /home/peni/<a href="http://olr.day.anomalies.2018.nc" target="_blank">olr.day.anomalies.2018.nc</a><br> file global attributes:<br> creation_date : Sel Des 17 08:24:41 WIB 2019<br> source_file : <a href="http://olr.day.mean.nc" target="_blank">olr.day.mean.nc</a><br> title : OLR: Daily Anomalies: 2018-2018<br> dimensions:<br> time = 365 // unlimited<br> lat = 73<br> lon = 144<br> variables:<br> double time ( time )<br> units : hours since 1800-01-01 00:00:0.0<br> long_name : Time<br> actual_range : ( 1528872, 1919688 )<br> delta_t : 0000-00-01 00:00:00<br> standard_name : time<br> axis : T<br><br> float lat ( lat )<br> units : degrees_north<br> actual_range : ( 90, -90 )<br> long_name : Latitude<br> standard_name : latitude<br> axis : Y<br><br> float lon ( lon )<br> units : degrees_east<br> long_name : Longitude<br> actual_range : ( 0, 360 )<br> standard_name : longitude<br> axis : X<br><br> float OLR_anom ( time, lat, lon )<br> _FillValue : 32766<br> long_name : Anomalies: Daily OLR<br> units : W/m^2<br><br> float OLR_anom_sm ( time, lat, lon )<br> _FillValue : 32766<br> long_name : Anomalies from Smooth Daily Climatology<br> units : W/m^2</div><div><br></div><div><br></div><div><b>this is what I got for zonal wind anomaly:</b></div><div>Variable: f3<br>Type: file<br>filename: uwnd.anomalies<br>path: /home/peni/<a href="http://uwnd.anomalies.nc" target="_blank">uwnd.anomalies.nc</a><br> file global attributes:<br> creation_date : Sel Des 17 12:56:50 WIB 2019<br> source_file : <a href="http://uwnd.day.mean.200.nc" target="_blank">uwnd.day.mean.200.nc</a><br> dimensions:<br> time = 365 // unlimited<br> lat = 73<br> lon = 145<br> variables:<br> double time ( time )<br> long_name : Time<br> units : minutes since 2018-01-01 00:00<br><br> double lat ( lat )<br> units : degrees_north<br> long_name : Latitude<br><br> double lon ( lon )<br> units : degrees_east<br> long_name : Longitude<br><br> float U_anom ( time, lat, lon )<br> long_name : Anomalies from Daily Climatology<br> _FillValue : -9.99e+08<br><br> float U_anom_sm ( time, lat, lon )<br> long_name : Anomalies from Smooth Daily Climatology<br> _FillValue : -9.99e+08</div><div><br></div><div><br></div><div><b>This is my script:</b></div><div>;-------------------------------------------------------------<br>; mjoclivar_2.ncl<br>;-------------------------------------------------------------<br>; These files are loaded by default in NCL V6.2.0 and newer<br>; load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl" <br>; load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl" <br>; load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl" <br>;******************************************************<br><br>begin<br>;******************************************************<br>; User specifications<br>;******************************************************<br> NC = True ; create netCDF? <br> PLOT = False ; sample plots?<br><br> ymdStrt = 20180101 ; start yyyymmdd<br> ymdLast = 20181231 ; last <br> yrStrt = ymdStrt/10000<br> yrLast = ymdLast/10000<br> <br> nhar = 4 ; number of fourier comp <br> <br> var = "uwnd2" ; name of file <br> vName = "U" ; name for plots<br><br> ;diri = "/project/cas/shea/CDC/" ; input directory <br> ;diri = "/Users/shea/Data/AMWG/" ; input directory <br> diri = "/home/peni/" ; new input directory<br> fili = var+".<a href="http://day.mean.850.nc" target="_blank">day.mean.850.nc</a>" ; input file <br><br> if (NC) then<br> diro= "/home/peni/" ; output dir<br> ;filo= var+".<a href="http://day.anomalies.nc" target="_blank">day.anomalies.nc</a>" ; output file <br> filo= var+".<a href="http://anomalies.nc" target="_blank">anomalies.nc</a>" ; output file <br> end if<br><br> if (PLOT) then<br> wksType = "png" ; send graphics to PNG file<br> wksName = "mjoclivar" ; "mjo."+yrStrt+"_"+yrLast<br> end if<br><br>;***********************************************************<br>; Read user specified time and create required yyyyddd <br>;***********************************************************<br> f = addfile (diri+fili , "r") <br><br> time = f->time ; all times on file<br> ymd = cd_calendar(time, -2) ; yyyymmdd<br> iStrt = ind(ymd.eq.ymdStrt) ; index start<br> iLast = ind(ymd.eq.ymdLast) ; index last <br> delete(time)<br> delete(ymd)<br><br>;***********************************************************<br>; Read user specified time and create required yyyyddd <br>;***********************************************************<br> time = f->time(iStrt:iLast) ; time:units = "hours since"<br> TIME = cd_calendar(time, 0) ; type float <br> year = floattointeger( TIME(:,0) )<br> month = floattointeger( TIME(:,1) )<br> day = floattointeger( TIME(:,2) ) <br> ddd = day_of_year(year, month, day) <br> yyyyddd = year*1000 + ddd ; needed for input<br><br>;***********************************************************<br>; Read data: dble2flt <br>;*********************************************************** <br> x = dble2flt( f->$var$(iStrt:iLast,:,:) ) ; convert to float <br> printVarSummary( x ) <br><br>;***********************************************************<br>; Compute daily climatology: raw and then 'smoothed' <br>;***********************************************************<br> xClmDay = clmDayTLL(x, yyyyddd) ; daily climatology at each grid point<br> printVarSummary(xClmDay) <br><br>;***********************************************************<br>; Compute smoothed daily climatology using 'nhar' harmonics<br>;***********************************************************<br> xClmDay_sm = smthClmDayTLL(xClmDay, nhar)<br> printVarSummary(xClmDay_sm)<br><br>;***********************************************************<br>; Compute daily anomalies using raw and smoothed climatologies<br>;***********************************************************<br> xAnom = calcDayAnomTLL (x, yyyyddd, xClmDay) <br> printVarSummary(xAnom)<br> printMinMax(xAnom, True)<br><br> xAnom_sm = calcDayAnomTLL (x, yyyyddd, xClmDay_sm) <br> xAnom_sm@long_name = "Anomalies from Smooth Daily Climatology"<br> printVarSummary(xAnom_sm)<br> printMinMax(xAnom_sm, True)<br><br> delete( x ) ; no longer needed<br><br><br>;***********************************************************<br>; Create netCDF: convenience use 'simple' method<br>;***********************************************************<br><br> dimx = dimsizes(xAnom)<br> ntim = dimx(0)<br> nlat = dimx(1)<br> mlon = dimx(2)<br><br> if (NC) then<br><br> lat = f->lat<br> lon = f->lon<br><br> system("/bin/rm -f "+diro+filo) ; rm any pre-exist file, if any<br> fnc = addfile (diro+filo, "c")<br><br> filAtt = 0<br> filAtt@title = vName+": Daily Anomalies: "+yrStrt+" <br> filAtt@source_file = fili<br> filAtt@creation_date = systemfunc("date")<br> fileattdef( fnc, filAtt ) ; copy file attributes <br><br> setfileoption(fnc,"DefineMode",True)<br> <br> varNC = vName+"_anom"<br> varNC_sm = vName+"_anom_sm"<br><br> dimNames = (/"time", "lat", "lon"/) <br> dimSizes = (/ -1 , nlat, mlon/) <br> dimUnlim = (/ True , False, False/) <br> filedimdef(fnc,dimNames,dimSizes,dimUnlim)<br><br> filevardef(fnc, "time" ,typeof(time),getvardims(time)) <br> filevardef(fnc, "lat" ,typeof(lat) ,getvardims(lat)) <br> filevardef(fnc, "lon" ,typeof(lon) ,getvardims(lon))<br> filevardef(fnc, varNC ,typeof(xAnom) ,getvardims(xAnom)) <br> filevardef(fnc, varNC_sm,typeof(xAnom) ,getvardims(xAnom)) <br><br> filevarattdef(fnc,"time" ,time) ; copy time attributes<br> filevarattdef(fnc,"lat" ,lat) ; copy lat attributes<br> filevarattdef(fnc,"lon" ,lon) ; copy lon attributes<br> filevarattdef(fnc,varNC ,xAnom) <br> filevarattdef(fnc,varNC_sm,xAnom_sm) <br><br> fnc->time = (/time/) <br> fnc->lat = (/lat/)<br> fnc->lon = (/lon/)<br> fnc->$varNC$ = (/xAnom /)<br> fnc->$varNC_sm$ = (/xAnom_sm/)<br> end if</div><div>end</div><div><br></div><div><br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Thu, Dec 19, 2019 at 9:59 AM Rahpeni Fajarianti <<a href="mailto:rahpenifajarianti@gmail.com" target="_blank">rahpenifajarianti@gmail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="auto"></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">---------- Forwarded message ---------<br>Dari: <strong class="gmail_sendername" dir="auto">Rahpeni Fajarianti</strong> <span dir="auto"><<a href="mailto:rahpenifajarianti@gmail.com" target="_blank">rahpenifajarianti@gmail.com</a>></span><br>Date: Kam, 19 Des 2019 9:10 AM<br>Subject: Fwd: [ncl-talk] fatal:divide: Division by 0, Can't continuefatal:Div: operator failed, can't continue fatal:["Execute.c":8635]:Execute: Error occurred at or near line 129<br>To: reza tisa <<a href="mailto:rezatisa1@gmail.com" target="_blank">rezatisa1@gmail.com</a>><br></div><br><br><div dir="auto"></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">---------- Forwarded message ---------<br>Dari: <strong class="gmail_sendername" dir="auto">Dennis Shea</strong> <span dir="auto"><<a href="mailto:shea@ucar.edu" rel="noreferrer" target="_blank">shea@ucar.edu</a>></span><br>Date: Kam, 19 Des 2019 7:57 AM<br>Subject: Re: [ncl-talk] fatal:divide: Division by 0, Can't continuefatal:Div: operator failed, can't continue fatal:["Execute.c":8635]:Execute: Error occurred at or near line 129<br>To: Rahpeni Fajarianti <<a href="mailto:rahpenifajarianti@gmail.com" rel="noreferrer" target="_blank">rahpenifajarianti@gmail.com</a>><br>Cc: Ncl-talk <<a href="mailto:ncl-talk@ucar.edu" rel="noreferrer" target="_blank">ncl-talk@ucar.edu</a>><br></div><br><br><div dir="ltr"><div>The GOLDEN RULE of data processing is: <b><span style="color:rgb(255,0,0)">***LOOK AT YOUR DATA***</span></b></div><div><br></div><div>This is *user* responsibility. <b>*ALL* the data are 0.0</b></div><div><br></div><div><b>%></b> ncl rdData.ncl</div><div>=================</div><div>Variable: olr<br>Dimensions and sizes: [time | 365] x [lat | 13] x [lon | 144]<br>Coordinates: <br> time: [1910952..1919688]<br> lat: [-15..15]<br> lon: [ 0..357.5]<br>(0) A<span style="color:rgb(0,0,255)"><b>nomalies: Daily OLR (W/m^2) : min=0 max=0</b></span><br>(0) ------------------------------<br><br>Variable: u200<br>Dimensions and sizes: [time | 365] x [lat | 13] x [lon | 145]<br>Coordinates: <br> time: [ 0..524160]<br> lat: [ -15.. 15]<br> lon: [ 0.. 360]<br>Number Of Attributes: 2<br> long_name : Anomalies from Daily Climatology<br> _FillValue : -9.99e+08<br>(0) <span style="color:rgb(0,0,255)"><b>Anomalies from Daily Climatology : min=0 max=0</b></span><br>(0) ------------------------------<br><br>Variable: u850<br>Dimensions and sizes: [time | 365] x [lat | 13] x [lon | 145]<br>Coordinates: <br> time: [ 0..524160]<br> lat: [ -15.. 15]<br> lon: [ 0.. 360]<br>Number Of Attributes: 2<br> _FillValue : -9.99e+08<br> <br>(0) <b><span style="color:rgb(0,0,255)">Anomalies from Daily Climatology : min=0 max=0</span></b><br></div><div><br></div><br><div><br></div><div><br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Wed, Dec 18, 2019 at 4:58 PM Rahpeni Fajarianti <<a href="mailto:rahpenifajarianti@gmail.com" rel="noreferrer noreferrer" target="_blank">rahpenifajarianti@gmail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div>Thank you so much for responding. I attach some files here. The data and script that I used. <br></div><div>The data : <a href="http://olr.day.anomalies.2018.nc" rel="noreferrer noreferrer" target="_blank">olr.day.anomalies.2018.nc</a> , <a href="http://uwnd.anomalies.nc" rel="noreferrer noreferrer" target="_blank">uwnd.anomalies.nc</a> (for 200hPa u wind anomalies) , <a href="http://uwnd2.anomalies.nc" rel="noreferrer noreferrer" target="_blank">uwnd2.anomalies.nc</a> (for 850hPa u wind anomalies)</div><div>my script: mjoclivar_14.ncl</div><div><br></div><div>I am very new in NCL and very pleased for help. Thank you.<br></div><div><div class="gmail_chip gmail_drive_chip" style="width:396px;height:18px;max-height:18px;background-color:rgb(245,245,245);padding:5px;color:rgb(34,34,34);font-family:arial;font-style:normal;font-weight:bold;font-size:13px;border:1px solid rgb(221,221,221);line-height:1"><a href="https://drive.google.com/file/d/1nstwpHwEVEB9nbuWQPnB9nntN2YVBz79/view?usp=drive_web" style="display:inline-block;max-width:366px;overflow:hidden;text-overflow:ellipsis;white-space:nowrap;text-decoration:none;padding:1px 0px;border:medium none" rel="noreferrer noreferrer" target="_blank"><img style="vertical-align: bottom; border: medium none;" src="https://ssl.gstatic.com/docs/doclist/images/icon_10_generic_list.png"> <span dir="ltr" style="color:rgb(17,85,204);text-decoration:none;vertical-align:bottom">mjoclivar_14.ncl</span></a><img style="opacity: 0.55; float: right; display: none;"></div><div class="gmail_chip gmail_drive_chip" style="width:396px;height:18px;max-height:18px;background-color:rgb(245,245,245);padding:5px;color:rgb(34,34,34);font-family:arial;font-style:normal;font-weight:bold;font-size:13px;border:1px solid rgb(221,221,221);line-height:1"><a href="https://drive.google.com/file/d/1RQIk85UqimEzDzqWZ9V97ZpbqRG8h-cz/view?usp=drive_web" style="display:inline-block;max-width:366px;overflow:hidden;text-overflow:ellipsis;white-space:nowrap;text-decoration:none;padding:1px 0px;border:medium none" rel="noreferrer noreferrer" target="_blank"><img style="vertical-align: bottom; border: medium none;" src="https://ssl.gstatic.com/docs/doclist/images/icon_10_generic_list.png"> <span dir="ltr" style="color:rgb(17,85,204);text-decoration:none;vertical-align:bottom">olr.day.anomalies.2018.nc</span></a><img style="display: none; opacity: 0.55; float: right;"></div><div class="gmail_chip gmail_drive_chip" style="width:396px;height:18px;max-height:18px;background-color:rgb(245,245,245);padding:5px;color:rgb(34,34,34);font-family:arial;font-style:normal;font-weight:bold;font-size:13px;border:1px solid rgb(221,221,221);line-height:1"><a href="https://drive.google.com/file/d/1ZxF2tZnCr7q1LFQWtGhMOLvaIE3mwK-p/view?usp=drive_web" style="display:inline-block;max-width:366px;overflow:hidden;text-overflow:ellipsis;white-space:nowrap;text-decoration:none;padding:1px 0px;border:medium none" rel="noreferrer noreferrer" target="_blank"><img style="vertical-align: bottom; border: medium none;" src="https://ssl.gstatic.com/docs/doclist/images/icon_10_generic_list.png"> <span dir="ltr" style="color:rgb(17,85,204);text-decoration:none;vertical-align:bottom">uwnd2.anomalies.nc</span></a><img style="display: none; opacity: 0.55; float: right;"></div><div class="gmail_chip gmail_drive_chip" style="width:396px;height:18px;max-height:18px;background-color:rgb(245,245,245);padding:5px;color:rgb(34,34,34);font-family:arial;font-style:normal;font-weight:bold;font-size:13px;border:1px solid rgb(221,221,221);line-height:1"><a href="https://drive.google.com/file/d/1rm-f84yy1PkKlA8ZBG87nRGDbnD8auqx/view?usp=drive_web" style="display:inline-block;max-width:366px;overflow:hidden;text-overflow:ellipsis;white-space:nowrap;text-decoration:none;padding:1px 0px;border:medium none" rel="noreferrer noreferrer" target="_blank"><img style="vertical-align: bottom; border: medium none;" src="https://ssl.gstatic.com/docs/doclist/images/icon_10_generic_list.png"> <span dir="ltr" style="color:rgb(17,85,204);text-decoration:none;vertical-align:bottom">uwnd.anomalies.nc</span></a><img style="display: none; opacity: 0.55; float: right;"></div></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Thu, Dec 19, 2019 at 2:22 AM Dennis Shea <<a href="mailto:shea@ucar.edu" rel="noreferrer noreferrer" target="_blank">shea@ucar.edu</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div>This is offlineI</div><div><br></div><div><div>If you do not know NCL well, please read the excellent tutorial at <br></div><div><a href="http://www.ncl.ucar.edu/Document/Manuals/NCL_User_Guide/" rel="noreferrer noreferrer" target="_blank"><b>http://www.ncl.ucar.edu/Document/Manuals/NCL_User_Guide/</b></a></div><div>---</div></div><div> We [ ncl-talk ] is frustrated by this thread. Please understand, nobody is paid to answer ncl-talk question. We do it on our own time to help people. Debugging is VERY time consuming. In particular, somebody else's code.<br></div><div>e<br></div><div>The NCL example scripts indicate an approach to use but are not necessarily plug in a file or files and use them. Often, users must make some changes.</div><div>---<br></div><div><div>This is the result:</div><div>Variable: olr<br>Type: float<br></div><div>[SNIP]<br></div><div>Number Of Attributes: 5<br> _FillValue : 32766<br><br></div><div>(0,0) 32766 ,==== All values are missing [_FillValue ]<br></div><div>[SNIP] <br></div><div><br></div><div>===</div><div>Well .... at what point in your script did things go 'wrong'?</div><div>Use printVarSummary(), printMinMax(), etc</div><div>===</div><div><br></div><div>I will look at the issue(s) if you make your data available via (say) google drive or ftp.</div><div><br></div><div>ftp <a href="http://ftp.cgd.ucar.edu" rel="noreferrer noreferrer" target="_blank">ftp.cgd.ucar.edu</a></div><div>anonymous</div><div>your_email</div><div>cd incoming</div><div>put ...</div><div>put ...</div><div>put ...</div><div>quit<br></div><div><br></div><div><br></div><div><br></div><div><br></div><div><br></div><div><br></div></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Wed, Dec 18, 2019 at 8:14 AM Rahpeni Fajarianti via ncl-talk <<a href="mailto:ncl-talk@ucar.edu" rel="noreferrer noreferrer" target="_blank">ncl-talk@ucar.edu</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div>Thanks for your help. I had use print and printVarSummary statement. <br></div><div><br></div><div>This is the result:</div><div>Variable: olr<br>Type: float<br>Total Size: 210240 bytes<br> 52560 values<br>Number of Dimensions: 2<br>Dimensions and sizes: [time | 365] x [lon | 144]<br>Coordinates: <br> time: [1910952..1919688]<br> lon: [ 0..357.5]<br>Number Of Attributes: 5<br> _FillValue : 32766<br> units : W/m^2<br> long_name : Anomalies: Daily OLR<br> average_op_ncl : dim_avg_n over dimension(s): lat<br> wgt_runave_op_ncl : wgt_runave_n<br>(0,0) 32766<br>(0,1) 32766<br>(0,2) 32766<br>(0,3) 32766<br>(0,4) 32766<br>(0,5) 32766<br>(0,6) 32766<br>(0,7) 32766<br>(0,8) 32766<br>(0,9) 32766<br>(0,10) 32766<br>(0,11) 32766<br>(0,12) 32766<br>(0,13) 32766<br>(0,14) 32766<br>(0,15) 32766<br>(0,16) 32766<br>(0,17) 32766<br>(0,18) 32766<br>(0,19) 32766<br>(0,20) 32766<br>(0,21) 32766<br>(0,22) 32766<br>(0,23) 32766<br>(0,24) 32766<br>(0,25) 32766<br>(0,26) 32766<br>(0,27) 32766<br>(0,28) 32766<br>(0,29) 32766<br>(0,30) 32766<br>(0,31) 32766<br>(0,32) 32766<br>(0,33) 32766<br>(0,34) 32766<br>(0,35) 32766<br>(0,36) 32766<br>(0,37) 32766<br>(0,38) 32766<br>(0,39) 32766<br>(0,40) 32766<br>(0,41) 32766<br>(0,42) 32766<br>(0,43) 32766<br>(0,44) 32766<br>(0,45) 32766<br>(0,46) 32766<br>(0,47) 32766<br>(0,48) 32766<br>(0,49) 32766<br>(0,50) 32766<br>(0,51) 32766<br>(0,52) 32766<br>(0,53) 32766<br>(0,54) 32766<br>(0,55) 32766<br>(0,56) 32766<br>(0,57) 32766<br>(0,58) 32766<br>(0,59) 32766<br>(0,60) 32766<br>(0,61) 32766<br>(0,62) 32766<br>(0,63) 32766<br>(0,64) 32766<br>(0,65) 32766<br>(0,66) 32766<br>(0,67) 32766<br>(0,68) 32766<br>(0,69) 32766<br>(0,70) 32766<br>(0,71) 32766<br>(0,72) 32766<br>(0,73) 32766<br>(0,74) 32766<br>(0,75) 32766<br>(0,76) 32766<br>(0,77) 32766<br>(0,78) 32766<br>(0,79) 32766<br>(0,80) 32766<br>(0,81) 32766<br>(0,82) 32766<br>(0,83) 32766<br>(0,84) 32766<br>(0,85) 32766<br>(0,86) 32766<br>(0,87) 32766<br>(0,88) 32766<br>(0,89) 32766<br>(0,90) 32766<br>(0,91) 32766<br>(0,92) 32766<br>(0,93) 32766<br>(0,94) 32766<br>(0,95) 32766<br>(0,96) 32766<br>(0,97) 32766<br>(0,98) 32766<br>(0,99) 32766<br>(0,100) 32766<br>(0,101) 32766<br>(0,102) 32766<br>(0,103) 32766<br>(0,104) 32766<br>(0,105) 32766<br>(0,106) 32766<br>(0,107) 32766<br>(0,108) 32766<br>(0,109) 32766<br>(0,110) 32766<br>(0,111) 32766<br>(0,112) 32766<br>(0,113) 32766<br>(0,114) 32766<br>(0,115) 32766<br>(0,116) 32766<br>(0,117) 32766<br>(0,118) 32766<br>(0,119) 32766<br>(0,120) 32766<br>(0,121) 32766<br>(0,122) 32766<br>(0,123) 32766<br>(0,124) 32766<br>(0,125) 32766<br>(0,126) 32766<br>(0,127) 32766<br>(0,128) 32766<br>(0,129) 32766<br>(0,130) 32766<br>(0,131) 32766<br>(0,132) 32766<br>(0,133) 32766<br>(0,134) 32766<br>(0,135) 32766<br>(0,136) 32766<br>(0,137) 32766<br></div><div>(0,138) 32766 ect..</div><div><br></div><div>What should I fix to my script? <br></div><div>Thanks. <br></div><div><br></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Wed, Dec 18, 2019 at 8:50 AM Adam Phillips <<a href="mailto:asphilli@ucar.edu" rel="noreferrer noreferrer" target="_blank">asphilli@ucar.edu</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr">Hi Rahpeni,<div>Please always include ncl-talk on your replies; that way others can assist and see the correspondence. </div><div><br></div><div>You stated:</div><div><div>Thank you so much, it helps ! </div><div><br></div><div>I am still getting some errors.<br></div><div>I got :</div><div>fatal:Subscript out of range, error in subscript #1<br>fatal:An error occurred reading olr<br>fatal:["Execute.c":8635]:Execute: Error occurred at or near line 141</div></div><div><br></div><div><br></div><div>Note that it is for you to attempt to diagnose the issue before responding to ncl-talk. The error message is telling you what the problem is: You are attempting to access an index of a dimension that is bigger than the dimension itself. The lines around 141 are (from the email you sent me):</div><div><br></div><div>ncl 140> do ml=0,mlon-1<br>ncl 141> cdata(ml ,:) = (/ olr(:,ml) /)<br>ncl 142> cdata(ml+ mlon,:) = (/ u850(:,ml) /)<br>ncl 143> cdata(ml+2*mlon,:) = (/ u200(:,ml) /)<br>ncl 144> end do<br></div><div><br></div><div>NCL is zero based, so error in subscript #1 is referring to the 2nd index, line 141 shows that you are accessing the 2nd index of olr like this: (/ olr(:,ml) /)</div><div><br></div><div>ml at some point is bigger than the 2nd dimension size of olr. You will have to figure out why. Use print and printVarSummary statements to help to diagnose the issue. </div><div><br></div><div>Hope that helps!</div><div>Adam</div><div><br></div><div><br></div><div><br></div><div><br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Tue, Dec 17, 2019 at 1:40 AM Rahpeni Fajarianti via ncl-talk <<a href="mailto:ncl-talk@ucar.edu" rel="noreferrer noreferrer" target="_blank">ncl-talk@ucar.edu</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div>Hi all ! <br></div><div>I want to make MJO clivar 14, but I have some problem.<br></div><div>This is :</div><div><div>fatal:divide: Division by 0, Can't continue</div>fatal:Div: operator failed, can't continue<br>fatal:["Execute.c":8635]:Execute: Error occurred at or near line 129</div><div><br></div><div>This is my script:<br></div><div>ncl 0> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/diagnostics_cam.ncl"</div>ncl 1> ;******************************************************<br> u850 = dim_rmvmean_n(u850, 0)<br> u200 = dim_rmvmean_nncl 2> ;<br>(uncl 3> ; mjoclivar_14.ncl<br>ncl 4> ;<br>ncl 5> ;***********************************************************<br>ncl 6> ; Combined EOFs<br>hncl 7> ; Latest Update: July, 2016: Eun-Pa Lim; Bureau of Meteorology, Australia<br>ncl 8> ;***********************************************************<br>iancncl 9> ;;<br>nncl 10> ;; The following are automatically loaded from 6.2.0 onward<br>ncl 11> ;;load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl" <br>ncl 12> ;;load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl" <br>ncl 13> ;;load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl" <br>ncl 14> <br>ncl 15> undef("read_rename")<br>ncl 16> function read_rename(f[1]:file, varName[1]:string \<br>ncl 16> ,iStrt[1]:integer, iLast[1]:integer \<br>ncl 16> ,latS[1]:numeric , latN[1]:numeric )<br>ncl 17> ; Utility to force specific named dimensions<br>ncl 18> ; This is done for historical reasons (convenience) <br>ncl 19> begin<br>ncl 20> work = f->$varName$(iStrt:iLast,{latS:latN},:) ; (time,lat,lon)<br>ncl 21> work!0 = "time" ; CAM model names<br>ncl 22> work!1 = "lat"<br>ncl 23> work!2 = "lon"<br>ncl 24> return(work)<br>ncl 25> end<br>ncl 26> ; =========================> MAIN <==============================<br>ncl 27> begin<br>ncl 28> neof = 2<br>ncl 29> <br>ncl 30> latS = -15<br>ncl 31> latN = 15<br>ncl 32> <br> ncl 33> ymdStrt = 20180101 ; start yyyymmdd<br>ncl 34> ymdLast = 20181231 ; last <br>ncl 35> <br>ncl 36> yrStrt = ymdStrt/10000<br>ncl 37> yrLast = ymdLast/10000<br>ncl 38> <br>ncl 39> pltDir = "/home/peni/" ; plot directory<br>*ncl 40> pltType = "png" <br>ncl 41> pltName = "mjoclivar" <br>ncl 42> <br>,ncl 43> diri = "/home/peni/" ; input directory<br>ncl 44> <br>ncl 45> filolr = "<a href="http://olr.day.anomalies.2018.nc" rel="noreferrer noreferrer" target="_blank">olr.day.anomalies.2018.nc</a>"<br>ncl 46> filu200 = "<a href="http://uwnd.anomalies.nc" rel="noreferrer noreferrer" target="_blank">uwnd.anomalies.nc</a>"<br>ncl 47> filu850 = "<a href="http://uwnd2.anomalies.nc.nc" rel="noreferrer noreferrer" target="_blank">uwnd2.anomalies.nc.nc</a>"<br>ncl 48> <br>ncl 49> ;************************************************<br>ncl 50> ; create BandPass Filter<br>ncl 51> ;************************************************<br>ncl 52> ihp = 2 ; bpf=>band pass filter<br>ncl 53> nWgt = 201<br>ncl 54> sigma = 1.0 ; Lanczos sigma<br>ncl 55> fca = 1./100.<br>ncl 56> fcb = 1./20.<br>ncl 57> wgt = filwgts_lanczos (nWgt, ihp, fca, fcb, sigma )<br>ncl 58> <br>)ncl 59> ;***********************************************************<br>ncl 60> ; Find the indices corresponding to the start/end times<br>ncl 61> ;***********************************************************<br>ncl 62> f = addfile (diri+filolr , "r") <br>ncl 63> TIME = f->time ; days since ...<br>ncl 64> <br> ncl 65> YMD = cd_calendar(TIME, -2) ; entire (time,6)<br>ncl 66> <br>ncl 67> iStrt = ind(YMD.eq.ymdStrt) ; index start<br>ncl 68> iLast = ind(YMD.eq.ymdLast) ; index last <br>ncl 69> delete([/ TIME, YMD /])<br>ncl 70> <br>ncl 71> ;***********************************************************<br>ncl 72> ; Read anomalies<br>ncl 73> ;***********************************************************<br>ncl 74> <br>ncl 75> work = read_rename(f,"OLR_anom",iStrt,iLast,latS,latN) ; (time,lat,lon) <br>ncl 76> OLR = dim_avg_n_Wrap(work, 1) ; (time,lon)<br>ncl 77> delete(work)<br>ncl 78> <br>ncl 79> f = addfile (diri+filu850 , "r") <br>ncl 80> work = read_rename(f,"U_anom",iStrt,iLast,latS,latN) ; (time,lat,lon) <br>ncl 81> U850 = dim_avg_n_Wrap(work, 1) ; (time,lon)<br>ncl 82> delete(work)<br>ncl 83> <br>ncl 84> f = addfile (diri+filu200 , "r") <br>ncl 85> work = read_rename(f,"U_anom",iStrt,iLast,latS,latN) ; (time,lat,lon) <br>ncl 86> U200 = dim_avg_n_Wrap(work, 1) ; (time,lon)<br>ncl 87> <br>ncl 88> dimw = dimsizes( work )<br>ncl 89> ntim = dimw(0)<br>ncl 90> nlat = dimw(1)<br>ncl 91> mlon = dimw(2)<br>ncl 92> delete(work)<br>ncl 93> <br>ncl 94> lon = OLR&lon <br>ncl 95> time = OLR&time <br>ncl 96> date = cd_calendar(time, -2) ; yyyymmdd<br>ncl 97> <br>ncl 98> ;************************************************<br>ncl 99> ; Apply the band pass filter to the original anomalies<br>ncl 100> ;************************************************<br>ncl 101> olr = wgt_runave_n_Wrap ( OLR, wgt, 0, 0) ; (time,lon)<br>) <br>;************************************************<br> <br> imax_olr_eof1 = maxind(ceof(0,0,:)) <br> imax_olr_eof2 ncl 102> u850 = wgt_runave_n_Wrap (U850, wgt, 0, 0)<br>ncl 103> u200 = wgt_runave_n_Wrap (U200, wgt, 0, 0)<br>ncl 104> <br>ncl 105> ;************************************************<br>ncl 106> ; remove temporal means of band pass series: *not* necessary <br>ncl 107> ;************************************************<br>ncl 108> olr = dim_rmvmean_n( olr, 0) ; (time,lon)<br>ncl 109> u850 = dim_rmvmean_n(u850, 0)<br>ncl 110> u200 = dim_rmvmean_n(u200, 0)<br>ncl 111> <br>ncl 112> ;************************************************<br>ncl 113> ; Compute the temporal variance at each lon<br>ncl 114> ;************************************************<br>ncl 115> var_olr = dim_variance_n_Wrap( olr, 0) ; (lon)<br>ncl 116> var_u850 = dim_variance_n_Wrap(u850, 0)<br>ncl 117> var_u200 = dim_variance_n_Wrap(u200, 0)<br>ncl 118> <br>ncl 119> ;************************************************<br>ncl 120> ; Compute the zonal mean of the temporal variance<br> Oncl 121> ;************************************************<br>ncl 122> zavg_var_olr = dim_avg_n_Wrap( var_olr , 0) <br>ncl 123> zavg_var_u850 = dim_avg_n_Wrap( var_u850, 0)<br>ncl 124> zavg_var_u200 = dim_avg_n_Wrap( var_u200, 0)<br>ncl 125> <br>ncl 126> ;************************************************<br>ncl 127> ; Normalize by sqrt(avg_var*)<br>ncl 128> ;************************************************<br>ncl 129> olr = olr/sqrt(zavg_var_olr ) ; (time,lon)<br>ncl 130> u850 = u850/sqrt(zavg_var_u850)<br>ncl 131> u200 = u200/sqrt(zavg_var_u200)<br>ncl 132> <br>cncl 133> ;************************************************<br>ncl 134> ; Combine the normalized data into one variable<br>ncl 135> ;************************************************<br>ncl 136> cdata = new ( (/3*mlon,ntim/), typeof(olr), getFillValue(olr))<br>ncl 137> do ml=0,mlon-1<br>ncl 138> cdata(ml ,:) = (/ olr(:,ml) /)<br>ncl 139> cdata(ml+ mlon,:) = (/ u850(:,ml) /)<br>ncl 140> cdata(ml+2*mlon,:) = (/ u200(:,ml) /)<br>ncl 141> end do<br>ncl 142> <br>=ncl 143> ;************************************************<br>ncl 144> ; Compute **combined** EOF; Sign of EOF is arbitrary<br>ncl 145> ;************************************************<br>ncl 146> eof_cdata = eofunc(cdata , neof, False) ; (neof,3*mlon)<br>ncl 147> print("==============")<br>ncl 148> printVarSummary(eof_cdata)<br>ncl 149> printMinMax(eof_cdata, True)<br>ncl 150> <br>ncl 151> eof_ts_cdata = eofunc_ts(cdata,eof_cdata,False) ; (neof,3*ntim)<br>ncl 152> print("==============") <br>ncl 153> printVarSummary(eof_ts_cdata)<br>ncl 154> printMinMax(eof_ts_cdata, True)<br>ncl 155> <br>incl 156> ;************************************************<br>ncl 157> ; For clarity, explicitly extract each variable. Create time series <br>ncl 158> ;************************************************<br>ncl 159> <br>tncl 160> nvar = 3 ; "olr", "u850", "u200"<br>ncl 161> ceof = new( (/nvar,neof,mlon/), typeof(cdata), getFillValue(cdata))<br>ncl 162> <br>dncl 163> do n=0,neof-1<br>ncl 164> ceof(0,n,:) = eof_cdata(n,0:mlon-1) ; olr<br>ncl 165> ceof(1,n,:) = eof_cdata(n,mlon:2*mlon-1) ; u850<br>ncl 166> ceof(2,n,:) = eof_cdata(n,2*mlon:) ; u200<br>ncl 167> end do<br>ncl 168> <br>(ncl 169> ceof!0 = "var"<br>ncl 170> ceof!1 = "eof"<br>ncl 171> ceof!2 = "lon" <br>ncl 172> ceof&lon = olr&lon<br>ncl 173> <br>incl 174> ceof_ts = new( (/nvar,neof,ntim/), typeof(cdata), getFillValue(cdata))<br>ncl 175> ceof_ts(0,:,:) = eofunc_ts_Wrap( olr(lon|:,time|:),ceof(0,:,:),False) ; (0,neof,ntim)<br>ncl 176> ceof_ts(1,:,:) = eofunc_ts_Wrap(u850(lon|:,time|:),ceof(1,:,:),False) ; (1,neof,ntim)<br>ncl 177> ceof_ts(2,:,:) = eofunc_ts_Wrap(u200(lon|:,time|:),ceof(2,:,:),False) ; (2,neof,ntim)<br>ncl 178> <br>ncl 179> ;**********************************************t*<br>ncl 180> ; Add code contributed by Marcus N. Morgan, Florida Institute of Technology; Feb 2015<br>ncl 181> ; Calculate % variance (pcv_ )accounted for by OLR, U850 and U200<br>ncl 182> ;************************************************<br>ncl 183> <br>ncl 184> pcv_eof_olr = new(neof,typeof(ceof))<br>ncl 185> pcv_eof_u850 = new(neof,typeof(ceof))<br>ncl 186> pcv_eof_u200 = new(neof,typeof(ceof))<br>ncl 187> <br>ncl 188> do n=0,neof-1<br>ncl 189> pcv_eof_olr(n) = avg((ceof(0,n,:)*sqrt(ceof@eval(n)))^2)*100<br>ncl 190> pcv_eof_u850(n) = avg((ceof(1,n,:)*sqrt(ceof@eval(n)))^2)*100<br>ncl 191> pcv_eof_u200(n) = avg((ceof(2,n,:)*sqrt(ceof@eval(n)))^2)*100<br>ncl 192> ;;print("pcv: neof="+(n+1)+": "+pcv_eof_olr(n)+" "+pcv_eof_u850(n)+" "+pcv_eof_u200(n))<br>ncl 193> end do<br>ncl 194> <br> ncl 195> ;************************************************<br>ncl 196> ; Change sign of EOFs for spatial structures of PC1 and PC2 <br>ncl 197> ; to represent convection over the tropical Indian Ocean and the tropical western Pacific Ocean, respectively <br>ncl 198> ; (Ad hoc approach) <br>ncl 199> ;************************************************<br>ncl 200> <br>ncl 201> imax_olr_eof1 = maxind(ceof(0,0,:)) <br>ncl 202> imax_olr_eof2 = maxind(ceof(0,1,:)) <br>ncl 203> <br>)ncl 204> lonmax_eof1 = ceof&lon(imax_olr_eof1) ; longitude of max value (i.e. suppressed convection)<br>ncl 205> lonmax_eof2 = ceof&lon(imax_olr_eof2)<br>ncl 206> <br>sncl 207> if (lonmax_eof1.ge.60 .and. lonmax_eof1.lt.180) then ; Change the sign of EOF1 <br>ncl 208> ceof(:,0,:) = -ceof(:,0,:) ; if OLR is positive<br>ncl 209> ceof_ts(:,0,:) = -ceof_ts(:,0,:) ; over the tropical Indian Ocean<br>ncl 210> eof_cdata(0,:) = -eof_cdata(0,:)<br>ncl 211> eof_ts_cdata(0,:) = -eof_ts_cdata(0,:)<br>ncl 212> end if<br>ncl 213> <br>encl 214> if (lonmax_eof2.ge.120 .and. lonmax_eof2.lt.180) then ; Change the sign of EOF2<br>ncl 215> ceof(:,1,:) = -ceof(:,1,:) ; if OLR is positive <br>ncl 216> ceof_ts(:,1,:) = -ceof_ts(:,1,:) ; over the tropical western Pacific Ocean<br>ncl 217> eof_cdata(1,:) = -eof_cdata(1,:)<br>ncl 218> eof_ts_cdata(1,:) = -eof_ts_cdata(1,:)<br>ncl 219> end if<br>ncl 220> <br>ncl 221> print("==============")<br>rncl 222> printVarSummary(eof_cdata)<br>ncl 223> printMinMax(eof_cdata, True)<br>ncl 224> <br>ncl 225> ;************************************************<br>ncl 226> ; Compute cross correlation of each variable's EOF time series at zero-lag<br>ncl 227> ;************************************************<br>ncl 228> r_olr_u850 = escorc(ceof_ts(0,:,:) , ceof_ts(1,:,:) ) ; (neof)<br>ncl 229> r_olr_u200 = escorc(ceof_ts(0,:,:) , ceof_ts(2,:,:) )<br>ncl 230> r_u850_u200 = escorc(ceof_ts(1,:,:) , ceof_ts(2,:,:) )<br>ncl 231> <br>sncl 232> print("==============")<br>ncl 233> do n=0,neof-1<br> ncl 234> print("neof="+n \<br>ncl 234> +" r_olr_u850=" +sprintf("%4.3f",r_olr_u850(n)) \<br>ncl 234> +" r_olr_u200=" +sprintf("%4.3f",r_olr_u200(n)) \<br>ncl 234> +" r_u850_u200="+sprintf("%4.3f",r_u850_u200(n)) )<br>ncl 235> end do<br>ncl 236> print("==============")<br>ncl 237> <br>ncl 238> ;************************************************<br>5N: "+yrStrt+"-"+yrLast <br> <br> do ncl 239> ; Compute cross correlation of the multivariate EOF; EOF 1 vs EOF 2<br>ncl 240> ;************************************************<br>ncl 241> <br> ncl 242> mxlag = 25<br>ncl 243> rlag_01 = esccr(eof_ts_cdata(0,:),eof_ts_cdata(1,:), mxlag) ; (N,mxlag+1)<br>ncl 244> rlag_10 = esccr(eof_ts_cdata(1,:),eof_ts_cdata(0,:), mxlag) ; (N,mxlag+1)<br>ncl 245> ccr_12 = new ( (/2*mxlag+1/), float) <br>ncl 246> <br>cncl 247> ccr_12(mxlag:) = rlag_10(0:mxlag) <br>ncl 248> ccr_12(0:mxlag) = rlag_01(::-1) ; reverse order<br>ncl 249> ;;print(ccr_12)<br>ncl 250> <br>sncl 251> <br>Pncl 252> ;************************************************<br>ncl 253> ; Normalize the multivariate EOF 1&2 component time series<br>ncl 254> ; Compute (PC1^2+PC2^2): values > 1 indicate "strong" periods<br>ncl 255> ;************************************************<br>ncl 256> eof_ts_cdata(0,:) = eof_ts_cdata(0,:)/stddev(eof_ts_cdata(0,:))<br>ncl 257> eof_ts_cdata(1,:) = eof_ts_cdata(1,:)/stddev(eof_ts_cdata(1,:))<br>ncl 258> <br>ncl 259> mjo_ts_index = eof_ts_cdata(0,:)^2 + eof_ts_cdata(1,:)^2 <br>ncl 260> mjo_ts_index_smt = runave(mjo_ts_index, 91, 0) ; 91-day running mean<br>ncl 261> <br> ncl 262> nGood = num(.not.ismissing(mjo_ts_index)) ; # non-missing<br>ncl 263> nStrong = num(mjo_ts_index .ge. 1.0)<br>ncl 264> print("nGood="+nGood+" nStrong="+nStrong+" nOther="+(nGood-nStrong))<br>ncl 265> <br>rncl 266> ;************************************************<br>ncl 267> ; Write PC results to netCDF for use in another example.<br>ncl 268> ;************************************************<br>ncl 269> mjo_ts_index!0 = "time"<br>ncl 270> mjo_ts_index&time = time <br>ncl 271> mjo_ts_index@long_name = "MJO PC INDEX" <br>ncl 272> mjo_ts_index@info = "(PC1^2 + PC2^2)" <br>ncl 273> <br> ncl 274> PC1 = eof_ts_cdata(0,:)<br>ncl 275> PC1!0 = "time"<br>ncl 276> PC1&time = time<br>ncl 277> PC1@long_name = "PC1"<br>ncl 278> PC1@info = "PC1/stddev(PC1)"<br>ncl 279> <br>oncl 280> PC2 = eof_ts_cdata(1,:)<br>ncl 281> PC2!0 = "time"<br>ncl 282> PC2&time = time<br>ncl 283> PC2@long_name = "PC2"<br>ncl 284> PC2@info = "PC2/stddev(PC2)"<br>ncl 285> <br>ncl 286> diro = "./"<br>ncl 287> filo = "MJO_PC_INDEX.nc"<br>ncl 288> system("/bin/rm -f "+diro+filo) ; remove any pre-existing file<br>ncl 289> ncdf = addfile(diro+filo,"c") ; open output netCDF file<br>of legend.<br> ncl 290> ; make time an UNLIMITED dimension <br>ncl 291> filedimdef(ncdf,"time",-1,True) ; recommended for most applications<br>ncl 292> ; output variables directly<br>ncl 293> ncdf->MJO_INDEX = mjo_ts_index <br>ncl 294> ncdf->PC1 = PC1 <br>ncl 295> ncdf->PC2 = PC2 <br>ncl 296> <br>ncl 297> ;------------------------------------------------------------<br>ncl 298> ; PLOTS<br>ncl 299> ;------------------------------------------------------------<br>ncl 300> <br>ncl 301> yyyymmdd = cd_calendar(time, -2)<br>tncl 302> yrfrac = yyyymmdd_to_yyyyfrac(yyyymmdd, 0.0)<br>ncl 303> delete([/ yrfrac@long_name, lon@long_name /])<br>ncl 304> <br>ncl 305> day = ispan(-mxlag, mxlag, 1)<br>ncl 306> ;day@long_name = "lag (day)"<br>ncl 307> <br>ncl 308> pltPath = pltDir+pltName<br>ncl 309> <br>ncl 310> wks = gsn_open_wks(pltType,pltPath)<br>ncl 311> gsn_define_colormap(wks,"default")<br>ncl 312> plot = new(3,graphic) <br>ncl 313> <br>ncl 314> ;************************************************<br>ncl 315> ; Multivariate EOF plots<br>ncl 316> ;************************************************<br>ncl 317> rts = True<br>ncl 318> rts@gsnDraw = False ; don't draw yet<br>ncl 319> rts@gsnFrame = False ; don't advance frame yet<br>ncl 320> rts@gsnScale = True ; force text scaling <br>ncl 321> <br>ncl 322> rts@vpHeightF = 0.40 ; Changes the aspect ratio<br>ncl 323> rts@vpWidthF = 0.85<br>ncl 324> rts@vpXF = 0.10 ; change start locations<br>ncl 325> rts@vpYF = 0.75 ; the plot<br>ncl 326> rts@xyLineThicknesses = (/2, 2, 2/)<br>ncl 327> rts@xyLineColors = (/"black","red","blue"/)<br>ncl 328> rts@gsnYRefLine = 0. ; reference line <br>ncl 329> rts@trXMaxF = max(lon)<br>ncl 330> rts@trXMinF = min(lon)<br>ncl 331> <br>ncl 332> rts@pmLegendDisplayMode = "Always" ; turn on legend<br>ncl 333> rts@pmLegendSide = "Top" ; Change location of <br>ncl 334> rts@pmLegendParallelPosF = 1.16 ; move units right<br>ncl 335> rts@pmLegendOrthogonalPosF = -0.50 ; move units down<br>ncl 336> rts@pmLegendWidthF = 0.15 ; Change width and<br>ncl 337> rts@pmLegendHeightF = 0.15 ; height of legend.<br>ncl 338> rts@lgLabelFontHeightF = 0.0175<br>ncl 339> <br>ncl 340> <br>ncl 341> rtsP = True ; modify the panel plot<br>ncl 342> ; rtsP@gsnMaximize = True ; large format<br>ncl 343> rtsP@gsnPanelMainString = "Multivariate EOF: 15S-15N: "+yrStrt+"-"+yrLast <br>ncl 344> <br>ncl 345> do n=0,neof-1<br>ncl 346> rts@xyExplicitLegendLabels = (/"OLR: "+sprintf("%4.1f", pcv_eof_u200(n)) +"%" \<br>ncl 346> ,"U850: "+sprintf("%4.1f", pcv_eof_u850(n))+"%" \<br>ncl 346> ,"U200: "+sprintf("%4.1f", pcv_eof_olr(n))+"%" /)<br>ncl 347> rts@gsnLeftString = "EOF "+(n+1)<br>ncl 348> rts@gsnRightString = sprintf("%3.1f",ceof@pcvar(n)) +"%"<br>ncl 349> plot(n) = gsn_csm_xy (wks,lon,ceof(:,n,:),rts)<br>ncl 350> end do<br>ncl 351> gsn_panel(wks,plot(0:1),(/2,1/),rtsP) ; now draw as one plot<br>ncl 352> <br>ncl 353> ;-----------------------------------------<br>ncl 354> ; The following doesn't work with some older versions of NCL<br>ncl 355> ; With old versions, the user must delete each individually.<br>ncl 356> ;----------------------------------------- <br>ncl 357> delete([/ rts@xyExplicitLegendLabels, rts@pmLegendDisplayMode \<br>ncl 357> , rts@xyLineThicknesses , rts@gsnLeftString \<br>ncl 357> , rts@gsnRightString , rts@xyLineColors \<br>ncl 357> , rts@trXMaxF , rts@trXMinF /] )<br>ncl 358> <br>ncl 359> lag = ispan(-mxlag,mxlag,1)<br>ncl 360> lag@long_name = "lag (days)"<br>ncl 361> <br>ncl 362> plot(0) = gsn_csm_xy (wks, lag ,ccr_12,rts)<br>ncl 363> rtsP@gsnPanelMainString = "Cross Correlation: Multivariate EOF: 15S-15N: " \<br>ncl 363> + yrStrt+"-"+yrLast <br>ncl 364> rtsP@gsnPaperOrientation = "portrait" ; force portrait<br>ncl 365> gsn_panel(wks,plot(0),(/1,1/),rtsP) ; now draw as one plot<br>ncl 366> <br>ncl 367> ;************************************************<br>ncl 368> ; MJO "strong" index <br>ncl 369> ;************************************************<br>ncl 370> rts@gsnYRefLine = 1.0<br>ncl 371> rts@gsnYRefLineColor = "black"<br>ncl 372> rts@xyMonoDashPattern = True<br>ncl 373> rts@xyLineColors = (/"black", "blue"/)<br>ncl 374> rts@xyLineThicknesses = (/1, 2/)<br>ncl 375> rts@pmLegendDisplayMode = "Always" ; turn on legend<br>ncl 376> rts@pmLegendWidthF = 0.12 ; Change width and<br>ncl 377> rts@pmLegendHeightF = 0.10 ; height of legend.<br>ncl 378> rts@pmLegendParallelPosF = 0.86 ; move units right<br>ncl 379> rts@pmLegendOrthogonalPosF = -0.40 ; move units down<br>ncl 380> rts@xyExplicitLegendLabels = (/"daily", "91-day runavg" /)<br>ncl 381> <br>ncl 382> mjo_ind_plt = new ( (/2,ntim/), typeof(mjo_ts_index))<br>ncl 383> mjo_ind_plt(0,:) = mjo_ts_index<br>ncl 384> mjo_ind_plt(1,:) = (/ mjo_ts_index_smt /)<br>ncl 385> plot(0) = gsn_csm_xy(wks, yrfrac,mjo_ind_plt,rts)<br>ncl 386> <br>ncl 387> rtsP@gsnPanelMainString = "MJO Index: (PC1^2+ PC2^2) : 15S-15N: "+yrStrt+"-"+yrLast <br>ncl 388> gsn_panel(wks,plot(0),(/1,1/),rtsP) ; now draw as one plot<br>ncl 389> <br>ncl 390> end<br><div><br></div><div>Thank you<br></div><br></div>
_______________________________________________<br>
ncl-talk mailing list<br>
<a href="mailto:ncl-talk@ucar.edu" rel="noreferrer noreferrer" target="_blank">ncl-talk@ucar.edu</a><br>
List instructions, subscriber options, unsubscribe:<br>
<a href="http://mailman.ucar.edu/mailman/listinfo/ncl-talk" rel="noreferrer noreferrer noreferrer" target="_blank">http://mailman.ucar.edu/mailman/listinfo/ncl-talk</a></blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div><div><div><span><font color="#888888">Adam Phillips <br></font></span></div><span><font color="#888888">Associate Scientist, </font></span><span><font color="#888888">Climate and Global Dynamics Laboratory, NCAR<br></font></span></div></div><div><span><font color="#888888"><a href="http://www.cgd.ucar.edu/staff/asphilli/" rel="noreferrer noreferrer" target="_blank">www.cgd.ucar.edu/staff/asphilli/</a> </font></span><span><font color="#888888">303-497-1726 </font></span></div><span><font color="#888888"></font></span><div><div><span><font color="#888888"><br></font></span><div><span><font color="#888888"><a href="http://www.cgd.ucar.edu/staff/asphilli" rel="noreferrer noreferrer" target="_blank"></a></font></span></div></div></div></div></div></div></div></div></div></div></div>
</blockquote></div></div>
_______________________________________________<br>
ncl-talk mailing list<br>
<a href="mailto:ncl-talk@ucar.edu" rel="noreferrer noreferrer" target="_blank">ncl-talk@ucar.edu</a><br>
List instructions, subscriber options, unsubscribe:<br>
<a href="http://mailman.ucar.edu/mailman/listinfo/ncl-talk" rel="noreferrer noreferrer noreferrer" target="_blank">http://mailman.ucar.edu/mailman/listinfo/ncl-talk</a></blockquote></div>
</blockquote></div>
</blockquote></div>
</div>
</div>
</blockquote></div>
_______________________________________________<br>
ncl-talk mailing list<br>
<a href="mailto:ncl-talk@ucar.edu" target="_blank">ncl-talk@ucar.edu</a><br>
List instructions, subscriber options, unsubscribe:<br>
<a href="http://mailman.ucar.edu/mailman/listinfo/ncl-talk" rel="noreferrer" target="_blank">http://mailman.ucar.edu/mailman/listinfo/ncl-talk</a></blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr"><div dir="ltr"><div><div dir="ltr">Barry H. Lynn, Ph.D<div><div>Senior Associate Scientist, Lecturer,</div><div><div><span style="color:rgb(136,136,136)">The Institute of the Earth Science, </span><br style="color:rgb(136,136,136)"><span style="color:rgb(136,136,136)">The Hebrew University of Jerusalem, </span><br style="color:rgb(136,136,136)"><span style="color:rgb(136,136,136)">Givat Ram, Jerusalem 91904, Israel </span><br style="color:rgb(136,136,136)"></div><span style="color:rgb(136,136,136)">Tel: 972 547 231 170</span><br style="color:rgb(136,136,136)"><span style="color:rgb(136,136,136)">Fax: (972)-25662581</span></div></div><div><span style="color:rgb(136,136,136)"><br></span></div><div>C.E.O, Weather It Is, LTD<br>Weather and Climate Focus<br><a href="http://weather-it-is.com" target="_blank">http://weather-it-is.com</a><br>Jerusalem, Israel<br>Local: 02 930 9525<br>Cell: 054 7 231 170<br>Int-IS: x972 2 930 9525<br><br></div></div></div></div></div>
</blockquote></div>