Introducción

En el Caribe Sur de Costa Rica transcurre el límite entre la placa Caribe y la microplaca de Panamá, también denominado Cinturón Deformado del Norte de Panamá (CDNP, Adamek, Frohlich y Pennington, 1988). Este cinturón de fallas inversas tiene una notable continuidad a lo largo del Caribe costarricense, con un rumbo noroeste-sureste, desde Sixaola hasta la ciudad de Limón (Fig. 1). En esta región se han documentado al menos seis terremotos destructivos desde el siglo XVIII (i.e., 1746, 1798, 1822, 1916, 1953, 1991; Montero, 1989; Boschini, 1989; Camacho y Víquez, 1993; Boschini y Montero, 1994; Peraldo y Montero, 1994), siendo el terremoto de Limón del 22 de abril de 1991, con una magnitud momento (Mw) de 7,7 (e.g., Goes, Velasco, Schwartz y Lay, 1993), el último gran evento sísmico en afectar esta zona del istmo centroamericano (Fig. 1).

Fig. 1:
Marco tectónico del sureste de Costa Rica. Los elementos tectónicos rotulados son: Cinturón Deformado del Norte de Panamá (CDNP), Fosa Mesoamericana (FMA) y Zona de Fractura de Panamá (ZFP). La línea blanca punteada marca la cuenca de Limón Sur (Denyer y Alvarado, 2007) y las líneas anaranjadas la localización de los perfiles tomográficos mostrados en la Fig. 5. Los triángulos azules representan las estaciones sismológicas, correspondientes principalmente a la RSN en Costa Rica, la UPA y la red de Chiriquí en Panamá y al proyecto SFB574 de la Universidad de Kiel en el Pacífico de Costa Rica. Las estrellas amarillas denotan los epicentros de los sismos mencionados en el texto, incluyendo el terremoto de Limón de 1991 (Mw 7,7). Los cuadros blancos representan las localidades San José (SJ), Limón y Sixaola, y su tamaño es relativo a la cantidad de población. Los triángulos café muestran los volcanes del Cuaternario y las líneas grises los contornos del techo de la placa del Coco de acuerdo con Lücke y Arroyo (2015). La región gris oscuro representa la cordillera de Talamanca con alturas mayores a 500 m.s.n.m. y los números, la profundad batimétrica en las áreas oceánicas. La velocidad de convergencia de la placa del Coco proviene de DeMets et al. (1994).

El límite entre las placas Caribe y Panamá, en las cercanías de la ciudad de Limón, cambia de rumbo y se adentra en la parte central de Costa Rica como una zona de deformación ancha (~100 km) compuesta por fallas principalmente dextrales y sinestrales denominada Cinturón Deformado del Centro de Costa Rica (CDCCR; e.g., Marshall, Fisher y Gardner, 2000; Montero, 2001). Debajo de las placas Caribe y Panamá se subduce la placa Coco con una velocidad de 83-89 mm/año (DeMets, Gordon, Argus y Stein, 1994) en conjunto con montes submarinos en el Pacífico Central y el levantamiento del Coco en el sureste de Costa Rica (Fig. 1). Más hacia el sureste del país, se encuentra un punto triple donde se unen la Fosa Meosamericana (FMA) y la Zona de Fractura de Panamá (ZFP), esta última limitando las placas Coco y Nazca. Esta compleja interacción entre las cuatro placas tectónicas mencionadas, resulta en una sismicidad muy alta que hace que Costa Rica sea una zona ideal para los estudios estructurales a partir de la detección de sismos locales.

La Red Sismológica Nacional (RSN) es un programa de investigación inscrito en la Universidad de Costa Rica (UCR) que tiene sus raíces en una larga historia de colaboración entre la UCR y el Instituto Costarricense de Electricidad (ICE; Linkimer, Arroyo, Alvarado, Arroyo y Bakkar, 2018). El centro de registro de la RSN se ubica en la Sección de Sismología, Vulcanología y Exploración Geofísica de la Escuela Centroamericana de Geología de la UCR, hacia donde actualmente transmiten en tiempo real 168 estaciones sismológicas (Fig. 1), de las cuales 120 son administradas por la UCR y 48 por el ICE. Esta robusta red de instrumentos permite la detección de los sismos con gran precisión, facilitando el estudio de la estructura de la litósfera y de las fuentes sismo-generadoras mediante técnicas como la tomografía sísmica.

