1 INTRODUCTION

The Guangdong-Hong Kong-Macao Greater Bay Area (GHMGBA) is an important economic and cultural center in South China, with a large population and intensive urbanization (Liu et al., 2020). In recent years, with accelerated urbanization and rapid development of infrastructure, greater demands have been placed on understanding the basic geological environment of the GHMGBA (Huang et al., 2019). It is of great importance to investigate the 3-D geological structures of the GHMGBA underground space to scientifically plan the construction of urban agglomeration and prevent property losses caused by large-scale geological disasters in the GHMGBA (Liu et al., 1994; Zhou et al., 2018; Chen et al., 2021).

The GHMGBA is located in the interaction zone between the South China Block (SCB) and the South China Sea (SCS). Influenced by the subduction of the Paleo-Pacific plate and the Cenozoic rift-spreading of the SCS since the Mesozoic (Zhou and Li, 2000; Zhang and Wang, 2007; Xia and Zhao, 2014; Zhang et al., 2018; Guo et al., 2019), the GHMGBA has experienced strong alternating processes of compression and extension (Yang, 2019). Extensive tectonic compression formed prominent fold belts, and subsequent continental expansion resulted in the upwelling of mantle materials and multiple periods of granitic magmatism (Huang et al., 2020). These processes caused the formation of a series of extensional basins in the study area. In addition, large-scale magmatic activities have made the region have a higher ground heat flow background, creating favourable conditions for the storage of geothermal resources (Wang et al., 2018). Three groups of faults with NE-, NEE-, and NW-trending control and influence the development of the basins and hot springs in the study area (Zhang, 1999; Yuan et al., 2013). Under the influence of dextral strike-slip of NE- and NEE-faults (Li et al., 2017), the main basins in the region experienced multiple clockwise rotations during their formation and evolution, resulting in multiple changes in their depositional centres (Hou et al., 1994). The mechanical properties of the Enping-Xinfeng fault and the Southern piedmont fault of Wugui mountain near the hot springs west of the Zhujiang (Pearl) River estuary (ZRE) are transtensional faults since Neotectonic movements (Yan, 1986; Ren et al., 2005). This is critical for the subsequent upwelling of hot springs fluid.

Numerous studies have been conducted in the GHMGBA and surrounding areas. Previous wide-angle seismic explorations (Zhang and Wang, 2007; Xia et al., 2010; Guo et al., 2019) and receiver function results (Huang et al., 2014, 2020) found that crust thickness of the coastal area decreases towards the sea, and there are high Poisson' s ratio structure characteristics in the crust on the east coast of the ZRE. This suggests that there was an intrusion of late Mesozoic mafic magma on the crust. The 2-D structures of the P-wave velocity on the east and west sides of the ZRE (Zhang et al., 2018; Lv et al., 2020) proved the existence of multiple structural discontinuities and low-velocity anomalies in the crust. Various seismological methods in Guangdong Province have revealed the existence of a low-velocity layer in the mid-crust (Zhao et al., 2004; Zhang et al., 2018; Zhou et al., 2020), which was generally believed to play an important role in regulating the crustal extension pattern in South China. Xia et al. (2012, 2014) discovered the existence of local high-velocity anomalies in the lower crust between Hong Kong and Dangan Island using an onshore-offshore seismic experiment that reflects the intensive basaltic underplating in the island area of Hong Kong. Using ambient noise tomography, Lü et al. (2022) studied the spatial correlation between seismic activity and crustal velocity anomalies in the Zhujiang (Pearl) River Delta and suggested that seismic clusters occur mainly along the NE-trending faults. Lu et al. (2017) and Wang et al. (2018) used hydrochemistry and numerical simulation methods to study the water cycle pattern, depth, and recharge of water sources of hot springs in the west of the ZRE. They suggested that hot spring systems in the coastal areas of Guangdong are mostly dominated by the deep-water cycle model, while recharge water of hot springs is mostly meteoric water and the depth of water cycle is more than 4 km. However, these studies are either limited to 2-D profiles or P-wave structures or focus on the entire crust and large-scale areas, which result in a poor understanding of the detailed 3-D S-wave velocity structure, basin characteristics, and genesis of hot springs in GHMGBA.

