Three-Dimensional Fault Model and Activity in the Arc-Shaped Tectonic Belt in the Northeastern Margin of the Tibetan Plateau

The arc-shaped tectonic belt, located in the northeastern margin of the Tibetan Plateau, is one of the leading edges of the plateau’s outward growth and uplift expansion, with a large number of active faults and frequent seismic activity. Researchers have carried out numerous studies on active faults in this region, and a wealth of reliable basic data has been accumulated. However, integrating multidisciplinary data to establish a 3D geometrical structure model that is concerned about seismogenic tectonics and can be tested, has become the key to restricting the regional seismic hazard evaluation. Based on a series of published active tectonic research, we analyze in detail the surface and deep coupling relationships of the major active faults in this region and establish three sets of 3D fault structure models, which are built respectively by active fault mapping and dip angles (the V1 model), 7 magnetotelluric profiles and 7 auxiliary profiles (the V2 model), and multi-source data (the V3 model) and continuously close to the real geological facts. From the model perspective, it is suggested that the controversial Haiyuan fault is a crustal-scale left-slip fault and the shape of the Liupanshan fault reflects the absorption of the left-slip component of the Haiyuan fault. Comparing the same fault plane of these three models, we find that the V3 model is more consistent with geological facts, showing that by assisting the multi-source data 3D geological modeling technique we can establish a 3D geological model closest to the real regional structure. Finally, combining the V3 model, the fault segmentation, and the empirical formulas of the moment magnitude-rupture parameters, we segment the faults and calculate the potential moment magnitudes of the major active faults in the study region. The faults with relatively higher seismic hazards are the Liupanshan fault, the southeastern segment of the Xiangshan-Tianjingshan fault, the westernmost segment of the Haiyuan fault, and the West Qinling fault, of which the estimated potential moment magnitudes are generally more than 7.0. Our study provides a referenced 3D geological model for exploring the deep structures of the region, regional geological research, and earthquake disaster prevention.