Pese a la alta sismicidad de Costa Rica y al progreso reciente en la detección de sismos a nivel nacional, la zona Caribe de Costa Rica aún carece de estudios sismológicos que aborden su estructura con detalle. Se han realizado varios trabajos previos de tomografías con sismos locales para toda la región de Costa Rica. Por ejemplo, Protti, Schwartz y Zandt (1996) y Husen, Quintero, Kissling y Hacker (2003) llevaron a cabo tomografías de todo el territorio nacional, empleando las bases de datos y estaciones sismológicas existentes hasta esos años, respectivamente. Syracuse et al. (2008) iluminaron tomográficamente Nicaragua y la mitad noroeste de Costa Rica usando redes temporales. Otras tomografías locales de detalle se han obtenido para áreas menores, usando sismos registrados por redes combinadas de estaciones a lo largo de la costa Pacífica y de fondo marino, tales como DeShon et al. (2006) para la península de Nicoya, y Arroyo et al. (2009) y Dinc et al. (2010) para el Pacífico Central.

En el contexto del 30 aniversario del terremoto de Limón de 1991, y en el marco de un proyecto de investigación (número 830-B8-500) inscrito en la UCR para la elaboración de un modelo tridimensional de velocidades de la toda región costarricense, presentamos en este artículo un estudio focalizado en la región Caribe Sur, en donde a la fecha no existe una imagen tomográfica detallada de la zona en donde se generó el terremoto de 1991. De esta forma, mostramos un modelo de velocidades de onda P en una dimensión que abarca toda la región sureste del país, dos secciones tomográficas que atraviesan la zona del terremoto de 1991, y además aspectos generales sobre la resolución de los cálculos efectuados. Nuestros resultados, además de mostrar por primera vez una gran anomalía de bajas velocidades sísmicas inclinadas hacia el suroeste bajo el Caribe Sur de Costa Rica, proveen un marco de referencia adicional para mejorar las localizaciones de los sismos en esta región del país e interpretar el origen de los grandes terremotos del Caribe costarricense.

Metodología

La Tomografía con Sismos Locales (TSL) es una técnica robusta que permite obtener localización de sismos de alta calidad y simultáneamente determinar la estructura de velocidades de una región de manera tridimensional. La TSL consiste en invertir tiempos de viajes de las ondas sísmicas para resolver la velocidad de ondas P (Vp) y S (Vs), o bien la relación Vp/Vs, y para determinar hipocentros. Los parámetros del modelo (i.e., tiempo de origen, latitud, longitud y profundidad de los sismos y estructura de Vp, Vs y Vp/Vs) se perturban con el fin de reducir con mínimos cuadrados las diferencias entre los tiempos de viaje observados en las estaciones sismológicas y los tiempos de viaje calculados usando un determinado modelo de referencia del interior de la Tierra (e.g., Eberhart-Phillips 1990, Kissling, 1988).

El Modelo Mínimo de Velocidades en una Dimensión (M1D) se define como el modelo de velocidades con corrección para estaciones que refleja más cercanamente la información obtenida por estudios previos en el área y a la vez resulta en el promedio mínimo de la media cuadrática (RMS, por sus siglas en inglés) de los residuos de tiempo de viaje para todos los sismos en la base de datos (Kissling, 1988; Kissling, Ellsworth, Eberhart-Phillips y Kradolfer, 1994 y Kissling, Kradolfer y Maurer, 1995). El M1D se emplea como modelo de referencia para la TSL y además, con él se relocaliza la base de datos para seleccionar el subconjunto de datos que se invierte en la TSL. Todo esto con el fin de poder linealizar el problema de interdependencia hipocentros-velocidades con velocidades e hipocentros cercanos a sus valores reales (Kissling et al., 1994).

Para derivar el M1D de la onda P se utilizó el programa VELEST (Kissling et al., 1995), que resuelve el problema de interdependencia hipocentros-velocidades a través de varias inversiones simultáneas de los parámetros hipocentrales, el modelo unidimensional de velocidades y las correcciones para las estaciones. Al final, se obtiene un modelo de capas paralelas de Vp con las correcciones correspondientes para las estaciones. Las correcciones compensan por la alta heterogeneidad y variaciones laterales en los primeros kilómetros de profundidad de la geología en toda la región estudiada.

