* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * subdiv_pdn.gs * * Plots annual, seasonal, or monthly rainfall anomalies for * subdivisions of India. Requires a DODS-enabled version of GrADS * and three supplemental files: indiaregions, india_polygons.asc, * and readpolygons.gsf. See http://monsoondata.org/customize for * further information. * * Usage: subdiv_pdn * Where: type can be a (annual), s (seasonal), or m (monthly). * year is between 1871 and 2000 * month# is between 1-12 * If no arguments are given, the script will prompt the user. * * Example: subdiv_pdn s 1886 * * Written by J.M. Adams March 2003 * function main(args) * Make sure GrADS is in landscape mode 'q gxinfo' pagesize = sublin(result,2) xsize = subwrd(pagesize,4) ysize = subwrd(pagesize,6) if (xsize != 11 & ysize != 8.5) say 'You must be in LANDSCAPE MODE for this script to work properly.' return endif 'reinit' rc = gsfallow("on") * Parse the args if (args = '') prompt 'Annual, Seasonal, or Monthly Climatology? (a/s/m) -> ' pull type if (type!='a' & type!='s' & type!='m') say 'Invalid choice. Exiting ...' return endif prompt 'Enter a 4-digit year -> ' pull year if (type='m') prompt 'Enter the month (1-12) -> ' pull mon endif else type = subwrd(args,1) year = subwrd(args,2) mon = subwrd(args,3) if (year='') prompt 'Enter a 4-digit year -> ' pull year endif if (type='m' & mon='') prompt 'Enter the month (1-12) -> ' pull mon endif endif if (type='m') months = "jan feb mar apr may jun jul aug sep oct nov dec" month = subwrd(months,mon) endif * Set up the map 'set display color white' 'clear' * Load the polygon info into a script variable rc = readpolygons() * Set colors and levels dorgb() _pcols = '21 20 19 22 18 17 16' _plevs = ' -40 -20 -10 10 20 40 ' * Open a data file and set the dimension environment if (type='a'); 'sdfopen http://monsoondata.org:9090/dods/iitm/subdiv/ann_pdn'; endif if (type='s'); 'sdfopen http://monsoondata.org:9090/dods/iitm/subdiv/seas_pdn'; endif if (type='m'); 'sdfopen http://monsoondata.org:9090/dods/iitm/subdiv/mon_pdn'; endif 'set dfile 2' 'set x 1' 'set y 1' 'set z 1' if (type='m') 'set time 'month%year else 'set time jan'year endif * Draw the polygons div = 2 while (div <= 34) if (div = 2 | div = 12 | div = 15 | div = 16) 'set line 99' else 'd p'div if (subwrd(result,1) = "Result") val = subwrd(result,4) if (val < -40) ; 'set line 21' ; endif if (val >= -40 & val < -20) ; 'set line 20' ; endif if (val >= -20 & val < -10) ; 'set line 19' ; endif if (val >= -10 & val <= 10) ; 'set line 22' ; endif if (val > 10 & val <= 20) ; 'set line 18' ; endif if (val > 20 & val <= 40) ; 'set line 17' ; endif if (val > 40) ; 'set line 16' ; endif else 'set line 23' endif endif _drawpolyf.div div = div + 1 endwhile * Draw the continental outline 'set dfile 1' 'set lat '_minlat' '_maxlat 'set lon '_minlon' '_maxlon 'set map 0' 'set mpdset indiaregions' 'set poli on' 'draw map' * Cover the neighboring coastlines 'set line 0' 'draw recf 7.6 1.5 9.1 4.38' ;* Burma 'draw polyf 5 1.43 4.6 0.68 5.5 0.68' ;* Sri Lanka 'draw polyf 1.96 5.1 2.27213 5.1493 1.96 5.8' ;* Pakistan * Draw a colorbar polycolor() * Draw a title if (type='m') longmonths = "January February March April May June July August September October November December" label1 = subwrd(longmonths,mon) label1 = label1' 'year' Rainfall' else if (type='a') label1 = year' Annual Rainfall' else label1 = year' Seasonal Rainfall' endif endif label2 = '(% departure from normal)' 'set string 99 bc' 'set strsiz 0.13 0.15' ; 'draw string 6.7 '_midpt-0.00' 'label1 'set strsiz 0.12 0.14' ; 'draw string 6.7 '_midpt-0.25' 'label2 * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * function dorgb * These are the BLUE shades 'set rgb 16 10 10 255' 'set rgb 17 110 110 255' 'set rgb 18 180 180 255' * These are the RED shades 'set rgb 19 255 180 180' 'set rgb 20 255 110 110' 'set rgb 21 255 10 10' * These are the GRAY shades 'set rgb 22 220 220 220' 'set rgb 23 30 30 30' 'set rgb 98 150 150 150' 'set rgb 99 1 1 1' * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * function polycolor * Get plot size info 'query gxinfo' rec2 = sublin(result,2) rec3 = sublin(result,3) rec4 = sublin(result,4) xsiz = subwrd(rec2,4) ysiz = subwrd(rec2,6) ylo = subwrd(rec4,4) xhi = subwrd(rec3,6) xd = xsiz - xhi * Set up constants cnum = 7 xl = xhi + xd/2 - 1.6 xr = xl + 0.2 xwid = 0.2 xmid = xl + 0.1 ywid = 0.5 if (ywid*cnum > ysiz*0.88) ywid = ysiz*0.88/cnum endif ymid = ysiz/2 ymid = ymid - 1.5 _midpt = ymid yb = ymid - ywid*cnum/2 'set string 1 l 5' vert = 1 * Plot colorbar 'set strsiz 0.10 0.11' num = 1 while (num <= cnum) col = subwrd(_pcols,num) val = subwrd(_plevs,num) 'set line 'col yt = yb + ywid if (num=cnum) 'draw polyf 'xl' 'yb' 'xr' 'yb' 'xmid' 'yt else if (num=1) 'draw polyf 'xl' 'yt' 'xr' 'yt' 'xmid' 'yb else 'draw recf 'xl' 'yb' 'xr' 'yt endif endif if (num