Open access peer-reviewed chapter - ONLINE FIRST

Unveiling the Evolution Journey from Pangea to Present Himalayan Orogeny with Relation to Seismic Hazard Assessment

Written By

Tandrila Sarkar, Abhishek Kumar Yadav, Suresh Kannaujiya, Paresh N.S. Roy and Charan Chaganti

Submitted: December 31st, 2021 Reviewed: January 14th, 2022 Published: March 29th, 2022

DOI: 10.5772/intechopen.102683

Earth's Crust and its Evolution - From Pangea to the Present Continents Edited by Mualla Cengiz

From the Edited Volume

Earth's Crust and its Evolution - From Pangea to the Present Continents [Working Title]

Dr. Mualla Cengiz and Dr. Savas Karabulut

Chapter metrics overview

17 Chapter Downloads

View Full Metrics


The objective is to understand incessant seismic activities in Northwest and Central Himalayan regions. GPS data acquired (2017–2020, Nepal; 2015–2019, Uttarakhand) from 65 GNSS stations are used to generate velocity solutions with respect to International Terrestrial Reference Frame 2014 & Indian fixed reference frame to determine the site’s precise position. These velocities are further used to calculate the strain rate and prevailing convergence rate by the respective Triangulation method and Okada’s formulation. The estimated mean maximum and minimum principal strain rate are 12.19 nano strain/yr. and − 102.94 nano strain/yr. respectively. And the respective mean shear strain and dilatation are 115.13 nano strain/yr. −90.75 nano strain, which implies that Higher Himalaya observes high compression rate compared to Outer and Lesser Himalayan region. Estimations have also elucidated presence of extensional deformation in the Northwestern part of the Himalayan arc. Accordingly, in Central Himalaya, paleoliquefaction investigations have deciphered turbidites, confirming that the seismic ruptures did not reach the surface during the 2015 Gorkha earthquake. The best-fit locking depth of 14 km and convergence rate of 21 mm/yr. (Nepal) & 18 mm/yr. (Uttarakhand) are obtained. The strain budget analysis indicates that Northwest and Central Himalaya can beckon a megathrust earthquake in the future.


  • triangulation method
  • Okada formulation
  • palaeoliquefaction
  • best-fit estimation
  • Himalaya

1. Introduction

The collision between the Indian and Eurasian Plate gave birth to the world’s youngest orogeny or fold mountain known as Himalaya. The Himalayan orogeny is very young compared to Aravallis and the Appalachian Mountains. During the mid-Permian, a shallow sea known as Tethys was present in the latitudinal area now covered by the Himalaya. In this era, the Pangea had also started splitting up into small landmasses. In this course, the Northern Eurasian landmass or Angara and the Southern Indian landmass or Gondwana deposited a reasonable amount of sediments into the Tethys, reducing the gap between the Eurasian and Indian landmass. The mountain building process took ∼55 Ma, which started in the Upper Cretaceous period and continued until the Middle Miocene. In the last stage, the formation of Siwaliks also took place that continued from Early Miocene to Neogene age [1, 2].

The convergence activity in the Himalayan arc has made it very vulnerable or important to high seismic activity [3], leading to megathrust earthquakes in the brittle part of the crust that has accumulated elastic strain for a prolonged period. Also, the complex rheology of Himalayan wedges and crustal structure is another reason for the seismic hazard potential. Past studies have shown that Uttarakhand and Western Nepal had been the most potential zone for megathrust earthquakes for 200–500 years [4, 10, 13, 38]. Therefore, in the past half-century, the Northwest Himalaya was ruptured by various small and large earthquakes. Similarly, in the Central Himalaya, the recent Gorkha earthquake of ∼Mw >7 that took place in Kathmandu, Nepal on 25 April 2015 resonated with the 1833 megathrust earthquake, thereby indicating an ongoing convergence rate of ∼14–20 mm/yr. (also consistent with the estimated value) as shown in Section no. 6 & 7 in the Himalayan region [5].

Figure 1.

A schematic diagram illustrating the study area along the northwest-central Himalayan stretch and the locations of the 65 GNSS stations. The square on the left figure is enlarged on the right highlighting the northwest Himalaya or Uttarakhand region along with the GNSS stations that are located in this realm. Note the white vectors indicate the velocity with respective to ITRF2014 and the white circles on it defines the error ellipse. While the red vectors represent the velocities derived with respect to India fixed reference frame (Jade et al. [16]) and the red circles on it define the seismicity of varying magnitude. The yellow lines feature the three Main thrusts of Himalaya; Himalayan frontal thrust (HFT), Main boundary thrust (MBT) and Main central thrust (MCT).

Figure 2.

The strain rate analysis using the triangulation approach. The red and blue vectors indicate maximum and minimum principal strain rate respectively. Note that the blue vectors are highly dominant towards the south eastern part of Nepal region indicating a high rate of compression. While in the northwest part of Uttarakhand region the red vectors are much dominant thus indicating an extension. The maximum and minimum values estimated are tabulated inTable A1.

