* This script creates 2 panels for a global DEM
* First panel is a global visualization.
* Second panel is a slide along a parallel.

* Created by: H. Sala - Instituto Antártico Argentino


*Before its execution you should download the file "elev.1-deg.nc" from:
* http://www.jisao.washington.edu/data_sets/elevation/                             

reinit
clear
'sdfopen elev.1-deg.nc'

'set grid off'

***Upper pannel****************************************************
*set vpage xmin xmax ymin ymax
'set vpage 0 11 2.5 8.5'
'set grads off'
'set gxout shaded'
*set frame on|off|circle
'set frame on'
*Defines an arbitrary gray scale to be used for topography.
'set rgb 18 80 80 80'
'set rgb 20 100 100 100'
'set rgb 22 110 110 110'
'set rgb 24 120 120 120'
'set rgb 25 130 130 130'
'set rgb 26 140 140 140'
'set rgb 28 150 150 150'
'set rgb 30 160 160 160'
'set rgb 32 170 170 170'
'set rgb 34 185 185 185'
'set rgb 35 195 195 195'
'set rgb 36 205 205 205'
'set rgb 38 210 210 210'
'set rgb 40 220 220 220'
'set rgb 42 235 235 235'
'set rgb 44 255 255 255'

'set csmooth on'

*Sets graphics type
'set gxout shaded'

*Asigns gray scale values for different hights
'set clevs 0 250 500 1000 1500 2000 2500 3000 3500 4000 4500 5000 6000 7000 8000 8500'
'set ccols 0 18 20 22 24 25 26 28 30 32 34 35 36 38 40 42 44'

*set xaxis start end <increment>
'set xaxis 0 0'

*Displays terrain topography as a gray scale
'd data'

*Sets the size for the strings
'set strsiz 0.17'
*Sets color of the strings
'set string 3'

*set line color# <style> <thickness>
'set line 11 1 2'

*Sets the extrems of the line
lon1= 0
lat1= -34.5
'q w2xy 'lon1' 'lat1
say result
*Loads the XY coords in x1 and y1
x1 = subwrd(result,3)
y1 = subwrd(result,6)
*Establece el otro extremo de la recta
lon2= 360
lat2= -34.5
'q w2xy 'lon2' 'lat2
say result
*Loads the XY coords in x2 and y2
x2 = subwrd(result,3)
y2 = subwrd(result,6)
* Dibuja la recta
*draw line x1 y1 x2 y2
'draw line 'x1' 'y1' 'x2' 'y2

'set vpage off'


***Lower pannel***************************************************
*set vpage xmin xmax ymin ymax
'set vpage 0 9.7 0 3.8'
'set grads off'
'set gxout line'

*set vrange y1 y2 - Specifies the range of y-axis scaling
'set vrange -7000 4000'

*set xlopts color <thickness <size>>
'set xlopts 1 4 0.13'

*Indica el rango de longitudes en q se va a graficar
'set lat -180 180'

*Para visualizar un corte a una determinada lat, darle un valor cte.
*La lat de BsAs es: -34.5 grados S
'set lat -34.5'

*set line color# <style> <thickness>
'set ccolor 11'

*set annot color <thickness> for the axis border, axis labels, and tickmarks
*'set annot 2'
'd data'

*draw xlab and ylab strings
'draw xlab Longitud'
'draw ylab Elev.'

*set line color# <style> <thickness>
'set line 1 2 1'
*draw line x1 y1 x2 y2
'draw line 2 2.54 10.5 2.54'

'set vpage off'


'printim DEM-2panels_ByN.png x1000 y800'