Por su parte, para realizar la TSL de la onda P se usó la solución iterativa amortiguada de mínimos cuadrados codificada en el programa SIMULPS14 (Thurber, 1983; Evans, Eberhart-Phillips y Thurber, 1994; Haslinger, 1998). Este código resuelve el problema no lineal de interdependencia hipocentros-velocidades tratando los hipocentros como incógnitas. Cada iteración consiste de una inversión para las velocidades y para las localizaciones hipocentrales. SIMULPS14 aplica además el algoritmo trazador de rayos de disparo tridimensional de Virieux y Farra (1991).

En SIMULPS14, el modelo de velocidades se parametriza como una rejilla de nodos. En este estudio, se eligió un espaciado entre nodos de 20 x 20 km, considerando la distribución de sismos y estaciones y con el objetivo de lograr la mayor uniformidad posible en cobertura radial y puntos de intersección entre rayos. Las velocidades se mantuvieron fijas en los nodos que no fueron cruzados por ningún rayo y, además, para aquellos cuya suma ponderada de derivadas (derivative weighted sum, DWS) fue menor a 8.

El parámetro de amortiguamiento para la inversión fue elegido inspeccionado la curva de varianza de datos contra varianza del modelo después de una iteración, probando un rango de valores de amortiguamiento entre 1 y 10000 (Eberhart-Phillips, 1986). Se escogió un amortiguamiento de 300, por ser el valor que minimiza la varianza de datos sin incrementar demasiado la varianza del modelo.

Para evaluar la calidad de la solución de la TSL en este trabajo, se usaron una combinación de herramientas clásicas tales como conteo de rayos que cruzan una celda (hit count), DWS, función de esparcimiento (spread function) y los elementos diagonales de la matriz de resolución (e.g., Reyners, Eberhart-Phillips y Stuart, 1999), y pruebas sintéticas de tablero de ajedrez (Haslinger y Kissling, 2001). De acuerdo con Husen, Smith y Waite (2004), se diseñaron dos modelos de tablero de ajedrez, cada uno consistente en capas con anomalías alternas de alta y baja velocidad (±10%), con un nodo de la rejilla abierto para probar si hay fugas de velocidad horizontales; al mismo tiempo, las capas con anomalías se alternan con capas libres para corroborar la presencia de fugas verticales. Se repitió el ejercicio cambiando el patrón de tablero de ajedrez para las capas inicialmente libres, para explorar así las fugas en todas las capas. Las pruebas se llevaron a cabo con anomalías 2 x 2 nodos y con anomalías de un solo nodo.

Datos

La base de datos usada en este estudio proviene del catálogo sísmico de la RSN. Actualmente, la red sísmica de la RSN incluye 168 estaciones, de las cuales 38 son de banda ancha, 122 de periodo corto y 8 son acelerógrafos. En las localizaciones rutinarias se incorporan además lecturas de estaciones administradas por otras redes sismológicas ubicadas en Costa Rica, Panamá y Nicaragua, compartidas directamente o disponibles a través del consorcio IRIS (Incorporated Research Institutions for Seismology). La detección y la localización automática de los sismos se realiza rutinariamente usando SeisComP3 (Gempa, 2019). Posteriormente, se realiza la lectura manual de los sismos y la localización y el cálculo de la Mw usando los programas Mulplt y Hyp (Lienert y Havskov, 1995) integrados en el software sismológico SeisAn (Havskov, Voss y Ottemöller, 2020). Para esta localización de sismos se utiliza un modelo de Vp basado en Matumoto, Ohtake, Latham y Umaña (1977) y una Vp/Vs de 1,75.

De los sismos localizados rutinariamente contenidos en el catálogo de la RSN, se seleccionaron 1208 eventos entre 1998 y 2020, que conforman la base de datos inicial (Fig. 2). Estos eventos fueron localizados inicialmente por al menos 15 estaciones y presentan en su localización inicial un RMS menor a 1,0 s y un gap de cobertura de estaciones igual o menor a 180°. Se realizó una revisión de la calidad de las lecturas de las ondas sísmicas de todos los sismos seleccionados, aplicando un esquema de pesos con un valor de 0 para las incertidumbres mínimas (< 0,1 s) y de hasta 4 para lecturas muy inciertas (> 0,3 s) que no se pueden usar para modelar.