INTRODUCTION
A large number of studies have shown that faults have non-single planes and complex geometric structures in 3D space, which control the rupture behavior of strong and large earthquakes (e.g., Evans et al., 2020;Ross et al., 2020). How to visually describe the fine 3D geometry of the fault underground has always been an unavoidable problem for structural geologists (e.g., Fischer and Wilkerson, 2000;Plesch et al., 2007;Guan et al., 2010). Since the concept of 3D geological modeling was proposed in the 1990's (Houlding, 1994), with the rapid development of computer technology, the establishment of 3D geological models that can reveal information about rock mass structure relationships and strata spatial distribution in computers has also become more and more sophisticated. 3D geological models have been the product of the mutual intersection of geography, geology, computer graphics, mathematics, computational geometry, environmental science, and other disciplines (Li et al., 2016b).
3D geological modeling methods are classified according to model data sources and can be divided into drilling-based, fieldbased geological data, geological section, and multi-source data modeling (eg, Carena and Suppe, 2002;Riesner et al., 2017;Yue et al., 2005;Wu et al., 2015, Wu et al., 2016. If it is only explained by the geometrically constrained fault structure, the results are generally uncertain and diversified (Evans et al., 2020;Li et al., 2007;Li et al., 2019b;Ramsay and Huber, 1987;Ross et al., 2020;Tavani et al., 2005). At present, the technology of building 3D geological structure models by comprehensively using a series of multi-source data such as shallow geology, artificial earthquake, and drilling data has become increasingly mature, and it has great application value whether it is for a single geological structure or complex regional active faults (e.g., Hubbard et al., 2016;Jia et al., 2010;Li et al., 2010, Li et al., 2013Li et al., 2019a;Li et al., 2021;Plesch et al., 2007;Shaw and Plesch, 2012;Shaw et al., 2015;Wang et al., 2020), which provides a new idea for the study of geometric forms and distribution in the 3D underground space of faults.
The arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau is one of the frontal parts of the entire plateau growing outward and uplifting expansion (e.g., Tapponnier et al., 2001;Zhang et al., 2003), located at the northernmost of the North-South seismic zone and has bred numerous active faults (Deng, 2007). Researchers have carried out a large number of active fault surveys, geophysical explorations, and seismological research in this area, so abundant and reliable original data has been accumulated. But at present, no researcher has used these data to carry out detailed 3D geological modeling work of this region. In addition, due to the boundaries of multiple active blocks in this area, the intracontinental tectonic movement and fault activity is strong, and a series of historical strong earthquakes have occurred, such as the 734 A.D. Tianshui M 7 1/2 earthquake, the 1125 A.D. Lanzhou M 7 earthquake, the 1622 A.D. North Guyuan M 7 earthquake and the 1920 A.D. Haiyuan M 8 1/2 earthquake Mei et al., 2012;Yuan et al., 2013).
At present, some of the structures of the arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau are still somewhat controversial. For example, for the extended depth of the Haiyuan fault, different researchers have used different methods to reveal whether it cuts through the Moho, and the deep structure of the Haiyuan fault needs further research. Zhan et al. (2005) and Jin et al. (2012) revealed that the Haiyuan fault is a deep fault that cuts the Moho by magnetotelluric profiles, the comprehensive geophysical profile shows that the Moho beneath the Haiyuan fault is a complex transition zone (Li et al., 2002), and the offset of the Moho beneath the Haiyuan fault can be directly observed by receiver function profile (Ye et al., 2015). Based on field geological investigation, Burchfiel et al. (1991) conclude that the Haiyuan fault is a fault terminating at the base of the upper crust, and some seismic reflection profiles (Gao et al., 2013;Guo et al., 2016) and magnetotelluric profiles (Min, 2012;Xia, 2019) also indicate that the Haiyuan fault does not cut through the Moho. Through 3D geological modeling, we can effectively clear the structure of the arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau and provide the most credible 3D geological model of this region, providing a basis for subsequent regional research work. Similar to the above problem, although there is a large amount of data accumulated in the region, a paradox of uncoupled data in the modeling process. How to efficiently integrate multiple sources of data to make the model close to reality is also a key challenge in the multi-source data 3D geological modeling technology. Ultimately, our model will serve earthquake prevention and disaster reduction in the study area with frequent earthquakes. Use the 3D geological model to study the seismic activity of major active faults, then focus on the prevention in advance in cities with large earthquake risks to reduce the loss of life and property. Therefore, the arcshaped tectonic belt in the northeastern margin of the Tibetan Plateau has become an ideal area for studying 3D geological modeling and fault activity of active faults.
To establish the 3D fault model of the main active faults in the arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau, we will collect 2D profiles revealing the coupling relationship between the surface and deep in the study area, geological surveys of active faults, and regional Digital Elevation Model (DEM) data. Then we will use the above data to extract the 3D information of faults and detachment layers, assist the SKUA-GOCAD software to generate models, and explore the reasonable method of integrating multi-source data to generate the 3D geological models. Based on the models, we will update the geological understanding of the relevant active faults, analyze the deep segmental properties of the main active faults in the study area, and estimate the seismic moment magnitude of each segment of the fault, to provide a theoretical basis for regional seismic risk assessment. TJSF), the Yantongshan fault (YTSF), the Liupanshan fault (LPSF), and the Xiaoguanshan fault (XGSF) and the Qishan-Mazhao fault (QS-MZF) ( Figure 1). Due to the special geographical location and strong Late Quaternary tectonic deformation in this area, it has always been a hotspot for the study of neotectonics and active tectonics . Researchers have studied the geometric structure, motion properties, activity periods, and historical earthquakes of the above-mentioned faults to a considerable extent.
The WQLF is located between the HYF and the eastern Kunlun fault, which is an ultra-crust fault that cuts through the Moho in the southern part of the study area, with a depth of about 50 km Huang et al., 2017;Ye et al., 2015;Zhan et al., 2014). The fault extends from the vicinity of Baoji in the west of the Weihe Basin to Tongren, Qinghai Province in the west of Linxia, with a total length of about 420 km, a general strike of 280°-310°, and dip angles of 50°-85° (Shao et al., 2011). Since the Quaternary, it has been mainly dominated by left-slip with a relatively small dip-slip component. According to historical records, there have been eight historical strong earthquakes along the WQLF, including five earthquakes of magnitude 6 or higher, which are the 128 A.D. Gangu M 6 1/2 earthquake, the 143 A.D. Gangu M 7 1/4 earthquake, the 734 A.D. Tianshui M 7 1/2 earthquake, the 1765 A.D. Gangu-Wushan M 6 1/2 earthquake, and the 1936 A.D. Kangle M 6 3/4 earthquake Yuan et al., 2007;Zhang, 2020b).
The MXSF zone is located in the southern mountainous area of Lanzhou and consists of four secondary faults, of which the northern margin fault of the Maxianshan is the main fault, with the largest scale and the strongest activity (Yuan et al., 2003;Yuan et al., 2013). The northern margin fault of the Maxianshan has a total length of about 115 km, a general strike of N60°W, a tendency of SW, and a steep dip, which is a lower crustalscale left-slip fault extending to about 40 km with a thrust component (Gao et al., 2013;Guo et al., 2016;Huang, 2017;Xin et al., 2021;Yuan et al., 2003;Yuan et al., 2013;Zhao, 2014). According to the textual research of historical seismic data, many earthquakes of moderate strength or above have occurred in this area, such as the 1125 A.D. Lanzhou M 7 earthquake, the 1590 A.D. Southeast Yongjing M 6 1/2 earthquake, the 1524 A.D. Maxianshan M 5 earthquake, the 1629 A.D. Lanzhou M 5 1/2 earthquake, the 1901 A.D. West Lanzhou M 5 earthquake, and others .
The HYF starts from Xiaokou, Guyuan in the eastern Liupanshan piedmont, passes through Haiyuan, crosses the Yellow River, and extends northwest to the Xijishui leftstepped pull-apart basin in Jingtai, with a total length of about 240 km. The HYF is mainly left-slip motion with thrust component and dips to SSW with a large dip angle, of which the strike gradually changes from NNW to NWW from east to west and generally presents a slight arc protruding NE. The 1920 A.D. Haiyuan M 8 1/2 earthquake occurred on this fault, forming a maximum left-lateral displacement of 10 m and a 230 km long surface rupture zone.
The XS-TJSF is also known as the Zhongwei-Tongxin fault with a total length of about 200 km, whose overall shape is similar to the HYF, presenting an arc slightly protruding NE in the middle. In this study, taking the compression stage area of two secondary faults intersecting at an angle of 17°at Shuangjingzi, Zhongwei as the boundary, the fault is divided into the northwest segment and the southeast segment (Li, 2014). The northwestern segment of the XS-TJSF strikes near E-W, dipping to S with steep dip angles, and is dominated by left-slip motion with either normal or thrust component, which also was the seismogenic fault of the 1709 A.D. Zhongwei M 7 1/2 earthquake. The southeastern segment is a thrust fault dipping to SW with leftslip property, and the general strike is near NW (Min et al., 2001;Li, 2014;Li et al., 2016a).
The YTSF is the third arc-shaped fault protruding northeastward in the northeastern margin of the Tibetan Plateau. The fault is mainly dominated by thrusting motion with a total length of over 150 km, a tendency of NW-NNW, a dip direction of SW, and dip angles of 40°-70° (Li and Shen, 2001;Qin et al., 2017). The linear contour of the fault is not obvious, and there are no active profiles and various geomorphological manifestations in the Holocene, indicating that the activity has been weak since the Late Quaternary (Min, 1998).
The LPSF is NNW-near SN as a whole with a length of about 175 km. The northernmost end of the fault is connected to the leftslip HYF near Xiaokou, Guyuan, and extends southward through Huating, Gansu Province to the vicinity of Longxian, Shaanxi Province, and is connected to the Longxian-Baoji fault zone. The LPSF is the main active fault in the Liupanshan area, dipping to SW and mainly dominated by thrust motion. The measured dip angles in the north are about 50°-55°and those in the south are 40°-60° (Shi, 2011b;Wang, 2018). The Liupanshan area has experienced the 600 A.D. Qinlong M 6 3/4 earthquake, the 1219 A.D. Guyuan M 7 earthquake, the 1306 A.D. Guyuan M 7 earthquake, and the 1622 A.D. North Guyuan M 7 earthquake, which are mainly concentrated in the northern segment of the fault. There are few historical records of strong earthquakes in the middle and south segments, which are unbroken active faults with the potential risk of strong earthquakes (Du et al., 2018;Lei et al., 2007;Mei et al., 2012;Wang et al., 2018).
The XGSF is located on the eastern side of the Liupanshan, nearly parallel to the eastern Liupanshan piedmont fault. It is a W-dipping thrust fault with a length of about 60 km and dip angles of 40°-50°. The Cretaceous Liupanshan group thrusts onto the Tertiary red bed (Zhang, 1990).
Located on the southwestern edge of Ordos, the QS-MZF is the largest, longest, and most active fault in the Longxian-Baoji fault zone (Shi, 2011a;Li, 2017). The fault has a total length of about 180 km, a strike 300°-320°, and dip angles of 70°-85°. According to field surveys, the features of various geomorphological bodies along the QS-MZF trace show that the fault movement mode is mainly left-slip, with a normal fault motion component (Li, 2017).

