<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">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" target="_blank">olr.day.anomalies.2018.nc</a> , <a href="http://uwnd.anomalies.nc" target="_blank">uwnd.anomalies.nc</a> (for 200hPa u wind anomalies) , <a href="http://uwnd2.anomalies.nc" 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" 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" 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" 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" 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" 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/" 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" 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" 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" 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" 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" target="_blank">olr.day.anomalies.2018.nc</a>"<br>ncl 46>    filu200 = "<a href="http://uwnd.anomalies.nc" target="_blank">uwnd.anomalies.nc</a>"<br>ncl 47>    filu850 = "<a href="http://uwnd2.anomalies.nc.nc" 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" 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"><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/" 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" 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" 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>
</blockquote></div>
</blockquote></div>