Fig. 2:
Distribución geográfica de los sismos preseleccionados para la inversión tomográfica. Esta base de datos consta de 1208 sismos (círculos rojos) registrados entre 1998 y 2020 por al menos 15 estaciones (triángulos) y localizados con un vacío de cobertura menor a 180°. Además, se muestran como referencia los sismos rotulados en la figura 1 (estrellas amarillas). A. Vista de planta mostrando los rayos sísmicos (líneas grises) entre los epicentros y las estaciones. Se muestra además la Fosa Mesoamericana (FMA) y el Cinturón Deformado del Norte de Panamá (CDNP). La estrella amarilla representa el terremoto de Limón de 1991 (Mw 7,7). B. Perfil en dirección norte-sur mostrando todos los hipocentros. C. Perfil en dirección este-oeste mostrando todos los hipocentros. D. Histograma de la cantidad de sismos por profundidad.

Resultados

Los sismos que componen la base de datos inicial se ubican principalmente en las regiones sureste y central de Costa Rica (Fig. 2A). Existe un vacío de sismos hacia el oeste y el sur del terremoto de Limón, ya que durante el periodo analizado, entre 1998 y 2020, no se localizaron sismos en esa zona con las características definidas en este estudio para ser usados en la inversión tomográfica (i.e., registrados por más de 15 estaciones con un vacío de cobertura menor a 180°). A pesar de esto, en el conjunto de datos seleccionados sí se tiene cobertura suficiente de rayos sísmicos que atraviesan toda la región del Caribe Sur de Costa Rica (Fig. 2A).

La distribución de los sismos usados por profundidad varía entre ~1 y 110 km (Fig. 2B, 2C y 2D). Los sismos más profundos (70-110 km) se ubican principalmente en la parte central de Costa Rica, entre las latitudes 9,0° y 10,5° norte (Fig. 2B) y longitudes 83,5° y 84,0° oeste (Fig. 2C) y corresponden con el ~4% de la totalidad de la base datos. El ~79 % de los sismos empleados en la inversión tienen profundidades entre 0 y 40 km y el ~17 %, entre 40 y 70 km.

Modelo mínimo de velocidades en una dimensión (M1D)

Para encontrar el M1D de la onda P, se invirtieron 18493 observaciones de onda P de 760 sismos bien ubicados en una región de 260 x 180 km. Estos sismos fueron localizados con un mínimo de 12 estaciones, incluyendo al menos una estación a una distancia de 1,5 veces la profundidad del evento y con un gap de cobertura menor a 180°. Se incluyeron solamente observaciones con la menor incertidumbre de lectura (< 0,10 s).

Se realizaron inversiones utilizando un modelo inicial basado en el de Quintero y Kissling (2001) para Costa Rica, y luego se examinó una variedad de geometrías de capa siguiendo un esquema de prueba y error. Durante estas pruebas se obtuvo una buena convergencia entre los 3 y los 70 km de profundidad, coincidiendo con la distribución de los hipocentros del set de datos (Fig. 2D). El modelo que resultó en un RMS menor (0,189 s) tiene 10 interfaces a profundidades desde 0 a 100 km con velocidades desde 3,00 hasta 8,29 km/s (Fig. 3A). El RMS final y la varianza para el conjunto de datos fueron 0,189 s y 0,043 s2 respectivamente, en contraste con los valores iniciales de 0,424 s y 0,217 s2.

Para demostrar la independencia de la solución final con respecto de los modelos de iniciales, se realizaron inversiones aumentando y disminuyendo las velocidades iniciales hasta en ~1,5 km/s para cada capa, conservando el gradiente original. Todos los resultados muestran una fuerte convergencia con el modelo de velocidad mencionado anteriormente, que se considera como el M1D para este conjunto de datos (Fig. 3A, línea roja).

Fig. 3:
Resultados del Modelo Mínimo de Velocidades en una Dimensión (M1D) para la onda P. A. Modelos de velocidad inicial y final. El área gris muestra el espacio de velocidades explorado a partir de diferentes velocidades iniciales (líneas grises). Todos los modelos iniciales convergen en un M1D (línea roja) o en modelos similares (líneas magenta punteadas). La línea azul muestra el gradiente de velocidad usado como modelo inicial para la TSL. B. Correcciones para las estaciones. El rombo negro marca la estación de referencia, EDP2. Las cruces rojas y azules denotan correcciones positivas y negativas, respectivamente, y su tamaño es relativo al valor de la corrección. Las correcciones de estación cercanas a cero (entre -0,05 y 0,05 s) se indican con círculos y las correcciones máximas y mínimas se indican en la esquina inferior izquierda del mapa como referencia. La estrella amarilla representa el terremoto de Limón de 1991 (Mw 7,7).

