Generic Mapping Tools Eduardo Suarez Mauricio Gende * [email protected] 11 de octubre de 2011

´Indice 1. ¿Qu´ e es Generic Mapping Tools?

2

2. Graficas X vs. Y 2.1. La base es lo primero. . . . . . . . . . . . . . . . . . . . . . . 2.2. Mejorando el sistema de referencia. . . . . . . . . . . . . . . .

4 5 7

3. Poniendo los datos en el gr´ afico.

10

4. Cambiando la representaci´ on de los datos.

12

5. Empezamos a graficar mapas

16

6. Mallas de datos

20

7. Curvas de nivel 7.1. Vecino m´as pr´oximo . . . . . . . . . . 7.2. Triangulaci´on de Delaunay . . . . . . . 7.3. Curvatura continua con m´ınima tensi´on 7.4. M´ascaras y cortes . . . . . . . . . . . .

22 23 24 25 27

*

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

Facultad de Ciencias Astron´ omicas y Geof´ısicas, Paseo del Bosque s/n, 1900 La Plata, Argentina

1

GMT

Referenciaci´ on

8. Pongamos un poco de color 27 8.1. Necesitamos una paleta . . . . . . . . . . . . . . . . . . . . . . 28 8.2. Armando un mapa con colores . . . . . . . . . . . . . . . . . 28 8.3. La iluminaci´on para marcar el relieve . . . . . . . . . . . . . . 29

1.

¿Qu´ e es Generic Mapping Tools?

Generic Mapping Tools (GMT) es una herramienta inform´atica que permite representar datos en forma gr´afica; debido a que naci´o en una ´epoca donde las computadoras en general no pose´ıan una interfaz gr´afica o ´esta era muy primitiva presenta las siguientes caracter´ısticas: El programa no posee un sistema de men´ ues gr´aficos que se maneja pulsando el rat´on sobre ellos sino que est´a pensando como una serie de comandos Cada uno de los comandos realiza una tarea simple y espec´ıfica cuyo resultado gr´afico se superpone o “apila” al anterior (ver figura 1). Cada comando tiene modificadores, algunos de estos son obligatorios y otros son optativos. Adem´as cada modificador tiene uno o m´as par´ametros, cada uno con la sintaxis particular del modificador. Los comandos ejecutados en forma secuencial terminan por armar una figura. El gr´afico se crea en un archivo en lenguaje postscript y no en la pantalla. Aunque el punto fuerte de GMT est´a en la representaci´on gr´afica de datos geoespaciales, tambi´en es posible generar gr´aficos m´as sencillos, como series temporales. Eventualmente el programa sirve para procesar informaci´on num´erica sin necesidad de visualizarla, por ejemplo generar una malla de valores 2

GMT

Referenciaci´ on

regulares a partir de datos que no presentan una distribuci´on equiespacial.

Figura 1: El gr´afico se construye solapando capas de soluciones parciales. En t´erminos generales GMT tiene Archivos de entrada relacionados con los datos Una lista con o´rdenes a realizar mediante diferentes comandos Una salida, generalmente gr´afica Un esquema completo aparece en la figura 2.

3

GMT

Referenciaci´ on

Figura 2: Entradas y salidas de GMT (extra´ıdo de Smith, 1995.) Aprendiendo mediante el estudio de ejemplos: Vamos a aprender a usar GMT mediante la explicaci´ on detallada de algunos ejemplos. Luego de haberlos vistos y comprendido el alumno deber´ıa poder recurrir al manual para aprender por s´ı solo detalles el programa. Advertimos que seg´ un dice su creador “GMT no es para flojos”, pero citando a una geof´ısica de nuestra Facultad “Al principio es duro pero al final pod´ es hacer hasta que te cebe mates”.

2.

Graficas X vs. Y

Empezamos con un ejemplo que usted podr´ıa resolver mediante el uso de Octave, o GNUplot o Excel. Lo hacemos con GMT porque creemos que es la manera m´as amena, si es que la hubiera, de familiarizarse con objetos ya conocidos por usted. La idea es graficar un per´ıodo completo de la funci´on f(x) = seno(x).

4

GMT

2.1.

Referenciaci´ on

La base es lo primero.

Comenzamos generando el sistema de referencia para el gr´afico. La salida de GMT se muestra en la figura 3.

Figura 3: Los ejes de un sistema cartesiano

5

GMT

Referenciaci´ on

