; Concepts illustrated:
;   - Calculating a two-dimensional correlation in time
;   - Reordering an array
;   - Copying attributes from one variable to another
;   - Copying coordinate arrays from one variable to another
;   - Removing the end boxes in a labelbar

begin

  in1  = addfile("b003_TS_200-299.nc","r")
  in2  = addfile("b006_TS_035-134.nc","r")
  
  tmp1 = in1->TS           ;  Surface temperature (K)  [time | 100] x [lat | 64] x [lon | 128]
  tmp2 = in2->TS           ;  Surface temperature (K)  [time | 100] x [lat | 64] x [lon | 128]
  ;printVarSummary(tmp1)
  ;printVarSummary(tmp2)


; reorder to get time as right most dimension

  ts1 = tmp1(lat|:,lon|:,time|:)       ;[lat | 64] x [lon | 128] x [time | 100]
  ts2 = tmp2(lat|:,lon|:,time|:)       ; [lat | 64] x [lon | 128] x [time | 100]
  ;printVarSummary(ts1)
  ;printVarSummary(ts2)
  

; calculate cross correlations

  maxlag = 2                             ; note: the max lag should not be more than N/4

  ccr = esccr(ts1,ts2,maxlag)            ; calc cross correlations
  copy_VarAtts(ts1, ccr)                    
  copy_VarCoords_1(ts2,ccr)              ; copy meta data and coordinate variables 
  
  printVarSummary(ccr)
  printMinMax(ccr,0)

; plot the correlations

  wks  = gsn_open_wks("png","corel")             ; send graphics to PNG file
 
  res                      = True                ; make plot mods

  res@cnFillOn             = True                ; turn on color
  res@cnFillPalette        = "BlWhRe"            ; set color map

  res@cnLinesOn            = False               ; turn off contour lines
 
  
  res@cnLevelSelectionMode = "ManualLevels"      ; manually set cn levels
  res@cnMinLevelValF       = -1.                 ; min level
  res@cnMaxLevelValF       =  1.                 ; max level
  res@cnLevelSpacingF      = .1                  ; contour level spacing

  res@cnLabelBarEndStyle   = "ExcludeOuterBoxes" ; remove the two end boxes from the labelbar
  res@pmLabelBarWidthF     = 0.8

  lag                      = 0
  res@tiMainString         = "Correlations at lag "+lag
  plot = gsn_csm_contour_map(wks,ccr(:,:,lag),res) 

end