The method of ambient noise tomography has unique advantages in urban underground detection. It can obtain underground structures in urban areas but does not depend on earthquakes or active seismic sources (Shapiro et al., 2005). Moreover, observation with a dense array can effectively improve ambient noise tomography resolution (Yao et al., 2006; Lin et al., 2007, 2013; Huang et al., 2010; Li et al., 2016; Yang et al., 2020; Zhu et al., 2021; Xia et al., 2022). In this study, we performed ambient noise tomography using 40-day continuous waveform data from 130 seismic stations to obtain 3-D S-wave velocity structures of the upper crust in GHMGBA. The S-wave velocity structures and overall morphological structures of the main basins in GHMGBA are revealed, and the genetic and deep source structures of hot springs are probed. Furthermore, we briefly discussed the relationship between the development of urban agglomerations and the basic geological environment in the GHMGBA from the perspective of crustal structures.

2 GEOLOGICAL SETTING

The depth of the Moho in the GHMGBA is between 26 and 32 km and gradually thins seaward. The crustal Poisson' s ratio is less than 0.26 overall, indicating that the GHMGBA crust is dominated by felsic composition (Huang et al., 2020). The study area underwent Caledonian folding orogeny and was strongly influenced by the Indosinian and Yanshanian movements (Yang, 2019). The current tectonic pattern was basically established during the Yanshanian movement (Qiu, 1992). From the early Cretaceous, the SCB underwent a tectonic transition from the Tethys domain to the Pacific domain, and then its tectonic pattern changed from EW to NE (Shu, 2012). The deformation characteristics of the GHMGBA are mainly manifested as fold deformation and fault block differentiation in the adjacent of fault zone (Qiu, 1992; Yao et al., 2013). Orogenies of large fault blocks have intensively changed the geomorphology on the surface and are accompanied by extensive granitic magmatic activities (Zhou and Li, 2000; Shu et al., 2009; Xia and Zhao, 2014; Deng et al., 2019).

Under the influence of complex geological evolution, three groups of faults with NE-, NEE-, and NW-trending are widely developed in the GHMGBA (Lu et al., 1991; Zhang and Wu, 1994; Sun et al., 2012; Yao et al., 2013; Cao et al., 2014). They are dominated by faults with a NE direction, which were mainly formed in the Mesozoic era and activated during the Yanshanian movement, along which vast amounts of Mesozoic granites were intruded. NEE-trending faults were also formed during the Mesozoic and were reactivated as tensional faults during seafloor spreading in the SCS (Cao et al., 2018). The NW-striking faults are thought to have originated later, mainly in the Late Yanshanian and Himalayan orogeny, with smaller scales and cutting deeper than those of the NE- and NEE-trending faults (Chen et al., 2021). These three groups of faults cut deep into the crust and controlled the development of sedimentary basins by affecting regional subsidence. The stratigraphic facies of the basins are dominated by lacustrine sediments from the Cretaceous to the Early Neogene (Zhu et al., 1989). The topographic uplift caused by the first Himalayan movement led to the absence of the Middle Neogene to Pleistocene strata in the region (Yao et al., 2013). Moreover, volcanic rocks are widespread in the basins. Lithologic analyses show that a set of tholei-calc-alkaline series of rocks with island arc characteristics are extensively developed in SSB (Zhu et al., 1989), which reflects lithosphere extension (Chen and Xie, 2011; Chen, 2018).

3 DATA AND METHOD 3.1 Data process

We installed 110 mobile seismometers in the GHMGBA from November 30, 2019, to January 20, 2020. After data processing and quality control, 92 short-period (down to 5 s) three-component portable seismic instruments and 38 existing fixed stations (down to at least 60 s) with good data quality were selected (average station spacing of about 15 km) (Fig. 1). The stations cover most areas of the GHMGBA and can well obtain the 3-D S-wave velocity structures in GHMGBA.