El c´odigo de GMT ser´a: 1. #!/bin/bash 2. # Ejemplo 1 de GMT 3. # Prop´ osito: Mostrar un gr´ afico de serie temporal 4. gmtset PAGE_ORIENTATION=landscape MEASURE_UNIT=cm 5. psbasemap -R0/360/-1.25/1.25 -JX14 -Ba90/a0.5WS -V >ejemplo_1a.ps La explicaci´on es la siguiente: La l´ınea 1 indica que usaremos el int´erprete de comandos Bash. Las l´ıneas 2 y 3 son comentarios y como tales empiezan con #. La l´ınea 4 define algunas variables de GMT, esto no es necesario hacerlo si uno est´a conforme con todas las variables que GMT toma por omisi´on. Estas variables est´an en el archivo .gmtdefaults4 que se encuentra en el directorio “home” de cada usuario. La l´ınea 5 arma finalmente el sistema de referencia usando el comando psbasemap. Analicemos cada modificador del mismo: • -R nos indica la regi´on de trabajo. Debemos indicar la regi´on en la forma x m´ınimo/x m´aximo/y m´ınimo/y m´aximo. Notemos que no hemos elegido para el eje de las Y el intervalo [-1, 1] solo para visualizar mejor la funci´on. • -J nos indica la manera de proyectar los datos al plano y la escala (o tama˜ no) de nuestro gr´afico. En este caso no necesitamos proyectar los datos ya que est´an en el plano, por lo tanto elegimos el par´ametro X y luego indicamos el tama˜ no (en cent´ımetros) del ancho del gr´afico. Si en lugar de indicar el tama˜ no, quisi´eramos usar una escala, en lugar de X (letra equis may´ uscula), usar´ıamos x (letra equis min´ uscula) e indicar´ıamos la escala como 1:valor.

6

GMT

Referenciaci´ on

• -B nos indica la ubicaci´on de las marcas, de las anotaciones en los ejes y de la malla o ret´ıcula. Con el par´ametro a indicamos el paso de las anotaciones, que en este caso son cada 90 grados en la horizontal y 0.5 unidades en la vertical. Como indicamos solamente el paso de las anotaciones, las marcas en los ejes corresponden solamente a ´estas. Las anotaciones estar´an al oeste o izquierda, por eso la letra W y abajo o al sur, por eso la letra S. • -V nos informa por pantalla que est´a haciendo el programa. Finalmente redirigimos la salida del programa con el s´ımbolo >hacia un archivo, cuyo nombre es ejemplo 1a.ps.

2.2.

Mejorando el sistema de referencia.

Si ahora queremos agregar marcas secundarias en los ejes, utilizamos en el modificador -B el par´ametro f, y para generar una malla o ret´ıcula utilizamos el par´ametro g. Una posible salida de GMT resultar´a como muestra la figura 4. El c´odigo de GMT ser´a: 1. gmtset PAGE_ORIENTATION=landscape MEASURE_UNIT=cm 2. psbasemap -R0/360/-1.25/1.25 -JX14 -Ba90f30g60/a0.5f0.25g0.25WS -V >ejemplo_1b.ps La u ´nica diferencia es en el modificador -B. Indicamos que las marcas secundarias (f) est´en cada 30 grados en la horizontal y cada 0.25 unidades en la vertical, y que la ret´ıcula (g) sea cada 60 grados y 0.25 unidades. Por u ´ltimo veamos como generar una escala horizontal en radianes y de paso “cerrar” la figura por derecha. La salida de GMT resultar´a en la figura 5. El c´odigo de GMT ser´a: 1. gmtset PAGE_ORIENTATION=landscape MEASURE_UNIT=cm 2. psbasemap -R0/360/-1.25/1.25 -JX14 -Ba90f30g60/a0.5f0.25g0.25WS -K -V >ejemplo_1c.ps 7

GMT

Referenciaci´ on

Figura 4: Una malla sobre los ejes 3. psbasemap -R0/6.28/-1.25/1.25 -J -Ba1.57/a0NE -O -V >> ejemplo_1c.ps Aparecen ahora algunos elementos muy interesantes que ejemplifican la forma de trabajo en GMT. Por un lado en la l´ınea 3 aparece el modificador -K con el cual indicamos que el gr´afico no debe cerrarse ya que luego vendr´a otro comando GMT que agregar´a elementos a la imagen. Por otro en la l´ınea siguiente aparece el modificador -O con el que indicamos que la salida de este comando debe superponerse a la informaci´on anterior. En esta misma l´ınea no redireccionamos la salida con un > que es un redireccionamiento destructivo, sino con un >> que es un redireccionamiento que adjunta la salida actual a las anteriores. Yendo al objetivo inicial, con la u ´ltima l´ınea se agregaron los ejes arriba (N) y a la derecha (E) en el modificador -B. En el mismo modificador indicamos que no hubiera anotaciones a la derecha a0 y que arriba las mismas deber´ıan tener un paso π/2. 8

