The geological structures of the Pyrenees and their peripheral basins examined through EMAG2v2 magnetic data

The major goal of this work is to provide insight into the structural and geological anatomy of the Pyrenees based on the magnetic anomalies and inferred lineaments from Earth Magnetic Anomaly Grid 2 arcmin resolution (EMAG2v2) magnetic data. We focus on providing qualitative and semi-quantitative evidence of the magnetic signature of the Pyrenees mountain range domains and structures. The integration of reduced-to-the-pole and processed maps, as well as the Bouguer anomaly map with geological data, has proved to be significantly useful in order to shed light on the main anomaly sources. Considering their magnetic response and texture, several anomalies can be linked to buried geological bodies or changes in the magnetic character of the basement. We have estimated their source bodies' depth through Euler and power spectrum calculations.

We have identified and confirmed eight magnetic zones with different features and interpreted them in terms of the geological and structural setting of the area: one on the Aquitaine Basin with the highest magnetic response in the area, linked to the Sub-Pyrenean Thrust; three along the North Pyrenean Zone linked to mantle materials; one to the east of the Axial Zone representing the boundary between the Pyrenean mountain range and the Mediterranean Sea; two more located in the Catalan Coastal Range related to the volcanic fields in the area; and one in the Gulf of Lyon caused by a thicker crustal block. The result is an overall interpretation of the Pyrenees main magnetic domains.

1 Introduction

The characterization of the Earth's crustal structure and composition through geophysical data has remained a significant challenge in geoscience over the last century. Among various geophysical exploration methods, aeromagnetic surveying has emerged during the last century, initially performed during World War II (Nabighian et al., 2005). Magnetic surveying has proven to be highly effective to investigate the Earth's crust at different scales, from mapping basement structures to constraining the extent of saline diapirs (Kimbell et al., 2002; Ayala et al., 2000; Bertrand et al., 2020; Al-Bahadily et al., 2024). Aeromagnetic data surveying allows large areas to be covered in a short period of time without being invasive (e.g., Vine and Matthews, 1963; Vine, 1966; Gibson and Millegan, 1998; Langel and Hinze, 1998; Golynsky et al., 2002; Purucker and Whaler, 2007; Maus et al., 2009). The processing of the acquired data allows the generation of magnetic-anomaly maps that provide valuable information concerning both continental terrains and Earth's oceanic floors (e.g., Müller et al., 2008). In particular, large-scale magnetic-anomaly maps are of great interest for analyzing major lineaments and structural trends (e.g., Maus et al., 2009) that can be used to improve geological mapping (e.g., Ayala et al., 2000; García-Lobón and Ayala, 2007).

Concerning the Pyrenees, numerous geological and geophysical studies have been conducted to infer their crustal structure at depth (e.g., Choukroune and ECORS team, 1989; Casas et al., 1997; Chevrot et al., 2014, 2018). However, their overall study based on magnetic data has received less attention (Zeyen and Banda, 1989; Zeyen et al., 1991).

In this study we utilize data from the Earth Magnetic Anomaly Grid (EMAG2v2) (Maus et al., 2009) to interpret the magnetic patterns related to the geological structure of the Pyrenees. This study highlights the effectiveness of magnetic-anomaly maps in interpreting the structural characteristics of the Pyrenees and gives an overall view of its magnetic complexity.

2 Data and methodology

The magnetic data used in this analysis came from EMAG2v2 (2009) instead of the EMAG2v3 (Meyer et al., 2017a) due to the lack of data of the latter in the study area. The homogenized grid of the EMAG2v2 incorporates the aeromagnetic survey data from France (Institut du Physique du Globe de Paris – Bureau de Recherches Géologiques et Minières, 1960) and Spain (Socías et al., 1991 – Instituto Geográfico Nacional, 1986–1987) as well as international satellite and ship data. It consists of a homogenized grid displaying the magnetic anomalies 4000 m above the geoid with a spatial resolution of 2 arcmin (at this latitude equivalent to ∼ 3700 m) (Maus et al., 2009). The magnetic-anomaly map of the studied area has been integrated into a digital database to compare it with the main Pyrenean structures, gravimetric anomalies, and the main outcrops of intrusive and volcanic rocks.

To enhance the interpretation of the magnetic signal of the Pyrenees, derivative maps and estimated calculations of the location and depth of the magnetic sources have been generated from the total magnetic field. This enables us to improve and complete the correlation between geophysical and geological information. These calculations were made with the Oasis Montaj© software by Seequent (Hinze et al., 2013), whose formula is primarily based on the works of Blakely (1996), Thompson (1982), and Spector and Grant (1970) and further developed considering the references included in the description of each methodology.

2.1 Reduction to pole (RTP)

This magnetic transformation aims to eliminate the inherent asymmetry present in total magnetic-field anomalies, locating the anomalies above the causative bodies. This is achieved by converting the observed anomaly to the equivalent anomaly measured at the north magnetic pole (Baranov, 1957; Baranov and Naudy, 1964). This transformation assumes that the remnant magnetization is negligible compared to the induced magnetization. The advantage of this transformation is that the resulting map is easier to correlate with the surface geology and other geophysical data, like gravity, which helps with the interpretation.

2.2 Magnetic derivatives

The magnetic derivatives have proved effective for locating the diverse anomaly sources since they accentuate the outline of their boundaries and sometimes even estimate their depth, the inclination, and susceptibility contrast between rocks (e.g., Verduzco et al., 2004). In this study, we have applied these calculations on the reduced-to-the-pole magnetic field.