Las correcciones para estaciones asociadas con el modelo de velocidades final (Fig. 3B) muestran un patrón atribuible a la geometría de la red, la distribución de la sismicidad y la geología bajo las estaciones y su elevación. En este trabajo se escogió como referencia la estación EDP2, instalada en Potrero Grande de Buenos Aires, que se encuentra lo más cerca posible al centro del arreglo de estaciones que se emplea en este estudio y a una altura media (201 m.s.n.m.) entre la estación más profunda (3488 m bajo el nivel del mar) y la más elevada (3506 m.s.n.m.). El empleo de estaciones de fondo marino (proyecto SFB574 de la Universidad de Kiel, Alemania) y la baja cantidad de sismos de profundidades inferiores a 3 km, hizo necesario fijar la primera capa con velocidades similares a las de los sedimentos de los primeros kilómetros de la plataforma y el talud continental (3 km/s). Así, las correcciones de estaciones son en general negativas. Sin embargo, las estaciones ubicadas sobre formaciones sedimentarias muy potentes, tales como la Cuenca de Limón Sur, retienen correcciones positivas notables (+1,06 a +1,16 s). Los retrasos negativos más altos (-2,34 a -1,05 s) se observan, por ejemplo, para las estaciones de fondo marino más profundas, asociados con las altas velocidades sísmicas esperables en la corteza oceánica.

Tomografía con sismos locales (TSL)

Para derivar la TSL se utilizó un set de 1208 sismos con 35964 arribos de onda P (Fig. 2) que fue seleccionado al relocalizar la base de datos inicial usando el M1D. Este set de datos de alta calidad incluye sismos con un gap de cobertura menor a 180°, un mínimo de 15 arribos de onda P y al menos una estación a una distancia de 1,5 veces la profundidad del evento.

Los resultados de las pruebas con tableros de ajedrez explicados en la metodología son muy similares usando ya sea cuatro nodos o un nodo. Como su objetivo principal es explorar las regiones del modelo donde es posible resolver estructuras en pequeña escala (~20 km), se presentan ejemplos de los resultados de tablero de ajedrez con anomalías de un solo nodo (Fig. 4). En general, los resultados de estas pruebas indican que la fuga horizontal de velocidades es mínima, con excepción de una región orientada noroeste-sureste bajo el flanco pacífico de Talamanca, a los 60 km de profundidad. Las fugas verticales son muy bajas a profundidades superiores a los 3 km, con excepción de un área cerca de la trinchera frente al Pacífico Central y en el oeste de Panamá, observables en las pruebas con anomalías de 2 x 2 nodos. La restitución de la amplitud de las anomalías es de buena a excelente (40-100%) entre los 10 y los 60 km de profundidad. En resumen, las pruebas de tablero de ajedrez permiten describir con alto grado de certidumbre los resultados tomográficos para la región Caribe Sur de Costa Rica entre los 6 y los 60 km de profundidad (Fig. 4).

Fig. 4:
Secciones horizontales con resultados de dos pruebas con modelos sintéticos de tablero de ajedrez para anomalías del tamaño de un nodo de la rejilla y con magnitud ±10 % en Vp. La ubicación de las anomalías iniciales de alta (+ 10%) y baja (-10%) velocidad se muestran como cuadrados azules y rojos, respectivamente. Los nodos de la rejilla usada en la tomografía se muestran como equis grises y las estaciones, como triángulos. A. Ejemplos a 10, 20 y 40 km de profundidad del modelo con anomalías a 3, 10, 40 y 80 km de profundidad. B. Ejemplos a 10, 20 y 40 km de profundidad del modelo con anomalías a 0, 6, 20, 60 y 100 km de profundidad.

El resultado final del TSL tuvo un RMS ponderado de todos los residuos de tiempos de viaje de 0,21 s, un valor relativamente alto. Esto podría deberse a diferencias entre analistas en el uso de la escala de pesos de las lecturas y/o a variaciones de velocidad de menor escala aún no modeladas, factible en vista especialmente de la alta complejidad tectónica del volumen modelado.

