Spanish (pdf)
Article in xml format
Article references
Automatic translation
Send this article by e-mail
Cited by SciELO 
Cited by Google 
Similars in
SciELO 
Similars in Google 
Permalink
INTRODUCCIÓN
En la Cuenca de Urabá han sido realizados diferentes estudios geológicos y geofísicos de tipo regional (Duque, 1980); algunos hacen referencia a la acreción de rocas Mesozoicas de los Andes (Duque, 1980; Meissnar et al., 1976; Feininger y Bristow, 1980). Actualmente, se tienen dos bloques adyacentes (Bloque Panamá Chocó y Bloque norte de los Andes), divididos por un conjunto de fallas que pueden ser la fuente sísmica del área (Nivia, 1989). El conocimiento de la geología del área aún es deficiente, y la falta de información limita la unicidad en la solución de los algoritmos de inversión, independientemente del método usado (Tarantola, 2005). Por lo tanto, establecer un modelo (parámetros) que se ajuste lo mejor posible a los datos medidos en superficie se convierte en una tarea cuya solución está en las condiciones o parámetros iniciales que se asignen en el modelo, debido a la información geológica que se tiene. Este trabajo es parte de la investigación, donde se desarrolla una metodología para encontrar la estructura geológica del subsuelo para la Cuenca de Urabá, a partir de la distribución de velocidades y la aplicación de la inversión de los eventos sísmicos naturales de la región (González, 2012). En este tipo de análisis, la relación costo/beneficio es baja, para profundidades considerablemente grandes, al ser comparado con otras herramientas o técnicas. Así, se hace posible determinar la distribución de velocidades, la cual está asociada a la litología del subsuelo. Los tiempos de arribo ofrecen información acerca de los parámetros hipocentrales y la distribución de velocidades del subsuelo (Suter et al., 2008). Asumiendo un tipo de dependencia funcional entre los datos observados, el hipocentro y la distribución de velocidades puede generarse un modelo que represente la estructura interna del subsuelo, y sea capaz de reproducir los tiempos de arribo, observados en cada estación sismológica ubicada en la región.
GENERALIDADES DEL ÁREA DE ESTUDIO
La Cuenca de Urabá se localiza en el extremo noroccidental, en la Costa Atlántica de Colombia, en el límite con Panamá, entre los departamentos de Antioquia, Córdoba y Chocó (FIGURA 1). Entre las cuencas de Urabá y Sinú - San Jacinto se dispuso de la Red Sismológica de Urabá (RSU), donde se adquirió la información en 22 estaciones (de 30 instaladas), con mejor relación señal-ruido. Las estaciones comprenden un sensor triaxial de banda ancha y digitalizadores CME, Quanterra y K2. Para complementar la información alrededor de la región se utilizaron 4 estaciones de la Red Sismológica Nacional de Colombia (RSNC) (TABLA 1).
FIGURA 1 Localización del área de estudio (Cuenca de Urabá - cuadro inferior: rojo -) y disposición de las estaciones en la Red Sismológica RSN (puntos en color naranja).
TABLA 1 Estaciones sismológicas utilizadas para el estudio, información obtenida de RSNC (Iyer y Hirahara, 1993; Hernández, 2009; RSNC, 2009a).
* Estaciones de la RSNC
Estaciones RSU
El comportamiento de las estaciones enmarca tres periodos de operación, donde el número de estaciones aumenta con cada periodo (TABLA 2 ). Los datos durante la adquisición presentan un retardo, resultado de la no sincronía entre las estaciones (FIGURA 2). Este error es debido a la diferencia entre los tiempos de arribo de la onda P para un mismo sismo, que excede la diferencia esperada, teniendo en cuenta la separación entre estaciones de la RSU. Los tiempos de arribo de la onda P para estaciones 1 y 2, son t P1 y t P2 , respectivamente.
TABLA 2 Funcionamiento de las estaciones RSU durante el tiempo de adquisición. Abreviatura de funcionamiento durante la semana: c = completa, p = parcial. -= menos de un día (modificado e RSNC 2009a, RSNC 2009b).
METODOLOGÍA
La metodología se basa en el registro de la perturbación de la superficie en función del tiempo. Al ocurrir un sismo, la perturbación se propaga desde la fuente sísmica hasta la estación receptora. Si la diferencia de amplitud del evento respecto al ruido es alta, se puede decidir con una incertidumbre relativamente baja, el tiempo de arribo de la onda P y de la onda S.
La información de los tiempos de arribo de las ondas P y S, por cada estación de la red sismológica ubicada en la cuenca, permiten encontrar hipocentros (localizaciones de los sismos), tiempos de origen y distribución de velocidades que validen la información obtenida, y está relacionada con la geología de la región (FIGURA 3). Para ello han sido adquiridas señales de amplitud vs tiempo, de manera continua por 9 meses. En la red sismológica de Urabá se determinaron los mismos eventos registrados en la Red Sismológica Nacional de Colombia (RSNC, 2009a, 2009b), determinando el arribo de las ondas P y S para cada estación, donde el arribo es claro para ser recolectado. Después se hizo una lista de los sismos con el nombre de las estaciones y su correspondiente tiempo para las ondas P y S (diferencia del tiempo de arribo entre las ondas S y P). Con la información de los tiempos de arribo y ubicación de las estaciones, se utilizó el método de inversión simultánea para diferentes casos: con un solo bloque homogéneo, con variación en profundidad de la velocidad (1D) y con variación de la velocidad en profundidad y lateralidad (2.5D) (Anders et al., 2016), (Castillo et al., 2010).
FIGURA 3 Localización epicentral de sismos registrados y localizados por la RSNC entre 1993-2006. Las estrellas representan sismos localizados por USGS en el año 1992-1994. Las líneas negras representan las fallas principales; 1 Murindó, 2 Murrí, 3 Urrao, 4 San Ruperto, 5 Abriaqui, 6 Cañasgordas, 7 Mutatá (Franco et al., 2006).
Durante el año 2009 (entre enero hasta septiembre) fueron detectados 75 sismos provenientes de esta región, en el cuadrante limitado por las coordenadas 7º a 9º en dirección sur - norte, y -78º a -76º en dirección oeste - este; generalmente con profundidades menores de 30 km, como se muestra en los boletines de sismicidad de la Red Sismológica Nacional de Colombia (Iyer y Hirahara, 1993). Con el fin de establecer la disposición y rasgos geológicos del subsuelo en la Cuenca de Urabá, se propuso un modelo de velocidades inicial (Franco et al., 2006), (Xiong et al., 2013), (Lou et al., 2016), mediante el registro de tiempo del primer arribo de las ondas P y S de los sismos de la región. Para ello, se tienen las siguientes consideraciones:
Depurar los registros sísmicos para obtener información con la mínima incertidumbre posible.
Selección de un modelo geológico inicial, representativo de la región, a partir de la información superficial y datos disponibles de la zona.
Solucionar la ecuación de onda, y encontrar los tiempos de arribo simulados de las ondas P y S, asumiendo el modelo de velocidad y localizaciones hipocentrales como parte del modelo final a encontrar. Se consideró el medio conformado por bloques y caracterizados únicamente por su velocidad, además se estableció que cada bloque se comporta como un sólido de Poisson, homogéneo e isotrópico. La simulación por medio de la solución de la ecuación de onda permite realizar análisis espectral y de amplitudes de los eventos sísmicos, sin embargo, en el presente trabajo solo se tiene en cuenta los tiempos de arribo de las ondas símicas.
Aplicar un algoritmo de inversión para encontrar hipocentros y el modelo de velocidades óptimo que ajuste los tiempos simulados de arribo de las ondas P y S, y sea plausible con la geología de la región. Para evitar fuertes cambios de velocidad entre bloques adyacentes se implementó una técnica conocida como suavizado (Stein y Wysession, 2003).
Modelo final de velocidades pseudo 3D y localizaciones hipocentrales.
Finalmente se aclara que todos los algoritmos utilizados en el estudio fueron desarrollados y validados por los autores del presente trabajo.
Clasificación de datos
Debido al error de sincronía de las estaciones, no se tiene confianza en los tiempos de arribo de la onda P, por tratarse de un valor absoluto. Sin embargo, la diferencia entre tiempo de arribo de la onda S y P es confiable, por ser un valor relativo. Con las estaciones CAP, DBB, HEL, SOL, pertenecientes a la RSNC (TABLA 1) se obtienen datos de tiempos absolutos confiables y el tiempo de arribo de la onda P puede ser adquirido. Así, los tiempos observados se componen de:
Diferencias de tiempo de arribo TSP = TS - TP para las estaciones RSU y RSNC.
Tiempo de arribo TP para las estaciones RSNC que sean confiables para ser recolectados.
Con estas características se tiene información de 54 sismos. Los tiempos observados serán la base para encontrar el mejor modelo que los represente. El modelo final será encontrado a partir del más sencillo. Así se supone inicialmente una velocidad constante y los hipocentros ubicados en el mismo punto. Mediante el programa de inversión se encontrará el primer modelo aproximado. Este modelo ubicará los hipocentros en sus posiciones aproximadas, considerando las dimensiones espaciales de la región de estudio. Por lo general la velocidad tiene mayor variación en profundidad que lateralmente, entonces se propone un primer modelo en profundidad (2.5D) y restringiendo los cambios laterales, con el fin de obtener un nuevo acercamiento al modelo real (q3D) (Bahrani y Adler, 2012; Pavlis y Mason, 2017).
Inversión
El programa de inversión se divide en dos etapas: La primera corresponde al algoritmo directo y la segunda al algoritmo inverso. En el algoritmo directo se realizó una simulación en C++ de la ecuación de onda (Ecuación 1) utilizando el método de diferencias finitas:
La simulación permite obtener los tiempos simulados de arribo de la onda P, y con ese valor se calcula el tiempo de arribo de la onda S, para diferentes sismos. Estos tiempos simulados dependen a su vez de la localización (hipocentro) de los sismos y de la velocidad de la onda P en el medio.
En el algoritmo inverso se comparan los resultados que se obtienen de la simulación (tiempos simulados) con los datos medidos en superficie. Así, se puede establecer la diferencia entre lo simulado y lo medido en campo. La diferencia modifica los valores iniciales del modelo, es decir, la velocidad de la onda P en el medio y las coordenadas de los hipocentros de los sismos. Los nuevos valores iniciales simulan los tiempos de arribo de la onda P para cada uno de los sismos, y se comparan los nuevos tiempos simulados con los tiempos registrados en superficie, obteniendo una nueva diferencia o error, la cual permitirá recalcular los nuevos parámetros iniciales para la siguiente simulación. Esta metodología se realiza iterativamente hasta que la diferencia entre los datos simulados y los registrados en superficie sean lo más parecidos posibles, es decir, que el error cuadrático medio sea el más pequeño. Para invertir simultáneamente la estructura de velocidad, los hipocentros y el tiempo de origen (Ramadan, 2016), se trabajó el problema en una dimensión, e independientemente se ha publicado el método de inversión en tres dimensiones (Crosson, 1976). De esta manera se generalizó el método de inversión de hipocentros a inversión simultánea de velocidades e hipocentros (Franco et al., 2006).
En el modelo se suponen N sismos cuyos hipocentros están dados por longitud x j , latitud y j y profundidad z j y sus tiempos de arribo t j con j = {1,2,3,…,N}. Además el espacio se divide en L bloques y cada bloque tiene velocidad v l , donde l={1,2,3,…,L}. Estas variables conforman el modelo m 0:
siendo 4N+L, el número de elementos.
Para encontrar m 0 , se tienen M estaciones (k=1,2,…M), con la posición (x k , y k , z k , t k ), y el tiempo de arribo observado en la k-ésima posición por el j-ésimo sismo. Existe una relación funcional T(m 0 ) que permite reproducir los tiempos observados de acuerdo con el modelo m 0 . Suponiendo que esta relación es lineal, se define la matriz sensibilidad de la siguiente forma:
donde Aj es una matriz Mx4 dada por:
y C es una matriz NxMxL; la i-ésima fila de C está dada por:
con i(j,k)=k+(j-1)M, donde δjki, está dado por:
De esta, la relación lineal entre los tiempos observados y el modelado es:
Para hallar m de la ecuación 7 es necesario invertir G, pero esta matriz no es necesariamente cuadrada. Para determinar m se utiliza el método de mínimos cuadrados, el cual consiste en multiplicar primero por la transpuesta y luego invertir la matriz. Conforme a esto la ecuación 7, se transforma en:
donde ΔT, es la variación entre el tiempo observado (medido en superficie) y el simulado,
y el error entre ambos tiempos se calcula de la siguiente manera:
y Δm es la variación del modelo, entonces:
donde I es la matriz identidad, y ϵ se conoce como parámetro de amortiguación, el cual se obtiene empíricamente. Este ajuste se realiza para evitar problemas de estabilidad en el cálculo de la matriz inversa.
De tal manera que se propone un modelo inicial m 0 , con el cual se halla el tiempo teórico y la diferencia con el tiempo observado, es decir ΔT 0. Ahora se halla un nuevo modelo m 1 por medio de la ecuación 11, se itera hasta que el error sea mínimo. Entonces el problema se reduce a encontrar la matriz de sensibilidad G. Esta matriz contiene elementos que representan la variación del tiempo con cada uno de los elementos del modelo. Usando la aproximación de primer orden de la serie de Taylor, se tiene:
donde T j,k (m), se puede obtener solucionando la ecuación de onda por diferencias finitas.
RESULTADOS
Modelo con velocidad constante
Dado que las dimensiones del estudio no superan los 300 km de longitud, y su ubicación se encuentra en la zona oeste colombiana, se decide utilizar coordenadas planas Gauss - Krügercon, referencia en Datum Bogotá, origen Oeste. Sin embargo, con el ánimo de facilitar la lectura numérica de las coordenadas, se redefine el origen de coordenadas en latitud, de tal manera que el origen del presente estudio se sitúa a 3000 km dirección Norte con respecto al origen Oeste con Datum Bogotá. De esta manera los límites de la zona de estudio se encuentran en los siguientes intervalos:
El modelo inicial se considera de la manera más sencilla posible, con velocidad constante v p = 5 km/s y los 54 sismos ubicados en el centro de la región, de manera análoga, como se hizo el modelo inicial para validar el programa de inversión. Además, el valor inicial para el tiempo de origen es el menor de aquellos sismos que han sido registrados, recolectados y etiquetados como onda P por alguna estación de RSNC. Después de 70 iteraciones, el error se mantuvo constante (FIGURA 4), y los hipocentros fueron ubicados.
Se tiene una mayor densidad de sismos en la región parte más profunda tiene muy baja densidad de sismos, suroeste (FIGURA 5A), con tendencia dominante hacia entonces se descartó del estudio, de tal manera que la superficie (FIGURA 5B). A profundidad se divide la para los modelos siguientes se considera hasta 80 km zona en 5 partes, cada una con 19 km de espesor. La de profundidad.
Este modelo permitió estimar los hipocentros de los eventos sísmicos. En seguida se da libertad para que el modelo estime la mejor distribución de velocidades y recalcule los orígenes hipocentrales.
Modelo con variación en profundidad
A partir del modelo inicial con velocidad constante del medio, se aplica el algoritmo de inversión con libertad de seleccionar las mejores distribuciones de velocidades a profundidad y localizaciones hipocentrales. Después de 111 iteraciones (FIGURA 6 ), se obtuvo que la velocidad aumenta con la profundidad y un fuerte cambio a 20 km de profundidad (FIGURA 7 y TABLA 3 ).
Modelo Pseudo-tridimensional (q3D)
Partiendo del modelo de velocidades con variación en profundidad (2.5D) modelo inicial (TABLA 3) se aplica el algoritmo de inversión con libertad de variación tanto en profundidad como en lateralidad (q3D), compuesto por 4 x 4 x 4 bloques. A partir de las secciones en profundidad se pudieron determinar lineamientos hasta una profundidad aproximada de 20 km y en el centro de la zona de estudio (FIGURA 8). Estos corresponden a rasgos estructurales, especialmente lineamentos y fallas, debidos a la respuesta de esfuerzos y la litología, los cuales pueden asociase a un material de mayor velocidad, por lo tanto, mayor competencia y que puede considerarse como el basamento. Además, el perfil b-b’ (cercano a 8°N) permitió delimitar la Cuenca de Urabá, de tal manera que ésta posiblemente alcanza una profundidad máxima de 15 km en esta zona (FIGURAS 8 y 9). Se presenta únicamente el contraste de velocidad en el plano superficial dado que la densidad de sismos disminuye rápidamente en superficie y planos más
FIGURA 8 Localización de sismos y distribución de velocidades A. Proyección hacia la superficie (vista en planta). B. Proyección lateral sobre el perfil a-a’.
CONCLUSIONES
El algoritmo directo e inverso permite la inversión de hipocentros, tiempos de origen y distribución de velocidades con variación 3D. Para el algoritmo directo se consideró propagación de onda por diferencias finitas con fronteras no reflectivas.
Tiempos de llegada teóricos son hallados con errores menores a 1% para el cálculo de velocidad y profundidad de hipocentros para modelos sintéticos. En el algoritmo inverso se tuvo en cuenta la matriz de suavidad para el modelo de distribución de velocidad en 3D.
Mediante el modelo de velocidad (constante inicialmente) se obtiene que los sismos son producidos en la región suroeste, cerca de la Falla de Murindó, a una profundidad máxima de 80 km, con velocidad representativa de 5-92 km/s.
El modelo final, evidencia características geológicas como la Falla de Murindó, y otras fallas satélites (en los límites, forma, geometría y definición). Se proponen dos fallas de acomodación debido a la subsidencia de la Placa Caribe y el movimiento del Bloque ChocóPanamá.
Se sugiere realizar la solución de la ecuación de onda por diferencias finitas con fronteras no reflectivas, y error de truncamiento de grado 4 o mayor, con el fin de disminuir el error del tiempo de llegada teórico sin perjudicar el tiempo de cómputo.
El desarrollo e implementación de la metodología descrita en el presente trabajo puede utilizarse para estaciones no sincronizadas. De tal manera que datos tomados por una red portátil de estaciones sismológicas e independientes pueden ser procesados.
