;*************************************************
; corel_2.ncl
;************************************************
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"  
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
begin
;************************************************
; variable and file handling
;************************************************
   in  = addfile("b003_TS_200-299.nc","r")
   ts = in->TS
   ts@_FillValue = getFillValue(ts)              ; assign _FillValue
;************************************************
; extract time series from 3d data
;************************************************
   ts1=ts(:,45,64)                               ; extract time series
   ts2=ts(:,23,117)
;************************************************
; calculate cross correlations
;************************************************
   maxlag=25                                     ; set lag
   totlag=2*maxlag                               ; set total lag
   ccr1 = esccr(ts1,ts2,maxlag)                  ; calc pos lag cross cor
   ccr2 = esccr(ts2,ts1,maxlag)                  ; calc neg lag cross cor
;************************************************
;  combine pos and neg into one time series
;************************************************
   ccrtot = new( (/totlag/),float)               ; allocate memory 
   ccrtot(0:maxlag-1) = ccr2(0:maxlag-1:-1)      
   ccrtot(maxlag:)    = ccr1(0:maxlag-1)
 
   x = ispan(-maxlag,maxlag,1)                   ; define x axis
;************************************************
; plot the correlations
;************************************************
   wks  = gsn_open_wks("ps","corel")                  ; open a ps plot
   res                 = True                         ; make plot mods
   res@tiMainString    = "37.7N 180E vs 23.72S 149W"  ; title
   res@tiXAxisString   = "LAG"                        ; x-axis label

   plot = gsn_xy(wks,x,ccrtot,res)                       ; plot correlation
;************************************************
end