Los resultados tomográficos finales se muestran en dos secciones transversales (Fig. 5) con valores de Vp absolutos y relativos. Estos últimos se presentan como porcentaje de variación con respecto del M1D, con una escala de colores clásica, donde el rojo representa velocidades menores al M1D y el azul, velocidades mayores. Las secciones tomográficas se extienden desde la trinchera hasta la proyección en superficie del CDNP y discurren a través de los límites noroeste y sureste de la cuenca Limón Sur (Fig. 1). Del lado Pacífico, se observa una banda de bajas velocidades que desciende hasta profundidades de ~70 km bajo el flanco suroeste de la cordillera de Talamanca. Las velocidades absolutas van desde 6,2 km/s a los 10 km de profundidad, hasta 8,0 km/s a los 70 km de profundidad. Su coincidencia con una franja de sismicidad descendiente permite suponer que se trata de la corteza oceánica engrosada del levantamiento del Coco. Velocidades similares han sido observadas para la corteza oceánica anómala en trabajos anteriores (e.g., Husen et al., 2003). Además, el contorno del techo de la losa, según el modelo de Lücke y Arroyo (2015), coincide muy bien con el límite superior de la anomalía y con la ubicación de los terremotos de Osa del 2002 (Mw 6,4; Fig. 5A; Arroyo, Grevemeyer, Ranero y von Huene, 2014) y Golfito de 1983 (Mw 7,4; Fig. 5B; Adamek, Tajima y Wiens, 1987).

Fig. 5:
Perfiles tomográficos a lo largo del sureste de Costa Rica (ver ubicación en Fig. 1). La estructura de velocidades se muestra como el porcentaje de cambio relativo al M1D según la escala de colores en la base de la figura, omitiendo las áreas con menor resolución. La línea blanca encierra el área de buena resolución. Los valores absolutos de Vp están interpolados en intervalos de 0,25 km/s. Se muestran los hipocentros (círculos rojos) ubicados a 20 km a ambos lados de cada perfil. Los nodos de la rejilla de la tomografía se muestran como equis grises. La línea punteada es el techo de la placa del Coco de acuerdo con Lücke y Arroyo (2015). Sobre cada perfil se indica la topografía como referencia. Las estrellas representan los sismos mencionados en el texto. T y C marcan la posición en superficie de la trinchera y la costa, respectivamente, CT la cordillera de Talamanca y CDNP la posición del Cinturón Deformado del Norte de Panamá. A. Perfil X. B. Perfil Y. C. Interpretación tectónica esquemática de los perfiles tomográficos y del origen del terremoto de Limón (estrella).

Del lado Caribe, se detectaron zonas conspicuas de baja velocidad. Entre los kilómetros 200 y 230 a lo largo del perfil X (Fig. 5A) y 150 y 200 en el perfil Y (Fig. 5B), entre las estribaciones noreste de Talamanca y la proyección en superficie del CDNP, se observa una reducción de velocidades de hasta 20% entre los 0 y los 20 km de profundidad. Las velocidades absolutas van de 2,8 a 6,0 km/s. Estas bajas velocidades están delimitadas lateralmente hacia el suroeste por una zona de alta velocidad bajo la cordillera de Talamanca, pero continúan hacia mayores profundidades (20 a 50 km) debajo de Talamanca, entre los kilómetros 120 y 170 a lo largo de los perfiles, donde se observan reducciones de velocidad de hasta 5 %. Las velocidades absolutas en ese sector varían entre 6,5 y 7,5 km/s. Esta zona tiene una forma elongada y buzante hacia el suroeste, más claramente observable en el perfil Y (Fig. 5B).

Discusión

El terremoto de Limón en el contexto de la imagen tomográfica

El terremoto de Limón del 22 de abril de 1991 (Mw 7,7, Fig. 1) ocurrió en el extremo oeste del CDNP, el cual marca el límite norte de la microplaca de Panamá (i.e., Adamek, et al., 1988). De acuerdo con Kobayashi et al. (2014), la microplaca de Panamá se mueve hacia el noreste como un bloque tectónico único, en respuesta a la colisión del levantamiento del Coco y a un redireccionamiento por colisión con el Bloque Norte de los Andes. En este contexto, el borde norte de la microplaca de Panamá está sometido a esfuerzos compresivos, que se expresan como un notable cinturón de fallas inversas orientadas predominantemente noroeste-sureste y pliegues asociados de rumbos similares interpretados en los primeros 10 km de profundidad (e.g., Silver, Reed, Tagudin y Heil, 1990; Fernández, Bottazzi, Barboza y Astorga, 1994; Silver, Galewsky y McIntosh, 1995).