The vertical derivative (VDR_Z) can highlight the shallowest anomaly sources because it enhances the high-wavenumber component of the spectrum and exposes positive values above the sources. The horizontal derivative (HDR_X and HDR_Y) outlines the anomaly source bodies since it produces a phase transformation and enhances the higher frequencies (Fanton et al., 2014). The direction of the attributes to highlight plays a crucial role in the latter since it exposes the directional variation in the magnetic field with regard to the magnetic causative bodies, providing a distinctive anomaly texture (Hayatudeen et al., 2021). The horizontal derivative in x (N90° E) direction will enhance N–S discontinuities, and the y (N180° E) direction will enhance E–W discontinuities. Additionally, directional derivatives at N110° E and N200° E would enhance the structures following the main Pyrenean directions. As a combination of both vertical and horizontal derivatives of the magnetic anomalies, the analytic signal shows symmetric maximum values above the anomaly sources (Ansari and Alamdar, 2009).

Initially described by Miller and Singh (1994) as an innovative method to distinguish the anomalies as well as their edges and peaks, the tilt derivative (TDR) combines the horizontal and vertical derivatives, later refined by Verduzco et al. (2004). Consequently, the tilt angle depends on both the horizontal and vertical derivatives of magnetization, so weak magnetic bodies are weighted the same as strongly magnetized ones (Blakely et al., 2016).

The positive vertical-derivative value's peak above the source has a value of zero along its edge, which turns negative outside of the source. The horizontal gradient is positive along the edges and is zero above the source (Miller and Singh, 1994). The tilt angle derivative is zero on the edge of the anomalies. It is well suited to locate anomalies despite their depth due to the equalization of the amplitudes of short- and long-wavelength anomalies (Lahti and Karinen, 2010). However, the horizontal derivative of the tilt derivative (TDR_THDR), which is invariant of geomagnetic inclination, is capable of providing sharp and well-defined peaks at the source edges. Despite being accurate in defining shallow sources, it may be ambiguous regarding deep bodies (Ma and Li, 2012).

2.3 Euler deconvolution

Euler deconvolution is a technique to estimate the depth of the center of the anomalous sources (Thompson, 1982; Reid et al., 1990) based on the x, y, and z derivatives of the data and a parameter called the structural index (SI) (Cooper, 2008). The SI is related to the geometry of the anomaly sources (Cooper, 2008; Reid et al., 1990). Each Euler solution is calculated based on the data within a window (WS) (Reid et al., 1990) that is recommended to be at least half of the size of the anomaly to compute (FitzGerald et al., 2004). According to the cited literature, SI = 0 is adequate to define contacts or steps of the causative bodies, whereas SI = 0.5 and SI = 1 are rather suitable to determine sills and dikes.

2.4 Power spectrum

The spectral analysis of the magnetic data provides information about the depth of the magnetic sources attending to magnetic contrasts (Maus and Dimri, 1995; Peredo et al., 2021). The spectrum pattern is a good approach to estimate the depth of the causative bodies. When plotting the natural logarithm of the spectrum against the radial frequency (or the wavenumber), the different linear segments of the spectrum showed in the graph can be associated with different geological discontinuities. The depth (z) is determined by the slope of each segment (m) according to the formula $m/\mathrm{4}\mathit{\pi }$ when using power spectrum and lineal wavenumbers (Spector and Grant, 1970; Peredo et al., 2021).

3 Geological setting of the study area

The study area comprises 110 184 km² within the coordinates −2.19–43.78°, 3.46–41,35° (Fig. 1). In addition to the main Pyrenees' zones, some areas at the eastern termination of the mountain range will also be described due to the occurrence of important magnetic anomalies like Olot, La Garrotxa volcanic field (e.g., Zeyen and Banda, 1988), and the western part of the Gulf of Lion passive margin (e.g., Canva et al., 2020) (see Fig. 1).

Figure 1 Major geological units differentiated in the Pyrenees and surrounding areas. Relief map of the study area. Data obtained from the IGN (Instituto Geográfico Nacional, Spain), GEOSERVICE (France), and EMODNET (European Marine Observation and Data Network, offshore). The location of the study area is shown in the upper right corner. Below the main map is the ECORS Pyrenees cross-section. The white lines on the upper map (modified from Muñoz, 1992) represent the approximate southern and northern boundaries of the profile. NPFT: North Pyrenean Frontal Thrust; NPF: North Pyrenean Fault; LT: Larra Thrust; MB: Morreres Back Thrust; LTF: La Têt Fault; SPT: Sub-Pyrenean Thrust; MT: Mazamet Thrust; TF: Toulouse Fault; SPFT: South Pyrenean Frontal Thrust; VPF: Vallès–Penedès Fault. UTM coordinates in meters, Zone 31N, ETRS89 datum.

The Pyrenees represent a WNW–ESE double-vergence mountain range with a dominant southward vergence resulting from the convergence between the Iberian and European plates (e.g., Muñoz, 1992). It developed as an asymmetric chain from the Late Cretaceous to the Early Miocene (e.g., Muñoz, 1992) (Fig. 2a). Deformation along the convergent plates' boundaries was accommodated by the underthrusting of the Iberian lithosphere underneath Europe and the thrusting and stacking of the upper crust (e.g., Choukroune and ECORS team, 1989; Roest and Srivastava, 1991; Rosenbaum et al., 2002; Wehr et al., 2018). The general structure, stratigraphy, and geodynamics of the Pyrenees have been extensively studied during the last 5 decades (e.g., Séguret, 1972; Garrido-Megías, 1973; Muñoz, 1992; Barnolas et al., 1996; Ford et al., 2022). In parallel, many geophysical studies have been carried out to infer the structure at the depth of the Pyrenees, mainly based on seismic, gravimetric, and magnetotelluric methods (e.g., Choukroune and ECORS team, 1989; Casas et al., 1997; Ledo et al., 2000; Chevrot et al., 2014, 2018; Wang et al., 2016; Campanyà et al., 2011).

The formation of the Pyrenees involved a complex polyphasic evolution, including the reactivation and/or inversion of previous Variscan and Mesozoic structures (e.g., Bond and McClay, 1995; Poblet, 1991; Muñoz, 1992; Clerc and Lagabrielle, 2014; Clerc et al., 2012; Ledo et al., 2000; Wehr et al., 2018; Le Maire et al., 2021; Cochelin et al., 2017). The Variscan orogeny took place during the Carboniferous associated with magmatic activity of a different nature (e.g., Zwart, 1986; García-Sansegundo, 1996; García-Sansegundo et al., 2011; Casas et al., 2019), while different extensional stages occurred from the Permian to Cretaceous, the latter being associated with the opening of the North Atlantic and Bay of Biscay (e.g., Roest and Srivastava, 1991; Ziegler, 1988).

