Caracterización resistiva 3D del área geotérmica de Paipa, Colombia

Introducción

La fase exploratoria de un proyecto geotérmico se encuentra caracterizada por los estudios de reconocimiento, prefactibilidad y factibilidad (OLADE y BID, 1994); donde se conjugan diversas técnicas geológicas, geoquímicas y geofísicas en la identificación de las estructuras características de un sistema geotérmico (fuente de calor, reservorio permeable, capa sello, zonas de recarga, etc.). Los métodos geofísicos se implementan para la identificación indirecta de las condiciones estructurales y la geometría del sistema. Particularmente, el método magnetotelúrico (MT) ha sido ampliamente utilizado en la caracterización resistiva de diversos ambientes geotermales (Volpi et al., 2003; Abdul-Azeez y Harinarayana, 2007; Arango et al., 2009; Cumming y Mackie 2010; Peacock et al., 2012); ya que la resistividad eléctrica del subsuelo está en gran medida controlada por la porosidad, permeabilidad, saturación del medio, salinidad, presión y temperatura del medio (Chave y Jones, 2012; Simpson y Bahr, 2005).

Trabajos como los de Spichak y Manzella (2009), Muñoz (2014) y Van Leeuwen (2016) muestran una recopilación de estudios magnetotelúricos realizados en diferentes ambientes geotérmicos, así como las posibles firmas resistivas asociadas a estos. Por ejemplo, los sistemas hidrotermales de alta temperatura (Cumming y Mackie 2010; Volpi et al., 2003), presentan una capa sello conductora dada por la presencia de arcillas sobre un reservorio geotermal más resistivo, reflejando un aumento de la temperatura con la profundidad (Ussher et al., 2000). Sin embargo, en algunos sistemas geotérmicos asociados a ambientes volcánicos la capa sello puede estar ausente (Muñoz, 2014). Por otro lado, en sistemas geotérmicos nomagmáticos (Abdul-Azeez y Harinarayana, 2007; Arango et al., 2009; Peacock et al., 2012) existe una amplia gama de configuraciones geométricas, lo que hace difícil la identificación de patrones resistivos a seguir (Muñoz, 2014). Cada sistema debe ser analizado rigurosamente con el apoyo de las diferentes técnicas de estudio.

El área geotérmica de Paipa, ubicada al noreste de Colombia (Figura 1), fue catalogada con una prioridad media-alta a partir de estudios de geología de superficie y de la geoquímica de la fase acuosa de las fuentes termales (OLADE, GEOTÉRMICA ITALIANA, ICEL y CONTECOL, 1982). Por otro lado, un estudio de prefactibilidad sugirió la existencia de un reservorio geotermal con potencial para la generación de energía eléctrica (Japan Consulting Institute, 1983). Dado lo anterior, el Servicio Geológico Colombiano ha realizado diversos estudios geológicos, geoquímicos y geofísicos, encaminados a la caracterización del sistema geotérmico de Paipa.

El área de estudio se encuentra en un ambiente sedimentario con presencia de depósitos piroclásticos y cuerpos efusivos (domos riolíticos) (Cepeda y Pardo, 2004; Rueda-Gutiérrez, 2017). Los manantiales termales presentes en la zona han registrado temperaturas máximas entre 70-76°C (Alfaro, 2002a; Ortiz y Alfaro, 2010). De acuerdo a la localización de las manifestaciones termales (Figura 1), estas estarían controladas estructuralmente por cruces de fallas en el área de influencia de formaciones sedimentarias (Velandia, 2003; Alfaro et al., 2012), permitiendo el ascenso de fluidos. Adicionalmente, Alfaro (2002a; 2002b) establece la presencia de sulfato sódico en el análisis de las muestras de los manantiales, proveniente de evaporitas al oeste de la secuencia sedimentaria, el cual se mezclaría con los fluidos geotérmicos, enmascarando la química del fluido hidrotermal y dificultando la correcta estimación de la temperatura en el reservorio. Sin embargo, a partir de geotermómetros de sílice, se ha estimado un reservorio con temperaturas entre 120°C y 230°C (Alfaro, 2002b).

Llanos et al. (2015), a partir de modelación geológicogeofísica con información de métodos potenciales (Vásquez, 2012), identificaron anomalías de alta densidad correspondientes a basamento profundo e intrusiones en la zona central de trabajo (domos Alto Los Volcanes y Alto Los Godos, en la Figura 1), así como anomalías de baja densidad asociadas a la secuencia sedimentaria al centro-oeste del área.