GMT

Referenciaci´ on

Figura 5: Gr´afico con doble escala Algo a destacar es que este cambio en la escala, de grados sexagesimales a radianes, solo fue posible gracias a dos cosas: la redefinici´on del ´area de trabajo con -R, donde ahora el eje X va de 0 a 2π y el hecho que se utiliz´o la misma escala. Esto u ´ltimo est´a impl´ıcito en el modificador -J sin ning´ un par´ametro, lo que se˜ nala que siguen activos los par´ametros usados con anterioridad 1 . Creemos que este u ´ltimo ejemplo es, para quien da una primera lectura al texto, un ejemplo que debe repasar con cuidado y tratar de entender opci´on a opci´on, realizando una comparaci´on con el segundo ejemplo, marcando claramente las diferencias y entendiendo por qu´e ellas ocurren. Quiz´as usted ya est´e pensando en tirar la toalla luego de todo lo engorroso que ha sido apenas generar el sistema de referencia para un simple gr´afico, no 1

GMT tiene memoria, es decir toma los valores de los par´ametros que tiene cada modificador cuando estos han sido definidos con anterioridad y s´olo cambian cuando se los redefine

9

GMT

Referenciaci´ on

se desaliente, la mayor parte de lo aprendido se repite en todos los comandos y es la parte esencial de GMT. Dejamos ahora un conjunto de ejercicios para que pueda repasar lo aprendido 1. Modifique el primer ejemplo y genere un gr´afico que tenga la mitad de tama˜ no que el primero. 2. Agregue al nuevo gr´afico un cambio en la extensi´on de las variables de manera de graficar solo los valores positivos en el rango que va de los 0 a los 180 grados 3. Al gr´afico saliente del punto anterior sobreimponga otro sistema de coordenadas cuya escala sea un 25 % menor y cuya regi´on sea igual a la del primer gr´afico de esta secci´on [0, 360] y [-1.25, 1.25]. Este grafico sobreimpreso debe estar corrido 2 cent´ımetros hacia arriba y hacia la derecha. Ayuda: utilice los par´ametros -X y -Y para correr un gr´afico sobreimpuesto a uno anterior. Tenga en cuenta que el corrimiento es relativo al primer gr´afico y que el origen est´a en la esquina inferior izquierda. La salida que debe poder generar debe ser como la de la figura 6. ¿Qu´e se aprendi´o hasta ahora? 1. C´omo ajustar variables al entorno GMT. 2. C´omo definir la regi´on de un mapa. 3. C´omo definir la proyecci´on de un mapa. 4. C´omo poner marcas en los ejes. 5. La manera de trabajar sobreimponiendo una salida de GMT a otra, opciones K, O y redireccionamiento >.

3.

Poniendo los datos en el gr´ afico.

Vayamos pues por la opci´on de graficar la funci´on trigonom´etrica seno. Obtenemos el gr´afico de la figura 7. 10

GMT

Referenciaci´ on

Figura 6: La soluci´on al ejercicio planteado Nuestro c´odigo para esto habr´a sido 1. gmtset PAGE_ORIENTATION=landscape MEASURE_UNIT=cm 2. psbasemap -R0/360/-1.25/1.25 -JX14 -Ba90f30g60/a0.5f0.25g0.25WS -K -V >ejemplo_2a.ps 3. psbasemap -R0/6.28/-1.25/1.25 -J -Ba1.57/a0NE -O -K -V >> ejemplo_2a.ps 4. psxy seno.prn -R0/360/-1.25/1.25 -J -O -V >> ejemplo_2a.ps Ahora en la u ´ltima l´ınea ejecutamos el comando psxy, el mismo toma del archivo seno.prn los datos en formato de tabla ASCII a dos columnas y grafica los puntos X vs. Y uni´endolos con una recta. Observemos que en la pen´ ultima l´ınea tuvimos que agregar el modificador -K ya que no es la u ´ltima l´ınea y que tuvimos que volver a cambiar la regi´on dado que nuestra tabla tiene grados sexagesimales en la columna que corresponde a las abscisas. 11

GMT

Referenciaci´ on

Figura 7: Gr´afico de una funci´on trigonom´etrica Como sucede siempre la u ´ltima l´ınea que produce salida gr´afica no tiene -K y s´ı tiene -O.

4.

Cambiando la representaci´ on de los datos.