Data Source and Filtering
The data used in this study can be mainly divided into two parts: 1) surface traces of active faults, and 2) deep structures of active faults. The surface traces are derived from the accurate and reliable geological mapping of active faults (Deng et al., 2002;Zheng et al., 2016, Zheng et al., 2019. The deep structures of active faults come from published research of seismic reflection profiles, geological profiles, magnetotelluric profiles, tomographic imaging profiles, receiver function profiles, accurate location of small earthquakes, and focal mechanism solutions in the arcshaped tectonic belt in the northeastern margin of the Tibetan Plateau. The geological mapping of active faults provides accurate fault surface traces for this study, which is not only the highest priority data for limiting fault surface locations but also the prerequisite for constraining the deep structural profile data of faults. With current technology, it is impossible to directly observe the deep structure of the earth. The indirect observation of the deep geometric structures of the faults through numerous 2D profiles published by researchers is the main way to constrain the deep structures of active faults in this study. Before the establishment of the 3D geological models, this study fully collected a series of data revealing the deep geometric structures of active faults in the study area, mainly including digital image data of the profile locations and the 2D profiles. Digital image data involves multiple disciplines, technical equipment updates, and different processing institutions, making the collected data have the properties of multiple sources, multiple resolutions, and multiple standards. Therefore, through the evaluation and weighting of multisource data, this research takes the data closer to the true deep geometric structure of the fault as the main body, extracts control points at a density of about 2 km horizontally, and extracts control points of low-priority data at a horizontal interval of about 4 km, to obtain the 3D space coordinate sets corresponding to the faults and detachment layers.

Vectorization of Active Fault Surface Traces
According to the most accurate geological mapping of active faults, the surface traces of each active fault are outlined in the geological map. Considering that each active fault has segmental properties, and the current resolution of the 2D profiles cannot reach the accuracy of revealing the segmental structure of deep fault sections, the fault trace of each fault is fitted into a continuous trace distribution in this study ( Figure 1). The geological map is registered through ArcGIS software and corrected to the CGCS2000 GK CM 105E projection coordinate system. Add control points according to the surface trace, add latitude and longitude coordinates for the control points and export the coordinate point set. Combined with DEM data, the elevation information is added to the control points of the corresponding surface trace to obtain complete 3D data.

Vectorization of Digital Images
Firstly, the profile locations are usually indicated in regional geological maps of the published literature, by which we can locate the specific distribution of each profile. Through ArcGIS software, register the corresponding regional geological map, project the longitude and latitude coordinate point set of active fault surface trace control points extracted in step (1), correct the profile location, and extract the profile distribution coordinates. Secondly, import 2D profile images into Adobe Illustrator software, add active fault control points with the corresponding horizontal point density according to the proportional scale in the figures and count the depth information ( Figure 2). Then, the 2D profile images are registered with the profile distribution coordinates, and the latitude and longitude coordinates of the control points are extracted via ArcGIS software. Data processing in ArcGIS software is unified to CGCS2000 GK CM 105E projected coordinate system. Finally, based on the above data, 3D coordinate point set files of the control points of each corresponding fault and detachment layer are sorted out.

Build 3D Models
The extracted 3D data point set was imported into SKUA-GOCAD software ( Figure 3). Then assist the function of the "Structure and Stratigraphy" module to draw strata and fault planes (Li et al., 2016b). Due to the possible sparsity of the modeling data sources, the constructed 3D geological fault plane model is relatively rough. The SKUA-GOCAD software uses a combination of discrete smoothing and Kriging interpolation methods to interpolate the sparse and empty area of the data to smooth the fault plane model (Franke, 1982;Mallet, 1992).

Theory of Earthquake Magnitude Estimation
In the past decades, to analyze the potential seismic hazards, a series of formulas of the moment magnitude-rupture parameters have been proposed by researchers, among which the most widely used is the empirical formulas proposed by Wells and Coppersmith (1994) that the earthquake magnitude corresponds to the surface rupture length (SRL), subsurface rupture length (RLD), rupture area (RA), and average displacement (AD) respectively. Cheng et al. (2020) collected seismic data from 91 earthquakes in mainland China and adjacent areas and used general orthogonal regression to derive the magnitude-rupture scaling relations. Their results were close to the statistical results given by global data without considering the five major earthquakes in eastern China before 1900 (Blaser et al., 2010;Hanks and Bakun, 2002;Wells and Coppersmith, 1994). Therefore, in this study, Eqs 1-5, which are more representative of the empirical seismic magnitude-rupture scale relationship in China, were used to estimate the seismic moment magnitudes (Cheng et al., 2020). For strike-slip earthquakes Mw 3.77 + 0.98Log(RA)(Mw < 6.7) (1) Mw 1.86 + 1.64Log(RA)(Mw ≥ 6.7) (2) Log(RLD) (−2.45 ± 0.17) + (0.61 ± 0.03)Mw (3) For thrust earthquakes At the early stage of this study, the initial 3D fault model V1 of active faults in the arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau was established (Figure 4, hereinafter referred to as the "V1 model"). The V1 model involves the deep structures of four large-scale active faults, namely, the HYF, the XS-TJSF, the LPSF, and the WQLF, and a set of detachment layers with continuously changing depths is also set up. The surface part of the initial V1 model is controlled by the main 3D control points extracted from the fault surface trace. The deep structure setting of the model is as follows. For the HYF, the WQLF, and the northwestern segment of the XS-TJSF, which are left-slip in nature, the dip angle is uniformly set to 90°. For the LPSF and the southeastern segment of the XS-TJSF with a thrust nature, the dip angle variation along the strike of the faults is determined according to the actual geological surveys. In this model, it is considered that the detachment layer with a depth of about 24 km controls the development of all faults in this region.
In the initial V1 model, firstly, the density of control points extracted from the surface traces of each fault is insufficient, which cannot accurately reflect the surface strike change of each fault, and secondly, the deep structures of each fault are only constrained by the geological dip angles. Therefore, the constraints at all levels of the V1 model are not accurate enough, and only the structure of active faults in the region can be briefly displayed. The model is still far from reaching our modeling goal, but it verifies for us the feasibility of establishing the 3D fault model of major active faults in the arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau.

3D Fault V2 Model
Given the shortcomings of the V1 model mentioned in Section 3.1 above, rational optimization of this model is the main direction of this research. The most effective way to optimize the fault surface trace is to encrypt the control points, to improve the accuracy of the fault surface trace. Therefore, we encrypt and extract the latitude and longitude coordinates of the control points at the turning points and curvature changes of the trace, then project them onto DEM to obtain the complete 3D FIGURE 5 | Distribution of modeling profiles and small earthquakes in the V2 model. Small earthquakes occurred from 1 January 1990 to 30 June 2018, with magnitudes less than 6.7 .
Frontiers in Earth Science | www.frontiersin.org May 2022 | Volume 10 | Article 893558 spatial coordinates with depth information. For the deep structure information of each fault, this study combines published geological, geophysical, and seismic data to interpret the deep structures of main active faults and detachment layers in 7 magnetotelluric profiles and 7 auxiliary profiles ( Figure 5) in the study area. Then extract 3D spatial coordinates of the deep structures, and obtain the accurate location of small earthquakes in the study area simultaneously . Finally, with the aid of the SKUA-GOCAD software modeling function, the 3D fault model V2 of the active faults in the arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau was obtained ( Figure 6, hereinafter referred to as the "V2 model"). The modeling scope of the V2 model has been extended to more active faults in the study area, such as the MXSF, the YTSF, the XGSF, and the QS-MZF, which are relatively small-scale faults. According to the published research, two sets of detachment layers with depths of 12 and 24 km were installed in the V2 model to replace the detachment layer with a depth of 24 km in the V1 model (Du et al., 2018;Tian et al., 2021;Xia, 2019). In addition, the accurate locations of small earthquakes are added to the model to visualize the relationship between the deep geometry of active faults and the distribution of small earthquakes, which provides model support for studying the relationship between the fault structure and earthquakes and evaluating the regional seismic risk.

3D Fault V3 Model
The optimization direction of the third model is mainly about the deep structure of the relevant active faults. To achieve this purpose, we searched a large number of published literature and obtained a series of 2D profiles passing through the arcshaped tectonic belt in the northeastern margin of the Tibetan Plateau. The surface traces, accurate location of small earthquakes, and detachment layers distribution of the third model follow the extracted information of the V2 model. For the deep structures of the faults, this study adopts 7 seismic reflection profiles, 3 geological profiles, 22 magnetotelluric profiles, 9 accurate location profiles of small earthquakes, 4 tomographic imaging profiles, and 1 receiver function profile from many collected profiles (Figure 7, Table 1), aiming to use multidisciplinary data to finely characterize the 3D subsurface structures of active faults. Based on the above data, we built the 3D fault model V3 of the active faults in the arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau (Figure 8, hereinafter referred to as the "V3 model").

Model Optimization and Comparison
Through three versions of the models established successively, it is not difficult to see that our results have been gradually transformed from simple modeling to a 3D fine model under multiple constraints. First, from the view of the whole model, the major active faults in the arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau have been all involved in the modeling process from the V1 model to the V3 model, whether it is a large-scale active block boundary fault across the study area or a small-scale fault that is only under the control of local neotectonic movement. Secondly, starting from the V2 model, referring to the accurate geological mapping of active faults, the surface traces of main active faults in the model are constrained by more 3D control points along the fault strike. Meanwhile, the accurate locations of small earthquakes are considered, and the relative positions are visualized in the model. Finally, the most important deep geometric structure, their constraints are gradually strengthened. The V1 model only depends on the geological profiles, the dip angles of left-slip faults are set to 90°uniformly, and thrust faults refer to the surface dip angles to infer the deep extension of faults. The V2 model mainly refers to a series of magnetotelluric profiles and auxiliary profiles to reveal the deep geometry of active faults. The deep structure of each fault in the V3 model is constrained by multi-source 2D profiles. In general, our current V3 model has meticulously depicted the distribution along the strike and extension along the longitudinal direction of the main active faults in the arcshaped tectonic belt in the northeastern margin of the Tibetan Plateau from surface to subsurface.
During our model optimization process, the interpolated formed fault planes have also changed dramatically with more detailed changes as the constraint data increased. Given the fact that in the V1 model only the dip angles from the geological field survey can hardly outline the deep geometric structures of each fault, the following fault plane comparisons only involve the V2 model and the V3 model.
The WQLF with left-slip nature is a deep boundary fault that divides the Tibetan Plateau block from the West Qinling block in the southern part of the study area. The fault cuts through the Moho and extends much deeper than our modeling scale. The top of the WQLF plane constructed by fault surface traces, geological profiles, and DEM data shows the elevation variation of the fault    Figures 9A,B). By analyzing and comparing the normal distribution of dip angles on the two fault planes, the dip angles are mainly distributed in the highangle area ( Figure 9C). For the dip angles of the WQLF plane in the V2 model, the average value is 88.8°, the standard deviation is 1.15°, and 99.7% of the dip angles are distributed in the range of 85.4°-90.0°. It shows that the distribution of dip angles is concentrated the nearly vertical. For the dip angles of the WQLF plane in the V3 model, the average value is 77.4°, the standard deviation is 6.99°, and 99.7% of the dip angles are distributed in the range of 60.0°-87.3°. Although the dip angles Frontiers in Earth Science | www.frontiersin.org May 2022 | Volume 10 | Article 893558 11 of the WQLF plane in the V2 model have a more concentrated distribution, the dip angles of the WQLF plane in the V3 model constrained by multi-source data are more in line with the geological surveys (Shao et al., 2011).
The LPSF with a thrust nature is the boundary fault, which divides the Tibetan Plateau block from the Ordos stable block in the east of the study area. Similar to the analysis of the WQLF plane above, the LPSF surface traces of the V2 model and V3 model show consistent elevation relief, and the shape of fault planes is significantly different (Figures 10A,B). The normal distribution of dip angles shows that the distribution center of dip angles in the V3 model is about 16°higher than that of the V2 model, but the distribution widths of dip angles are extremely similar ( Figure 10C). For the dip angles of the LPSF plane in the V2 model, the average value is 37.3°, the standard deviation is 12.17°, and 99.7% of the fault dip angles are distributed in the range of 27.9°-70.0°. For the dip angles of the LPSF plane in the V3 model, the average value is 53.6°, the standard deviation is 12.05°, and 99.7% of the fault dip angles are distributed in the range of 19.1°-72.3°.
In the V2 model, the near-surface dip angles of northern LPSF are about 40°, and those of southern LPSF are about 60°. The dip angles of the LPSF gradually increase southward along the strike. Meanwhile, as the LPSF extends deeper, the dip angles decrease, and the LPSF plane gradually merges into the detachment layer with a depth of about 12 km ( Figure 10A). However, there is a completely different situation in the V3 model. The dip angles of LPSF gradually decrease southward along the strike. The nearsurface dip angles of northern LPSF are about 60°and the northern LPSF plane gradually merges into the detachment layer with a depth of about 22 km. The near-surface dip angles of southern LPSF are about 40°, and the southern LPSF plane gradually merges into the detachment layer at a depth of about 12 km. According to geological investigations (Shi, 2011b;Wang, 2018), the shallow dip angles in the V3 model are more in line with geological facts. We suggest that the major active faults in the arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau are largely controlled by two sets of detachment layers at depths of about 12 and 22 km. Except for the WQLF in the scale of super crust and the HYF and the MXSF in the scale of at least lower crust, thrust faults terminate between 12 and 22 km in depth, indicating that they are mainly controlled by the above two sets of detachment layers. And the QS-MZF and the northwestern segment of the XS-TJSF with left-slip nature extend to about 12 and 22 km, respectively. In this paper, our research is limited to the HYF (narrow sense), through which most of the 2D profiles show that it has not yet cut through the Moho. However, some researchers have suggested that the western segment of the HYF (broad sense) may be a deep fault cutting through the Moho Zhan et al., 2008;Zhao et al., 2019, Zhao et al.,2020. The LPSF continues southeastward along the tail of the HYF, of which the northwestern segment has a similar extended depth to the tail of the HYF and surface dip angles of about 60°, while the southeastern segment has an extended depth of only about 12 km and gentler dip angles of about 40°. The change of the deep structure of the LPSF from northwest to southeast probably reflects the fact that the strike-slip component of the near E-W trending HYF is absorbed by the near NW-SE trending LPSF, and the nature changes from predominantly left-slip to predominantly thrust Zheng et al., 2013.