El objetivo de este trabajo es realizar la caracterización resistiva del área geotérmica de Paipa, a partir de la inversión 3D de datos magnetotelúricos, utilizando el código ModEM (Egbert y Kelbert, 2012; Kelbert et al., 2014), e identificar las principales estructuras resistivas del sistema geotérmico.

Marco geológico

El área de estudio se encuentra principalmente en un ambiente sedimentario al sur del municipio de Paipa (Boyacá, Colombia). De acuerdo a la caracterización geológica y estructural realizada por Velandia (2003), al NO y SO del área de estudio (Figura 1) se encuentran sedimentos del Cuaternario y del Neógeno (Formación Tilatá), representados por secuencias de arenas, limos arcillas, conglomerados correspondientes a la actividad más reciente de tipo aluvial, lacustre y fluvio-lacustre. En el sector central y oriental, se encuentran rocas sedimentarias del Paleógeno (Formación Bogotá) y Cretácico (Formaciones Guaduas, Labor y Tierna, Los Pinos, Plaeners, Conejo, Churuvita Une y Tibasosa), constituidas por intercalaciones de areniscas cuarzosas, limolitas, arcillolitas, niveles esporádicos de calizas e impregnaciones de sales y lutitas. El espesor de esta secuencia sedimentaria se estima entre 1 – 1,5 km aproximadamente al NO del área (Velandia, 2003). Bajo la secuencia sedimentaria se encuentra el basamento formado por rocas sedimentarias (niveles arenosos fracturados, intercalados con arcillolitas y limolitas) y metamórficas (filitas, esquistos y gnéises) del Paleozoico. Dicho basamento correspondería a la continuación en profundidad del Macizo de Floresta, el cual aflora fuera del área de estudio al NE, formando el núcleo del Anticlinal Tibasosa-Toledo (Figura 1). Adicionalmente, en la Tabla 1 se puede observar una breve descripción litológica de cada una de las formaciones geológicas, así como su importancia hidrogeológica.

Sobre la cubierta sedimentaria, Cepeda y Pardo (2004) observaron una secuencia de unidades piroclásticas correspondientes a dos etapas eruptivas posteriores al Eoceno tardío, con intercalaciones de sedimentitas de ambientes fluvial, lacustre y coluvial (Figura 1). Por otro lado, al oeste del área, Velandia (2003) describe la presencia de una brecha hidrotermal (Figura 1) caracterizada por fragmentos muy alterados de posible origen ígneo, posterior a los eventos volcánicos, con un alto grado de alteración argílica avanzada (González et al., 2008; Rodríguez y Alfaro, 2015), sugiriendo un posible aporte de calor radiogénico al sistema hidrotermal (Alfaro et al., 2017), producto de una intrusión granítica en la zona.

Rueda-Gutiérrez (2017) caracterizó y cartografió tres sectores con cuerpos dómicos de origen de arco volcánico, denominados Alto Los Volcanes, Alto Los Godos y Quebrada Honda; con edades estimadas de 1,76 y 1,8 Ma para el Alto Los Volcanes y Quebrada Honda, respectivamente; y 2,8, 2,7 y 2,6 Ma para el sector del Alto Los Godos. Sin embargo, refiere la edad más reciente de estas intrusiones en 1 millón de años a partir de dataciones complementarias por trazas de fisión.

El área se encuentra deformada por fallas regionales con dirección NE-SO, asociadas a la apertura de una cuenca marina en el Cretácico y la inversión tectónica consecuencia de un régimen compresivo que dio lugar al levantamiento de la Cordillera Oriental en el Neógeno (Ortiz, 2014). Velandia (2003) observó una serie de estructuras longitudinales y transversales a los lineamientos regionales (Figura 1). Para las fallas longitudinales con orientación NNE-SSO (El Bizcocho, El Batán, Rancho Grande, Buena Vista, Agua Tibia, Lanceros) describe sus trazos como sinuosos y localmente discontinuos por estar cubiertos por depósitos neógenos y cuaternarios, tanto sedimentarios como vulcano-clásticos, dificultando la interpretación de relaciones estructurales entre las unidades preneógenas (Velandia, 2003). Las fallas transversales con orientación NO-SE (Cerro Plateado, Paipa-Iza, Las Peñas) controlan morfología y drenaje, posiblemente facilitando el paso de fluidos hidrotermales (Alfaro et al., 2012), y permitiendo el ascenso de los magmas que dieron origen al volcanismo reconocido en la zona. Entre tanto, las fallas transversales con orientación NESO (El Hornito, Canocas, Santa Rita, El Tuno) afectan sólo la secuencia sedimentaria (Velandia, 2003).