Forty-day continuous waveform data are deployed to conduct this tomography. We first cut the original data into one-day segments, and then decimated the sampling rate to 10 Hz. Before downsampling, the instrument response, mean, and trend have been removed. Spectrum whitening and temporal normalization were adopted to broaden the frequency band of ambient noise signals and reduce the impact of earthquake signals, instrument irregularity, and non-stationary noise sources on subsequent cross-correlations (Bensen et al., 2007). The final cross-correlation function (CCF) of each station pair was obtained (Fig. 2) by calculating the CCF between the two stations and stacking the CCF of each day. To improve the signal-to-noise ratio of CCFs, we stacked the positive and negative lags of CCFs. After that, Rayleigh wave group velocity dispersion curves were measured by frequency-time analysis (FTAN) from the cross correlograms (Herrmann, 2013): we first converted the original data into the waveform with the frequency band of 0.1–2 Hz, then obtained the group velocity of each period by tracking the maximum amplitude. Dispersion data with a station spacing less than 1.5 times wavelength were eliminated. The average group velocity in each period was calculated on the basis of the extracted dispersion curves, and the dispersion data with a difference of more than ±0.2 km/s from the average group velocity in each period were removed. Finally, 5 100 group velocity dispersion curves were obtained (Fig. 3).

3.2 Ambient noise tomography method

Based on group velocities dispersion data, the Fast Marching Method (FMM) (Herrmann, 2013) was used for ray tracing in order to consider the impact of underground medium heterogeneity on the surface wave propagation path. During the inversion, the parameters of smoothing and damping are considered and tested, which were set as 10 and 5. Then we inverted the group velocity structures of the GHMGBA for periods of 0.6–8 s (Supplementary Figs.S1–S2). The travel-time residuals before and after Rayleigh wave velocity inversion are shown in

Supplementary Fig. S3, and the fittings between the measured and calculated dispersion curves are shown in Supplementary Figs.S4–S7. Based on this, we first inverted the average group velocities to obtain a 1-D initial S-wave velocity model (Fig. 4a). Then using this 1-D S-wave velocity model as the initial model and the Rayleigh wave group velocities map of 0.6–8-s periods as inputs to invert the 3-D S-wave velocity model following the linearized inversion method of Herrmann (2013). In the inversion processes, the model is divided into 5 layers at depth with each layer of 1.5 km (5 layers are set at 1.5, 3, 4.5, 6, 7.5 km), and the velocity ratio of S-wave to P-wave V_p/V_s is kept at 1.75. Density parameters were obtained by V_p conversion using the Nafe-Drake relation (Brocher, 2005), as follow:

ρ=1.6612V_p–0.4721V_p²+0.0671V_p³–0.0043V_p⁴+ 0.000106V_p⁵.

Sensitive kernel curves at depths in different periods show that our data can constrain the velocity structures shallower than 8 km (Fig. 4b).

4 RESULT

The lateral resolution of ambient noise tomography mainly depends on the conditions of ray coverage and crossover (Li et al., 2009). The denser the coverage with the better the crossing, the higher the imaging resolution. Figure 5 shows the ray coverage for periods of 1–8 s. We can conclude that the ray coverage of 3–6-s periods is the best. Checkerboard resolution tests are also performed to check the robustness and resolution of the tomography. Figure 6 shows the results of the checkerboard test for periods of 1–8 s with the grid size of 0.2°×0.2°, from which we can determine that most areas are well resolved except for the margin of the seismic array. Due to the uneven distribution of stations in the east of the ZRE, the NE-striking rays are more numerous than those in other directions, which leads to slight NE-trending smearing in the checkerboard test results. We also conducted comparative tests with different grid sizes (Supplementary Figs.S8–S9). We found that the grid spacing of 0.2°×0.2° has a well recovery in the whole study area. Considering the resolution and our research objects, finally, we chose the grid spacing of 0.2°×0.2°.

