*Generates 2 panels that display 1-a global DEM elevation and 2-values along a parallel.
*Created by Hernán Sala, Instituto Antartico Argentino

reinit
clear

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


'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 0 0'
'set rgb 20 100 0 0'
'set rgb 22 110 0 0'
'set rgb 24 120 0 0'
'set rgb 25 130 0 0'
'set rgb 26 140 0 0'
'set rgb 28 150 0 0'
'set rgb 30 160 0 0'
'set rgb 32 170 0 0'
'set rgb 34 185 0 0'
'set rgb 35 195 0 0'
'set rgb 36 205 0 0'
'set rgb 38 210 0 0'
'set rgb 40 220 0 0'
'set rgb 42 235 0 0'
'set rgb 44 255 0 0'

'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 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 extremes 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)
*Sets the other extreme of the line
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)
* Draws the line
*Draws line x1 y1 x2 y2
'draw line 'x1' 'y1' 'x2' 'y2

'set vpage off'


***Lower panel**************************************************************
*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'

*Range of longs to be plotted
'set lat -180 180'

*To visualize a cut to a given lat, set a constant value.
*Lat of BsAs es: -34.5 (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 output-2panels.png x1000 y800'