Metodología

El modelo resistivo del área geotérmica se estimó por medio de la inversión de datos magnetotelúricos. El método magnetotelúrico es un método geofísico de fuente pasiva, el cual registra en superficie las variaciones temporales del campo electromagnético natural, por medio del cual se puede estimar la distribución de la resistividad del subsuelo (Vozoff, 1991) desde decenas de metros hasta centenas de kilómetros. A partir de la relación entre las componentes horizontales del campo magnético y eléctrico, en el dominio de las frecuencias, se halla un tensor de impedancia complejo (Z) de 2x2, cuya expresión matemática está dada por:

Cada una de las componentes de . es transformada en resistividad aparente (ρa) y fase (φ), que posteriormente son modeladas.

donde μ0 es la permeabilidad magnética del vacío y . es la frecuencia angular.

Adicionalmente, la medida de la componente vertical del campo magnético y su relación con las componentes horizontales del mismo campo, genera la función de transferencia geomagnética (VTF en inglés). La cual permite estimar las variaciones laterales de la resistividad en el medio (Vozoff, 1972).

En el área de estudio se adquirieron 78 sondeos audiomagnetotelúricos y magnetotelúricos (AMT+MT) (Figura 2) a lo largo de diferentes campañas entre los años 2012 y 2016 (Moyano y Vallejo, 2015; González-Idárraga y Rodríguez-Rodríguez, 2017). La información se adquirió con equipos receptores de la marca Phoenix Geophysics y sensores magnéticos con capacidad de registro hasta 10000 s. Los dipolos eléctricos y bobinas fueron orientadas al norte para el eje x y al este para el eje y. Los tiempos de medición fueron de una hora para las mediciones AMT y no menos de 12 horas para las mediciones MT. A partir de un procesamiento robusto de las señales temporales adquiridas, bajo un esquema de referencia remota local (Gamble et al., 1979), se obtuvieron las respectivas curvas de resistividad aparente, fase y VTF.

En la Figura 3 se presentan las curvas de resistividad aparente, fase y VTF de 6 sondeos representativos de las adquisiciones realizadas en el área de estudio. Estos muestran buenas estimaciones hasta periodos de 10 s para todas las componentes del tensor de impedancia, mientras que en 17 sondeos, la componente Zyx alcanza estimaciones iguales o mayores a 100 s (PA0064, PA0077, PA0130, PA0077, PA0145, PA0150, PA0160, PA0162, entre otros. Ubicaciones en la Figura 2). Entretanto, las curvas de VTF presentan buenas estimaciones hasta periodos entre 1 y 10 s. También se identificó en algunas mediciones una baja estimación en el rango de banda muerta del AMT, como se observa en las componentes de VTF de los sondeos PA0020, PA0073 y PA0098. Para periodos mayores a 10 s las curvas presentaron mayores errores y dispersión en los datos, por lo que no se tuvieron en cuenta al momento de realizar la inversión.

Finalmente, sondeos como PA0020, PA0027 y PA0077 presentaron corrimiento estático (static shift, en inglés) en las curvas de resistividad aparente. Sin embargo, la inversión 3D puede corregir la distorsión galvánica que genera dicho corrimiento, así como los efectos inductivos debido a estructuras cercanas a la superficie (Siripunvaraporn, 2012).

Análisis de dimensionalidad

Con ánimo de estimar la dimensionalidad geoeléctrica del medio, se realizó un análisis a partir del Tensor de Fase (Caldwell et al., 2004). El tensor de fase está definido a partir de las relaciones de fase entre las componentes horizontales de los campos eléctrico y magnético. La representación gráfica de las cantidades escalares valor máximo del tensor (φmax), valor mínimo del tensor (φmin) y el ángulo de inclinación skew β (medida de la asimetría del tensor) reflejan la dimensionalidad del medio, así como la dirección de strike del contraste de conductividad (Caldwell et al., 2004). En la Figura 4 se puede observar la representación gráfica del tensor de fase, donde las circunferencias definen un medio 1D con un ángulo skew de cero, un medio 2D es representado por una elipse con valores de skew de cero o cercanos (±3°) y valores diferentes para el skew reflejan un medio 3D (Caldwell et al., 2004; Bibby et al., 2005; Booker, 2014).

