[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

Rahpeni Fajarianti rahpenifajarianti at gmail.com
Fri Dec 20 17:22:32 MST 2019


Thank you.

Now, I use long term means for daily values for 1981-2010 as base
climatology data. It works.
I got min max values for :
- OLR anomaly (min=-185.5106410148043   max=77.65257890949087)
- 200hPa zonal wind anomaly  (min=-48.11269687686956
max=57.21568789939967)
- 850hPa zonal wind anomaly (min=-31.22578327383305
max=41.24650545117836).

then, I use mjoclivar_14.ncl to get multivariate eof for data 2018. I
attach my script.

Some error shown like this:
sub drveof: ier,jopt=  10  0 returned from vcmssm/crmssm
warning:eofunc: 10 eigenvectors failed to converge
(0) ==============
Variable: eof_cdata
Type: double
Total Size: 6960 bytes
            870 values
Number of Dimensions: 2
Dimensions and sizes: [2] x [435]
Coordinates:
Number Of Attributes: 5
  _FillValue : -999000000
  method : no transpose
  matrix : covariance
  pcvar : ( -9.99e+08, -9.99e+08 )
  eval : ( -999000000, -999000000 )
(0)
(0) min=-999000000   max=-999000000
warning:eofunc_ts: 2 eigenvectors failed to converge
(0) ==============

Variable: eof_ts_cdata
Type: double
Total Size: 5840 bytes
            730 values
Number of Dimensions: 2
Dimensions and sizes: [2] x [365]
Coordinates:
Number Of Attributes: 3
  _FillValue : -999000000
  matrix : covariance
  ts_mean : ( -999000000, -999000000 )
(0)
(0) min=-999000000   max=-999000000
warning:eofunc_ts: 2 eigenvectors failed to converge
warning:eofunc_ts: 2 eigenvectors failed to converge
warning:eofunc_ts: 2 eigenvectors failed to converge
warning:avg: Entire input array contains missing values, can't compute
average

fatal:Subscript out of range, error in subscript #0
fatal:An error occurred reading unknown
fatal:["Execute.c":8635]:Execute: Error occurred at or near line 207

this is the error line :
ncl 207>   lonmax_eof1 = ceof&lon(imax_olr_eof1)      ; longitude of max
value (i.e. suppressed convection)
ncl 208>   lonmax_eof2 = ceof&lon(imax_olr_eof2)

The data I used:
 olranom.nc
<https://drive.google.com/file/d/1xbmRiRO9VrTQ6qqKTfTq_rraiRSZFOXg/view?usp=drive_web>
 uanom.200.nc
<https://drive.google.com/file/d/1u1IYOGmjLuWAR4izwLkAgEMzWDxrkZlX/view?usp=drive_web>
 uanom.850.nc
<https://drive.google.com/file/d/1ZTtEukGralmUuwBWLGJ9ei6-N0yf199b/view?usp=drive_web>
The script I used:

On Sat, Dec 21, 2019 at 2:21 AM Dennis Shea <shea at ucar.edu> wrote:

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


More information about the ncl-talk mailing list