A dense network of GNSS stations with ∼ > 5 years long period also accomplices in crustal deformation study and help in detecting the slip rate asperity on decollement to unravel the zones/areas with a high risk of seismicity. A detailed study is needed to understand the factors and precursors that trigger megathrust earthquakes. This research work estimated the surface deformation to assess the future seismic hazard potential in Northwest-Central Himalayan region using a close-knit network of 65 GNSS stations (both permanent and campaign).


2. Tectonic setting

Himalaya is located along the southern fringe of Tibetan plateau with a stretch of ∼2500 km length and 250–300 km width, is bounded by Nanga Parbhat (Indus gorge) in northwest and Namcha Barwa (Tsangpo gorge) in the Northeast. The morphologic and structural framework of Himalaya is classified as Kashmir-Punjab of length ∼ 550 km, Kumaun-Garhwal of length ∼ 320 km, Assam segments of length ∼ 400 km and Nepal of length ∼ 800 km [6]. At the east and west syntaxes, the Tsangpo River separates Himalaya from the Indo-Burmese range in the east, and the Indus River separates Hindukush in the west, respectively. The Himalaya with peaks high resides on compressed, thrusted, folded and 70 km thick (twice the thickness from the normal continental crust) continental crust, witnessed through geophysical surveys (seismic reflection and gravity) [41, 42, 43]. MFT lies adjacent to Indo-Gangetic Plain a foreland formed due to convergence between the Indian and Eurasian Plate. The Main Himalayan Thrust acts as a décollement or detachment that dips at a shallow depth. The South Tibetan Detachment, the Main Central Thrust, the Main Boundary Thrust, the Himalayan Frontal Thrust coalesce and meet [7, 8]. The Formation of the Himalaya took place in five distinct phases, which are described in Table 1 relative to their geology and age modified after [9, 10, 11].

Evolution PhaseAge of formationTectonic/geological processGeology (metamorphism grade, Fossils)
Trans-Himalayan Uplift55–35 MaUpliftment due to thrusting on South-directed thrust such as Gangdese Fault, Transportation of sediments by Indus and Tsangpo river, rapid cooled and eroded igneous rocks, formation of major basins such as Kargil, Kailash and LhasaGranitic and volcanic rock intruded by metamorphic and sedimentary rock of the South Tibetan block
Tethyan Himalaya Uplift45–35 MaExtensive folding, uplift of granite gneiss domes, formation of North Himalayan thrust, ‘antecedent’ pattern of Himalayan rivers such as Indus, Satluj, Kali, Ganga, Gandaki, and BrahmaputraCambrian-Eocene sandstone, shale and limestone (contains fossils which includes cretaceous age ammonites), Miocene Leucogranite, Post early middle Eocene-Biotite grade and lower; Early Miocene and younger- Kyanite-sillimanite grade
Higher Himalaya Uplift24–17 MaUplift of the metamorphic rocks buried at 20–25 km depth along the MCT, separation of Tethyan Himalaya from Higher Himalaya along STD that is characterized by tectonic extension, gravitational gliding and back-folding of the sediments of Tethyan Himalaya and upliftment and exhumation of the Higher Himalaya, Klippe formation, deposition of sediments in the foothills of Himalaya10–20 km thick schist and Gneiss, Miocene leucogranite, Central Crystalline Zone, Post-Early middle Eocene-garnet grade
Lesser Himalaya uplift11–7 MaThrusting across MBT, began of Monsoon, rapid erosion of Himalaya, excessive rate of sedimentation in the Siwalik basin,Quartzite, marble, slate, schist and gneiss, Eocene-Chlorite grade and lower
Neotectonics2.6–0 MaUplift of Siwalik range along HFT, rapid erosion, deep and narrow gorges, abundant of coarse-grained fluvial sedimentsSandstone, mudstone (abundant of mammal fossils), Unmetamorphosed

Table 1.

Evolution phase of the Himalaya with associated tectonic processes and geology.


3. GPS velocity measurements within close-knitted GNSS stations