A niveles más profundos, el CDNP ha sido descrito en la parte central de Panamá por Camacho, Hutton y Pacheco (2010) como una zona de subducción, en la cual la placa Caribe se hunde debajo de la microplaca de Panamá. Estos autores evidenciaron que a través del centro de Panamá, los hipocentros aumentan de profundidad hacia el sur hasta los 80 km, dibujando una zona de Wadati-Benioff bien definida en ese sector. Un aspecto muy relevante del CDNP es que los sismos de magnitudes altas (M > 7) son frecuentes, llegando incluso a desencadenar tsunamis en 1822, 1882, 1916 y 1991 (e.g., Camacho, 1994). Aparte del terremoto de Limón de 1991, otros grandes terremotos en el Caribe costarricense ocurrieron en 1746, 1798 y 1822 (Montero, 1989; Boschini, 1989; Camacho y Víquez, 1993; Boschini y Montero, 1994; Peraldo y Montero, 1994) y en Panamá en 1882, 1904 y 1916 (Camacho y Víquez, 1993).

Las bajas velocidades debajo del Caribe Sur de Costa Rica observadas en los perfiles tomográficos obtenidos en este estudio (Fig. 5) se pueden correlacionar con dos zonas a diferentes profundidades. Entre la superficie y los ~10 km de profundidad, las velocidades de 2,8 a ~4,0 km/s observadas se asocian con las rocas sedimentarias del Cretácico al Cuaternario de la cuenca de Limón Sur, cuyo espesor se ha estimado en ~10 km (Fernández et al., 1994; Denyer y Alvarado, 2007). En este trabajo se interpreta que esta es la zona altamente fallada y plegada que ha sido tradicionalmente descrita como el CDNP en Costa Rica, y además que su estructura se asemeja a un prisma acrecional de una zona de subducción.

Por debajo de esa profundidad (> 10 km), la sismicidad aumenta y las bajas velocidades forman una gruesa franja (~10-20 km) inclinada hacia el suroeste, visible hasta los 50 km de profundidad (Fig. 5). La sismicidad de 1998 a 2020 relocalizada en este trabajo ocurre entre los 10 y los 30 km de profundidad, y se distribuye sugiriendo también un plano buzante hacia el suroeste. Esta configuración de bajas velocidades y sismicidad se interpreta aquí como una zona de convergencia donde la placa Caribe se hunde bajo la microplaca de Panamá con un ángulo inicial de ~20° que parece incrementar a ~50° hasta alcanzar profundidades máximas de 50 km bajo las elevaciones mayores de la cordillera de Talamanca (Fig. 5).

El terremoto de Limón tuvo una solución de tensor de momento inversa con un plano nodal de rumbo 103° (noroeste-sureste), inclinado 25° hacia el suroeste (Dziewonski, Ekström y Salganik, 1992). Diferentes estudios y redes sismológicas ubican su hipocentro a una profundidad de 10-20 km (Dziewonski et al., 1992; Fan, Beck y Wallace, 1993; Goes et al., 1993; Protti y Schwartz, 1994; Suárez et al., 1995). En los perfiles tomográficos, el hipocentro del terremoto se ubica siguiendo la tendencia general de los hipocentros relocalizados y dentro de la franja inclinada de bajas velocidades.

La interpretación tectónica conjunta de la imagen tomográfica, la sismicidad aquí presentada y el hipocentro del terremoto de Limón es consistente con las conclusiones de trabajos anteriores (Plafker y Ward, 1992; Fan et al., 1993; Goes et al., 1993; Ponce et al., 1994; Protti y Schwartz, 1994; Suárez et al., 1995; Tajima y Kikuchi, 1995), que describen el origen del terremoto en un régimen altamente compresivo, al mismo tiempo que también es coherente con el trabajo de Camacho et al. (2010), que ilumina una zona de Wadati-Benioff en el centro de Panamá. Esta es la primera imagen tomográfica con resolución suficiente en la que se puede interpretar el CDNP (< 10 km de profundidad) en el Caribe costarricense como un sobrecorrimiento (¿prisma acrecional?) perteneciente a la microplaca de Panamá, bajo la cual se hunde la placa Caribe hasta profundidades de 50 km. El terremoto de Limón es por consiguiente un evento originado en una zona sismogénica interplacas Caribe-Panamá (Fig. 5C).

