******* Phan loai cap TC cua JMA:
* 2 : Tropical Depression (TD)
* 3 : Tropical Storm (TS)
* 4 : Severe Tropical Storm (STS)
* 5 : Typhoon (TY)
* 6 : Extra-tropical Cyclone (L)
* 7 : Just entering into the responsible area of Japan Meteorological Agency (JMA)
* 8 : Not used
* 9 : Tropical Cyclone of TS intensity or higher
***********************

   "reinit"
   "set display white"; "c"
   "set mpdraw on"; "set frame on"; "set grid on";"set xlab on";"set ylab on";

    mark.2=0; mark.3=0; mark.4=0; mark.5=0; mark.6=0; mark.7=0; mark.9=0
*    colo.2=4; colo.3=10; colo.4=12; colo.5=2; colo.6=1; colo.7=12; colo.9=12
    colo.2=5; colo.3=7; colo.4=8; colo.5=6; colo.6=0; colo.7=12; colo.9=12
    thic.2=4; thic.3=4; thic.4=4; thic.5=4; thic.6=4; thic.7=4; thic.9=4
    styl.2=1; styl.3=1; styl.4=1; styl.5=1; styl.6=1; styl.7=1; styl.9=1;

    NInt_tc=4
    inten.1="2"; inten.2="3"; inten.3="4"; inten.4="5"; inten.5="6"; inten.6="7"; inten.7="9"
    cdn.1="(TD)"; cdn.2="(TS)"; cdn.3="(STS)"; cdn.4="(TY)"

    x00=1.0; x01=8.3; y00=1; y01=10
*    Main_title ((x00+x01)/2-0.5,y01-0.35, "TC tracks during period 1975-2014 (Based on RSMC data)")
  
   "sdfopen 4TC_paper.nc"
   lat00=0; lat01=35; lon00=90; lon01=130; dlat=lat01-lat00; dlon=lon01-lon00
   "set lon "lon00" "lon01; "set lat "lat00" "lat01

   "set ccolor 0"; "set clab off"; "set grads off"
   "set mpdset hires"; "set xlint 5"; "set ylint 5"
   "d tc"
*   "hatch "tc" 1.1 3 -angle 45,135 -int 0.1 -color 15"

*  Define something needed:
   'query gxinfo'
   rec3 = sublin(result,3) 
   rec4 = sublin(result,4) 
   xlim1 = subwrd(rec3,4)
   xlim2 = subwrd(rec3,6)
   ylim1 = subwrd(rec4,4)
   ylim2 = subwrd(rec4,6)
   xd = xlim2 - xlim1
   yd = ylim2 - ylim1
*   say 'xd = 'xd'  yd = 'yd
*   say 'xlim2 = 'xlim2'  ylim2 = 'ylim2
   'query dims'
   rec3 = sublin(result,2) 
   rec4 = sublin(result,3) 
   xl1 = subwrd(rec3,6)
   xl2 = subwrd(rec3,8)
      if (xl1 < 0)
      xl1 = xl1  + 360.0
      endif
      if (xl2 < 0)
      xl2 = xl2  + 360.0
      endif
   yl1 = subwrd(rec4,6)
   yl2 = subwrd(rec4,8)
   lonw = xl1
   lats = yl1
   xld = xl2 - xl1
   yld = yl2 - yl1
*   say 'xld = 'xld'  yld = 'yld
*   say 'lats = 'lats'  lonw = 'lonw

*  Legend first:
   tit.1="TC Genesis locations during period 1975-1994"; tit.2="TC Genesis locations during period 1995-2014"
   tit.3="TC Decay locations during period 1975-1994"  ; tit.4="TC Decay locations during period 1995-2014"
   M.1=3; M.2=9; M.3=1; M.4=6