Figures 7 and 8 show the inversion results of S-wave velocity structures at five different depths (1.5–7.5 km) and three vertical profiles crossing major basins, respectively. On the plane, the main basins in the GHMGBA show obvious low-velocity anomalies at shallow depths (0–3 km). S-wave velocities of sedimentary formation in basins are about 2.8–3.1 km/s, and faults at the basin boundaries correspond well to the boundaries of the low-velocity anomaly. As the depth increases, the scale of low-velocity anomalies decreases significantly. If we take the contour of 3.1 km/s as the lower boundary of the basins, the depth and morphological characteristics of the basins can be roughly delineated. The range of low-velocity anomalies in SSB is the largest, while the velocity value is the smallest. The low-velocity anomalies in SSB characterized by 3.1 km/s extend to the depth of 4 km, with a funnel shape in the NW direction (Fig. 8a). However, it is irregular in the EW direction (Fig. 8b), as a velocity contour in the east side of SSB is gentler than the west side. The SSB velocity gradient changes abruptly in the Xijiang fault (F1) in the west and the Baini-Shawan fault (F2) in the east. The DGB is divided by the Luogang-Taiping fault (F3), and its velocity structures differ significantly on east and west sides of the fault. On the west side of F3, the DGB has obvious low-velocity anomalies extending over 2 km in depth, while the east side does not show obvious low-velocity characteristics (Fig. 8b). The HZB is located between the Heyuan-Shaowu fault (F8) and the Lianhuashan fault (F9). Its S-wave velocity structures also show clear low-velocity characteristics (Fig. 8b–c). Whereas it should be noted that the HZB is located on the edge of the study area, and our data cover only the southern part of the basin, the description of detailed characteristics of imaging results may be unclear, so a detailed analysis was not conducted. The low-velocity anomalies of the EPB extend for more than 2 km (Fig. 8c), with a general V-shape distribution. The Enping-Xinfeng (F4) fault apparently controls the north-eastern boundary of the EPB. The overall S-wave velocities of the XHB are slightly higher than that of the EPB, and its depth is shallower than that of the EPB (Fig. 8c).

We projected the locations of hot springs on maps showing the S-wave velocity structures at different depths (Fig. 7). The hot springs are mostly located on the edge of shallow low-velocity anomalies in a 1.5-km depth slice (Fig. 7a). Vertical profiles over areas with intensive hot springs (Fig. 9) show that at a depth of about 4.5 km below the hot springs there are pronounced low-velocity anomalies. According to the shape characteristics of low-velocity anomalies, they can be roughly divided into two types: (1) areas of intensive hot springs around the EPB (Fig. 9a–b), where the S-wave velocity of the low-velocity anomalies about 3.1 km/s, extending northward, and with a thickness of almost 2 km. The 4.5-km low-velocity anomalies have a certain connectivity with the shallow low-velocity anomalies of the EPB; and (2) areas of intensive hot springs along the coast of Zhuhai (Fig. 9c), where the S-wave velocity of the low-velocity anomalies at 4.5-km depth is slightly less than 3.1 km/s with a banded distribution. Its thickness is less than the thickness of the low-velocity layer below the EPB.

5 DISCUSSION 5.1 Development characteristics and tectonic genesis of the basins

Our results show that the S-wave velocity structures of the upper crust in the study area have obvious lateral heterogeneity, and the major low-velocity anomalies correspond well to the basin locations outlined by previous studies (Zhang, 1999; Lü et al., 2022). Our model revealed the deep 3-D structure of major basins and their coupling relationships with surrounding faults, which is difficult to achieve in previous studies. The ranges of low-velocity anomalies are in line with the basin locations outlined by previous studies (Zhang, 1999; Lü et al., 2022). As the largest basin in the study area, outfield surveys and borehole data revealed that sediment thickness in the SSB was about 4 km (Tang, 1994; Zhang, 1999; Wei et al., 2018), which agrees well with the range of low-velocity anomalies outlined by a 3.1-km/s velocity contour in our study (Fig. 8a–b). This study, combined with the velocity interval of the Mesozoic basin basements (Özalaybey et al., 1997), selected a 3.1-km/s velocity contour to roughly constrain the deep morphological characteristics of the SSB. The S-wave velocity structures (Fig. 8a) indicate that the SSB is funnel-shaped with almost symmetrical velocity contours in the NW direction. Whereas the velocity structures of the SSB have obvious asymmetric morphological characteristics in the EW direction (Fig. 8b). If the deepest sedimentary point of the SSB is considered as the axis (black thick solid line in Fig. 8b), the low-velocity anomaly outlined by the 3.1-km/s velocity contour is 15 km wider on the east side (L2) than that on the west side (L1). Similarly, based on the 3.0-km/s velocity contour, the east side of the SSB (S2) is 25 km wider than the west side (S1). According to regional sedimentary stratigraphy studies, during the continuous expansion of the SSB from the Late Cretaceous to the Tertiary, magmatic activities in the SSB centre showed a linear north-south distribution, which indicated that the SSB spreading direction was mainly east-west (Zhu et al., 2004). Spatial and temporal characteristics of volcanic rock distribution show that volcanic rock distribution in the SSB is consistent with sedimentary centre migration both in space and time, and volcanic rocks gradually migrate west of the SSB with basin expansion (Dong et al., 2006; Zhou et al., 2009). Combined with the deep morphological characteristics of the SSB revealed by the S-wave velocity structures, we proposed that the SSB has an east-west asymmetric extension pattern. In other words, the expansion on the east side is more intense than on the west side, and the asymmetric pattern becomes more severe with the basin expansion.