The early stages of Atlantic rifting resulted in the emplacement of subvolcanic tholeiitic dolerites, known as ophites, derived from a subcontinental lithospheric mantle, which are now found as isolated outcrops within the Upper Triassic evaporites in the Pyrenees and other regions of Iberia (e.g., Beziat, 1983; Alibert, 1985). The Mesozoic extension culminated in the Pyrenees with the formation of a hyperextended margin during the Albian–Cenomanian, characterized by an extreme crustal thinning and subcontinental lithospheric mantle exhumation within the northern Pyrenees (e.g., Clerc and Lagabrielle, 2014; Jammes et al., 2009; Lagabrielle et al., 2010; Pedrera et al., 2017).

Figure 2 (a) Chronological and structural map (modified from IGME-BRGM, 2009). (b) Lithological and structural map (modified from IGME-BRGM, 2009) of the Pyrenees study area. A simplified version of Fig. 1 is shown in the upper right corner. NPFT: North Pyrenean Frontal Thrust; NPF: North Pyrenean Fault; ECT: Eaux-Chaudes Thrust; GT: Gavarnie Thrust; CF: Couflens Fault; MF: Mérens Fault; BT: Bono Thrust; BxT: Bóixols Thrust; MsT: Montsec Thrust; SMT: Sierra Marginales Thrust; LT: Larra Thrust; MB: Morreres Back Thrust; LTF: La Têt Fault; SPT: Sub-Pyrenean Thrust; MT: Mazamet Thrust; TF: Toulouse Fault; SPFT: South Pyrenean Frontal Thrust; VPF: Vallès–Penedès Fault. UTM coordinates in meters, Zone 31N, ETRS89 datum.