Si queremos aumentar el grosor, el color o, el tipo de del trazo usamos el modificador -W. As´ı nuestro c´odigo queda: 1. gmtset PAGE_ORIENTATION=landscape MEASURE_UNIT=cm 2. psbasemap -R0/360/-1.25/1.25 -JX14 -Ba90f30g60/a0.5f0.25g0.25WS -K -V >ejemplo_2b.ps 3. psbasemap -R0/6.28/-1.25/1.25 -J -Ba1.57/a0NE -O -K -V >> ejemplo_2b.ps 4. psxy seno.prn -R0/360/-1.25/1.25 -W4 -J -O -V >> ejemplo_2b.ps 12

GMT

Referenciaci´ on

Donde el grosor de trazo est´a dado en puntos. La salida gr´afica se ver´a como la figura 8.

2

Figura 8: Aumentando el grosor de la l´ınea. Otro cambio que podemos introducir usando el modificador -W es cambiar el tipo de trazo: en lugar de una l´ınea continua graficamos una l´ınea con rayas (ver figura 9). Esto lo logramos agregando al modificador -W el par´ametro (signo menos) pero antes del mismo hay que elegir el grosor del trazo y luego separar con , (coma) los par´ametros. As´ı nuestro programa para un trazo de 5 puntos de grosor y con rayas queda: 1. gmtset PAGE_ORIENTATION=landscape MEASURE_UNIT=cm 2. psbasemap -R0/360/-1.25/1.25 -JX14 -Ba90f30g60/a0.5f0.25g0.25WS -K -V >ejemplo_2c.ps 3. psbasemap -R0/6.28/-1.25/1.25 -J -Ba1.57/a0NE -O -K -V >> ejemplo_2c.ps 2

En tipograf´ıa el punto es la menor unidad de medida. Hoy en d´ıa, el punto PostScript corresponde a 127/360 mm = 0,3527 mm

13

GMT

Referenciaci´ on

Figura 9: Trazo punteado 4. psxy seno.prn -R0/360/-1.25/1.25 -J -W5,- -O -V >> ejemplo_2c.ps Por u ´ltimo veamos como poner un s´ımbolo en cada dato (modificador -S en psxy) y como poner un texto explicativo al gr´afico. El resultado queda as´ı: Y el c´odigo del programa es: 1. gmtset PAGE_ORIENTATION=landscape MEASURE_UNIT=cm 2. psbasemap -R0/360/-1.25/1.25 -JX14 -Ba90f30g60/a0.5f0.25g0.25WS -K -V >ejemplo_2d.ps 3. psbasemap -R0/6.28/-1.25/1.25 -J -Ba1.57/a0NE -O -K -V >> ejemplo_2d.ps 4. psxy seno.prn -R0/360/-1.25/1.25 -Sc0.02 -J -K -O -V >> ejemplo_2d.ps 5. pstext titulo.txt -R -J -X7 -Y9 -N -O -V >> ejemplo_2d.ps

14

GMT

Referenciaci´ on

El modificador -Sc0.02 nos dice que el s´ımbolo es un c´ırculo y que las dimensiones del mismo son 0.02 puntos. Con el comando pstext titulo.txt escribimos sobre la misma regi´on y para la misma proyecci´on un texto. Al usar los modificadores -X7 -Y9 desplazamos el origen del gr´afico 7 cm en la direcci´on X y, 9 cm en la direcci´on Y. El cero esta en 2,54 cm, 2,54 cm as´ı que con este desplazamiento, ahora se ubicar´a en 9,54 cm, 11,54 cm. Cuando usamos el modificador -N le indicamos al comando pstext que NO restrinja el texto al lugar que ocupa el gr´afico, o sea, puede imprimir por fuera de los bordes del mismo. Para entender qu´e y c´omo se imprime hay que leer el archivo ASCII titulo.txt En este caso este tiene una sola l´ınea que indica: 0 0 15 0 1 CM Funci\360n trigonom\345trica seno Donde los n´ umeros indican de izquierda a derecha: 1. Posici´on del texto en X. 2. Posici´on del texto en Y. 3. Tama˜ no de letra. ´ 4. Angulo respecto de la horizontal. 5. Tipo de letra. En n´ umero 1 significa Helvetica. 6. Justificaci´on: (L)eft, (C)enter, (R)ight, (T)op, (M)edium, (B)ottom Notar que los acentos se ponen siguiendo un n´ umero octal cuya tabla de conversi´on est´a en el manual de GMT. Tarea 1. Arme una tabla que muestre los valores num´ericos de un a´ngulo, en radianes, versus la funci´on coseno. Los l´ımites del ´angulo ser´an [0 pi/2] y el paso de muestreo de la tabla un d´ecimo de radi´an. 2. Grafique la tabla para los l´ımites o´ptimos uniendo los datos con una l´ınea de rayas y puntos (-.-.-.-. ) . Coloque una estrella cuando el coseno tome el valor π/4. El resultado al que debe llegar esta representado en la figura 10. 15