By analyzing the correlation between the faults and the velocity structures of the basins, it can be determined that the fault systems in the GHMGBA have obvious control and impact on the basins. Both the EPB and the XHB are enclosed by the NE- and NEE-trending faults. The Enping-Xinfeng fault (F4) separates the EPB from the high-velocity anomalies of Tianlu mountain on its northwest side, where Indosinian and Yanshanian granites (Long, 2004) are widely exposed (Fig. 7). The Luogang-Taiping fault (F3) crosses the DGB and shows significant differences in S-wave velocity structures between the east and west sides of the fault. On the west side, the basin shows obvious low-velocity anomalies, while this is not the case for the east side (Fig. 8b), which reveals that uneven uplift exists on the east and west sides of the DGB under the impact of the Luogang-Taiping fault (F3). The S-wave velocity gradient beneath the SSB shows obvious and rapid changes on the east, west, and south sides of the Baini-Shawan fault (F2), the Xijiang fault (F1), and the Enping-Xinfeng fault (F4) (Fig. 8a–b), which indicate that these faults significantly limited the boundary of the SSB. Similarly, faults around the DGB and the XHB also show strong control on the basin boundaries.

The deep velocity structures of many small and medium-sized basins do not show distinct connectivity but show relatively independent geometric characteristics. Their velocity structures are controlled by the NE- and NEE-trending faults (Fig. 8), which is in line with the NE- and NEE-trending basin and range tectonics of South China (Shu et al., 2004; Wang and Shu, 2012). The results from this study show that the S-wave velocity structures in the GHMGBA represent the obvious characteristics of anomalies interphase distribution of high and low velocities along the NE-SW direction (Fig. 8d–e), and the boundaries of high and low velocity anomalies are limited by the main NE-trending faults. This may reflect the extension and compression stress filed in the NW direction. These results are in line with the deep dynamics of strong crustal extension and thinning in South China caused by the rupture and delamination of the NW-trending Paleo-Pacific subduction plate and the asthenosphere upwelling since the early Cretaceous (Li and Li, 2007). Intensive crustal extension enabled the creation of a series of extensional fault-depression basins, followed by massive intrusion of granitic magma and accumulation of volcanic eruptions, which formed the NE-trending granitic volcanic and intrusive complex zone (Shu et al., 2004; Shu and Wang, 2006; Zhou et al., 2009; Shu, 2012). Our results show the high-velocity anomalies beneath the high topography in the GHMGBA, which is assumed to be a Late Mesozoic solidified magma channel system in the crust. In addition, the large granite intrusion zone in Zhuhai and Macao also shows relatively high velocity characteristics. The 3-D S-wave velocity structures depicted shallow and deep basin-range structures in the GHMGBA, reflecting the reshaping effects of the thermal material upwelling on the crust, and providing a clear perspective for understanding the regional deep basin-range structures.

5.2 Genesis and types of hot springs