3D Fault Model Uncertainty
To a certain extent, the V3 model has visualized the geometric structure of main active faults in the study area from the overall model. However, there are still a series of uncertainties in the V3 model, which affect the accuracy of the V3 model and determine whether it has practical application value. First of all, the multisource 2D profiles, which are the basis for establishing the V3 model, have different geological interpretations for the same fault, and the resolution of profiles obtained by different methods is different. Through unified evaluation and screening of these data in the early stage, we can slightly reduce the impact of these uncertainties. Secondly, as the ordinate of some seismic profiles is the two-way travel time, we uniformly adopt the empirical seismic wave velocity of 5 km/s to make the time-depth conversion, instead of using the actual seismic wave velocity. This approach will inevitably distort the converted fault plane in the model and cannot make it closest to the real structure of the fault plane. Finally, during the modeling process, we found that the density of control profiles in different parts of the study area is different, such as 2D profiles through the convergence area of the LPSF, HYF, and southeastern segment of the XS-TJSF are extremely dense. In contrast, the control profiles through the northwestern segment of the XS-TJSF and the western HYF are relatively sparse. As a result, the control points in different areas of the V3 model are unevenly distributed, and different parts of the fault plane constructed by interpolation also have different accuracy.

3D Fault V3 Model Application
Through the visual V3 model, we can easily obtain the cutting profiles of the arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau in any direction, which reflects the deep geometry of regional active faults (Figures 11, 12). Therefore, the V3 model can provide pre-study results for subsequent researchers who plan to carry out relevant geological research in this area, such as collecting long-distance 2D profiles in a certain direction, helping them to arrange and adjust subsequent fieldwork.
In addition, the arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau is the front of the northeastward expansion of the plateau with active tectonic movements, developing a series of active faults. The evolutionary history of this region since the late Cenozoic is of great significance for studying the expansion of the Tibetan Plateau. At present, the geological research of this area is only limited to the regional surface geological evolution and the deep structures of several active faults Shi et al., 2013;Wang et al., 2013;Wang et al., 2014;Zheng et al., 2013), and have not combined the deep structures of the overall active fault and surface geology to jointly explain the regional geological evolution history. The 3D fault model may serve as a bridge to combine the two, effectively linking regional surface geology, deep structures, most of the active faults, and so on. However, if we want to study the smallscale geological structures and the finer structure of a single fault in this area, our model is far from reaching the corresponding resolution and accuracy. The actual distribution of active faults in this region is far more complex than the V3 model. To reveal the deep segmental and geometric features of active faults more precisely in the 3D fault model depends on the improvement of the regional 2D profiles coverage density.

