begin
   
    ;Open input file.
    ;************************
     DATADir = "../DATA/"
     FILES = systemfunc (" ls -1 " + DATADir + "wrfout_d01* ")
     numFILES = dimsizes(FILES)

     a = addfiles(DATADir+FILES+".nc","r")

     T = wrf_user_getvar(a,"T2",-1)
     T = T-273.15
     wrf_smooth_2d(T, 3)             ; smooth T

     lat    = wrf_user_getvar(a,"XLAT",0)
     lon    = wrf_user_getvar(a,"XLONG",0)

     lat@units= "degrees_north"
     lon@units= "degrees_east"

     xMaxTime = dim_max_n_Wrap(T, 0)

     ;Start the graphics
     ;*******************
     wks = gsn_open_wks("png","Tmax-da_31")
   
   
     ;************************************************
     fill_specs = (/"water","land"/)
     ;mask_specs = (/"zambia"/)
     mask_specs = (/"sudan","ethiopia", "kenya","tanzania","uganda", "rwanda",\
	 "burundi","central african republic","chad","congo","congo democratic republic"/)
     ;mask_specs = (/\
     ;"algeria","angola","angola-exclave-called-cabinda","benin","botswana",\
     ;"burundi","cameroon","central-african-republic","chad","congo","djibouti",\
     ;"egypt","equatorial-guinea","ethiopia","gabon","gambia","ghana","guinea",\
     ;"guinea-bissau","ivory-coast","kenya","lesotho","liberia","libya",\
     ;"madagascar","malawi","mali","mauritania","mauritius","morocco",\
     ;"mozambique","namibia","niger","nigeria","rwanda","senegal","sierra-leone",\
     ;"somalia","south-africa","sudan","swaziland","tanzania","togo","tunisia",\
     ;"uganda","upper-volta","western-sahara","zaire","zambia","zimbabwe"/)
   
     ;Define colormap.
     ;************************************************
      color_map = (/(/1.00,1.00,1.00/),(/0.00,0.00,0.00/),(/0.70,0.70,0.70/), \
                  (/0.75,0.50,1.00/),(/0.50,0.00,1.00/),(/0.00,0.00,1.00/), \
                  (/0.00,0.50,1.00/),(/0.00,1.00,1.00/),(/0.00,1.00,0.60/), \
                  (/0.00,1.00,0.00/),(/0.70,1.00,0.00/),(/1.00,1.00,0.00/), \
                  (/1.00,0.75,0.00/),(/1.00,0.38,0.38/),(/1.00,0.00,0.38/), \
                  (/1.00,0.00,0.00/)/)
   
      gsn_define_colormap(wks,color_map)
      
     ;Set some graphical resources
     ;****************************
   
      res                 = True
   
      res@gsnMaximize     = True
      res@cnFillOn        = True
      res@cnLinesOn       = False
      res@tfDoNDCOverlay  = True          ; This is necessary if you don't set sfXArray/sfYArray
      res@cnFillDrawOrder = "Predraw"
      res@cnLineDrawOrder = "Predraw"
      res@gsnDraw         = False         ; don't draw yet
      res@gsnFrame        = False         ; don't advance yet
   
      res@tiMainString    = "WRF 2m Temperature"           ; draw a title
      res@gsnLeftString   = ""
      res@gsnRightString  = ""
      res@gsnCenterString = "Maximum Temperature with da ~C~D+02 Valid on 31/07/2022 From 00z on 29/07/2022"
   
     ;Add additional map resources
     ;****************************
      ;res = wrf_map_resources(a,res)
      ;print(res)
      ;exit

   
     ;Change some of the resources that were set (these were set to "gray")
     ;*********************************************************************
      res@mpDataBaseVersion           = "MediumRes"      ; use finer database
      res@mpDataSetName               = "Earth..4"
      res@mpOutlineOn                 = True
      res@mpOutlineBoundarySets       = "AllBoundaries"
      res@mpGeophysicalLineColor      = "blue"
      res@mpGeophysicalLineThicknessF = 3.0
      res@mpNationalLineColor         = "Sienna"
      res@mpNationalLineThicknessF    = 3.5
      res@mpLimbLineThicknessF        = 0.5
      res@mpLimbLineColor             = "Gray"
      res@mpOceanFillColor            = "Light blue"
      res@mpInlandWaterFillColor      = "blue"
      res@mpGridLatSpacingF           = 2.
      res@mpGridLonSpacingF           = 2.
   
     ;---Zoom in on map and plot again
      res@cnLineThicknessF = 2            ; make lines thicker, Default: 1
      res@mpMinLatF      = min(lat)   ; zoom in on lat/lon area
      res@mpMaxLatF      = max(lat)
      res@mpMinLonF      = min(lon)
      res@mpMaxLonF      = max(lon)
   
   ;   res@mpMinLatF      = -5   ; zoom in on lat/lon area
   ;   res@mpMaxLatF      =  6
   ;   res@mpMinLonF      =  33
   ;   res@mpMaxLonF      =  43
   
       res@cnLevelSelectionMode = "ManualLevels"      ; set contour levels, Options: ManualLevels, AutomaticLevels, EqualSpacedLevels, ExplicitLevels

       res@cnMinLevelValF       =  -9.       ; minimum contour level
       res@cnMaxLevelValF       =  39.       ; maximum contour level
       res@cnLevelSpacingF      =   3.       ; contour level spacing
 
      ;Turn off the fill boundaries, since we are going to explicitly set which boundaries we want to fill.
      ;****************************************************************************************************
       ;res@mpFillBoundarySets    = "NoBoundaries"
       res@mpFillAreaSpecifiers  = fill_specs
       res@mpSpecifiedFillColors = (/7,2/)
   
      ;Turn on masking, and set resource indicating areas to mask.
      ;************************************************************
       res@mpAreaMaskingOn       = 1
       res@mpMaskAreaSpecifiers  = mask_specs
   
       plot = gsn_csm_contour_map(wks, xMaxTime ,res)
   
      ;Define shapefiles
      ;*****************
       shp_filename = ("$HOME/DA/TMINTMAX/Shapefiles/Zambia/ZMB_ADM1.shp")
   
      ;Set shapefile resources
      ;***********************
       shpres                    =  True
       shpres@gsLineThicknessF   =  1             ; increase line thickness
       shpres@gsLineColor        =  "Black"       ; line colorgsLineThicknessF
       shpres@gsLineDashPattern  =  1
   
       shp= gsn_add_shapefile_polylines(wks,plot,shp_filename,shpres)
   
       draw(plot)
       frame(wks)
   
   end