Early Cretaceous paleomagnetic and geochronologic results from the Tethyan Himalaya: Insights into the Neotethyan paleogeography and the India–Asia collision

To better understand the Neotethyan paleogeography, a paleomagnetic and geochronological study has been performed on the Early Cretaceous Sangxiu Formation lava flows, which were dated from ~135.1 Ma to ~124.4 Ma, in the Tethyan Himalaya. The tilt-corrected site-mean characteristic remanent magnetization (ChRM) direction for 26 sites is Ds = 296.1°, Is = −65.7°, ks = 51.7, α95 = 4.0°, corresponding to a paleopole at 5.9°S, 308.0°E with A95 = 6.1°. Positive fold and reversal tests prove that the ChRM directions are prefolding primary magnetizations. These results, together with reliable Cretaceous-Paleocene paleomagnetic data observed from the Tethyan Himalaya and the Lhasa terrane, as well as the paleolatitude evolution indicated by the apparent polar wander paths (APWPs) of India, reveal that the Tethyan Himalaya was a part of Greater India during the Early Cretaceous (135.1–124.4 Ma) when the Neotethyan Ocean was up to ~6900 km, it rifted from India sometime after ~130 Ma, and that the India-Asia collision should be a dual-collision process including the first Tethyan Himalaya-Lhasa terrane collision at ~54.9 Ma and the final India-Tethyan Himalaya collision at ~36.7 Ma.