Se representaron los tensores de fase para los periodos 0,01, 0,1, 1 y 10 s (Figura 4). En el periodo 0,01 s se observa, de forma general, circunferencias con valores de skew cercanos a cero indicando una dimensionalidad 1D en las zonas más someras del área de estudio. Este comportamiento es esperado, al encontrarse la zona de estudio en un ambiente sedimentario. A periodos más altos (0,1, 1 y 10 s) se observan dimensionalidades 2D y 3D. Por otra parte, las flechas negras representan los vectores de inducción (Parkinson, 1959) que proporcionan información de contrastes en la conductividad del medio. Este parámetro refleja un gradiente conductivo hacia el NO del área de estudio a un periodo de 1 s, apoyado por la orientación del eje mayor de las elipses del tensor de fase en la misma dirección.

Inversión 3D

El algoritmo de inversión magnetotelúrica 3D ModEM (Egbert y Kelbert, 2012; Kelbert et al., 2014) fue empleado para generar el modelo de resistividad 3D del área de estudio. Se realizó la inversión conjunta del tensor de impedancia y de la función de transferencia geomagnética de 78 sondeos, previamente depurados en los periodos con una baja calidad. Se estableció un valor de 0,3 de covarianza. Se implementó un modelo inicial con una resistividad de 100 ohm.m, un mallado horizontal de 33 celdas en el eje N-S, 30 celdas en el eje E-O, dimensiones de 400 m x 400 m y 10 celdas de atenuación a cada lado con un incremento logarítmico de 1,3. El eje vertical se discretizó en 79 celdas. Entre el sondeo realizado a mayor altura y el de menor altura, las celdas poseen una resolución de 50 m y posteriormente se realizó un incremento logarítmico de 1,2. Adicionalmente, los sondeos se proyectaron a la superficie acorde a la topografía y se centraron en cada celda.

La inversión fue realizada con 18 periodos del tensor de impedancia (Zxx, Zxy, Zyx, Zyy) y el VTF, en un rango desde 0,001 s hasta 10 s. Se realizaron pruebas de inversión de sólo Z, sólo VTF y Z+VTF. La secuencia se inició invirtiendo el tensor de impedancia con errores de piso de hasta 50%, posteriormente se utilizó el modelo resultante como modelo inicial de una nueva inversión con errores de piso menores, esto con el ánimo de ir evaluando los ajustes de las anomalías en el modelo (Romero et al., 2019). Finalmente se definieron errores de piso en 10% para la antidiagonal de Z, 20% para la diagonal de Z y un valor de 3% para el VTF, obteniendo un modelo resistivo del área con un valor RMS de 1,68 después de 118 iteraciones.

La Figura 5 muestra curvas de ajuste del modelo resistivo calculado. De manera general, se detalla un ajuste coherente para las resistividades aparentes de las componentes Zxy y Zyx del tensor de impedancia y para las curvas de VTF, en el rango de periodo modelado. Las componentes de la diagonal principal del tensor de impedancia (Zxx y Zyy) presentan desajustes superiores. Estas componentes están más sujetas a errores y se les ha asignado un error de piso superior. Sin embargo, se recomienda la inversión de todas las componentes del tensor de impedancia, ya que las estructuras pueden ser mejor cartografiadas (Kiyan et al., 2014).

Resultados

Los resultados de la modelación 3D realizada se muestran en las Figuras 6 y 7. Al oeste del área se observa la presencia de una serie de anomalías conductoras someras (C1, C2 y C3) que no superan los 500 m de profundidad. Dichas anomalías son reflejo del análisis de los vectores de inducción que mostraron un gradiente de conductividad al NO de la zona (Figura 4).

Al centro-este y sur del área modelada (Figuras 6 y 7) se encuentra una respuesta resistiva (R1) de grandes dimensiones, por debajo de 500 metros de profundidad para las secciones al sur. Sobre dicha anomalía resistiva se encuentra una capa de resistividad intermedia (entre 10 y 100 ohm.m), aparentemente con una mayor profundidad al oeste del área. Sin embargo, como se observa a partir del análisis de sensibilidad presentado en la siguiente sección, dado el rango de periodos modelado, al norte y noroeste del área no es posible asegurar la caracterización de las estructuras por debajo de 1470 msnm.