*   C.1=4; C.2=2; C.3=4; C.4=2
    C.1=1; C.2=1; C.3=1; C.4=1
   N=4
   x00=2; x01=8; dx=(x01-x00)/NInt_tc
   y00=1.1; icd=0
   While (icd<N)
      icd=icd+1; xx0=x00+0.2; yy0=y00-(icd-1)*0.2
      draw_titl(xx0, yy0, tit.icd)
      "set strsiz 0.12 0.12"; "set line "C.icd
      "draw mark "M.icd" "xx0-0.3" "yy0" "0.12
   Endwhile

*   "drawrec 105 120 5 25 1"
*   "drawrec1 110 120 5 15 5"

********************************************************
*  Years to be plotted: Case=3

   yy0=1975; yy1=2014
   Nyear=yy1-yy0+1; Nyear2=Nyear/2

*  ENSO years: Case=1 for El Nino, Case=2 for La Nina
   NEL=15; NLA=18
   E.1.1=1951; E.1.2=1957; E.1.3=1963; E.1.4=1965; E.1.5=1972; E.1.6=1977; E.1.7=1982; E.1.8=1987; 
   E.1.9=1991; E.1.10=1994; E.1.11=1997; E.1.12=2002; E.1.13=2004; E.1.14=2009; E.1.15=2015

   E.2.1=1954; E.2.2=1955; E.2.3=1956; E.2.4=1964; E.2.5=1970; E.2.6=1971; E.2.7=1973; E.2.8=1974; 
   E.2.9=1975; E.2.10=1984; E.2.11=1985; E.2.12=1988; E.2.13=1998; E.2.14=1999; E.2.15=2000; 
   E.2.16=2007; E.2.17=2010; E.2.18=2011

*  ------ Chose Case for plotting: -------
   Case=3
   if (Case=1); Ncase=NEL; ica=1; Fout="EL"; endif
   if (Case=2); Ncase=NLA; ica=2; Fout="LA"; endif
   if (Case=3); Ncase=Nyear; Fout=yy0"-"yy1; endif

*  Reading and plotting:
*   filename = TC_tr.txt

* --------------------------------------
   filename = "/home/tanpv/Progs/TCs/TC_Dat/TC_tr_for_GrADS.txt"
* --------------------------------------
   code=0
   res=read (filename)
   code=sublin(res,1)
   if (code !=0); say Error !!!; pull dd; endif
   Ntc=sublin(res,2)

   While (code=0)
      res=read (filename)
      code=sublin(res,1)
      if (code=2); "draw map"; Rec_Hat (); "gxprint Tracks_"Fout".pdf pdf"; "close 1"; exit; endif
      if (code=9); say "something wrong !! "code; "close 1"; exit; endif
      xx=sublin(res,2)
      Nlf=subwrd(xx,1)
      res=read (filename)
      code=sublin(res,1)
      xx = sublin(res,2)
      id = subwrd(xx,1); yr = subwrd(xx,2); mm = subwrd(xx,3); dd = subwrd(xx,4); hh = subwrd(xx,5);
      lon1 = subwrd(xx,6); lat1 = subwrd(xx,7); CD1 = subwrd(xx,8) 
      if (CD1=9); CD1=3; endif
      lon1 = lon1*0.10; lat1 = lat1*0.10
      x1 = (lon1 - lonw)/ xld * xd; y1 = (lat1 - lats)/ yld * yd
      x1 = x1 + xlim1; y1 = y1 + ylim1

      ilf=1
      While (ilf<Nlf)
         ilf=ilf+1
         res=read (filename)
         code=sublin(res,1)
         xx = sublin(res,2)
         id = subwrd(xx,1)
         yr = subwrd(xx,2); mm = subwrd(xx,3); dd = subwrd(xx,4); hh = subwrd(xx,5);
         lon2 = subwrd(xx,6); lat2 = subwrd(xx,7)
         lon2 = lon2/10; lat2 = lat2/10; CD2 = subwrd(xx,8); if (CD2=9); CD2=3; endif

         x2 = (lon2 - lonw)/ xld * xd; y2 = (lat2 - lats)/ yld * yd
         x2 = x2 + xlim1; y2 = y2 + ylim1

         yyyy=0
         While (yyyy<Ncase)
           yyyy=yyyy+1; 
           if (Case=3);