GMT

Referenciaci´ on

Figura 10: Soluci´on al segundo problema. ¿Qu´e se aprendi´o es esta secci´on? 1. Como graficar una tabla bidimensional con l´ıneas y s´ımbolos. 2. Como agregar texto a un gr´afico.

5.

Empezamos a graficar mapas

Comenzamos ahora a realizar una tarea m´as espec´ıfica y por la que GMT es frecuentemente usado: la representaci´on cartogr´afica. En el primer ejemplo escribimos pscoast -R-180/180/-90/90 -JI24 -W -Bg30/g30 >ejemplo_3a.ps Y obtenemos un mapa de toda la Tierra (ver figura 11) utilizando la proyecci´on sinusoidal (JI). La opci´on W indica que se debe graficar el borde entre el agua y la tierra. Como la proyecci´on deforma mucho a la Tierra podemos partirla y generar todo el globo a partir de tres secciones del mismo (ver figura 12). 16

GMT

Referenciaci´ on

Figura 11: Mapa de la Tierra. 1. gmtset PAGE_ORIENTATION=landscape MEASURE_UNIT=cm 2. pscoast -R200/340/-90/90 -Ji0.07c -Bg30/g15 -W -K >ejemplo_3b.ps 3. pscoast -R-20/60/-90/90 -J -B -W -K -O -X9.85 >> ejemplo_3b.ps 4. pscoast -R60/200/-90/90 -J -B -W -O -X5.6 >> ejemplo_3b.ps Cambiemos la regi´on y solo grafiquemos la zona de Tierra del Fuego como se ve en la figura 13 pscoast -R-76/-63/-57/-51 -Jn-71/2 -Ba5 -W -Dc -V >ejemplo_3c.ps ¡Caramba, el mapa ha sido dibujado por un pintor cubista!. No, es que hemos elegido la resoluci´on de la costa (-D), que antes la tomaba por omisi´on, y la elecci´on ha sido la resoluci´on “cruda” (-Dc). En la figura 14 podemos ver que hay otras resoluciones disponibles: l (low), i (intermedia) y h (high). Tambi´en existe otra que es f (full) y es la m´as detallada. Cuanto mejor es la resoluci´on m´as pesado es el gr´afico, por lo que no siempre la m´as densa es la ideal. 17

GMT

Referenciaci´ on

Figura 12: Mapa con discontinuidades. 1. gmtset PAGE_ORIENTATION=landscape 2. gmtset MEASURE_UNIT=cm 3. pscoast -R-76/-63/-57/-51 -JI-71/11 -Ba5 -Dc -W -Y12 -K -V >ejemplo_3d.ps 4. pscoast -R-76/-63/-57/-51 -JI-71/11 -Ba5 -Dl -W -X14 -K -O V >> ejemplo_3d.ps 5. pscoast -R-76/-63/-57/-51 -JI-71/11 -Ba5 -Di -W -X-14 -Y-11 -K -O -V >> ejemplo_3d.ps 6. pscoast -R-76/-63/-57/-51 -JI-71/11 -Ba5 -Dh -W -X14 -O -V >> ejemplo_3d.ps Finalizamos poniendo en Tierra del Fuego el nombre de dos ciudades que se indican con un s´ımbolo (ver figura 15).

18

GMT

Referenciaci´ on

Figura 13: Tierra del Fuego. 1. pscoast -R-76/-63/-57/-51 -Jn-71/2 -Ba5 -W -Di -Y4 -K -V >ejemplo_3e.ps 2. psxy ciudades.txt -R-76/-63/-57/-51 -Jn-71/2 -: -Sn.4c0/0/0 -Gblack -K -O -V >> ejemplo_3e.ps 3. pstext leyendas_ciudades.txt -R-76/-63/-57/-51 -Jn-71/2 -: D0/0 -G0/0/0 -O -V >> ejemplo_3e.ps Notar que el s´ımbolo fue rellenado de negro (-Gblack). Le preguntamos ahora: ¿Qu´e proyecci´on y que s´ımbolo se ha elegido? Tarea: Arme un mapa de toda la Argentina en donde aparezcan los l´ımites internos entre provincias, adem´as de sus r´ıos y coloque para alguna provincia tanto un s´ımbolo como un cartel indicando la localizaci´on de su capital. Ayuda: preste atenci´on a los modificadores -N y -I para llevar adelante esta tarea. ¿Qu´e se aprendi´o es esta secci´on? 19

GMT

Referenciaci´ on