Análisis de sensibilidad

Al modificar los valores de resistividad a diferentes profundidades, hasta la base del modelo, y correr modelos directos para observar el desajuste en los periodos altos respecto al de las curvas calculadas en el modelo final, se puede estimar la profundidad de estudio alcanzada por el modelo.

El área modelada tiene una altura media de 2780 msnm, tomando este valor como cero, se modificó la resistividad en un valor de 1000 ohm.m a profundidades de -3150 msnm, -2230 msnm, -1470 msnm, -100 msnm y +1010 msnm. Adicionalmente, se realizaron análisis de sensibilidad con modelos escalonados (Figura 8A), donde el resistivo se encontraba más cercano a la superficie al norte y oeste del área modelada, dado que en dicha zona se encuentran cuerpos conductores superficiales y, al noroeste, sondeos MT con la menor longitud en términos de periodo, como se describió anteriormente. Al observar los diferentes desajustes en los sondeos se concluyó que, de manera general los ajustes a las diferentes profundidades, el modelo posee una profundidad de estudio de hasta la cota de -1470 msnm al centro-este y toda la zona sur del área (sondeos PA0020, PA0077 y PA0098, en la Figura 8B), mientras que al norte y noroeste solo tiene buena precisión hasta la cota +1010 msnm (sondeo PA0073 en Figura 8B).

Para el análisis de sensibilidad de los cuerpos conductores C1, C2 y C3, se modificó su resistividad a valores de 30 ohm.m (Figura 8C) y se realizó una inversión utilizando este como modelo inicial. En la Figura 8C, en la sección iii, se observa como los conductores se recuperan, confirmando la presencia de estos dentro del modelo.

Discusión

Las anomalías conductoras al oeste del área, como se observa en la Figura 6, parecen tener un control estructural dado por las fallas de bajo ángulo de buzamiento. El cuerpo conductor C1, formaría un aparente conducto de fluidos a través de las formaciones Guaduas, Los Pinos y particularmente Plaeners (al ser acuífero), con un posible control a lo largo de las fallas El Bizcocho y El Batán, terminando al norte con la descarga en los manantiales del sector de La Playa, los cuales han registrado temperaturas de hasta 76°C (Alfaro, 2002a; Ortiz y Alfaro, 2010). Por su parte, el cuerpo conductivo C2 representaría tránsito de fluidos termales a través de los depósitos Cuaternarios, a lo largo de las fallas Canocas y El Hornito, hasta su descarga en los manantiales del sector de ITP-Lanceros (temperatura máxima de descarga de 71°C). Dicha descarga se da en el cruce de fallas de El Hornito y El Bizcocho al norte del área de estudio (Figura 1). Estas dos anomalías (C1 y C2) parecieran estar conectadas, como se observa en las vistas en planta a 2484 y 2258 msnm (Figura 6). Por otro lado, de acuerdo a la extensión de la anomalía C2, se esperaría que el manantial de SALPA tuviera aportación de fluidos hidrotermales como en el sector de ITP-Lanceros. Sin embargo, este manantial posee altas concentraciones de sulfato de sodio y una temperatura baja (Alfaro, 2002a; Ortiz y Alfaro, 2010). Alfaro et al. (2017) proponen que este sector seria la descarga de fluidos provenientes de niveles inferiores de la secuencia sedimentaria, los cuales originarían disolución de evaporitas al oeste del área. Adicionalmente, dados los cortes superficiales del modelo (Figura 1 y Figura 6), no se observa una diferenciación resistiva entre las formaciones geológicas, ya que estas presentan grandes similitudes litológicas e hidrogeológicas.

La anomalía resistiva R1 hace referencia a la respuesta del basamento metamórfico al sureste del área, donde muestra una mayor presencia. Por otro lado, las respuestas resistivas R2 y R3 (Figura 6) se asocian a las respuestas superficiales de los sectores dómicos Alto Los Volcanes y Alto Los Godos, respectivamente.

La capa resistividad intermedia (entre 10 y 100 ohm.m), sobre la anomalía resistiva del basamento, se ha interpretado como el bloque de formaciones sedimentarias (González-Idárraga y Rodríguez- Rodríguez, 2017). Sin embargo, dado el análisis de sensibilidad, al noroeste del área no permitiría caracterizar las estructuras por debajo de las anomalías conductoras. Por lo anterior, esto no permite verificar la respuesta resistiva de la fuente salada del manantial de SALPA.

