;
; This script tests the problem reported in NCL-2759,
; where setting cnSmoothingOn = True causes one of the
; contour areas to not be filled.
;
; The script works in these cases:
;  - You don't set cnSmoothinOn to True (NCL V6.4.0)
;  - If you set cnSmoothingDistanceF to 0.003
;  - If you use raster fill and set cnRasterSmoothingOn = True
;
lt1 = 5
lt2 = 25
ln1 = 60
ln2 = 90

fin  = addfile("saji_rain_jun-sep.nc","r")
fmk  = addfile("saji_landsea_mask.nc","r")
rain = fin->var(:,{lt1:lt2},{ln1:ln2})
msk  = fmk->ls_mask({lt1:lt2},{ln1:ln2})
lmsk =  conform(rain,msk,(/1,2/))

rain_masked = rain
rain_masked = (/where(lmsk.eq.1, rain@_FillValue,rain) /)
rjja_masked = dim_avg_n_Wrap(rain_masked,0)

wks = gsn_open_wks("png","contour_smoothing_bug")

res               = True

res@gsnDraw       = False       ; will panel later
res@gsnFrame      = False
res@gsnAddCyclic  = False

res@cnFillOn      = True
res@cnLinesOn     = False

res@mpMinLatF     = lt1
res@mpMaxLatF     = lt2
res@mpMinLonF     = ln1
res@mpMaxLonF     = ln2

res@gsnLeftString = ""

;---Create two plots without smoothing
res@tiMainString  = "rjja: no contour smoothing"
plt_rjja_no_smth = gsn_csm_contour_map(wks,rjja_masked,res)

res@tiMainString  = "rain: no contour smoothing"
plt_rain_no_smth = gsn_csm_contour_map(wks,rain_masked(3,:,:),res)

;---Create two plots with smoothing turned on; here's where you see the issue in NCL V6.4.0
res@cnSmoothingOn = True
res@tiMainString = "rjja: contour smoothing"
plt_rjja_smth = gsn_csm_contour_map(wks,rjja_masked,res)

res@tiMainString = "rain: contour smoothing"
plt_rain_smth = gsn_csm_contour_map(wks,rain_masked(3,:,:),res)

;---Panel all four plots
pres = True
pres@gsnMaximize = True
nclv = get_ncl_version()
if(nclv.eq."6.3.0") then
  pres@txString = "NCL version " + nclv
else
  pres@gsnPanelMainString = "NCL version " + nclv
end if
pres@gsnPanelMainString = "NCL version " + nclv
gsn_panel(wks,(/plt_rjja_no_smth,plt_rain_no_smth,plt_rjja_smth,plt_rain_smth/),(/2,2/),pres)
 
;----------------------------------------------------------------------
; Side detour: try raster fill with smoothing
;----------------------------------------------------------------------
rres = res
delete(rres@cnSmoothingOn)
rres@cnRasterSmoothingOn = True 
rres@cnFillMode          = "RasterFill"
rres@tiMainString        = "rjja: raster contour smoothing"
plt_rjja_rst_smth = gsn_csm_contour_map(wks,rjja_masked,rres)

rres@tiMainString        = "rain: raster contour smoothing"
plt_rain_rst_smth = gsn_csm_contour_map(wks,rain_masked(3,:,:),rres)

gsn_panel(wks,(/plt_rjja_rst_smth,plt_rain_rst_smth/),(/2,1/),pres)

;----------------------------------------------------------------------
; Setting cnSmoothingDistanceF to a smaller value than the default 
; (0.01) seems to help.
;
; A value of 0.003 works for both rain and rjja plots. 
; 0.004 and 0.005 causes a problem with the rjja plot. 
; 0.007 has a problem for both plots.
;----------------------------------------------------------------------
distances = (/0.007,0.005,0.004,0.003,0.002/)
ndist     = dimsizes(distances)
plts_dist = new(2*ndist,graphic)

do nd=0,ndist-1
  res@cnSmoothingDistanceF = distances(nd)
  res@tiMainString = "rjja: contour smooth distance = " + res@cnSmoothingDistanceF
  plts_dist(nd*2) = gsn_csm_contour_map(wks,rjja_masked,res)
  res@tiMainString = "rain: contour smooth distance = " + res@cnSmoothingDistanceF
  plts_dist((nd*2)+1) = gsn_csm_contour_map(wks,rain_masked(3,:,:),res)
end do
gsn_panel(wks,plts_dist,(/5,2/),pres)