* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * 
* subdiv_timeseries.gs
* 
* Plots a variety of timeseries of rainfall data for 
* subdivisions of India. Requires a DODS-enabled version of GrADS.
* See http://monsoondata.org/customize for further information. 
*
* Usage: subdiv_timeseries <type> <var> <div> <byr> <eyr>
* Where: type can be a (annual), s (seasonal), or m (monthly).
*        var  can be t (total) or p (% departure from normal)
*        var  can be also be c (climatology) if type is 'm'
*        div  is the subdivision number (3-34)
*        byr  is the beginning year
*        eyr  is the ending year 
* If no arguments are given, the script will prompt the user.
*
* Examples: subdiv_timeseries m t 33 1900 1910
*           subdiv_timeseries m c 8
* 
* Written by J.M. Adams March 2003
*
function main(args)
'reinit'  
* Parse arguments
if (args = '')
  prompt 'Do you want Monthly, Seasonal, or Annual timeseries? (m/s/a) -> '
  pull type
  if (type='m')
    say 'With monthly data, you can plot the following:'
    say '  A 12-month Climatology  (c)'
    say '  Monthly Rainfall Totals (t)'
    say '  Monthly Percentage Different from Normal (p)'
    prompt 'Enter your choice (c/t/p) -> '
    pull var
  endif
  if (type='s')
    say 'With seasonal data, you can plot the following:'
    say '  Rainfall Totals (t)'
    say '  Percentage Different from Normal (p)'
    prompt 'Enter your choice (t/p) -> '
    pull var
  endif
  if (type='a')
    say 'With annual data, you can plot the following:'
    say '  Rainfall Totals (t)'
    say '  Percentage Different from Normal (p)
    prompt Enter your choice (t/p) ->
    pull var
  endif
  if (var='p' | var='t')
    prompt 'Enter the beginning year [1871 by default] -> '
    pull byear
    if (byr=''); byr=1871; endif
    prompt 'Enter the ending year [2000 by default] -> '
    pull eyr
    if (eyr=''); eyr=2000; endif
  endif
  prompt 'Enter the subdivision number (3-34) -> '
  pull div
  if (div<3 | div=12 | div=15 | div=16 | div>34)
    say 'No data for that Subdivision'
    return
  endif
else
  type = subwrd(args,1)
  var  = subwrd(args,2)
  div  = subwrd(args,3)
  byr  = subwrd(args,4)
  eyr  = subwrd(args,5)
endif  


if (type = 'a')
  if (var = 'p'); 'sdfopen http://monsoondata.org:9090/dods/iitm/subdiv/ann_pdn'; endif
  if (var = 't'); 'sdfopen http://monsoondata.org:9090/dods/iitm/subdiv/ann_total'; endif
  title = "Annual "
endif

if (type = 's')
  if (var = 'p'); 'sdfopen http://monsoondata.org:9090/dods/iitm/subdiv/seas_pdn'; endif
  if (var = 't'); 'sdfopen http://monsoondata.org:9090/dods/iitm/subdiv/seas_total'; endif
  title = "Seasonal (Jun-Sep) "
endif

if (type = 'm')
  if (var = 'c'); 'sdfopen http://monsoondata.org:9090/dods/iitm/subdiv/mon_clim'; endif
  if (var = 'p'); 'sdfopen http://monsoondata.org:9090/dods/iitm/subdiv/mon_pdn'; endif
  if (var = 't'); 'sdfopen http://monsoondata.org:9090/dods/iitm/subdiv/mon_total'; endif
  title = "Monthly "
  if (var='c'); title=title%'Climatological '; endif
endif

if (var = 'c'); ylabel='mm/month'; endif
if (var = 'p'); ylabel='% different from normal'; endif
if (var = 't'); ylabel='mm/month'; endif
 

* Set window characteristics
'set display color white'
'clear'