Alfaro et al. (2017) proponen un modelo conceptual para el área geotérmica a partir de estudios de geología, geoquímica y geofísica. La Figura 9 muestra una representación de dicho modelo conceptual sobre el modelo resistivo 3D obtenido en este trabajo. La recarga del sistema (flecha a en la Figura 9) sucedería en formaciones sedimentarias permeables (Tibasosa y Une) que afloran en el anticlinal de Tibasosa-Toledo al sureste del área, donde el terreno basculado al NNO permite la circulación de los fluidos en la mismadirección (Alfaro et al., 2017). La fuente de calor estaría asociada al calor remanente de las intrusiones del Altos los Volcanes y Alto Los Godos (Alfaro et al., 2017), cuya edad de evento más reciente es del orden de 1 millón de años (1,0±0,25 y 1,0±0,35 Ma) (Rueda- Gutiérrez, 2017). Los fluidos, dada la permeabilidad de las capas sedimentarias, así como la presencia de fallas de basamento (Paipa-Iza, Las Peñas), infiltrarían la zona intermedia entre los sectores dómicos (fuentes de calor) (Alfaro et al., 2017). Al encontrarse el sistema geotermal en un ambiente sedimentario es difícil identificar la respuesta resistiva de la capa sello y del reservorio geotermal (Van Leeuwen, 2016). En el modelo presentado (Figura 9), y bajo la suposición de la ubicación del reservorio entre los sectores dómicos (Alfaro et al., 2017), es posible observar de una forma más clara la presencia de una capa conductora somera (C3 en Figura 6 y 7) que se asociaría al sello del sistema y bajo esta se ubicaría el reservorio geotermal (Figura 7 y 8) con un aumento en la resistividad, como lo menciona Muñoz (2014).

Posteriormente, los fluidos termales continuarían su tránsito al norte por medio de la permeabilidad dispuesta por las fallas de Cerro Plateado y El Bizcocho hasta su descarga en los sectores de La Playa e ITP- Lanceros (flechas b y c en la Figura 9). Alfaro et al. (2017) establece que a lo largo de este tránsito, los fluidos recibirían aportaciones de aguas salobres más frías. Adicionalmente, proponen un segundo sistema de circulación, donde la fuente correspondería al calor radiogénico de la intrusión El Durazno (Figura 1) y que tendría como conducto el corredor entre las fallas Canocas y El Hornito (flecha d en Figura 9).

Conclusiones

El modelo resistivo del área geotérmica de Paipa obtenido a partir de la inversión 3D de datos magnetotelúricos, mostró la presencia de cuerpos conductores al centro-oeste del área, aparentemente controlados por las fallas de bajo ángulo (Canocas, El Hornito y El Bizcocho), formando conductos para el tránsito de fluidos hidrotermales hacia las zonas de descarga de La Playa e ITP-Lanceros. El basamento del área es identificado como un cuerpo resistivo al este y sur del área de estudio, sobre el que se identificó, en inmediaciones de los domos Alto Los Volcanes y Alto Los Godos, las respuestas eléctricas asociadas a la capa sello (cuerpo conductivo) y del reservorio geotermal.

Agradecimientos

El presente trabajo se realizó dentro de los objetivos misionales del Servicio Geológico Colombiano (SGC). El autor agradece a los integrantes del Grupo de Exploración de Recursos Geotérmicos del SGC, especialmente a la Química Claudia Alfaro y al MSc. Juan Camilo Herrera, por sus aportaciones y sugerencias durante la realización de este trabajo. Así como a los revisores, que con sus observaciones se logró una mejora sustancial del artículo.

Referencias

Información adicional

Forma de citar: González-Idárraga, C.E. (2020). Caracterización resistiva 3D del área geotérmica de Paipa, Colombia. Boletín de Geología, 42(3), 81-97 DOI: https://doi.org/10.18273/revbol.v42n3-2020003

Boletín de Geología
ISSN: 0120-0283
Vol. 42
Num. 3
Año. 2020

Caracterización resistiva 3D del área geotérmica de Paipa, Colombia

Carlos Eduardo González-Idárraga
Servicio Geológico Colombiano,Colombia