Figura 14: Diferentes resoluciones de costas. 1. C´omo graficar un mapa. 2. C´omo cambiar la resoluci´on de las costas. 3. C´omo cambiar l´ımites interprovinciales o nacionales. 4. C´omo indicar si graficar o no las masas de aguas internas.

6.

Mallas de datos

En las Ciencias de la Tierra es muy com´ un trabajar con gr´aficos llamados de dos dimensiones y media (2D 1/2). En ellos se representa una variable dependiente en funci´on de dos variables espaciales independientes, por ejemplo: latitud y la longitud o X e Y. Una manera habitual de hacer esto es mediante un gr´afico de curvas de nivel donde cada curva une puntos de igual valor de la funci´on f(lon,lat) o f(x,y) representada.

20

GMT

Referenciaci´ on

Figura 15: Simbolos en un mapa. El primer paso que suele requerir esta representaci´on es la elaboraci´on de una malla regular de valores estimados a partir de los datos observados que casi nunca son tomados en forma regular. Las dos maneras de generar mallas son: 1. Dando valores a los nodos de una ret´ıcula regular, opci´on llamada malla con registro por nodo o en ingl´es gridline registrated. Esta es la opci´on por omisi´on en GMT. 2. Dando valores a cada uno de los centros de los cuadril´ateros definido por cuatro nodos, opci´on conocida como malla con registro por pixel o e ingl´es pixel registrated. GMT tiene 5 m´etodos de interpolaci´on para generar mallas equiespaciadas: 1. Vecino m´as pr´oximo (nearneighbor) 2. Triangulaci´on de Delaunay (triangulate) 21

GMT

Referenciaci´ on

3. Curvatura continua con m´ınima tensi´on (surface) 4. Mallado en tensi´on de datos esf´ericos (sphinterpolate) 5. Splines de Green (greenspline) Los tres primeros m´etodos realizan la interpolaci´on sobre una superficie plana, el cuarto la realiza sobre la superficie esf´erica y el quinto la puede realizar en ambas superficies. Algunos de los par´ametros son comunes a los cuatro m´etodos: 1. La regi´on en la cual se interpola. Mediante el modificador R. 2. El nombre de la grilla de salida. Mediante el modificador G donde el par´ametro es el nombre del archivo de salida. Por ejemplo: -Gnombre del archivo 3. El intervalo entre los puntos o la cantidad de puntos donde se interpola. Mediante el modificador I donde el par´ametro indica el espaciamiento en minutos de arco (m) o en kil´ometros(k) o metros (e) o millas (i) o millas n´auticas (n). Por ejemplo: -I4/2k, que indica una malla con puntos cada 4 kil´ometros en el eje x y cada 2 kil´ometros en el eje y. Si el espaciamiento no depende de la direcci´on se indica con un solo n´ umero. As´ı -I3n indica que la ret´ıcula la forman cuadrados de 3 millas n´auticas de lado. Es importante entender que no existe un m´etodo de interpolaci´on mejor que otro sino que cada uno tiene ventajas y desventajas y la elecci´on del mismo debe realizarse caso por caso. No es el objetivo de este texto discutir c´omo seleccionar un m´etodo de interpolaci´on por lo que simplemente diremos como hacerlo desde el punto de vista operativo.

7.

Curvas de nivel

Armaremos un ejemplo en base a las cotas topogr´aficas de Tierra del Fuego tomadas en tierra firme, generando curvas de nivel con varios de los m´etodos y analizando sus modificadores y par´ametros. Las curvas pueden realizarse solo luego de haber generado la malla regular a trav´es del comando grdcontour. El mismo se puede ejecutar as´ı: 22

GMT

Referenciaci´ on

1. gmtset BASEMAP_TYPE=plain ANNOT_FONT_SIZE_PRIMARY=7p 2. psbasemap -R-68.75/-65.0/-55.1/-52.6 -JM-66.875/-53.85/7.5c -B1 -P -Y5 -K -V >salida.ps 3. pscoast -R -J -Df -W -Ia/1/0/0/0 -N1/1/0/0/0 -O -K -V >> salida.ps 4. grdcontour grilla_t.grd -C100 -A200+s4.5p -Gd5c -R -J -W1,0/0/0,-O -V >> salida.ps Analicemos el c´odigo: 1. El modificador -C acompa˜ nado del par´ametro 100 nos dice que hay que dibujar las l´ıneas de igual altura cada 100 metros 2. El modificador -A acompa˜ nado del par´ametro 200 nos dice que hay que etiquetar las l´ıneas de igual altura cada 200 metros. s4.5p nos dice que grosor deben del texto deben tener 4.5 y p indica la unidad: puntos 3. El modificador -Gd5c controla la ubicaci´on de las etiquetas, no permitiendo que existan etiquetas a menos de 5 cent´ımetros. Miremos ahora las salidas de las diferentes grillas y la forma de hacerlas