. Sketches of geology and sampling location for this study. (a), Regional geologic map of the Himalayan belt and adjacent areas modified from Yin and Harrison 1 and Yin 60 . Solid squares and circles show sampling locations of previous Cretaceous and Paleocene paleomagnetic studies on volcanic and sedimentary rocks, respectively (for sample location abbreviations see Table 2 Cretaceous paleomagnetic data are still necessary. Considering that the sedimentary rocks often suffer from compaction-induced inclination shallowing 20,[30][31][32] , whereas the volcanic rocks are immune from its effect, we carried out a combined geochronologic and paleomagnetic study on the Early Cretaceous Sangxiu Formation (Fm) lava flows in the Tethyan Himalaya. These new high-quality and well-dated paleomagnetic data can significantly attribute to the two issues mentioned above.
The Sangxiu Fm, which is defined as an intercalated volcanic-sedimentary sequence, is only distributed in the southeast of the Yangzhuoyongcuo Lake along the eastern part of central Tethyan Himalaya (Fig. 1b). It conformably overlies the Upper Jurassic Weimei Fm and underlies the Lower Cretaceous Jiabula Fm [1:250,000 scale Luoza county regional geological survey report (H46C004001), 2002]. The age of the Sangxiu Fm volcanics is ~133 Ma indicated by SHRIMP U-Pb zircon dating 33 . The earliest folding of the Sangxiu Fm occurred in the latest Early Cretaceous [H46C004001, 2002]. The sedimentary strata include coarse to fine-grain quartz sandstone, quartz greywacke, siltstone and shale. The volcanic strata consist of pillowed and massive, sparsely amygdaloidal basalts and minor dacites.
A total of 32 paleomagnetic sites were sampled from two sections located in the southeast of the Yangzuoyongcuo Lake (Fig. 1b) and ~150 km northwest of the sampling area of the Lakang Fm lava flows 12 . Twenty-three sites were collected from section A located at (28.8°N, 91.3°E), and another 9 sites were sampled from monoclinal strata of section B at (28.8°N, 91.1°E). Each sampling site spans several meters of stratigraphic thickness and covers at least one lava flow. The bedding attitudes of both sampling sections are obvious and can be well determined by measuring the intercalated sedimentary rocks ( Supplementary Fig. S1). Furthermore, two block samples of fresh volcanic rocks from sites ZL1 and ZL23 of section A were collected for zircon U-Pb chronology.

Results
U-Pb Zircon Geochronology. The zircon grains are euhedral to subhedral prism (50-200 μ m in length) with an aspect ratio of ~1-3 ( Supplementary Fig. S2). This, together with clear oscillatory zonings in Cathodoluminescence images, indicates a magmatic origin. Oscillatory zonings where without inclusions and cracks were selected and analyzed for each sample. Zircon U-Pb analyses yielded diverse age groups, indicating different sources for the zircons (Fig. 2). The weight mean 206 Pb/ 238 U ages of the main population are interpreted as the time of emplacement of the studied volcanics. Other ages of minor populations, which include less than 5 zircon grains, may yield from inherited or contaminative zircons. Samples ZL1 and ZL23 yield weighted Scientific RepoRts | 6:21605 | DOI: 10.1038/srep21605 mean 206 Pb/ 238 U ages of 135.1 ± 0.7 Ma and 124.4 ± 0.7 Ma, respectively. These new ages are well consistent with the age of ~133 ± 3 Ma reported by Zhu et al. 33 , indicating that the Sangxiu volcanics erupted during the Early Cretaceous. We used ~135.1-124. 4 Ma for the Sangxiu Fm in further discussion.
Rock magnetic results. Isothermal remanent magnetization (IRM) acquisition curves of the representative specimens rise very quickly below 200 mT and the saturation is essentially reached at ~160-300 mT, revealing that low-coercivity magnetic carriers are dominant ( Supplementary Fig. S3a,b). Progressive demagnetization of the saturation IRM (SIRM) by applying reverse fields indicates that the maximum coercive force is less than 60 mT. This, together with a clear Curie temperature at ~580 °C, indicates that the low-coercivity magnetite is dominant in the Sangxiu Fm volcanic samples ( Supplementary Fig. S3c,d). Hysteresis loops close in a magnetic field at ~250 mT, also indicating relatively low coercivity ( Supplementary Fig. S4a,b). The Day plot 34 reveals that the magnetic grains are located within the region of pseudo-single-domain (PSD) (Supplementary Fig. S4c). First-order reversal curves (FORCs) 35 show closed peak structures with relatively open contours, suggesting the presence of PSD magnetite and minor single-domain (SD) magnetite ( Supplementary Fig. S4d,e). Both PSD and SD grains are very efficient carriers of remanent magnetization, so they should most probably carry a stable remanence when the Sangxiu Fm lava flows were erupted.
Paleomagnetic Results. About 80% of specimens underwent stepwise thermal demagnetization, and about 20% were performed stepwise alternating field (AF) demagnetization. Both demagnetization methods yield the same ChRM directions (Fig. 3). Some volcanic specimens give a low-temperature component (LTC) below ~250 °C or a low-coercivity component (LCC) below ~20 mT. Both LTC and LCC in geographic coordinates are closed to the present-day geomagnetic field direction. After removing the LTC or LCC, a high-temperature component (HTC) or a high-coercivity component (HCC) can be isolated from most specimens, which is defined as the ChRM directions, between ~400 °C and ~580 °C or between ~30 mT and ~100 mT. The ChRM directions decay toward the origin and include antipodal normal and reverse polarities (Fig. 3). However, some specimens from sites ZL6, ZL16, ZL17 and ZL18 display erratic demagnetization patterns (Fig. 3p), and no reliable ChRM directions can be isolated from them. Based on the following filtering criteria to ChRM and site-mean directions:  Table 1 and the corresponding ChRM directions of 230 specimens are presented in Supplementary Table S1. Notably, the tilt-corrected ChRM directions of sites ZL4 and ZL8 show >45° angular deviations from the overall-mean direction ( Table 1), indicating that they probably recorded a transitional or excursional paleomagnetic direction. Therefore, these two sites are discarded for final analysis. Finally, twenty-six paleomagnetic sites yield an overall-mean direction of Dg = 268.7°, Ig = − 47.4°, kg = 1.9, α 95 = 30.0° in situ and Ds = 296.1°, Is = − 65.7°, ks = 51.7, α 95 = 4.0° after tilt correction ( Table 1, Fig. 4). This overall-mean direction passes both McElhinny 38 and McFadden 39 fold tests at the 95% and 99% confident level, which, combined with that the reversal test is also positive at the 95% confidence level 40 (Table 1), indicates that the ChRMs have a pre-folding origin and are probably primary magnetization acquired during the eruption of the Sangxiu Fm lava flows. The Fisherian site-mean paleopole for 26 sites is located at 5.9°S, 308.0°E with A 95 = 6.1°, corresponding to a paleolatitude of 48.5° ± 6.1°S for the study area (28.8°N, 91.3°E).
Paleomagnetic data using for paleogeographic and tectonic reconstructions must average paleosecular variations. Based on the following evidences: (1) the sampling of lava flows spans a long time from ~135.1 Ma to ~124.4 Ma as determined by U-Pb zircon ages; (2) the sampling sections span many lava flows interbedded with sedimentary rocks; (3) the ChRM directions include antipodal dual polarities; (4) twenty-six lava sites provide a virtual geomagnetic pole (VGP) scatter of 17.1° at ~48.5°S, which matches with the paleosecular variation model at similar paleolatitude 41 ; (5) the A 95 , which is obtained from the VGPs of 26 lava sites, is 6.1°, which is well consistent with a N-dependent A 95 envelope with a 95% confidence interval (3.3°, 10.5°) proposed by Deenen et al. 42 , we confidently conclude that the Fisherian site-mean pole (5.9°S, 308.0°E with A 95 = 6.1°) obtained from the 26 Sangxiu Fm lava flows has averaged paleosecular variation. Therefore, it should be a reliable Early Cretaceous pole for the Tethyan Himalaya.

Discussion
The present geological boundary between India and Asia is along relatively east-west ITSZ that separates the Tethyan Himalaya from the Lhasa terrane (Fig. 1a), where the Neo-Tethyan Ocean opened in the Late Triassic and reached its greatest width in the Early Cretaceous 1 . Therefore, the paleolatitude evolutions of the Indian craton, the Tethyan Himalaya and the Lhasa terrane have been widely used to constrain the Cretaceous-Paleogene paleogeography of the Neotethyan Ocean and the India-Asia collision process 12,14,17,43 . Because the present-day ITSZ stretches a length of ~2000 km from east to west (Fig. 1a), a reference point (29°N, 87.5°E) located in its middle part has been used to calculate the expected and observed paleolatitudes ( Table 2).
For the Lhasa terrane, considering that 1) the possible inclination shallowing is still a critical and unresolved problem for the validity of paleomagnetic data from sedimentary rocks 20,30-32 ; 2) the Paleogene paleomagnetic data are mainly from the Linzizong Gp volcanic rocks which span a long time range from ~69 to ~40 Ma and might be erupted after the Tethyan Himalaya-Lhasa terrane collision age (such as the Panna Fm volcanic rocks) 31,44,45 ; and 3) the Lhasa terrane accreted onto the Qiangtang terrane by the Early Cretaceous and its southern margin maintained a stable paleolatitude during the whole Cretaceous 11,12,17,20,43 , in this study we use two reliable Cretaceous volcanic paleomagnetic poles, which come from a large number of volcanic sites with clear bedding attitudes and satisfy all the 7-point data quality criteria proposed by Van der Voo 46 , to constrain the paleolatitude of the Lhasa terrane as Huang et al. 47   Yanhu (QS) 17 and Cuoqin (CQ) 11,20 areas yield paleolatitudes of 18.3° ± 2.1°N and 14.9° ± 5.5°N for the reference point (29°N, 87.5°E), respectively (Table 2), and their mean paleolatitude of 16.6° ± 5.9°N should be a credible estimate for the precollisional southern margin of Asia.
For the Tethyan Himalaya, only nine Cretaceous-Paleocene poles were available due to serious remagnetization (Table 2). Three poles WL 29 , TD 25 and GB 15 from the Early Cretaceous sedimentary rocks don't provide a robust field test, and two Paleocene poles BS 26 and TY 27 are from too small sites (specimens) to pass the basic selection criterion of specimen (site) number. Therefore, these five poles are discarded for further discussion. The remaining two volcanic (LK 12 and SX (this study)) and two limestone (ZS 14,15 and ZP [14][15][16] ) poles satisfy the 7-point data quality criteria proposed by Van der Voo 46 . Noticeably, Dupont-Nivet et al. 48 applied the Elongation/ Inclination (E/I) correction method to these limestone paleomagnetic data from the southern Tethyan Himalaya, and suggested that they did not suffer from significant inclination shallowing. Therefore, in this study we use these four paleomagnetic poles to position the Tethyan Himalaya.   49 , we made a 1.6° N-S shortening correction to the LK dataset whose sampling area is located at the southern margin of the Tethyan Himalaya (Fig. 1a). Because our sampling area is located at the middle Tethyan Himalaya, a 1.1° (~120 km) N-S shortening was corrected for the SX dataset. Such the N-S shortening corrections would relocated the northern margin of the Tethyan Himalaya at a more northerly paleolatitude of 46.5° ± 5.7°S at ~134-131 Ma and 44.2° ± 6.1°S at ~135.1-124.4 Ma (Fig. 5). Considering that these two paleolatitudes are very consistent within the paleomagnetic confidence level, as well as that the age range of the SX pole entirely covers that of the LK pole, we use their mean paleolatitude of ~45.4°S as a credible estimate for the Tethyan Himalaya at the reference point (29°N, 87.5°E) (Fig. 5). Comparing it with ~16.6°N for the southern margin of the Lhasa terrane shows a paleolatitude difference of ~62.0°, implying that the Neotethyan Ocean between the Tethyan Himalaya and the Lhasa terrane opened a latitudinal width up to ~6900 km during ~135.1-124.4 Ma (Fig. 5). This width is also consistent with ~6700 km and ~7000 km estimated by Chen et al. 11 and Yang et al. 12 , respectively.
Comparing the Early Cretaceous (135.1-124. 4 Ma) paleolatitude of ~45.4°S observed from the Tethyan Himalaya with the ~49.6°S and ~50.2°S predicted by the most widely used 23 and most recent global synthetic 24 Indian APWPs at 130 Ma for the reference point (29.0°N, 87.5°E) indicates a paleolatitude difference of ~4.2° (~470 km) and ~4.8° (~530 km), respectively. Although the ~4.2°-4.8° (~470-530 km) crustal shortening deduced from high-quality Early Cretaceous paleomagnetic data is within the paleomagnetic resolution (~6° or ~670 km), the most possible crustal shortening amount (~470-530 km) deduced from these high quality paleomagnetic results is basically consistent with ~480-650 km, which includes ~176 km within the Tethyan Himalaya and the ITSZ 49 and ~302-476 km between the MFT and the STDS 50 , estimated by balanced cross-section analyses, supporting that the Tethyan Himalaya belonged to a contiguous Indian continental lithosphere at 130 Ma 12,14,29 .
Based on the Late Cretaceous and Paleocene paleomagnetic data 15,16 reanalyzed by van Hinsbergen et al. 14 , the Tethyan Himalaya was located at paleolatitudes of 5.1° ± 3.5°S at ~68 Ma and 8.7° ± 1.7°N at ~59 Ma for the reference point (29°N, 87.5°E), respectively ( Table 2). A ~176 km (~1.6° of latitude) N-S shortening correction 49 would relocated the northern margin of the Tethyan Himalaya at paleolatitudes of 3.5° ± 3.5°S at ~68 Ma and 10.3° ± 1.7°N at ~59 Ma (Fig. 5). Comparing them with the two coeval paleolatitudes of ~21.0°S and ~10.0°S calculated from the Indian APWPs 24 reveals a paleolatitude difference of ~17.5° (~1940 km) for ~68 Ma and ~20.3° (~2250 km) for ~59 Ma, respectively. Obviously, the paleolatitude difference observed from the high-quality Late Cretaceous and Paleocene data is much greater than ~4.8° (~530 km) determined by the high-quality Early Cretaceous (~135.1-124.4 Ma) volcanic data and ~480-650 km estimated by balanced cross-section analyses 49,50 , as well as far more than the Indian subcontinent extension of ~500-950 km deduced from the fitting of India in-Gondwana and an analysis of bathymetric features in the eastern Indian Ocean 8 . Therefore, such a great paleolatitude gap between the Indian craton and the Tethyan Himalaya cannot wholly be attributed to the continental crustal shortening within the Himalaya terrane after the India-Asia collision [15][16][17]51 , as well as to an ocean existed between the Indian craton and the Tethyan Himalaya even earlier in the Early Cretaceous 13 and it subducted beneath the Tethyan Himalaya after the Eocene. However, it can be interpreted as a small post-Neotethyan ocean (basin) extended between the Indian craton and the Tethyan Himalaya after the Early Cretaceous, and it subducted beneath the Tethyan Himalaya after the Eocene 12,14,20 . Such a Cretaceous extension model implies that the India-Asia collision should be a dual-collision process, which consists of a first collision occurred between the Tethyan Himalaya and the Lhasa terrane and a final collision occurred between the India craton and Tethyan Himalaya (Fig. 5). The collisional age of the Tethyan Himalaya and the Lhasa terrane can be determined by their paleolatitudes overlaps. Based on high-quality Late Cretaceous and Paleocene paleomagnetic data observed from the Tethyan Himalaya, it moved northward at a velocity of ~17.0 cm/yr during ~68-59 Ma. Extrapolating the constant northward velocity ~17.0 cm/yr implies that the Tethyan Himalaya would intersect with the Lhasa terrane at ~54. 9     Because the crustal shortening in Asia is ~600-750 km in the last 50 Ma 52 (i.e. ~12-15 km/Myr), the shortening amount for the time interval of ~50-38 Ma can be up to ~144-180 km (~1.3-1.6° of latitude). We conservatively add 1.6° of latitude to the leading edge of the Tethyan Himalaya-Lhasa terrane for ~38-36 Ma to estimate the effect of postcollisional shortening within Asia (Fig. 5). A final collision between the India Craton and the Tethyan Himalaya-Lhasa terrane was occurred at ~22.6°N at ~36.7 Ma (Fig. 5). It is important to note that dating the second collision largely depends on the paleolatitude observed from the upper Linzizong Fm volcanics and crustal shortening in the Asia and the Himalaya terrane prior to the second collision.
Critically, such a two-stage collision model also match with many geological evidences, such as: (1) the northward drift velocity of the Indian craton sharply decreases at ~55 Ma 3,4 and ~40 Ma 24 ; (2) a small post-Neotethyan Ocean opened during the Late Cretaceous and then subducted at the time interval of ~55-37 Ma, which corresponds to a two-stage collision model 12,20 , can explain the reason why the 'missing' convergence didn't be documented by geological record of Asia and the Himalaya 14 ; (3) the Asian assemblage arrived on the Tethyan Himalaya during the Early Eocene 53 , while the Tethyan Himalaya detritus reached the Indian foreland region during the Middle Eocene 6 ; (4) the subduction-related ultrahigh-pressure rocks in the Himalayan terrane dated at the Early Eocene 54 , while the earliest deep tectonic burial, high-grade metamorphism, and anatexis in the Greater Himalaya and middle-southern Tethyan Himalaya occurred during the middle Eocene 55,56 .
In summary, we have obtained a high quality and well-dated Early Cretaceous (135.1-124.4 Ma) paleopole (5.9°S, 308.0°E with A 95 = 6.1°) for the Tethyan Himalaya. This new paleomagnetic pole not only satisfies the secular variation model, but fulfills all the seven quality criteria proposed by Van der Voo 46 to appraise the reliability of paleomagnetic data. Based on our new paleomagnetic data and previous high-quality Cretaceous and Paleocene paleomagnetic results from the Tethyan Himalaya and the Lhasa terrane, as well as the APWPs of the India craton, several main conclusions can be summarized as follows: (1) the Gongbuxue area of the Tethyan Himalaya during ~135-124 Ma was located at ~48.5°S; (2) the Tethyan Himalaya belonged to a contiguous Indian continental lithosphere at ~130 Ma and it ever separated from the Indian Craton sometime after ~130 Ma; (3) in the present-day Indian craton coordinates the Neotethyan Ocean between the Indian craton and the Tethyan Himalaya was up to a latitudinal width of ~6900 km (~62.0°) during ~135-124 Ma; (4) the India-Asia collision should be a dual-collision process including the first Tethyan Himalaya-Lhasa terrane collision occurred at ~54. 9 Ma and the final India-Tethyan Himalaya collision occurred at ~36.7 Ma.