* Set the time range
if (type='m')
  if (var='p' | var='t')
    'set time jan'byr' dec'eyr
    'q dims'
    tdims = sublin(result,5)
    t1 = subwrd(tdims,11)
    t2 = subwrd(tdims,13)
    'set t 't1-0.5' 't2+0.5
  endif
  if (var='c')
    'set t 0.5 12.5'
    'set tlsupp year'
  endif
endif
if (type='s' | type='a')
  'set time jan'byr' jan'eyr
  'q dims'
  tdims = sublin(result,5)
  t1 = subwrd(tdims,11)
  t2 = subwrd(tdims,13)
  'set t 't1-0.5' 't2+0.5
endif

* Set plot area
yt = 8
yb = 4
xl = 0.75
xr = 10.8
'set parea 'xl' 'xr' 'yb' 'yt

* Get range of values to set yaxis limits
if (type='m')
  if (var='c' | var='t')
    'set gxout stat'
    'd p'div
    range = sublin(result,8)
    max = subwrd(range,5)
    max = max + 0.05*max
    max = math_nint(max)
    'set vrange -2 'max
  endif
endif

* Set bar chart characteristics
'set gxout bar'
if (var='p'); 'set barbase 0'; endif
if (var='c'); 'set barbase 0'; endif
if (var='t' & type='m'); 'set barbase 0'; endif
'set bargap 40'
'set grads off'
'set rgb 99 1 1 1'
'set annot 99'
'set xlopts 99'
'set ylopts 99'

* Draw the time series
if (var='c' | var='t')
  'set ccolor 3'
  'd p'div
endif
if (var='p')
  'd p'div
  'set ccolor 4'
  'd maskout(p'div',p'div')'
  'set ccolor 2'
  'd maskout(p'div',-p'div')'
endif

*   Annotate
if (div= 3) ; name = 'North Assam'               ; endif
if (div= 4) ; name = 'South Assam'               ; endif
if (div= 5) ; name = 'Sub-Himalayan West Bengal' ; endif
if (div= 6) ; name = 'Gangetic West Bengal'      ; endif
if (div= 7) ; name = 'Orissa'                    ; endif
if (div= 8) ; name = 'Bihar Plateau'             ; endif
if (div= 9) ; name = 'Bihar Plains'              ; endif
if (div=10) ; name = 'East Uttar Pradesh'        ; endif
if (div=11) ; name = 'West Uttar Pradesh'        ; endif
if (div=13) ; name = 'Haryana'                   ; endif
if (div=14) ; name = 'Punjab'                    ; endif
if (div=17) ; name = 'West Rajsthan'             ; endif
if (div=18) ; name = 'East Rajsthan'             ; endif
if (div=19) ; name = 'West Madhya Pradesh'       ; endif
if (div=20) ; name = 'East Madhya Pradesh'       ; endif
if (div=21) ; name = 'Gujrat'                    ; endif
if (div=22) ; name = 'Saurashtra & Kutch'        ; endif
if (div=23) ; name = 'Konkan & Goa'              ; endif
if (div=24) ; name = 'Madhya Maharashtra'        ; endif
if (div=25) ; name = 'Marathwada'                ; endif
if (div=26) ; name = 'Vidarbha'                  ; endif
if (div=27) ; name = 'Coastal Andhra Pradesh'    ; endif
if (div=28) ; name = 'Telangana'                 ; endif
if (div=29) ; name = 'Rayalaseema'               ; endif
if (div=30) ; name = 'Tamil Nadu'                ; endif
if (div=31) ; name = 'Coastal Karnataka'         ; endif
if (div=32) ; name = 'North Interior Karnataka'  ; endif
if (div=33) ; name = 'South Interior Karnataka'  ; endif
if (div=34) ; name = 'Kerala'                    ; endif

xlabel = title 'Rainfall for Subdivision 'div': 'name 
ymid = yb + (yt-yb)/2 
xmid = xl + (xr-xl)/2

'set string 99 c 1 0'
'set strsiz 0.16 0.18' ; 'draw string 'xmid' 'yb-0.5' 'xlabel
'set string 99 c 1 90'
'set strsiz 0.11 0.13' ; 'draw string 'xl-0.56' 'ymid' 'ylabel