7.1.

Vecino m´ as pr´ oximo

Ejecutamos la sentencia: nearneighbor datos.dat -R-68.75/-65.0/55.1/-52.6 -Ggrilla_n.grd -S0.5K -N4/2 -I3c -: -V Donde: 1. datos.dat es una tabla, que no tiene porque estar equiespaciada, de tres columnas: longitud, latitud, altura sobre el elipsoide. 2. -S indica el tama˜ nos del c´ırculo donde busca datos y el par´ametro K (may´ uscula hace que la distancia sea sobre la superficie elips´oidica o esf´erica y no en el plano)

23

GMT

Referenciaci´ on

3. -N indica la cantidad de sectores en los que debe dividir en c´ırculo y luego viene la cantidad de valores totales que debe haber. Si no se cumple la condici´on aparece un NaN (Not a Number) El resultado puede verse en la figura 16. Como se aprecia las curvas de nivel quedan segmentadas. Notar que a pesar de que no se tomaron datos en el mar las curvas de nivel aparece sobre ´el en las poximidades de la costa. El par´ametro -S funciona como una m´ascara natural para que no existan valores de topograf´ıa muy lejos de la costa.

Figura 16: Curvas de nivel utilizando nearneighbor.

7.2.

Triangulaci´ on de Delaunay

Ejecutamos la sentencia: 24

GMT

Referenciaci´ on

triangulate datos.dat -R-68.75/-65.0/-55.1/-52.6 -Ggrilla_t.grd -I3c -: -V El resultado puede verse en la figura 17. Como se aprecia las curvas de nivel quedan menos segmentadas pero son tramos de rectas. Una vez m´as existen curvas de nivel sobre el mar en las poximidades de la costa, algo que no deber´ıa existir.

Figura 17: Curvas de nivel utilizando triangulaci´on de Delaunay.

7.3.

Curvatura continua con m´ınima tensi´ on

Ejecutamos la sentencia: surface datos.dat -R-68.75/-65.0/-55.1/52.6 -Ggrilla_s.grd -T0.25 -C0.1 -I3c -: -V

25

GMT

Referenciaci´ on

El resultado puede verse en la figura 18. Como se aprecia las curvas de nivel quedan suaves y sin v´ertices. Donde 1. -T indica la tensi´on, el valor de este par´ametro suele indicarse entre 0.25 para datos asociados a campos potenciales y 0.35 para datos de topograf´ıa rugosa. 2. -C indica un l´ımite de convergencia debido a que el algoritmo es iterativo, el valor de este par´ametro suele ser 0.1.

Figura 18: Curvas de nivel utilizando curvatura continua con m´ınima tensi´on. Pero como se ve surface no tiene reparos en estimar la topograf´ıa sobre el mar, a´ un a distancias muy grandes de la costa donde es claro que no existen 26

GMT

Referenciaci´ on

datos. Este hecho muestra lo peligroso que es ajustar algunas funciones que pueden estimar valores sin un solo dato remotamente cercano.

7.4.

M´ ascaras y cortes

GMT tiene varias maneras eliminar la representaci´on de datos en una regi´on dada. Una posibilidad, compleja y flexible, es la de colocar una m´ascara sobre una regi´on para que sobre ella no queden representados datos utilizando grdmask. Otra, m´as sencilla, es recortando del gr´afico los datos sobre el mar o la tierra. El siguiente c´odigo recorta los datos sobre el oc´eano. Para ellos en la l´ınea 4 est´a el comando pscoast con el modificador -S y el par´ametro c. Luego ponemos la malla de datos y finalmente cerramos con otro pscoast y el modificador -Q. 1. gmtset BASEMAP_TYPE=plain ANNOT_FONT_SIZE_PRIMARY=7p 2. psbasemap -R-68.75/-65.0/-55.1/-52.6 -JM-66.875/-53.85/7.5c -B1 -P -Y5 -K -V >salida.ps 3. pscoast -R -J -Df -W -Ia/1/0/0/0 -N1/1/0/0/0 -O -K -V >> salida.ps 4. pscoast -R -J -Df -Sc -O -K -V >> salida.ps 5. grdcontour grilla_s.grd -C100 -A200+s4.5p -Gd5c -R -J -W1,0/0/0,-O -V >> salida.ps 6. pscoast -Q -O -K -V >> salida.ps El resultado del recorte se visualiza en la figura 19.

8.