Geological data of the study (age of units, structure, and lithology) (Fig. 2a and b) were obtained from the cartography base map of the IGME-BRGM (2009) through the Geological Survey of Spain (IGME-CSIC; http://info.igme.es/cartografiadigital, last access: 22 April 2025).

3.1 Pyrenean main tectonic units

From north to south, five main tectonic units can be distinguished in the Pyrenees (Mattauer, 1968; Barnolas and Pujalte, 2004 and references therein) (Fig. 1): the Aquitaine Basin (AB), the North Pyrenean Zone (NPZ), the Axial Zone (AZ), the South Pyrenean Zone (SPZ), and the Ebro Basin (EB).

The Aquitaine Basin (Fig. 2a and b) represents the northern foreland basin developed on the overriding European plate. It has been extensively explored over the last 60 years (e.g., Biteau et al., 2006; Serrano et al., 2006). Its mid-Eocene to Miocene sedimentary fill can be attributed to a foreland subsidence directly linked to the growth of the Pyrenees (e.g., Biteau et al., 2006; Rougier et al., 2016). Geophysical data, like seismic profiles (e.g., Rougier et al., 2016), have revealed the existence of highly deformed Paleozoic basement rocks at depth and Mesozoic sub-basins related to an important rift system highly variable along-strike (e.g., Canérot et al., 2005). The Aquitaine Basin is mainly deformed at depth between the North Pyrenean Frontal Thrust (NPFT) and the Sub-Pyrenean Thrust (SPT).

The North Pyrenean Zone is characterized by thick inverted sequences of syn-rift to postrift Mesozoic rocks (e.g., Souquet et al., 1977; Déramond et al., 1988). This zone represents a relatively narrow fold-and-thrust belt, measuring between 25 and 35 km, which has developed along the southern boundary of the European (upper) plate between the North Pyrenean Fault (NPF) and the North Pyrenean Frontal Thrust (NPFT) (Fig. 2a and b). Regional high-temperature, low-pressure (HT–LP) metamorphism and about 40 outcrops of subcontinental peridotites are hosted in the northern Pyrenees related to the Early Cretaceous hyper-extension and mantle exhumation (e.g., Jammes et al., 2009, 2010; Clerc et al., 2012, Vauchez et al., 2012).

The Axial Zone (Fig. 2a and b) is located at the core of the alpine chain and comprises an antiformal stack of south-vergent basement thrust sheets. It mainly consists of metasedimentary Paleozoic rocks deformed and metamorphosed previously during the Variscan orogeny, along with Variscan granitoids (e.g., Barnolas et al., 1996; Casas et al., 2019). The Axial Zone elevation varies along its strike and exhibits a direct correlation with the depth of the European basement top (e.g., Soto et al., 2006; Saspiturry et al., 2019), indicating that the control of the deep crustal configuration of the Pyrenees lies on their structural architecture since the subduction is deeper in the center of the mountain chain and decreases laterally (e.g., Beaumont et al., 2000; Jammes et al., 2009; Saspiturry et al., 2019; Muñoz, 2019).

The South Pyrenean Zone constitutes a south-verging, thin-skinned fold-and-thrust belt detached from the basement along the Middle–Upper Triassic evaporites (e.g., Séguret, 1972; Muñoz, 1992). It is characterized by a Mesozoic succession with an along-strike thickness variation and a thick syn-orogenic Cenozoic sequence (e.g., Séguret, 1972; Garrido-Megías, 1973; Cámara and Klimowitz, 1985; Muñoz, 1992; Soto et al., 2002) (Fig. 2a and b).

Finally, the Ebro Basin (Fig. 2a and b) represents the southern foreland basin of the Pyrenees and encompasses a continuous sequence up to 5500 m of the latest Eocene, Oligocene, and Miocene continental sediments (e.g., Alonso-Zarza et al., 2002). Its autochthonous sediments are deformed by several folds and thrusts detached along the Eocene evaporites close to the South Pyrenean Frontal Thrust.

4 Data analysis and results

4.1 Reduction to the pole (RTP)

The reduced-to-the-pole magnetic field was calculated according to the declination and inclination parameters corresponding to the average date of data acquisition, 15 June 1973.

The magnetic properties of the different rocks are crucial for understanding the RTP map. These could include anything from sedimentary rocks to granulites or lherzolites. Lherzolites, granulites, basalts, and meta-basalts represent the most important magnetic-mineral-rich lithologies outcropping in the study area.

This way, as a consequence of the rocks, the resulting map displays highly defined and located anomalies above the source bodies (Fig. 3a). Magnetic-anomaly values range between −50.2 and 83.24 nT, while the mean value is −5.79 nT. It is noticeable that the overlaying of the magnetic-anomaly maxima is above the Bouguer anomaly maxima (Pedrera et al., 2017) (Fig. 3b and c). Despite the correlation between some western North Pyrenean Zone anomalies and some gravimetric highs along the Mediterranean coast, the gravimetric map displays a more dulled and homogeneous distribution of anomalies in comparison to the RTP magnetic map.

Figure 3 (a) Reduced-to-the-pole (RTP) magnetic map of the study area. (b) Bouguer anomaly map (modified from Pedrera et al., 2017). The main structures are also represented. A simplified version of Fig. 1 is shown in the upper right corner. MT: Mazamet Thrust; TF: Toulouse Fault; SPT: Sub-Pyrenean Thrust; NPFT: North Pyrenean Frontal Thrust; NPF: North Pyrenean Fault; LTF: La Têt Fault; MB: Morreres Back Thrust; LT: Larra Thrust; SPFT: South Pyrenean Frontal Thrust; VPF: Vallès–Penedès Fault. UTM coordinates in meters, Zone 31N, ETRS89 datum.

On the RTP map (Fig. 3a), wavelength values range from 7 to 1000 km. The Pyrenean anomalies can be grouped into magnetic zones that roughly coincide with the main Pyrenean structures.

In the North Pyrenean Zone (Figs. 1, 2b, and 3a and b), next to the NPFT, the high-intensity, short- to intermediate-wavelength anomalies can be associated with the outcrops of lherzolites, granulites, and ophites.

The Axial Zone and South Pyrenean Zone (Figs. 1 and 3a and b) are characterized by an elongated medium- to long-wavelength minimum with some short-wavelength, high-amplitude positive and negative anomalies on top. These anomalies could be related to the basement geometry and/or the accumulation of Triassic evaporites. Given that these anomalies have a NNW–SSE direction, they can also be related to Axial Zone structures. In the South Pyrenean Zone it is worth noting the so-called Lérida anomaly, a medium-wavelength high that could have originated from a basement high within the Ebro Basin (Zeyen and Banda, 1988).

Towards the northeast, next to the Catalan coast, there is a zone characterized by several maxima that correlates with outcrops of the Olot basaltic volcanic rocks of La Garrotxa and the volcanic field of L'Empordà and La Selva. Since these magnetic anomalies are larger than the outcrops of volcanic rocks, this could point to a wider volcanic area than previously thought, with the anomalies having also originated from subsurface volcanic rocks that could be part of a feeder system. In the Catalan Coastal Range (CCR), the relative maxima are related to outcrops of the Paleozoic basement. Towards its northern termination, their anomalies are masked by the magnetic signature of the La Garrotxa volcanic field (Olot).

4.2 Magnetic derivatives

Since it focuses on shallower sources, the vertical-derivative (VDR_Z) map (Fig. 4a) delimits near-surface anomalous bodies quite well. Three different magnetic-texture zones are distinguished: a set of NNW–SSE lineaments to the north of the southern border of the Axial Zone; a zone with relatively intermediate anomalies on the western side of the Ebro Basin; and NE–SW lineaments towards the eastern side of the EB and the Catalan Coastal Range as well as L'Empordà, La Selva, and La Garrotxa volcanic field.

The map of the horizontal derivative in x (HDR_X; Fig. 4b) displays predominantly N–S-oriented lineaments along the North Pyrenean Zone, the Axial Zone, and the Catalan Coastal Range. These zones also exhibit a strong ESE–WNW component on the map of the horizontal derivative in y (HDR_Y; Fig. 4c), whereas, L'Empordà, La Selva, and La Garrotxa volcanic field present a more chaotic magnetic pattern.

Horizontal derivatives in the N110° E direction (HDR_110) (Fig. 4d) and the N200° E direction (HDR_200) (Fig. 4e), which are the Pyrenees and their perpendicular main directions, enhance structures in those directions such as La Têt Fault or the Morreres Back Thrust, respectively. HDR_110 seems to be more defined over the Axial Zone and North Pyrenean Zone and especially defined eastwards, towards Olot. HDR_200 outlines the southern border of the Axial Zone almost perfectly but has a less concordant texture around the other structures.

The analytic signal (AS) (Fig. 4f) shows two different texture zones divided by an imaginary NW–SE line that seems to follow the trace of the Larra and Morreres thrusts, splitting an upper part characterized by rounded shapes above the main anomalies along the North Pyrenean Zone, the volcanic field (La Empordà, La Selva, La Garrotxa), and the Gulf of Lyon as well as a southern part corresponding to the Axial Zone, the South Pyrenean Zone and the Ebro Basin, with a diffuse and uneven texture.

Lineaments in the tilt derivative (TDR) map (Fig. 4g) appear to be outlining most of the main anomaly bodies that generate the magnetic anomalies on the RTP map (Fig. 3a). The intricate morphology of the North Pyrenean Zone, the volcanic field, and the Catalan Coastal Range anomaly sources are very well defined. It is worth noting that from the southern border of the Axial Zone northwards (Axial Zone, North Pyrenean Zone, Aquitaine Basin) it displays a NW–SE texture with an approximate N110° E direction, while a N–S orientation pattern prevails to the south of that boundary (South Pyrenean Zone, Ebro Basin), Olot (La Garrotxa) denoting the eastern end of both textures. Its horizontal derivative (TDR_THDR) (Fig. 4h) enhances the edges of the anomaly source bodies, exceptionally so along the central Pyrenees (See Sect. 2.2 “Magnetic derivatives”).

Figure 4 (a) Vertical-derivative map (VDR_Z). (b) Map of horizontal derivative in x (HDR_X). (c) Map of horizontal derivative in y (HDR_Y). (d) Horizontal-derivative map along N110° E (HDR_110). (e) Horizontal derivative map along N200° E (HDR_200). (f) Analytic signal (AS) map. (g) Tilt derivative (TDR) map. (h) Horizontal derivative of the tilt derivative map (TDR_THDR). A simplified version of Fig. 1 is shown in the upper right corner. In addition, the main lineaments of each derivative and main structures of the study area are also represented. MT: Mazamet Thrust; TF: Toulouse Fault; SPT: Sub-Pyrenean Thrust; NPFT: North Pyrenean Frontal Thrust; NPF: North Pyrenean Fault; LTF: La Têt Fault; MB: Morreres Back Thrust; LT: Larra Thrust; SPFT: South Pyrenean Frontal Thrust; VPF: Vallès–Penedès Fault. UTM coordinates in meters, Zone 31N, ETRS89 datum.

4.3 Euler deconvolution

In order to constrain the depth estimation of the magnetic source bodies, we have used the Euler deconvolution calculation. As indicated in the methodology, we have chosen the parameters SI = 0 (contacts/steps) with a window of 5 km and SI = 1 (sills/dikes) with a window of 10 km. Each window size was selected to ensure that the whole extension of the target bodies fits within.

For the contact/step solutions (SI = 0, WS = 5) (Fig. 5a) the depth estimation reaches up to 5 km to the north of the North Pyrenean Fault in different clusters, surrounding some of the RTP high positive anomalies, usually near the outcrops of lherzolites and granulites. The deepest solutions (from 3 to more than 20 km) are located in the western sector of the South Pyrenean Zone and the volcanic field, while shallower solutions (from 0.1 to 5 km) are found along the Ebro Basin and the Catalan Coastal Range.

For the sill and dike solutions (SI = 1, WS = 10), as the window size increases, so does the number of solutions (Fig. 5b). Despite showing a similar distribution to that in the previous case, their estimations are deeper than in the previous case. North of the North Pyrenean Fault, the solutions display deeper values to the west than to the east, gathered around the higher-value RTP anomalies in the zone. Along the volcanic field, deep solutions are found in areas without volcanic outcrops. Depth solutions in the western termination of the Ebro Basin, South Pyrenean Zone, Axial Zone, and North Pyrenean Zone (from 10 to more than 20 km) could be related to a deepening of the top of the basement. Values of 10 km are calculated at the border of the great offshore RTP anomaly in the Gulf of Lyon. In contrast, the shallowest estimations (0.1 km or less) are found in the center of the RTP anomaly in the Aquitaine Basin, crossed by the Sub-Pyrenean Thrust (Fig. 5b). However, the outcropping materials do not have a relevant magnetic response on the RTP map (Fig. 2b).

Figure 5 In order to highlight the most relevant magnetic anomalies, the RTP map is only portrayed above 0 nT in addition to (a) the Euler solution for SI = 0, WS = 5, and (b) the Euler solution for SI = 1, WS = 10. A simplified version of Fig. 1 is shown in the upper right corner. The main structures are also represented. MT: Mazamet Thrust; TF: Toulouse Fault; SPT: Sub-Pyrenean Thrust; NPFT: North Pyrenean Frontal Thrust; NPF: North Pyrenean Fault; LTF: La Têt Fault; MB: Morreres Back Thrust; LT: Larra Thrust; SPFT: South Pyrenean Frontal Thrust; VPF: Vallès–Penedès Fault. Volcanic materials, lherzolites, Lagabrielle et al. (2010) lherzolite and granulite mapping, Gulf of Lyon, La Selva, L'Empordà, La Garrotxa, and Olot are highlighted. UTM coordinates in meters, Zone 31N, ETRS89 datum.

4.4 Power spectrum

The radially averaged spectrum consists mainly of three components (Fig. 6) that give the maximum depth of the top of the magnetic basement. The first is a very steep segment (green) that shows a low wavenumber (0.001–0.016 km⁻¹). The maximum depth displayed for this segment reaches up to 21–23 km. The next segment (pink) ranges between wavenumbers of 0.016 to 0.084 km⁻¹. This segment displays the shallower sources and noise, reaching up to 15 km for the maximum depth. As for the third segment (blue), comprising wavelength values from 0.084 to 0.14 km⁻¹, shallower maximum depths are reached despite the peaks associated with noise, around 12 km.

Figure 6 Radially averaged power spectrum of the RTP magnetic field. The wavelength range spans from 7.14 to 1000 km. The maximum depth values reach 35 km.

5 Structural interpretation of the magnetic data

5.1 Magnetic fabric

The magnetic texture of the RTP map within the study area is complex and displays a pronounced heterogeneity in terms of source type, distribution, and shape (Fig. 3a, b, c). In order to highlight these parameters through the relationship between the magnetic transformations and the geological structures and lithologies, we have created three merged maps, which are displayed in Fig. 7: (1) a ternary image combining the RTP (Fig. 3a) grid, the analytic signal (Fig. 4f) grid, and the VDR_Z (Fig. 4a) together with the main lineaments of the HDR_X, HDR_Y, HDR_110, and HDR_200 (Fig. 4b, c, d, e) shown in Fig. 7a since they provide the best approach to the area's structural lineaments; (2) the geological map with the RTP anomalies above 0 nT to enhance the anomalous bodies' location and the gravimetric maximums (Fig. 7b); (3) the TDR (Fig. 4g) with the main structures and the Axial Zone mapping (Fig. 7c).

Figure 7 (a) Ternary image (RTP-AS-VDR_Z) together with the horizontal-derivative lineament mapping (color lines) and the selected magnetic anomalies (white lines). In the upper right corner, the triangle represents the color scale. (b) Lithological map of the study area in addition to the reduced-to-the-magnetic-pole above 0 nT mapping, gravimetric maximums (blue lines), and selected magnetic anomalies (red lines). A simplified version of Fig. 1 is shown in the upper right corner. (c) Comparison between the tilt derivative mapping (blue lines) and the main structures (black lines) with the TDR grid as a template and the AZ as a pattern mapping. The main structures are represented in black lines for the three maps. In blue: LA (Labourd), L (Lourdes), SG (Saint-Gaudens). In red: MT (Mazamet Thrust), TF (Toulouse Fault), SPT (Sub-Pyrenean Thrust), NPFT (North Pyrenean Frontal Thrust), NPF (North Pyrenean Fault), LTF (La Têt Fault), MB (Morreres Back Thrust), LT (Larra Thrust), SPFT (South Pyrenean Frontal Thrust), VPF (Vallès–Penedès Fault). I: Labourd; II: Lourdes; III: Saint-Gaudens; IV: Saint-Gaudens-2; V: Olot; VI: L'Empordà–La Selva; VII: east coast; VIII: Rascasse Horst. UTM coordinates meters, Zone 31N, ETRS89 datum.

The most outstanding anomaly zones have been selected with the following criteria. Firstly, as shown in Fig. 7a, those areas where high and similar values of RTP, analytic signal, and VDR_Z coincide (displayed in black) are regarded as sound evidence of subsurface anomalous bodies. However, their shape is more accurately constrained by the TDR_THDR peak values (white lines in Fig. 7a), which nearly matches with the shape of the black zones (Fig. 7a). Then, high RTP values (from 20 to 50 nT, reaching 70 nT in zone VIII and near 85 nT in zone IV) and the concurrency with Bouguer anomalies (anomalies I, II, III) (Fig. 7b) were also considered for the selection. These anomalies, labeled I to VIII in Fig. 7a and b, correlate with outcrops of magnetic bodies but also allow the buried part of these source bodies to be delimited.

Mapping the main lineations of the horizontal derivatives (x, y, N200° E, and N110° E direction) (Fig. 7a) allows the magnetic and tectonic fabric to be illustrated. Lineaments with a strong N–S component such as HDR_X and HDR_110 are mostly grouped in the eastern area, while those with a predominant E–W component, HDR_Y and HDR_200, are assembled along the Axial Zone and prevail in the western sector. This could imply that the presence of predominantly E–W structures to the west could be linked to the compressive structures related to the Pyrenees formation (e.g., Muñoz, 1992), while the predominantly N–S structures could be associated with extensive structures generated as a consequence of the opening of the northwestern Mediterranean Sea (e.g., Maillard and Mauffret, 2001).

The N–S structures highlighted by the HDR_X lineaments outline the edges of some anomalies that could be related to evaporites to the southeast of the Ebro Basin (Zeyen and Banda, 1989), while the easternmost ones next to the Catalan Coastal Range seem to be affected by the volcanic field of L'Empordà, La Selva, and La Garrotxa. This N–S lineation may be related to the volcanic materials' disposition at depth since the most notorious anomalies linked to the volcanic materials appear to the east of the lineation, while the western side appears to be less magnetic (Fig. 7a and b), separating La Garrotxa to the west and L'Empordà and La Selva to the east. This layout may be related to a fault system setting that places shallow bodies to the east of the lineaments, whereas western magnetic source bodies may be located in structurally deeper positions. Additionally, another N–S lineation seems to outline the northwestern edge of the Catalan coast anomaly and could be linked to the intrusion of the magnetic source materials (anomaly VII). Regarding the N110° E direction lineations (Figs. 4d and 7a), the south–central ones appear to constrain the evaporite bodies alongside the HDR_X and HDR_Y lineaments. On both sides of La Têt Fault (LTF), lineations with a N110° E direction are mapped (Fig. 7a), which we interpret to imply the presence of similar structures along the eastern end of the Axial Zone.

On the other hand, the lineations with an E–W component (HDR_Y and HDR_200) constrain the configuration of the Axial Zone in its western and central parts, where the thickness of the Pyrenean crustal prism is greatest (e.g., Soto et al., 2006; Pedrera et al., 2017). HDR_Y lineations mimic the pathway of the North Pyrenean Fault to the north, as well as the Larra Thrust and the Morreres Back Thrust to the south. When combined with the HDR_200 lineaments, they outline the anomaly of the southern border of the Axial Zone as well as the morphology of the Saint-Gaudens anomaly (III). To the south, they define a parallel texture which may be related to the presence of evaporites (e. g. Santolaria et al., 2024).

Lastly, the tilt derivative peak mapping exposes the pathway of the structures such as the North Pyrenean Frontal Thrust, the Morreres Back Thrust, and the Larra Thrust, which delimit the Axial Zone (Fig. 7c). In addition to describing these main structural features, this method portrays structural fabric inside the Axial Zone. A N200° E lineament domain is displayed westward and turns into a N110° E domain eastward, adjacent to the La Têt Fault (Fig. 7a and b).

5.2 Geological interpretation of the positive magnetic anomalies

The correlation between the magnetic response of the most remarkable magnetic anomalies (Fig. 8a and b) and their geological context (Fig. 2a) can be interpreted as follows.

Figure 8 Reduced-to-the-pole magnetic map of the study area with the (a) the Euler solution for SI = 0, WS = 5, and (b) the Euler solution for SI = 1 WS = 10 on top. A simplified version of Fig. 1 is shown in the upper right corner. The main structures are also represented. MT: Mazamet Thrust; TF: Toulouse Fault; SPT: Sub-Pyrenean Thrust; NPFT: North Pyrenean Frontal Thrust; NPF: North Pyrenean Fault; LTF: La Têt Fault; MB: Morreres Back Thrust; LT: Larra Thrust; SPFT: South Pyrenean Frontal Thrust; VPF: Vallès–Penedès Fault. Volcanic materials, lherzolites, and Lagabrielle et al. (2010) lherzolite and granulite mapping are displayed. I: Labourd; II: Lourdes; III: Saint-Gaudens; IV: Saint-Gaudens-2; V: Olot; VI: L'Empordà–La Selva; VII: east coast, VIII: Rascasse Horst.

5.2.1 Aquitaine Basin

Located on the southern border of the Aquitaine Basin, anomaly IV has previously been interpreted by Le Maire et al. (2021) to be linked to a slice of subcontinental mantle. This magnetic anomaly does not coincide with lherzolite or granulite outcrops (Figs. 2b and 7b), and it is associated with a mild gravimetric response (Fig. 3b and c); it may correspond to a portion of serpentinized mantle. Nonetheless, this magnetic anomaly coincides with the inferred trace of the Sub-Pyrenean Thrust (SPT) (Angrand et al., 2018).

On the RTP map, anomaly IV displays a magnetic response up to +83 nT, the highest in the study area. According to the Euler solutions, it exhibits the shallowest depths for SI = 0 (0.1 to 3 km), suggesting the presence of a body at those depths (Figs. 5a and 8a). SI = 1 shows a well-shaped circular geometry with an average depth up to 7 km (Figs. 5b and 8b). According to subsurface data, coinciding with the trace of the Sub-Pyrenean Thrust (SPT), Rougier et al. (2016) already interpreted the presence of Upper Albian–Cenomanian volcanic rocks at depth, which may have used previous normal faults as a conduit. Therefore, we interpret magnetic anomaly IV to be related to a highly magnetic and low–medium-density, vertical, volcanic, dike-like body as the source.

5.2.2 North Pyrenean Zone

The North Pyrenean Zone is characterized by the presence of magnetic anomalies I, II, and III. Their traces match exceptionally well with the outline of the gravimetric anomalies (Fig. 7b), also located between the NPFT and the NPF. The shape of these anomalies displayed on the RTP map is a WNW–ESE-oriented elongated shape (Fig. 3a); however the ternary image (Fig. 7a) suggests a more uneven shape that could reveal the source extension at greater depths.

Regarding the Euler solutions, these anomalies are not well defined for SI = 0. For SI = 1, the deep solutions are found approximately at the edges of the anomalies (up to 15 km), commonly between 1 and 5 km, and reach more than 20 km for the northwestern areas of anomaly I (Fig. 8a and b). In contrast, for SI = 1 the shallow depths are registered at the center of the clusters (up to 1 km) (Figs. 5b and 8b).

The Labourd anomaly (I) appears to be continued towards the northwest, probably deepening, as depicted in the Euler solutions (SI = 0 and SI = 1) (Fig. 8a and b). The Lourdes anomaly (II) seems to be connected to a small, southern, WNW–ESE-oriented body that appears to be at an intermediate depth (1 to 3 km) from the Euler solutions (Fig. 8a and b). The Saint-Gaudens anomaly (III) has an elongated shape that, in the Euler solutions, is delineated by a string of solutions with depths ranging between 0.6 and 7 km.

From a geological point of view, the North Pyrenean Zone is characterized by the presence of several outcrops of lherzolites, granulites, and serpentinites (see the mapping compilation of Lagabrielle et al., 2010) that have already been linked to magnetic anomalies I, II, and III (e.g., Casas et al., 1997; Le Maire et al., 2021), as their outcrops cover the extension of these anomalies (Fig. 7b). These rocks consist of basic and ultrabasic mantle and lower-crustal rocks emplaced at upper-crustal levels related to an extreme crustal thinning, and subcontinental lithospheric mantle exhumation occurred within the northern Pyrenees during the Albian–Cenomanian (e.g., Clerc and Lagabrielle, 2014; Jammes et al., 2009; Lagabrielle et al., 2010).

5.2.3 Axial Zone

As shown in Figs. 1a and 3a, the Axial Zone does not exhibit any remarkable anomalies along the Pyrenees mountain range (anomalies range between −20 and 0 nT). The absence of highly magnetic bodies derived from the mantle is consistent with the structural architecture of the Axial Zone, formed by a thick alpine stack of basement thrust sheets. Therefore, in contrast to the North Pyrenean Zone, the Axial Zone's mild magnetic response may be on account of an increase in the upper-crustal thickness and the fact that the granitoids and other intrusive rocks have very low magnetic susceptibility, most of them falling within the paramagnetic domain (Pueyo et al., 2022; Porquet et al., 2017). Some of the WNW–ESE magnetic lineaments present on the western side of this zone could be related to the deformed and elongated metamorphic Variscan domes (Cochelin et al., 2017).

Nevertheless, anomaly VII, located on the eastern coastline, presents values of at least +20 nT. Its trace extends northward from the main RTP anomaly, as portrayed in the ternary image (Fig. 7b), and represents the transition zone between the Axial Zone, the Mediterranean coast, and the westernmost part of the Gulf of Lyon. The outcrops of Cambrian meta-basaltic rocks (Debon et al., 1996; Padel et al., 2018), as well as buried Pleistocene volcanic rocks, which show small outcrops in this area, may be the anomaly source (Carreras and Druguet, 2014; Carreras et al., 1988).

Zone VII is not well defined in the Euler solutions. The values of the solutions (Fig. 8a and b) are shallower for SI = 0 (up to 5–7 km) and deeper for SI = 1 (up to 7–10 km), yet in both cases the source of the anomalies deepens towards the north.

It is important to highlight the existence of low-intensity anomalies that have the same orientation as the LTF to its west (Fig. 3a) as well as the northern area of anomaly V, which reaches depths of 3–5 km, probably indicating the presence of a bigger spread of volcanic materials beneath the surface (Fig. 7b).

5.2.4 South Pyrenean Zone

The RTP map (Figs. 1a and 3a) reveals a rather low magnetic response for the SPZ. Its south–central area displays values up to 6 nT (Fig. 7b). In this area, small and dull responsive magnetic anomalies have been linked to the presence of underlying Triassic ophites by Zeyen and Banda (1989), but we do not have enough data to confirm this hypothesis. The Euler solutions of SI = 0 mainly follow the SPFT with depths that can reach up to 2 km westwards. By contrast, SI = 1 describes two distinct areas, a deeper area (7 to 20 km) in the western half and a shallower area (0.2 to 10 km) to the east (Fig. 8a and b). Magnetic-anomaly sources within the area located at depths below 2–6 km (SI = 1) might be linked to sources located in the footwall of the South Pyrenean Thrust according to the geometry interpreted as Pyrenean crustal cross-sections (see Teixell et al., 2018).

5.2.5 Catalan Coastal Range

Characterized by high positive anomalies (20–50 nT) on the RTP map (Figs. 1a and 3a), the CCR structure is also highlighted by the derivatives (e.g., Fig. 4g) and the ternary image (Fig. 7a). These are linked to anomalies V and VI in Fig. 7a and b.

Anomaly V, located around the town of Olot, corresponds to the outcrop of basalts and basanites of the La Garrotxa volcanic unit (Zeyen et al., 1991). These Neogene and Quaternary volcanic materials are genetically correlated with those of anomaly VI. The latter compromises the L'Empordà and La Selva volcanic units, which are mainly Neogene (Dèzes et al., 2004; Araña et al., 1983; Zeyen et al., 1991). These volcanic terrains that extend offshore involve few volcanic outcrops; however, there is a large number of granitoids and intermediate to ultramafic plutonic rocks outcropping (Enrique, 1990). We interpret the high RTP anomaly values (up to 50 nT) to be related to the underlying basalts and outcrops of ultramafic rocks.

For the Euler solutions, anomaly V seems to display shallow solutions towards the edges (<0.1–0.2 km for SI = 0 and up to 1–5 km for SI = 1; Fig. 8a and b). Anomaly VI also displays depth values from 1 to 5 km for SI = 0 (Fig. 8a) around the edges, while depths of 5–10 km (SI = 1) are calculated for the volcanic materials (Fig. 8b), reaching up to 20 km in the westernmost areas.

5.2.6 Ebro Basin

The RTP map displays anomaly values above 0 nT with a mainly N–S orientation for most of the area, mostly on its western and eastern edges and especially so in the center, extending from those previously described to the south of the South Pyrenean Zone. The central anomalies located in the South Pyrenean Zone seem to extend southwards, showing the lowest magnetic-anomaly values across the area. However, westward within the area, a high anomaly value (+10 nT) appears to portray the northern edge of the Iberian chain (Figs. 1a and 3a) that could be related to a shallower top of the basement.

5.2.7 Gulf of Lyon

Anomaly VII represents the abrupt eastern end of the Pyrenees mountain range as the magnetic fabric changes drastically in the northwestern Mediterranean Sea (Canva et al., 2020). In addition to the contribution of the Cambrian metabasalts (Debon et al., 1996; Padel et al., 2018), this anomaly is suggested to be related to the Oligocene–Miocene rifting that opened the western Mediterranean Sea, with a progressive thinning of the crustal thickness from the Pyrenees to the Mediterranean Sea (Mauffret et al., 2001) that also generated the Pleistocene volcanic rocks. The presence of the Catalan Transfer Zone as well as other NW–SE structural lineations, crucial during the opening of the Gulf of Lyon (Canva et al., 2020), may be genetically linked to the NW–SE magnetic anomalies, such as VIII.

This anomaly comprises the Rascasse Horst, part of a NE–SW horst–graben display of thinned continental crust (Gorini et al., 1994; Séranne et al., 1995; Gailler et al., 2009; Moulin et al., 2015; Canva et al., 2020). The absence of volcanic evidence and the presence of a deep Moho point to a relatively thicker crustal block with mafic intrusions (Canva et al., 2020) as the source for anomaly VIII. The RTP intensity reaches values up to +70 nT, which represents the second-highest anomaly values in the study area (Figs. 3a and 7b). The anomaly VIII outline displays a homogeneous depth around 1–3 km for the Euler solutions of SI = 0 (Fig. 8a), as previously described by Canva et al. (2020).

6 Conclusions

The EMAG2v2 magnetic data transformations for the Pyrenees and their foreland basin area have provided new insights into the geological interpretation of the magnetic pattern. The magnetic response for the different domains has portrayed textures that are consistent with the main Pyrenean lineaments. The qualitative interpretation of the magnetic-field data has resulted in a comprehensive delineation of eight main distinctive magnetic anomalies: one in the Aquitaine Basin (anomaly VI), three along the North Pyrenean Zone (anomalies I, II, and III), one in the Axial Zone (anomaly VII), two in the Catalan Coastal Range (anomalies V and VI), and one in the Gulf of Lyon (anomaly VIII). Through the Euler solutions and power spectrum calculations, depth estimations for each zone have been calculated according to the most relevant morphologies, contacts, and sills/dikes. The integration of the geological map and the magnetic transformations has shed light on the intricate pattern of structures and materials that comprise the selected anomalies. Furthermore, the utilization of the tilt derivative method has provided a delineation of the borders of the Axial Zone materials, even though geologically speaking, in the western part the contact between the Axial Zone and the SPZ is indicated by the presence or absence of the sedimentary cover, while the combination of the horizontal-derivative lineaments has exposed two different structural domains. Including more detailed studies such as petrophysical, structural, or geophysical surveys around the areas and anomalies depicted in this work can reduce the uncertainties in the geophysical interpretation provided in this work, and it is expected to be part of our future works. Hence, this study is put forward to lay the groundwork for future research aimed at defining the magnetic fabric and deep structure of the Pyrenees mountain range in specific areas that could be of interest, for instance, to investigate its geothermal potential or mineral resources.