;----------------------------------------------------------------------
; panel_33.ncl
;
; Concepts illustrated:
;   - Combining two sets of paneled plots on one page
;   - Maximizing paneled plots after they've been created
;   - Drawing two labelbars in a combined panel plot
;   - Using lbBoxEndCapStyle to draw triangles at the end of a labelbar
;----------------------------------------------------------------------
;
; 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
;---Open file
    fname = "TS.cam3.toga_ENS.1950-2000.nc"
    fname = "panel_data.nc"
    f     = addfile(fname,"r")

;---Convert "time" to an ntim x 6 array of year,mon,day,hour,min,sec
    newtime = cd_calendar(f->time,0)
    years   = newtime(:,0)
    months  = newtime(:,1)
    
;---Indicate start and end years of interest. Best to only do 1 or 2 here
    month      = 1                  ; January
    smonth     = "Jan"
    start_years = (/1951,1961/)
    mid_years   = (/1971,1981/)
    end_years   = (/1991,2001/)
    nyears_cols = dimsizes(start_years)
    nyears_rows = 3


;---Read in temperature and convert to degrees C.
    t       = f->TS                 
    t       = t - 273.15
    t@units = "degC"

;---Start the graphics
    wks  = gsn_open_wks("png","panel")         ; send graphics to PNG file
    
    res                     = True
    res@gsnDraw             = False
    res@gsnFrame            = False
    res@mpFillOn            = False            ; no need   

    res@cnLevelSelectionMode= "ManualLevels"   ; manual set levels
    res@cnMinLevelValF      = -3.0
    res@cnMaxLevelValF      = 27.0
    res@cnLevelSpacingF     = 1.5              ; 20 contour levels        
    res@cnFillOn            = True             ; color fill plot
    res@cnLinesOn           = False
    res@cnLineLabelsOn      = False
    res@cnInfoLabelOn       = False
    res@cnFillPalette       = "BlAqGrYeOrRe"

    res@lbLabelBarOn        = False            ; turn off individual label bars

    res@gsnStringFontHeightF= 0.02
    res@gsnLeftString       = "TS"
    res@gsnRightString      = t@units

;---Resources for diff plots
    dres                    = res    
    dres@cnMinLevelValF     = -4.
    dres@cnMaxLevelValF     =  4.
    dres@cnLevelSpacingF    =  1.
    dres@cnFillPalette      = "temp_diff_18lev"

;---Create arrays to hold series of plots
    plots      = new(nyears_rows*nyears_cols,graphic)    
    diff_plots = new((nyears_rows-1)*nyears_cols,graphic)    

    do i=0,nyears_cols-1
;---Get data for start and end year of interest
      sy = ind(years.eq.start_years(i).and.months.eq.month)
      my = ind(years.eq.mid_years(i)  .and.months.eq.month)
      ey = ind(years.eq.end_years(i)  .and.months.eq.month)

      diff_mid = t(sy,:,:)                ; trick to copy  metadata
      diff_mid = t(my,:,:) - t(sy,:,:)    ; overwrite with diff values
      diff_end = t(my,:,:)                ; trick to copy  metadata
      diff_end = t(ey,:,:) - t(my,:,:)    ; overwrite with diff values

      diff_mid@long_name = mid_years(i) + "-" + start_years(i) + \
                       " TS field differences"
      diff_end@long_name = end_years(i) + "-" + mid_years(i) + \
                       " TS field differences"
;---Debug prints
      print("========================================")
      printMinMax(t(sy,:,:),0)
      printMinMax(t(my,:,:),0)
      printMinMax(t(ey,:,:),0)
      printMinMax(diff_mid,0)
      printMinMax(diff_end,0)
    
;---Create the start/mid/end year plots and the difference plots
      res@gsnCenterString = smonth + " " + start_years(i)
      plots(i) = gsn_csm_contour_map_ce(wks,t(sy,:,:),res)

      res@gsnCenterString = smonth + " " + mid_years(i)
      plots(i+nyears_cols) = gsn_csm_contour_map_ce(wks,t(my,:,:),res)

      res@gsnCenterString = smonth + " " + end_years(i)
      plots(i+2*nyears_cols) = gsn_csm_contour_map_ce(wks,t(ey,:,:),res)

      dres@gsnCenterString = "Difference: " + smonth + " " + \
                             mid_years(i) + "-" + start_years(i)
      diff_plots(i) = gsn_csm_contour_map_ce(wks, diff_mid ,dres)
      dres@gsnCenterString = "Difference: " + smonth + " " + \
                             end_years(i) + "-" + mid_years(i)
      diff_plots(i+nyears_cols) = gsn_csm_contour_map_ce(wks, diff_end ,dres)
    end do

    total_rows = tofloat(nyears_rows + (nyears_rows-1))

;---Resources for paneling
    pres                   = True              ; modify the panel plot
    pres@gsnPanelMainString= fname
    pres@gsnFrame          = False             ; don't advance frame yet
    pres@gsnDraw           = False             ; will draw later, in maximize mode
    pres@gsnPanelBottom    = 1.-(nyears_rows/total_rows)

    pres@gsnPanelLabelBar  = True              ; add common colorbar
    pres@lbOrientation     = "vertical"
    pres@pmLabelBarWidthF  = 0.075             ; make thinner
    pres@pmLabelBarHeightF = (1.-pres@gsnPanelBottom) * 0.8

;
; It is necessary to return the id of the paneled plots,
; so they can "live" for the rest of this script, for
; when we maximize them later.
; 
    panelid1 = gsn_panel_return(wks,plots,(/nyears_rows,nyears_cols/),pres) 

;
; Grab the top and bottom locations of the all but the top plots.
; The height area of these plots will be the height area of the
; difference plots. And the bottom of these plots will be the top
; of the difference plots.
;
    bb     = NhlGetBB(plots(nyears_cols:))
    bottom = min(bb(:,1))
    top    = max(bb(:,0))
    height = (top-bottom)

;---Panel the difference plots at the bottom
    delete(pres@gsnPanelMainString)
    pres@gsnPanelTop       := bottom - 0.01   ; the 0.01 gives a little spaces b/w the two sets of plots
    pres@gsnPanelBottom    := pres@gsnPanelTop - height
    pres@lbBoxEndCapStyle  = "TriangleBothEnds"
    pres@pmLabelBarHeightF = (pres@gsnPanelTop - pres@gsnPanelBottom) * 0.8

    panelid2 = gsn_panel_return(wks,diff_plots,(/nyears_rows-1,nyears_cols/),pres)

;    drawNDCGrid(wks)
;    draw(wks)
;    frame(wks)

;---This will maximize the size of all the paneled stuff.
    maximize_output(wks,True)
end