The west ZRE is rich in geothermal resources, most of which are medium-low temperature hydro-geothermal systems (Zhang et al., 2002; Yuan, 2013). Our imaging results show prominent low-velocity anomalies at a depth of about 4.5 km below the intensive outcrops of hot springs (Fig. 9). According to the borehole data from the geothermal field, the terrestrial heat flow in the GHMGBA is about 70–80 mW/m2, and the maximum geothermal gradient can reach 42 ℃/km (Wang, 2018). This is significantly higher than the normal geothermal gradient of 10–30 ℃/km. Nevertheless, it is far lower than high-temperature geothermal systems with additional heat sources (e.g., the gradient in the Tengchong geothermal field is 100–120 ℃/km with magma heat sources) (Long, 1988). Thus, low-velocity anomalies beneath hot springs cannot be bodies of molten magma. In addition, previous theoretical calculations and analyses of water flow and temperature in hot springs have indicated that the mechanism of hot springs formation in the region was a type of deep-water circulation (Zhang et al., 2002). Numerical simulation has also suggested that high-temperature fluid reservoirs are located beneath the hot springs (Lu et al., 2017). Moreover, Wang et al. (2018) confirmed that the depth of circulation of hot spring fluid through the pipeline model was 3.9–4.7 km, based on the characteristics of the medium-low temperature hydro-geothermal system in deep fault zones. Combined with previous studies, we have suggested that the low-velocity anomalies at a depth of 4.5 km beneath the hot springs in the study area are the thermal reservoirs of the hot springs. Around the thermal reservoirs, the fluid of hot springs reaches a temperature-pressure equilibrium, and fluids migrate mainly through advection (Lu et al., 2017). When the fluid is disturbed, for example when it encounters large faults, its equilibrium state is destroyed, and the buoyancy generated at high temperatures causes the fluid to migrate upward along the faults, only to be exposed as hot springs on the surface (Qi et al., 2017).

The hot springs in the west of the study area can be roughly divided into two types (shown in Fig. 7c) by comparing the hot springs with the characteristics of the topography and S-wave velocity. The first type of hot springs is mainly located on the periphery of the EPB. The fluid from the hot springs is more likely to upwell along the faults around the basin and then be exposed at the edge of the basin since the heat transfer efficiency around the basin is higher than in the centre (Cheng et al., 2021). Moreover, the hot springs are scattered around the thermal reservoirs at a depth of 4.5 km below the EPB, and the thickness of the thermal reservoirs (shown in Fig. 7c) is slightly wider than 2 km. Moho uplift and crust extension and thinning are revealed by the high gravity anomalies in the EPB (Tian, 2014), which can not only lead to basin subsidence, creating space for overlying deposition, but also creating a large number of listric normal faults on the basin periphery. Thus, the sedimentary stratum of low thermal conductivity created favourable conditions for the formation of thermal reservoirs, while large faults around the basin provided ascending channels for the hot spring fluid (Xi et al., 2018). Therefore, we proposed that the hot springs around the EPB belong to a representative geothermal system of the sedimentation basin conduction type, and we speculated that this type of hot springs can easily form thick thermal reservoirs because of the superior thermal storage conditions.

The second type of hot springs is mainly exposed along the coast of Zhuhai, and generally has a linear distribution parallel to the NE-trending faults (shown in Fig. 7c). The outcropping conditions of hot springs here differ from those on the periphery of the EPB. In this region, Yanshanian granites are widely exposed, and granites near the fault zones have a high degree of fragmentation and extreme hydrothermal alterations (Wang et al., 2018). We have speculated that the thermal reservoirs at a depth of 4.5 km are formed by high-temperature fluid hosting in granites with high porosity. And its overlaid is covered with a stratum of low permeability, which shows a relatively high velocity in the S-wave velocity image results, forming the thermal reservoir cap. In addition, hot springs are mainly exposed in areas with large topographic contrast (Huang et al., 1991), thus providing sufficient head pressure for deep water circulation. This results in a relatively shorter cycle of hot springs fluids, thus showing the characteristics of a thin thermal reservoir. We proposed that the hot springs here belong to the uplift mountain convection type due to the characteristics of regional thermal reservoirs and topography. Thermal storage conditions for this type of hot springs are not superior, so it is more likely that this type of hot springs will form thin thermal reservoirs with the mostly banded distribution.

5.3 Characteristics of the GHMGBA geological environment