Himalaya has a history of seismicity or earthquakes of very high magnitude that lead to devastation and casualties. To assess the seismic hazard potential, the interseismic strain that is getting accumulated for quite a long period is evaluated with the aid of a close-knitted network of continuous and campaign GNSS stations [10, 12, 13, 46]. These stations are located exquisitely along the Northwest and Central Himalayan stretch spanning the three major thrusts (HFT, MBT and MCT, Figure 1). The Continuous GNSS stations or CORS (Continuously Operating Reference Stations) are well equipped with robust & firmly constructed receivers from Trimble and Leica, Choke-ring geodetic antenna, <10° elevation angle and the time taken for data recording is at 15 s sampling interval. Consequently, the campaign or episodic GNSS stations are equipped with a Zephyr antenna. The data is recorded for 3–4 days repeatedly at regular intervals so that the hydrological mass distribution and seasonal variation effect in this region diminish accordingly. The data is acquired from 65 stations for 2015–2019 (Uttarakhand region) and 2017–2020 (Nepal region). GPS data acquired is processed in GAMIT/GLOBK (vs. 10.71) software along with nearby IGS (International GNSS Services) stations to calculate the precise position of these 65 stations [14]. It takes two steps to complete the GPS data processing. In the first step, the daily relative position of an individual station/site is estimated with respect to clock errors and orbital parameters of the satellite, along with error predicted due to ocean-tidal effect, ionospheric electron content and atmospheric water vapor. The Finite Element Solution (FES) 2012 and International Earth Rotation system (IERS) 2010 models are employed to reduce the respective ocean tidal and Earth tidal effects. While the Global Pressure (GPT) model and Global Mapping Function (GMF) correct the tropospheric delay, which results due to wet and dry water mass [15]. The final or the second step of this process uses GLOBK modules to analyze the loosely constrained solution for obtaining the site velocity with respect to International Terrestrial Reference Frame 2008 (ITRF08) and ITRF14 Uttarakhand and Nepal region, respectively. The velocities obtained from ITRF08 is converted to ITRF14 using the HTDP software ( Also, the velocities for each station is calculated with respect to fixed Indian reference Frame/Euler Pole as suggested by [16] (Table 2, Figure 1). The velocities estimated with respect to ITRF14 vary from 34.70 ± 1.5 (BADR) to 52.42 ± 0.07 (BRN2), thereby indicating crustal shortening rate in the Higher Himalaya is very high compared to the stations that are located in the Outer Himalayan region. Accordingly, the velocities obtained from India Plate fixed reference frame, vary from 23.3 (CHLM) to 0.1 (DEHR). Such variation signifies that the surface deformation in the Himalayan stretch occurs due to the thrusting of the Indian Plate below the Eurasian Plate.

Figure 3.

The colored patches show the dilatation at the Centre of the respective triangles as described inTable A1. The prominent red patch in northwest of the Uttarakhand region defines extension (with reference toFigures 2and4, here the rate of extension is high). And the prominent blue patch in the Southeastern part of Nepal region indicates high compression rate.

SITELONGLATEast_velError(east)North_velError(north)Fixed Indian (North)Fixed Indian (East)

Table 2.

Estimated velocities from both continuous and campaign stations with respect to ITRF 2014 and Indian fixed reference frame [16].

The starred (×) stations indicate the permanent ones.


4. Strain rate estimation on a localized scale using the triangulation method

The estimated site velocities depend on the deformation rate and the translation motion [17]. The translation motion does not depend on the position in a localized zone and has a similar direction and magnitude for a particular reference frame. Whereas the deformation rate indicates a change in the displacement (ui) or position of the site secularly, which is explained in 2D by the following equation.


where xjis the change in the site position and Lijis the displacement gradient tensor. The Lijor displacement gradient tensor is decomposed into asymmetric and symmetric parts (Eq. 2), thereby representing the rotation rate and the infinitesimal strain rate.


Within a close clustered network of GNSS stations that are distributed along the Himalayan stretch, the estimation of strain rate (on a local scale) and the rate of rotation can be carried with the aid of the Triangulation method [47]. In this method, the horizontal field velocity (u̇x, u̇y) determines the deformation rate at the centroid of the triangle at three nonlinear stations. Each of the horizontal velocity component is expressed in terms of deformation rate (ε̇xx,ε̇xy,ε̇yx,ε̇yy) and translation motion (tx,ty) with the recognized initial site location (xo,yo). So, the velocity field at three different stations gives six equations (i.e. 2 for each site) alongwith six unknown, and it is expressed by the following equation.


The above equation and its solution for three non-collinear GPS stations can be written in the form of


where u̇is a 1 × 6 column matrix of known instantaneous displacement/velocity vectors at three stations, mis a 1 × 6 column matrix of unknown translation vector and the deformation gradient tensor, and Gis a square matrix of 6 × 6 order, and G−1 is its inverse matrix which includes zeros, ones and the location vector coordinates of the three GNSS stations.

In a triangulation network, the translation gradient, rotational gradient, maximum shear strain rate and dilatation is estimated using the derived model parameters [18].

The translation gradient is calculated as:


The rotation rate is expressed as:


2D Lagrangian strain rate tensor described as;


The magnitude and orientation of maximum (e1) and minimum (e2) principal strain rate tensor. Also, compression and extension in a regime are determined by the respective negative and positive values of principal strain rates. And the maximum infinitesimal shear strain (γmax) is expressed by the following equation.


The first invariant of the 2D strain tensor is the areal strain or dilatation ((ellipse area-circle area)/circle area) and is explained by the following equation.

Area strain rate/Dilatation=e1+e2E10

5. Inferences from localized strain rate estimation

The strain rate and rotation rate are estimated at the centroid of 89 triangular zones of the 65 sites forming an angle of 30° or more (Table A1). The maximum strain rate varies from −47.90 to 195.22 nano strain/yr. with an azimuth variation from 23.21°N to 161.02° N. Similarly, the minimum strain rate varies from −333.606 to −8.9 nano strain/yr. with an azimuth variation from 178.55 °N to 0.07 °N. The value of maximum shear strain lies between 13.47 to 301.51 nano strain/yr., and the angle of rotation range from −9.9 E-0.7 to 9.12E-0.7 (°/yr). Along with this, the change in dilatation or areal strain differs between 136.07 to −381.405 nano strain. From these observations, it can be stated that the stations (mainly station number 71–86) located in the Higher Himalaya (South East) of the Nepal region experiences a high rate of compression (Figure 2, Table A1), clockwise motion with a mean shear value of 115.13 nano strain/yr. (Figure 4) and negative dilatation (highlighted by blue to yellow color, Figure 3). On proceeding towards the Northwest Himalaya, the stations located there observe a low compression rate (especially along the MBT boundary). At the same time, the cluster of stations (mainly station number 36–50) located in the Higher Himalaya of Uttarakhand region (Northwest) reflects positive dilatation (highlighted by dark red color, Figure 3), anticlockwise movement (Figure 4) and positive strain rate values (Figure 2) or extension.

S.NO.Rotation ± uncertainty (degrees/yr)Max horizontal extension (e1H) (nano-strain)Azimuth of S1H (degrees)Min horizontal extension (e2H) (nano-strain)Azimuth of S2H (degrees)Max shear strain (nano-strain)Area strain (nano-strain)Long (°)Lat (°)

Table A1.

The estimated maximum & minimum principal strain rate, maximum shear strain rate, rotation rate, dilatation with their corresponding azimuth and coordinated by Triangulation approach.

Figure 4.

The maximum shear strain rate is estimated by the triangulation approach. Note the green and blue lines represent right lateral (RL) and left lateral (LL) respectively. The respective estimates and direction of the angular rotation is described inTable A1.

Figure 5.

A 2D second invariant map estimated by the triangulation approach showing the spatial distribution of seismicity. Note the hollow red circles are the seismic events ranging in magnitude from 4 to >5.

Figure 6.

The rate of ongoing convergence in the Nepal region is estimated along the arc parallel or fault parallel direction. Note the green arrows are the estimated convergence rate while the red arrows define modeled convergence rate. The hollow blue circles indicate seismic events (2017–2021) of magnitude varying from 4 to >5.

The second invariant of the strain rate field, where the second invariant is defined as;


where ˙εφφ,˙εθθand ˙εφθare the horizontal components of the strain rate tensor. When the second invariant is estimated, the whole tectonic regime experiences moderate to high areal dilatation (color varies from yellow to orange to red, Figure 5), leaving out certain portions of the Eastern or NNE region (blue to yellow color, Figure 5). The 2D second invariant to the strain rate tensor has been estimated independently from the principal components of strain (Figure 5) to establish a relation between the rate of deformation, seismicity and the interseismic strain accumulating for a prolonged period [19, 44]. The map of the second invariant also focuses on the spatial distribution of earthquakes (Figure 5).

Figure 7.

The rate of ongoing convergence in the Uttarakhand region is estimated along the arc parallel or fault parallel direction. Note the green arrows are the estimated convergence rate while the red arrows define modeled convergence rate. The hollow blue circles indicate seismic events (1960–2021) of magnitude varying from 4 to >6.

In an orogen of convergence, observation of such paradox situation becomes quite inquisitive [20]. Past studies have shown that prominent extensional deformation exists in an orogen of convergence like Andean cordillera, Scandinavian Caledonides etc. [21, 22]. In the Himalayan realm, intra-continental tectonic features, prominent folds, and several other deformational geological structures indicate a compressional tectonic area. But in recent studies, it has been observed that normal faults and many other extensional geological structures are also distributed in different Himalayan sectors [23], but these are not consistent with ongoing seismicity and thrusting in the southward direction. Also, the presence of large scale lineament, shifting and tilting of river, subsidence of older rocks, upliftment of river terraces & active thrusts and normal faults overriding the young Holocene sediments demonstrates that extensional and compressional features co-exist with each other in Himalaya [24, 25].


6. Estimation of ongoing convergence rate and its future implications

The rate of convergence along the active plate boundaries determines the seismic potential and the rate of seismicity occurring today. Previous studies have stated that currently, MFT is the most active fault because here, the crustal shortening is taking place due to the plate movement and evidence of GPS vectors [40]. The site velocities (obtained with respect to ITRF) are distributed along the strike direction (fault parallel or strike slip fault) (Figures 6 and 7) and across the strike direction (fault normal, or oblique fault) (Figures 8 and 9), and the structural trend of MFT is taken as 303° and 288° in the Nepal & Uttarakhand region. It is assumed that the fault normal velocity as dip slip and fault parallel velocity as strike slip on the plate interface, ‘i.e.’ MHT (Main Himalayan Thrust). In fault normal velocity, the uncertainty is calculated by applying the following formula

Figure 8.

The rate of ongoing convergence in the Nepal region is estimated across the arc normal or fault normal direction. Note the green arrows are the estimated convergence rate while the red arrows define modeled convergence rate. The hollow blue circles indicate seismic events (2017–2021) of magnitude varying from 4 to >5.

Figure 9.

The rate of ongoing convergence in the Uttarakhand region is estimated across the arc normal or fault normal direction. Note the green arrows are the estimated convergence rate while the red arrows define modeled convergence rate. The hollow blue circles indicate seismic events (1960–2021) of magnitude varying from 4 to >6.


where σEis the standard deviation of the east velocity, σNis the standard deviation of the north velocity, and θis the bearing of the convergence vector, counterclockwise from east. At the site (GNSS stations) location, the subsurface deformation generally depends on several fault parameters: rake, dip, strike, fault length, width, etc. This is due to finite subsurface dislocation, and it can be explained by the back slip approach as expressed by [26] in the following expression


where u is the surface deformation vectors, s is the vector of parameters, and G is the Green’s functions matrix computed using the semi-analytical formulation published by [27].

The Okada dip slip dislocation approach [27] estimates the surface deformation at each GNSS site due to a finite rectangular dip slip dislocation on MHT for an elastic half-space. In this approach, it is presumed that the frontal portion of MHT is completely locked to some extent and aseismic creeping occurs with a uniform slip rate at downdip zone of MHT. After that, reduced chi-square uncertainty is calculated because a grid search analysis was required to estimate the best fit value of slip rate and locking depth. The expression used for calculating the chi-square uncertainty is as follows;


where nis the number of observations, Pis the number of free parameters, riis the residual between observed and calculated velocity, σiis the formal data uncertainty, and fis the data scaling factor. The scaling factors helps to avoid the inclusion of any additional uncertainty and prevent influencing various other types of data. Also, the formal GPS velocity uncertainty estimates are generally underestimated by a factor of 2 to 5 [28].


7. Inferences from ongoing plate convergence and seismicity assesment

The plate convergence rate (fault normal) or the elastic strain accumulated in the brittle frontal part of MHT is 21 mm/yr. and 18 mm/yr. for Nepal and the Uttarakhand region, respectively, with 14 km locked depth and dip angle 10° (Nepal) and 7° (Uttarakhand) (Figures 10 and 11). But the fault parallel motion is relatively small (3–4 mm/yr) and varies inconsistently (Figures 6 and 7), so the fault normal motion is only considered. The slip deficit in the Uttarakhand region is calculated as 3.9 m, indicating that a significant amount of strain had accumulated since the last magnitude earthquake (1803 earthquake) [45]. And the estimated slip deficit corresponds to seismic moment ∼4.9 × 1021 Nm and Mw ∼ 8.3 by taking fault length 412 km and assuming its width as 111 Km. Therefore, the strain budget analysis in the Uttarakhand region signifies that this area bears the potential to produce a megathrust earthquake, and the results are consistent with [11, 13, 29, 30]. While the Nepal region was struck by a recent high magnitude (Mw >7) Gorkha earthquake that immensely disturbed the tectonic activity of this region. It has been observed that the 2015 seismicity took place along the planar MHT and also gives good information about ground motion frequencies, the direction of rupture etc. It is also presumed that earthquakes with Mw > 8 can rupture the surface devastatingly. But the earthquakes with Mw between 7 and 8 rupture mainly within depth and rarely reach the surface. Also, Paleoliquefaction investigations have deciphered turbidites that formed due to 2015 seismic activity from the Rara Lake, which presents alternative evidence for ruptures that did not reach the surface [32, 34, 35, 36, 37]. Therefore, it can be concluded that interseismic strain is building up within the Nepal Himalayan region to beckon the next calamitous earthquake in the near future.

3D reference model (for both Nepal and the Uttarakhand region) illustrating the predominant occurrence of earthquakes ranging from low to high magnitude is developed using the parameters that have been applied in Okada’s formulation. This model discusses the evolution of slip rate spatiotemporally and the related seismic events that happened in and around the Main Central Thrust fault, thereby indicating that the locations of hypocentre fall in the Higher Himalayan region (Figures 10 and 11) [31, 33, 34]. But previous studies have shown that many historic megathrust earthquakes have ruptured partially and did not reach the surface [38, 39]. Therefore, through magnitude moment calculations and 3D reference model, it can be presumed that a megathrust event similar to the Gorkha earthquake or even more catastrophic than that will occur instantaneously.

Figure 10.

A 3D reference model showing the cross section of the Himalayan geometry in Nepal region. Note the locked width is 96 km, locked depth is 14 km, dip angle 10° and slip rate estimated is 21 mm/yr. the hollow blue circles represent seismic events from 2017 to 2020 ranging between 3 to >5. Also, the concentration of seismicity is high near the Main central thrust, thereby indicating hypocentre of these earthquakes falls in and around the higher Himalaya.

Figure 11.

A 3D reference model showing the cross section of the Himalayan geometry in Uttarakhand region. Note the locked width is 111 km, locked depth is 14 km, dip angle 7° and slip rate estimated is 18 mm/yr. the hollow blue circles represent seismic events from 2015 to 2021 ranging between 3 to >5. Also, the concentration of seismicity is high near the Main central thrust, thereby indicating hypocentre of these earthquakes falls in and around the higher Himalaya.


8. Conclusion

The surface deformation rate is measured along the Northwest -Central Himalayan belt using 65 GNSS stations located close to the wide spatial distribution. GPS data acquired (2015–2019, Uttarakhand and 2017–2020, Nepal) from these GNSS stations are used to estimate the velocity solution and determine the precise position of each site with respect to ITRF2014 and Indian fixed reference frame as suggested by Jade et al. [16]. These estimated velocities are then used to calculate the strain rate and prevailing convergence rate of the Himalayan stretch using the Triangulation method and Okada formulation. From the triangulation approach, the mean observed values of maximum principal strain are 12.19 nano strain/yr., minimum principal strain is −102.94 nano strain/yr., maximum shear strain is 115.13 nano strain/yr. and dilatation of −90.75 nano strain. Although the velocity solution from ITRF2014 implies an increase in convergence rate towards the Higher Himalaya, but results from triangulation have also featured comparatively high extensional deformation in the Northwest (Uttarakhand region) than in the Nepal region. Accordingly, using the grid search analysis, the best fit convergence calculated is 21 mm/yr. and 18 mm/yr. in Nepal and Uttarakhand regions. The locking depth is 14 km in both the regions with dip angle 10° (Nepal) and 7° (Uttarakhand). Considering the last megathrust earthquake in the Uttarakhand region (1803 earthquake), the slip deficit estimated is ∼3.9 m corresponding to ∼4.9 × 1021 Nm seismic moment and Mw > 8. Whereas in Nepal, after the 2015 Gorkha earthquake (Mw 7.3), the tectonic regime is disturbed greatly and has significantly impacted the infrastructure and the GNSS stations present here. Secondary evidences like turbidites from palaeoseismic investigation show that 2015 seismic activity has not ruptured the surface to a megathrust front. This implies building up interseismic stress waiting to release into a calamitous earthquake in the near future.



The authors are grateful to the Director of the Indian Institute of Remote Sensing for his consistent encouragement in carrying out this study. We acknowledge those scholars who have assisted us in the field surveys and those who have helped directly and indirectly in this work. We are thankful to Dr. Robert W. King (MIT) for providing GAMIT/GLOBK 10.70 software for processing GPS data. The authors heartily acknowledges the Editor-in-chief and Reviewer for helping in exemplifying their work exquisitely.


  1. 1. Jacquelyne KW, Tilling RI.“The Story of Plate Tectonics.” this Dynamic Earth. United States Geological Survey; 2016. DOI: 10.3133/7000097
  2. 2. Lovett RA.Texas and Antarctica Were Attached, Rocks Hint. National Geographic News, Magazine; 2011
  3. 3. Seeber L, Armbruster J. Great detachment earthquakes along the Himalayan arc and long-term forecasts. In: Simpson DW, Richards PG, editors. Earthquake Prediction: An International Review. Vol. 4. Washington DC: Am. Geophys. Union) Marice Ewing Series; 1981. pp. 259-277
  4. 4. Bilham R, England P. Plateau ‘pop-up’ in the great 1897 Assam earthquake. Nature. 2001;410(6830):807-809. DOI: 10.1038/35071057
  5. 5. Banerjee P, Bürgmann R, Nagarajan B, Apel E. Intraplate deformation of the Indian subcontinent. Geophysical Research Letters. 2008;35(18). DOI: 10.1029/2008GL035468
  6. 6. Le Fort P. Himalaya: The collided range. Present knowledge of the continental arc. American Journal of Science. 1975;275A:1-44
  7. 7. Lyon-Caen H, Molnar P. Gravity anomalies, flexure of the Indian plate, and the structure, support and evolution of the Himalaya and Ganga Basin. Tectonics. 1985;4(6):513-538
  8. 8. Lang KA, Huntington KW, Burmester R, Housen B. Rapid exhumation of the eastern Himalayan syntax is since the late Miocene. Geological Society of America Bulletin. 2016;128(9–10):1403-1422
  9. 9. Sorkhabi RB, Macfarlane A. Himalaya and Tibet: Mountain roots to mountain tops. Geological Society of America Special Paper. 1999;328:1-7
  10. 10. Kannaujiya S, Yadav RK, Champati ray PK, Sarkar T, Sharma G, Chauhan P. Unraveling seismic hazard by estimating prolonged crustal strain buildup in Kumaun-Garhwal, northwest Himalaya using GPS measurements. Journal of Asian Earth Sciences. 2022;223:104993. DOI: 10.1016/j.jseaes.2021.104993
  11. 11. Kannaujiya S, Philip G, Kannaujiya S, Yadav RK, Champati ray PK, Sarkar T, Sharma G, Chauhan P. Integrated geophysical techniques for subsurface imaging of active deformation across the Himalayan frontal thrust in Singhauli, kala Amb, India. Quaternary International. 2021;575–576:72-84. DOI: 10.1016/j.quaint.2020.05.003
  12. 12. Yadav RK, Gahalaut VK, Bansal AK. Tectonic and non-tectonic crustal deformation in Kumaun Garhwal Himalaya. Quaternary International. 2021;585:171-182. DOI: 10.1016/j.quaint.2020.10.011
  13. 13. Yadav RK, Gahalaut VK, Bansal AK, Sati SP, Catherine J, Gautam P, et al. Strong seismic coupling underneath Garhwal–Kumaun region, NW Himalaya, India. Earth and Planetary Science Letters. 2019;506:8-14. DOI: 10.1016/j.epsl.2018.10.023
  14. 14. Herring TA, King RW, McClusky SC. Introduction to GAMIT/GLOBK, Release 10.4. Cambridge, MA: Dept. of Earth Atmos. and Planet. Sci., Mass. Instit. of Technol.; 2010. p. 48
  15. 15. Boehm J, Heinkelmann R, Schuh H. Short note: A global model of pressure and temperature for geodetic applications. Journal of Geodesy. 2007;81:679-683
  16. 16. Jade S, Shrungeshwara TS, Kumar K, Choudhury P, Dumka RK, Bhu H. India plate angular velocity and contemporary deformation rates from continuous GPS measurements from 1996 to 2015. Scientific Reports. 2017;7(1). DOI: 10.1038/s41598-017-11697-w
  17. 17. Frank F. Deduction of earth strains from survey data. Bulletin of the Seismological Society of America. 1966;56(1):35-42
  18. 18. Cronin VS, Resor PG, Hammond WC, Kreemer C, Olds SE, Pratt-Sitaula B, et al. Developing a curricular module for introductory geophysics or structural geology courses to quantify crustal strain using earth scope PBO GPS velocities. In: Abstract ED41B-0681, Fall Meeting. San Francisco: AGU; 2012
  19. 19. Sreejith KM, Sunil PS, Agrawal R, Saji AP, Rajawat AS, Ramesh DS. Audit of stored strain energy and extent of future earthquake rupture in central Himalaya. Scientific Reports. 2018;8(1). DOI: 10.1038/s41598-018-35025-y
  20. 20. Joshi GR, Hayashi D. Development of extensional stresses in the compressional setting of the Himalayan thrust wedge: Inference from numerical modelling. Natural Science. 2010;02(07):667-680. DOI: 10.4236/ns.2010.27083
  21. 21. Norton M. Late Caladonian extension in western Norway: A response to extreme crustal thickening. Tectonics. 1986;5(2):192-204
  22. 22. Corredor F. Eastward extend of the late eocene early oligocene onset of deformation across the northern Andes: Constraints from the northern portion of the eastern cordillera fold belt, Colombia. Journal of South American Earth Sciences. 2006;16(6):445-457
  23. 23. Royden LH, Burchfiel BC. Thin skinned N-S extension within the convergence Himalayan region; gravitational collapse of a Miocene topographic front. In: Croward MP, Dewey JF, Hanback PL, eds. Continental Extensional Tectonics. Vol. 26(5-6) London: Geological Society; 1987. pp. 611-619
  24. 24. Lave J, Avouac JP. Active folding of fluvial terraces across the Siwalik Hills, Himalayas of Central Nepal. Journal of Geophysical Research. 2000;105(B3):5735-5770
  25. 25. Joshi GR, Hayashi D. Finite element modelling of the pull-apart formation: Implication for tectonics of Bengo Co pull-apart basin, Southern Tibet. Natural Science. 2010;2(6):654-666
  26. 26. Savage JC. Strain accumulation in Western United States. Annual Review of Earth and Planetary Sciences. 1983;11:11-43
  27. 27. Okada Y. Surface deformation due to shear and tensile faults in a half-space. Bulletin of the Seismological Society of America. 1985;75:1135-1154
  28. 28. Mao A, Christopher GA, Timothy HD. Noise in GPS coordinate time series. Journal of Geophysical Research. 1999;104:2792-2816
  29. 29. Yadav A, Kannaujiya S, Champati Ray PK, Yadav RK, Gautam PK. Estimation of crustal deformation parameters and strain build-up in northwest Himalaya using GNSS data measurements. Contributions to Geophysics and Geodesy. 2021;51(3):225-243. DOI: 10.31577/congeo.2021.51.3.2
  30. 30. Lindsey EO, Almeida R, Mallick R, Hubbard J, Bradley K, Tsang LLH, et al. Structural control on Downdip locking extent of the Himalayan megathrust. Journal of Geophysical Research - Solid Earth. 2018;123(6):5265-5278. DOI: 10.1029/2018JB015868
  31. 31. Dal Zilio L, van Dinther Y, Gerya T, Avouac JP. Bimodal seismicity in the Himalaya controlled by fault friction and geometry. Nature Communications. 2019;10(1). DOI: 10.1038/s41467-018-07874-8
  32. 32. Dal Zilio L, Hetényi G, Hubbard J, Bollinger L. Building the Himalaya from tectonic to earthquake scales. Nature Reviews Earth & Environment. 2021;2(4):251-268. DOI: 10.1038/s43017-021-00143-1
  33. 33. Rizza M et al. Post earthquake aggradation processes to hide surface ruptures in thrust systems: The M8.3, 1934, Bihar-Nepal earthquake ruptures at Charnath Khola (eastern Nepal). Journal of Geophysical Research - Solid Earth. 2019;124:9182-9207
  34. 34. Elliott J et al. Himalayan megathrust geometry and relation to topography revealed by the Gorkha earthquake. Nature Geoscience. 2016;9:174
  35. 35. Bilham R. Earthquakes in India and the Himalaya: Tectonics, geodesy and history. Annales de Geophysique. 2004;47:2-3
  36. 36. Stevens V, Avouac JP. Millenary mw > 9.0 earthquakes required by geodetic strain in the Himalaya. Geophysical Research Letters. 2016;43:1118-1123
  37. 37. Armijo R, Tapponnier P, Mercier J, Han TL. Quaternary extension in southern Tibet: Field observations and tectonic implications. Journal of Geophysical Research - Solid Earth. 1986;91:13803-13872
  38. 38. Rajendran CP, John B, Anandasabari K, Sanwal J, Rajendran K, Kumar P, et al. On the paleoseismic evidence of the 1803 earthquake rupture (or lack of it) along the frontal thrust of the Kumaun Himalaya. Tectonophysics. 2018;722:227-234
  39. 39. Rajendran CP, John B, Rajendran K. Medieval pulse of great earthquakes in the central Himalaya: Viewing past activities on the frontal thrust. Journal of Geophysical Research - Solid Earth. 2015;120:1623-1641
  40. 40. Thakur VC. Active tectonics of Himalayan frontal fault system. International Journal of Earth Sciences. 2013;102(7):1791-1810
  41. 41. Caldwell WB, Klemperer SL, Lawrence JF, Rai SS, Ashish. Characterizing the Main Himalayan thrust in the Garhwal Himalaya, India with receiver function CCP stacking. Earth and Planetary Science Letters. 2013;367:15-27
  42. 42. Rawat G, Arora BR, Gupta PK. Electrical resistivity cross-section across the Garhwal Himalaya: Proxy to fluid-seismicity linkage. Tectonophysics. 2014;637:68-79
  43. 43. Gilligan A, Priestley KF, Roecker SW, Levin V, Rai SS. The crustal structure of the Western Himalayas and Tibet. Journal of Geophysical Research - Solid Earth. 2014;120(5):3946-3964
  44. 44. Kreemer C, Holt WE, Haines AJ. An integrated global model of present-day plate motions and plate boundary deformation. Geophysical Journal International. 2003;154(1):8-34. DOI: 10.1046/j.1365-246x.2003.01917.x
  45. 45. Kumahara Y, Jayangondaperumal R. Paleoseismic evidence of a surface rupture along the Northwestern Himalayan frontal thrust (HFT). Geomorphology. 2013;180(181):47-56
  46. 46. Sharma G, Kannaujiya S, Gautam PKR, Taloor AK, Champatiray PK, Mohanty S. Crustal deformation analysis across Garhwal Himalaya: Part of western Himalaya using GPS observations. Quaternary International. 2021;575:153-159. DOI: 10.1046/j.1365-246X.2003.01917.x
  47. 47. Bibby HM. Crustal strain from triangulation in Marlborough. New Zealand. Tectonophysics. 1975;29:529–540

Written By

Tandrila Sarkar, Abhishek Kumar Yadav, Suresh Kannaujiya, Paresh N.S. Roy and Charan Chaganti

Submitted: December 31st, 2021 Reviewed: January 14th, 2022 Published: March 29th, 2022