Methods
Six paleomagnetic sites were collected with standard paleomagnetic cores from section A using a portable gasoline-powered drill and were oriented by both a magnetic compass and a Sun compass. A declination difference is less than 2° by comparing these two oriented results, indicating that the local magnetic disturbance can be neglected. Twenty-six paleomagnetic sites were collected with oriented blocks because the portable drills were broken during the field sampling. Cores (25 mm diameter) were further drilled from the block samples in the laboratory.
Standard 2.5-cm-diameter cores were cut into 2.2-cm-long specimens in the laboratory. Most specimens were subjected to either stepwise thermal demagnetization from 100 °C to ~580 °C in an ASC-TD 48 furnace with an internal residual field less than 10 nT or alternating field (AF) demagnetization from 5 mT to 110 mT using a D-2000 alternating field demagnetizer. Remanent magnetization measurements were carried out with 2G-755-4K cryogenic magnetometers. All stepwise demagnetization and remanent magnetization measurements were performed in a shielded room with residual fields less than 300 nT at the Paleomagnetic and Environmental Magnetism Laboratory (PEML) at the China University of Geosciences, Beijing (CUGB). ChRM directions of all the specimens were determined using principal component analysis 36 , and site-mean directions were calculated using Fisherian statistics 37 . Paleomagnetic data were analyzed using Enkin's 57 and Cogné's 58 computer program packages.
Acquisition of the IRM, backfield demagnetization of SIRM, thermal demagnetization of the three-axis IRM, hysteresis loops, and FORCs were performed on some representative standard specimens and corresponding powder specimens. The IRM, backfield demagnetization of SIRM and three-axis IRM were acquired using an IM10-30, and then were measured using a JR-6A spinner magnetometer at PEML of the CUGB. Hysteresis loops, FORCs and magnetic hysteresis parameters including saturation magnetization (Ms), saturation remanent magnetization (Mrs), coercivity (Bc), remanence coercivity (Bcr) were measured on representative powder samples at room temperature using a MicroMag Model 3900 Vibrating Sample Magnetometer at the Institute of Geophysics, China Earthquake Administration. FORCs' data were processed using the FORCinel software with a smoothing factor (SF) of 6 35 .
Zircons and images preparation are following the methods described in Ma et al. 17 . U-Pb zircon geochronology was performed using a laser-ablation-multicollector inductively coupled-plasma-mass (LA-ICP-MS) at the Institute of Tibetan Plateau Research, Chinese Academy of Sciences. More detailed analytical procedures and configuration of the LA-ICP-MS have been described in Wu et al. 59 .