Moment Magnitude Estimation of Main Active Faults
In the modeling process of this study, considering that the resolution and density of the modeling data cannot identify the segmental properties of active faults, the discontinuous surface traces of the faults are fitted into continuous traces. However, in the earthquake occurrence, not the entire fault is active at the same time, but only the energy release of the partial fault. Given this, it is necessary to clarify the segmentation of faults before quantifying the potential seismic hazard (Ding, 1993;Ding, 1995). For strike-slip faults, this study segments the fault planes according to the segmentation of fault surface traces (Howarth et al., 2021;Wesnouky, 1988;Wesnouky, 2006). For thrust faults, the dip angle contour of the fault planes was generated based on SKUA-GOCAD software, and the changing interfaces of dip angles are used to segment the fault planes (Qiu et al., 2016) (Figure 13). The accurate locations of small earthquakes indicate that there are still a considerable number of small-earthquake distributions within the region at least 30 km deep (Du et al., 2018;Sheng, 2015;Xin et al., 2010;Zhan et al., 2020;Zhang B., 2020), whereas our model sets a uniform cut-off of the fault plane extension depth to 30 km, for which the statistical fault plane area can be used to quantify the potential seismic moment magnitude. Assuming that the planes and surface traces of the seismogenic faults all rupture during an earthquake, we can calculate the potential seismic moment magnitude of each fault segment according to the empirical formulas proposed by Cheng et al. (2020) (Figure 14). If the major active faults in the arc-shaped tectonic belt in the northeastern margin of the Tibetan Plateau produce earthquakes, the estimated moment magnitudes are generally higher than magnitude 6.0.
(2) The MXSF is unsegmented as a whole. The moment magnitude is inferred to be 6.4 by the RAF, and 7.2 inferred by the RLF. (3) The HYF can be divided into five segments. From SWW to SE, the moment magnitudes of each segment are 7.3, 6.7, 6.4, 6.4, and 6.6 deduced from the RAF, and 7.0, 6.4, 6.5, 6.8, and 6.9 according to the RLF. (4) The LPSF can be divided into three segments. From NNW to SSE, the moment magnitudes of each segment are 6.7, 7.3, and 7.1 deduced from the RAF, and 6.5, 7.3, and 6.9 according to the RLF. (5) The QS-MZF can be divided into three segments. From NW to SE, the moment magnitudes of each segment are 6.1, 5.8, and 6.4 deduced from the RAF, and 6.8, 6.3, and 7.0 according to the RLF. (6) The northwestern segment of the XS-TJSF can be divided into five segments. From W to E, the moment magnitudes of FIGURE 14 | Estimated moment magnitude of major active faults in the northeastern margin of the Tibetan Plateau. Red numbers represent results from the rupture area-moment magnitude empirical formula, and black numbers represent results from the rupture length-moment magnitude empirical formula.
Frontiers in Earth Science | www.frontiersin.org May 2022 | Volume 10 | Article 893558 each segment are 5.9, 6.1, 5.9, 6.1, and 6.1 deduced from the RAF, and 6.4, 6.4, 6.2, 6.6, and 6.7 according to the RLF. (7) The southeastern segment of the XS-TJSF can be divided into three segments. From NW to SSE, the moment magnitudes of each segment are 6.8, 7.2, and 7.0 deduced from the RAF, and 6.7, 7.0, and 7.1 according to the RLF. (8) The YTSF can be divided into three segments. From NW to SSE, the moment magnitudes of each segment are 6.5, 7.0, and 6.7 deduced from the RAF, and 6.3, 6.9, and 6.6 according to the RLF. (9) The XGSF can be divided into two segments. From NW to SSE, the moment magnitudes of each segment are 6.8 and 6.8 deduced from the RAF, and 6.8 and 6.7 according to the RLF.
The moment magnitudes respectively estimated from the rupture area and rupture length of thrust faults are in good agreement, with the difference basically within ±0.3, while for left-slip faults the moment magnitudes calculated from the rupture length are mostly larger than those calculated from the rupture area, with significant deviation. In combination with the deep structures of the main active faults shown in the V3 model, the thrust faults are relatively shallow in extended depth. The potential seismic hazards of thrust faults in this region are greater than those of left-slip faults with high dip angles. Therefore, we pay more attention to the possible seismic hazards caused by thrust faults in this region.
For thrust faults, the estimated moment magnitudes are all less than or equal to 7.3. Except for the southernmost part of the southeastern segment of the XS-TJSF, all moment magnitudes estimated from the RAF are greater than or equal to those estimated from the RLF. These thrust faults, which have experienced many historical strong earthquakes (Figure 1), have the potential for large earthquakes of magnitude 7 or greater. Therefore, the estimations based on the RAF are closer to reality and more reliable. Among thrust faults, the LPSF and the southeastern segment of the XS-TJSF, where the maximum seismic moment magnitude can be respectively up to 7.3 and 7.2. Moreover, the historical records of strong earthquakes in the middle and south segments of the LPSF are few and the passing time is long, so there is a risk of strong earthquakes (Du et al., 2018;Lei et al., 2007;Mei et al., 2012;Wang et al., 2018).
At the same time, the seismic hazards caused by left-slip faults, which also have the potential for a strong earthquake of magnitude 7 or more, should not be ignored. For the WQLF, the QS-MZF, the MXSF, and the northwestern segment of the XS-TJSF, all moment magnitudes estimated from the RLF are greater than those estimated from the RAF. The eastern segment of the HYF is consistent with this phenomenon, but the western segment is opposite. We suggest that this may be because the western segment of the HYF, which is constrained by seismic reflection profile, has smaller dip angles of about 60°compared to other left-slip faults . Therefore, the estimations based on the RLF are more reliable for the left-slip faults with steep planes, and conversely, the estimations based on the RAF are more reliable. Among the left-slip faults, the westernmost segment of the HYF and the WQLF are at relatively higher seismic risks. The potential moment magnitudes of the westernmost segment of the HYF can reach up to 7.3. Considering that the WQLF is beyond the modeling range, it extends continuously to the west of the model, and the surface trace is fitted into a continuous trace. Therefore, the V3 model underestimates the length of the westernmost segment of the WQLF, and the moment magnitude calculated by the RLF is small, which is most likely larger than 7.0. All segments of the WQLF have potential seismic risks of moment magnitude 7.0 or greater, and one segment can reach a maximum of magnitude 7.6.
Compared with the historical earthquake records, we find that the moment magnitudes estimated by the V3 model are small. We believe that this error mainly comes from the following two reasons. Firstly, the fault traces are simplified in the model, and the rupture length and rupture area of the seismogenic fault are underestimated. Secondly, we calculate the moment magnitudes under the assumption that only a single fault segment ruptures, without considering the cascade rupture of active faults, and the seismogenic structure generally involves more than a single fault segment. At present, the accuracy of the V3 model is not yet sufficient to solve the problem. However, our result reveals the relative seismic hazard and the minimum seismogenic potential of each rupture in the region, which is of the reference value.