Our results show that the distribution of major urban agglomerations in the GHMGBA appears to have a good correlation with geological structures characterized by low velocities in shallow areas (less than 3 km). Moreover, the population density distribution (Fig. 10) revealed that the trend of GHMGBA population development in the past 20 years has mainly concentrated in areas with relatively low velocities in geological structures. This may indicate that there is a certain relationship between the urban agglomerations' development and the basic geological conditions in the GHMGBA.

To better investigate how basic geological conditions affect the development of urban agglomerations, we projected data of first-order rivers in the study area onto the S-wave velocity view map at a depth of 1.5 km (Fig. 11). The results show that the trend of first-order rivers in the GHMGBA is either along the periphery of the basins or directly traverse the basins, which may be related to faults on the periphery of the basins. As a fracture zone, faults are always under the impacts of weathering and erosion, which is why the channels of large rivers easily change along the direction of faults (Yao et al., 2013; Jin et al., 2015). It results in a subtle coupling between rivers and basins to some extent (Chang and Wang, 1992; Chen et al., 2008; Sidorchuk et al., 2009). Based on the coupling relationship, the basins are equipped with a set of conditions necessary for human survival, including abundant water and fertile soil, thus allowing a large number of people to settle there and form urban communities. While the cluster effects reflected on the crustal structures are the major urban agglomerations are mainly located in the structural province of relatively low velocity. It is worth noting that there are numerous factors that affect urban development, such as the political system, human capital, investment innovation, etc. (Lin et al., 2014; Zhou et al., 2015; Liang et al., 2018). The urban agglomerations we are discussing are more likely to refer to traditional urban agglomerations with a certain history of farming. Comparatively, "young" cities such as Shenzhen, Hong Kong, Macao, and Zhuhai benefit more from modern policies and the booming shipping industry. The basic geological conditions determined to some extent the urban planning and the environment of the human settlement. Although the correlations between the basic geological conditions and urban development need to be further explored, the obtained characteristics of the geological environment of major cities in the GHMGBA may provide a new reference for the future planning of large cities.

6 CONCLUSION

In this paper, we obtained the 3-D upper crust S-wave velocity structures of the GHMGBA using the method of ambient noise tomography based on 40-day continuous waveform data from 130 stations. Combined with regional geological and geophysical data, the characteristics of S-wave velocity structures in the main basins and beneath the hot springs are discussed. The main findings are as follows:

(1) The velocity structures and overall outlines of the main basins in the GHMGBA are well illustrated by the imaging results. Significant differences in velocity gradients exist between the east and west sides of the SSB, which suggests that the SSB has an asymmetric EW-trending extension pattern. The NE-, NEE-, and NW-trending faults represent intensive control effects of the basins;

(2) There are obvious low-velocity bodies at a depth of about 4.5 km beneath the hot springs, which are speculated to be thermal reservoirs of the hot springs. The hot springs in the west ZRE can be divided into two types according to their exposure conditions and morphological characteristics of thermal reservoirs. The first is the type of sedimentation basin conduction around the EPB and the second is the uplift mountain convection type along the Zhuhai coast. In addition, we have proposed that the hot springs of the sedimentation basin conduction type tend to form thicker thermal reservoirs, while the hot springs of the uplift mountain convection type tend to form thinner thermal reservoirs;

(3) We found that the distribution of traditional urban agglomerations in the GHMGBA shows a correlation with the geological structures of the low-velocity areas at shallow depths (less than 3 km) from a new perspective of crustal structures. It is speculated that the influence of faults leads to a certain coupling relationship between large rivers and basins, forming a favourable geological environment that is suitable for human settlement.

7 DATA AVAILABILITY STATEMENT

All data generated and/or analyzed during this studay are available on request from the corresponding author.

8 ACKNOWLEDGMENT

We appreciate two anonymous reviewers for their useful and constructive comments, which helped improve this manuscript. We are grateful to the Guangdong Earthquake Agency for providing the waveform data, and to all those who participated in data acquisition.

Electronic supplementary material

Supplementary material (Supplementary Figs. S1–S9) is available in the online version of this article at https://doi.org/10.1007/s00343-022-2087-0.