De los perfiles tomográficos mostrados en este trabajo se puede determinar que ~100 km de la placa Caribe se han subducido bajo la microplaca de Panamá. Considerando una velocidad de convergencia de 6,9 mm/año para el oeste del CDNP calculada por Kobayashi et al. (2014) y asumiendo que ha sido una tasa constante, se puede estimar que la subducción de la placa Caribe bajo la microplaca de Panamá inició hace al menos ~14,5 Ma. Por su parte, Kobayashi et al. (2014) estiman un inicio de la subducción bajo el centro de Panamá hace al menos 22-23 Ma. Así, esta subducción comenzó mucho antes del arribo del levantamiento del Coco a la Fosa Mesoamericana, un evento que los trabajos más recientes datan en el Pleistoceno Inferior (Vannucchi, Sak, Morgan, Ohkushi y Ujiie, 2013) o Plioceno Superior (Morell, 2015).

Estos hallazgos implican que el sureste de Costa Rica es una zona de doble subducción convergente, en la que dos losas subducidas (Coco y Caribe) se acercan una a la otra en profundidad. Este complejo escenario también ha sido identificado por ejemplo en el norte de Sulawesi, en Indonesia, en donde en el transcurso de una distancia de ~100 km, la losa de Célebes, que subduce la litósfera del mar de Célebes, se encuentra en profundidad (~50 km) con la losa de la placa Sangihe (Hall y Spakman, 2015). Una geometría tan compleja como esta podría tener implicaciones en el entendimiento de los singulares rasgos tectónicos de este sector del istmo centroamericano, que posee sólo 120 km de ancho y las mayores elevaciones (4 km) del sur de América Central y que carece además de vulcanismo activo y sismos de profundidad intermedia (> 70 km).

Conclusiones

Se obtuvo la primera imagen tomográfica con resolución suficiente para interpretar el límite entre la placa Caribe y la microplaca de Panamá como una zona de subducción en el Caribe Sur de Costa Rica. En esta imagen, el CDNP coincide con una zona de bajas velocidades hasta profundidades de 10 km, que corresponde con la región que ha sido descrita en trabajos previos como altamente fallada y plegada y que podría asemejarse a un prisma acrecional de una zona de subducción. A profundidades mayores, la placa Caribe se observa como una banda de bajas velocidades que se inclina hacia el suroeste y es visible hasta los 50 km de profundidad bajo Talamanca. Esta zona de bajas velocidades incluye el hipocentro del terremoto de Limón de 1991 (Mw 7,7) y la sismicidad relocalizada hasta los 30 km. En este contexto, se concluye que el terremoto de 1991 es un evento originado en una zona sismogénica interplacas.

La edad mínima del inicio de la subducción de la placa Caribe en el Caribe Sur de Costa Rica podría ser 14,5 Ma, mucho antes del arribo del levantamiento del Coco a la Fosa Mesoamericana. La interpretación aquí presentada implica además que el sureste de Costa Rica es una zona de doble subducción convergente, en la que las losas subducidas de las placas Coco y Caribe se acercan una a la otra en profundidad debajo de Talamanca. La geometría propuesta es un nuevo aporte para esclarecer la compleja configuración tectónica del sureste de Costa Rica y provee un marco de referencia adicional para interpretar el origen de los grandes terremotos y tsunamis del Caribe costarricense, así como para la estimación de la amenaza sísmica de la región.

Agradecimientos

Este trabajo es resultado de los proyectos de investigación 830-B8-500 ''Estructura de la zona de subducción en Costa Rica mediante tomografía de sismos locales'', 113-B5-704 ''Vigilancia Sísmica de Costa Rica'' y 113-A1-716 ''Apoyo de asistentes a la Sección de Sismología, Vulcanología y Exploración Geofísica'' y del programa de investigación 113-B9-911 ''Red Sismológica Nacional'', inscritos en la UCR y del proyecto GREAT de la Universidad de Rutgers. Además, ha sido posible gracias a los recursos económicos y de personal brindados por la UCR y el ICE, y al financiamiento que la RSN recibe de la Ley Nacional de Emergencias Nº 8488. Agradecemos al personal de la RSN (UCR-ICE) que permite mantener la red de estaciones y los sistemas de localización, especialmente a Luis Fernando Brenes y Jean Paul Calvo. Agradecemos las sugerencias y comentarios del Dr. Ronnie Quintero como revisor, que contribuyeron a mejorar este manuscrito. Las siguientes personas colaboraron en la localización de sismos: Magda Taylor, Brandon Acosta, Carolina Fallas y Roberto Masis. Gran parte de las figuras fueron realizadas con Generic Mapping Tool (GMT), de Wessel, Smith, Scharroo, Luis y Wobbe (2013).

Referencias bibliográficas

Fechas de Publicación

Histórico