CONCLUSION
In this paper, three versions of 3D fault models of the arcshaped tectonic belt in the northeastern margin of the Tibetan Plateau have been successfully established separately using geological data, magnetotelluric profiles, and multi-source data. Comparing these three models, it is suggested that the multi-source data 3D geological modeling technology is more advantageous in establishing 3D models close to the real regional structure than the single-source data 3D geological modeling technology. Meanwhile, the structures of the main active faults in the study area are reconnoitered from the model perspective. In particular, the V3 model shows that the Haiyuan fault does not cut through the Moho and is a crustal-scale left-slip fault. The shape of the Liupanshan fault plane, which is steep and deep in the north and opposite in the south, may reflect that the left-slip component of the Haiyuan fault with nearly E-W strike is absorbed by the near NW-SE trending Liupanshan thrust fault. Based on the V3 model, the fault segmentation, and the empirical formulas of moment magnitude-rupture parameters, the modeling faults are segmented and the potential moment magnitudes of each segment are calculated. The thrust faults with relatively higher seismic hazards are the Liupanshan fault and the southeastern segment of the Xiangshan-Tianjingshan fault, of which the maximum moment magnitudes are respectively estimated to be 7.3 and 7.2. Among the left-slip faults in the study region, the most westernmost segment of the Haiyuan fault has a potential moment magnitude of 7.3, and potential moment magnitudes of each segment of the West Qinling fault are at least 7.0, of which one segment can be up to 7.6.