*             genesis locations:
*              if (yr>=yy0 & yr<=yy1 & ilf=2);
              if (yr>=yy0 & yr<=yy1);
*                 say yr" "CD2
*                 "set line "colo.CD2" 1 4"
                 "set line 15 1 1"
                 if (x1>xlim1 & x1<xlim2 & x2>xlim1 & x2<xlim2 & y1>ylim1 & y1<ylim2 & y2<ylim2);
                    "draw line "x1" "y1" "x2" "y2;
                    if (ilf=2);
                       if (yr<=yy0+Nyear2-1); "set line 1"; "draw mark 3 "x1" "y1" "0.08; endif
                       if (yr> yy0+Nyear2-1); "set line 1"; "draw mark 9 "x1" "y1" "0.08; endif 
                    endif 
                    if (ilf=Nlf);
                       if (yr<=yy0+Nyear2-1); "set line 1"; "draw mark 1 "x2" "y2" "0.12; endif
                       if (yr> yy0+Nyear2-1); "set line 1"; "draw mark 6 "x2" "y2" "0.12; endif
                    endif
                 endif
                 x1 = x2; y1 = y2;
              endif
*             decay locations:
*              if (yr>=yy0 & yr<=yy1 & ilf=Nlf);
**                 say yr" "CD2
**                 "set line "colo.CD2" 1 4"
*                   "set line 15 1 1"
*                 if (x1>xlim1 & x1<xlim2 & x2>xlim1 & x2<xlim2 & y1>ylim1 & y1<ylim2 & y2<ylim2);
*                    "draw line "x1" "y1" "x2" "y2;
*                   if (yr<=yy0+Nyear2-1); "set line 1"; "draw mark 1 "x2" "y2" "0.12; endif
*                   if (yr> yy0+Nyear2-1); "set line 1"; "draw mark 6 "x2" "y2" "0.12; endif  
*                 endif
*                 x1 = x2; y1 = y2;
*              endif
              yyyy=Ncase
           endif

           if (Case=1 | Case=2);
              if (yr=E.ica.yyyy);
                 say yyyy" "yr" "E.ica.yyyy" "CD2
                 "set line "colo.CD2" 1 4"
                 if (x1>xlim1 & x1<xlim2 & x2>xlim1 & x2<xlim2 & y1>ylim1 & y1<ylim2 & y2<ylim2);
                    "draw line "x1" "y1" "x2" "y2;
                 endif
                 x1 = x2; y1 = y2;
              endif
           endif 
         Endwhile
         if (code=2); "draw map"; Rec_Hat (); "gxprint Tracks_"Fout".pdf pdf"; "close 1"; exit; endif
         if (code=9); say "something wrong !! "code; "close 1"; exit; endif
      Endwhile
   Endwhile
   "draw map"; Rec_Hat (); "gxprint Tracks_"Fout".pdf pdf"
   pull dd
*  ----------------------------------------------------
   function Rec_Hat ()
     lat00=0; lat01=35; lon00=90; lon01=130; dlat=lat01-lat00; dlon=lon01-lon00
     "set lon "lon00" "lon01; "set lat "lat00" "lat01

     "set ccolor 0"; "set clab off"; "set grads off"
     "set mpdset hires"; "set xlint 5"; "set ylint 5"
*     "d tc"
     "hatch "tc" 1.1 3 -angle 45,135 -int 0.15 -color 1"

     "drawrec 105 120 5 25 1"
    "drawrec1 110 120 5 15 5"
   return
*******************************************************
   function Main_title (x, y, title)
     "set string 1 c 5.5 0"
     "set strsiz 0.12 0.13"
     "draw string "x" "y" "title
   return
*******************************************************

   function draw_tit(x, y, title)
      "set string 1 c 5.5"; "set strsiz 0.12 0.12"
      "draw string "x" "y" "title
   return
*******************************************************

   function draw_titl(x, y, title)
      "set string 1 l 5.5"; "set strsiz 0.12 0.12"
      "draw string "x" "y" "title
   return
*******************************************************