PLOT = True

  if (PLOT) then
    pltdir       = "./"
    pltname      = "O3"
    plttype      = "x11"      ; send graphics to X11 or PNG file

  end if

;--------------------------------------------
; end user input
;--------------------------------------------
;
;--------------------------------------------
; set up
;--------------------------------------------

;-----read the data and define

; NCL script
;Najib Yusuf 19 May, 2022
;;**************************************

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" ; BRH

begin

; read model data

   spec_units = ("ppb")
   factor = (/1.e+09/)	
   wkdir = "/glade/scratch/najiby/"
   diri = "/glade/scratch/najiby/"
   fili = systemfunc("ls "+diri+"f.e20.FSCD_FullEmissNig.cam.h0.2013*.nc")
   f = addfiles(fili, "r")
  o3x = f[:]->O3
       o3x = o3x*1.e+09
  o3x_surf = o3x(0,55,:,:)
  xMon = f[:]->date
  
printVarSummary(o3x_surf)
printMinMax(o3x_surf,0)

printVarSummary(xMon)
print(xMon)

;-------------------------------------------------------
; seasons average
  ;xDJF = month_to_season(o3x_surf(0:N-1),:,:),"DJF)
  ;xDJF = month_to_seasonN (xMon, (/"DJF","MAM","JJA","SON"/))
  xDJF = month_to_seasonN(o3x,(/"DJF","MAM","JJA","SON"/))
   
   o3x_surf1 = "DJF"
   o3x_surf2 = "MAM"
   o3x_surf3 = "JJA"
   o3x_surf4 = "SON"
 
;define plot
 wks = gsn_open_wks("x11", "seasons_FullEmissNig")

; define plot generals
;********************************************************
  res                 = True
  gsn_define_colormap(wks,"WhiteBlueGreenYellowRed")
  res@gsnDraw                = False
  res@gsnFrame               = False
  res@gsnMaximize = False
  res@cnFillOn             = True                  ; turn on color fill
  res@cnLinesOn         = False                 ; no contour lines
  res@cnLineLabelsOn    = False
  res@tmXTOn    = False
  res@tmXBOn = False
  res@tmYLLabelFontHeightF = 0.015
  res@tmYRMode             = "Automatic"          ; turn off special labels on right axis
  res@lbLabelBarOn      = True
  res@lbOrientation     = "Vertical"
  res@cnInfoLabelOn        = False
  res@cnLevelSelectionMode      = "ExplicitLevels"       ; manual levels
  res@cnLevelSelectionMode  = "ManualLevels"   ;-- set contour levels manually
  res@cnMinLevelValF        =0.2            ;-- minimum contour level
  res@cnMaxLevelValF        =  160
  res@cnMissingValFillColor = 0
  res@txFontHeightF = 0.020
  res@txFont = "helvetica-bold"
  res@gsnPaperOrientation = "landscape"
  res@pmLegendWidthF = 0.15
  res@pmLegendHeightF = 0.15
  res@lgLabelFontHeightF = .015
  res@lgPerimOn = True
  res@txFontHeightF = 0.015
  res@mpDataSetName         = "Earth..4"   ; This new database contains
                                           ; divisions for other countries.
  res@mpDataBaseVersion     = "MediumRes"  ; Medium resolution database

  res@mpOutlineOn           = True         ; Turn on map outlines
  res@mpOutlineSpecifiers   = (/"Africa:states"/)
   res@mpLimitMode           = "LatLon"
  ;res@mpMaxLatF             =  15          ;  Africa limits
  ;res@mpMinLatF             =  0
  ;res@mpMaxLonF             = 17
  ;res@mpMinLonF             = 0
  res@mpMaxLatF             =  35          ;  Africa limits
  res@mpMinLatF             = -35
  res@mpMaxLonF             = 60
  res@mpMinLonF             = -30

;*********************************************************************
   res@cnLevelSelectionMode = "ExplicitLevels"
   res@cnFillOn             = True
   res@cnLinesOn            = False
   res@cnLineLabelsOn       = False
   res@lbTitleOn = True
   res@lbTitleString = spec_units
   res@lbLabelFontHeightF = 0.018
   res@lbTitleFontHeightF = 0.02
   res@lbLabelStride  = 2
   res@lbOrientation = "Horizontal"
; contour levels
 con = new(14,float)
 dcon = new (14,float)
 fa = (/5/)
 con=fa*(/0.5, 1.,2., 3.,4.,5.,6.,7.,8.,10.,15.,20.,25.,30./)

;**********************************************************************

; res@txString =month(mo)
  res@cnLevels = con(:)
; res@cnFillColors   = (/3,4,5,6,7,8,9,10,11,12,13,14,15,16/)
  res@cnLevels = con(:)
  res@gsnStringFontHeightF = 0.018
  
; if you want to plot the surface you do this

   res@gsnLeftString = "DJF_dry cold"
  plot1 = gsn_csm_contour_map_ce(wks,o3x_surf1({-35:35},{-30:60}),res)
  
   res@gsnLeftString = "MAM_dry warm"
  plot2 = gsn_csm_contour_map_ce(wks,o3x_surf2(:,:),res)
  
   res@gsnLeftString = "JJA_monsoon"
  plot3 = gsn_csm_contour_map_ce(wks,o3x_surf3(:,:),res)
  
   res@gsnLeftString = "SON_post monsoon"
  plot4 = gsn_csm_contour_map_ce(wks,o3x_surf4(:,:),res)
  
 ;--Panel all three plots in a m x n configuration
pres = True

;*********************************************
; add lat long info to plot
    txres               = True
    txres@txFontHeightF = 0.017

    panres = True
    panres@gsnFrame = False
    panres@gsnMaximize = True
    panres@gsnPanelTop = 0.96
     panres@gsnPanelBottom = 0.05
     gsn_panel(wks, (/plot1,plot2,plot3,plot4/),(/2,2/),pres)
    frame(wks)
;-------------------------------------------------------------------------
 delete(res)

end