;------LOAD THE LIBRARIES-----
load "$NCARG_ROOT/lib/ncarg/nclex/gsun/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
;-----------------------------


begin

;;******DIKKAT BU KOD SADECE IKI AYIN PANEL BICIMINDE CIZIMI ICIN*****
;******ilk iki sutun cizilecek degerler, 3-4 enlem boylam******
;********degiskenlerle oynamayiniz:)**************

homsiniffile="~/Desktop/DATA_MAKALESI/yagis_interpol/homogenized/hom_siniflama.txt"
dummy=asciiread(homsiniffile,-1,"float") ; -1 means read all rows into a one-dimensional variable
ncol=19 ; number of columns is 19, we know this
npts=dimsizes(dummy)/ncol ; get number of points
sinifdata=onedtond(dummy,(/npts,ncol/)) ; sinifdata will be npts x ncol 

ay1=1
ay2=7 ; which months will be shown?

z=sinifdata(:,ay1-1:ay2-ay1:ay2-1)    ; z values to be shown (one column for each month)
y=sinifdata(:,13) ; latitude values
x=sinifdata(:,14) ; longitude values

wks = gsn_open_wks("eps","sinif_map_hom") ; Open a workstation.
drawNDCGrid(wks) ; DRAWS NORMALIZED DEVICE COORDINATES (for ease)

res = True

xmin=25.50; vertices of map 
ymin=35.50
xmax=45.20
ymax=42.50

getvalues wks                 ; Retrieve default color map
    "wkColorMap" : cmap
end getvalues

res@gsnDraw = False   ; draw the graphics
res@gsnFrame= False
;res@gsnSpreadColors = True
;res@gsnMaximize= True ; maximize map 


cmap(0,:) = (/1.0,1.0,1.0/)      ; set the background to white and the foreground to black
cmap(1,:) = (/0.,0.,0./)
gsn_define_colormap(wks, cmap) 


;res@mpProjection = "CylindricalEquidistant" ; CHANGE MAP PROJECTION
res@mpProjection = "Mercator" ; CHANGE MAP PROJECTION
;res@mpCenterLonF = (xmin+xmax)/2   ; CENTER LONGITUDE
;res@mpCenterLatF = (ymin+ymax)/2   ; CENTER LATITUDE

res@mpDataBaseVersion = "MediumRes"

res@mpLimitMode = "LatLon" ; LIMIT THE MAP AREA BY LATITUDE LONGITUDE
res@mpMinLonF   = xmin
res@mpMaxLonF   = xmax
res@mpMinLatF   = ymin
res@mpMaxLatF   = ymax

res@pmTickMarkDisplayMode = "Always" ; enlemler-boylamlari gostermek icin
res@tmXBLabelsOn         = False              ; do not draw bottom labels
res@tmXBOn               = False
;res@tiMainString = ":F26:Precip Correlations with NCPI for Spring" ; SET A TITLE
;res@tiXAxisString = "Longitude" ; x ekseni baslik
;res@tiYAxisString = "Latitude"  ; y ekseni baslik
;res@tiXAxisFontHeightF = 0.012 ; x ekseni baslik font height
;res@tiYAxisFontHeightF = 0.012 ; y ekseni baslik font height

res@mpPerimOn   = False ; 	TURN MAP PERIMETER ON

;res@vpXF      = 0.2 ; CHANGE THE SIZE AND LOCATION OF THE MAP ON THE VIEWPORT
;res@vpYF      = 0.8
;res@vpWidthF  = 0.6
;res@vpHeightF = 0.6
;res@vpKeepAspect = True
;res@vpUseSegments = True

res@gsnScale = True

;res@pmLabelBarDisplayMode = "Always" ; Turn on/off a label bar.
;res@lbPerimOn             = False    ; Turn on/off label bar perim.

res@mpOutlineOn           = True
res@mpOutlineBoundarySets = "Geophysical"
;res@mpOutlineSpecifiers   = "Turkey"
res@mpGeophysicalLineThicknessF = 0.5
res@mpNationalLineThicknessF    = 0.5
res@mpFillOn              = True
res@mpFillAreaSpecifiers  = (/"Water","Land"/)
res@mpSpecifiedFillColors = (/"white","white"/) ; ilki deniz, oburu turkiye disi karalar icin
res@mpAreaMaskingOn       = True
res@mpMaskAreaSpecifiers  = "Turkey"

res@mpGridAndLimbOn = False
;res@mpGridLonSpacingF = 1.0
;res@mpGridAndLimbDrawOrder = "PreDraw"
res@mpGridMaskMode        = "MaskMaskArea"

res@mpFillDrawOrder       = "Draw"
res@mpOutlineDrawOrder    = "PostDraw"

map1 = gsn_csm_map(wks,res) ; Draw the first map.

res@tmXBLabelsOn         = True          ; do not draw bottom labels for the second map
res@tmXBOn               = True
map2 = gsn_csm_map(wks,res) ; Draw the second map.

polyres=True ; polymarker resources
stanumber = dimsizes(z(:,0))
dum1 = new(stanumber,graphic)
dum2 = new(stanumber,graphic)
markercodes = (/16,4,5,14/) ; the codes for the markers. CHANGE THEM HERE

do a=0, 1
do i=0, stanumber-1
        if (z(i,a).eq.1) then
           polyres@gsMarkerIndex = markercodes(0) ; non-testable
        else
        if (z(i,a).eq.100) then
           polyres@gsMarkerIndex = markercodes(1) ; homogeneous
        else
        if (z(i,a).gt.350) then       ; inhomogeneous after correction
           polyres@gsMarkerIndex = markercodes(2)  
        else polyres@gsMarkerIndex = markercodes(3)  ; homogeneous after correction
        end if
        end if
        end if
        if (a.eq.0) then
        dum1(i) = gsn_add_polymarker(wks,map1,x(i),y(i),polyres) ;Add the polymarkers
        else dum2(i) = gsn_add_polymarker(wks,map2,x(i),y(i),polyres) ;Add the polymarkers
        end if
end do
end do

;;LEGEND ITEMS AND  LOCATIONS;;
sx = 0.7                                  ; CHANGE HERE ONLY
symbol_y_values = (/0.55,0.53,0.49,0.51/) ; CHANGE HERE ONLY
symbol_x_values = (/sx,sx,sx,sx/)
labels = (/"Non-testable","Homogeneous","Homogeneous After Correction","Inhomogeneous"/)
lx = sx + 0.2
label_x_values = (/lx,lx,lx,lx/)
label_y_values = symbol_y_values

;;LEGEND RESOURCES;;
leg_symbol_res=True
leg_label_res=True
;leg_symbol_res@gsMarkerSizeF = 0.1
;leg_label_res@txFontHeightF = 0.08

;;DRAW THE LEGEND;;
do i=0, dimsizes(markercodes) - 1
        leg_symbol_res@gsMarkerIndex = markercodes(i)
        gsn_polymarker_ndc(wks,symbol_x_values(i),symbol_y_values(i),leg_symbol_res)
        gsn_text_ndc(wks,labels(i),label_x_values(i),label_y_values(i),leg_label_res)
end do

panelres = False
gsn_panel(wks,(/map1,map2/),(/2,1/),panelres)


end