Pongamos un poco de color

Una alternativa a la representaci´on utilizando curvas de nivel es el coloreado de las superficies entre isol´ıneas con diferentes colores. Llamaremos a este tipo de gr´aficos como curvas de nivel rellenas. 27

GMT

Referenciaci´ on

Figura 19: Curvas de nivel con recorte sobre el oc´eano.

8.1.

Necesitamos una paleta

Para armar estos gr´aficos habr´a que definir una paleta de colores asociada a una escala num´erica de la variable dependiente. Podremos proponer que la variaci´on de la paleta sea continua o discreta. El tipo de paleta es una elecci´on que en general se realiza respetando las convenciones hist´oricas que existen; as´ı por ejemplo hay paletas est´andares para representar la topograf´ıa, la concentraci´on de ciertos minerales, etc. GMT trae predefinida una buena cantidad de paletas (ver figura 20), para adaptarlas a la escala a nuestros datos utilizamos el comando makecpt. Por ejemplo la sentencia: makecpt -Cpaleta_predefinida.cpt -T0/500/50 -Z >paleta_propia.cpt Toma a una paleta que GMT ya posee (modificador -C ) y la adapta al rango que va de 0 a 500 con paso 50. El modificador -Z hace que los cambios en los colores sean suaves o “continuos”. Para ver otras posibilidades de paletas puede ir al sitio: cpt-city

8.2.

Armando un mapa con colores

Dada la paleta y la malla de datos mediante el comando grdimage se colorea el mapa seg´ un el valor de la variable dependiente. El resultado podr´a ser algo como lo que muestra la figura 21. Y el c´odigo podr´a tener una forma como esta: 28

GMT

Referenciaci´ on

Figura 20: Paletas predefinidas en GMT. 1. psbasemap -R-68.75/-65.0/-55/-52.6 -J M-66.875/-53.85/6.5c Ba1f30m -K -P -V >salida.ps 2. pscoast -R -J -Df -Gc -O -K -V >> salida.ps 3. grdimage malla_mia.grd -R -J -Cpaleta_propia.cpt -Sb -E300 O -K -V >> salida.ps 4. pscoast -Q -O -K -V >> salida.ps 5. psscale -D 7.5c/3.1c/5.5c/0.25c -B250 -Cpaleta_propia.cpt O -K -V >> salida.ps 6. pscoast -R -J -Df -S 208/229/242 -W 0.25p,black -N1/0.25p,black,.-O -V >> salida.ps Notar que en la l´ınea 5 aparece el comando psscale que se ocupa de poner la escala de colores y su equivalente num´erico.

8.3.

La iluminaci´ on para marcar el relieve

Cuando se ilumina un mapa de colores con luz proveniente de una direcci´on particular se genera un efecto de relieve que puede resaltar alg´ un aspecto particular de la figura. Para esto hace falta generar primero una grilla de iluminaci´on con el comando grdgradient, Por ejemplo 29

GMT

Referenciaci´ on

Figura 21: Mapa topogr´afico 2D 1/2 coloreado grdgradient mi_malla.grd -Ggradiente_mi_malla.grd -E45/75/0.55/0.6/0.4/10 -Ne0.6 -V Genera una malla de iluminaci´on llamada gradiente mi malla.grd Luego tomando el gr´afico de la figura 21 e ilumin´andolo queda la figura 22. EL c´odigo correspondiente es este 1. psbasemap -R-68.75/-65.0/-55/-52.6 -J M-66.875/-53.85/6.5c Ba1f30m -K -P -V >salida.ps 2. pscoast -R -J -Df -Gc -O -K -V >> salida.ps 3. grdimage datos.grd -R -J -Cmi_paleta.cpt -Igradiente.grd -Sb 30

GMT

Referenciaci´ on

Figura 22: Mapa topogr´afico 2D 1/2 con relieve -E300 -O -K -V >> salida.ps 4. pscoast -Q -O -K -V >> salida.ps 5. pscoast -R -J -Df -S 208/229/242 -W 0.25p,black -N1/0.25p,black,.-O -V >> salida.ps Si se observa en la figura 23 aparecen representados en una escala de grises el mismo juego de datos a derecha e izquierda pero con diferentes iluminaciones. A la derecha de la figura se representa el relieve mediante una iluminaci´on que supone que la luz viene del a´ngulo superior izquierdo, esta es la iluminaci´on cl´asica para dar relieve a las figuras. A la izquierda de la figura se representa el relieve mediante una iluminaci´on inusual donde la luz 31

GMT

Referenciaci´ on

viene del a´ngulo inferior derecho.

Figura 23: Diferentes iluminaciones

32