Method for analyzing target area under earth surface and associated device

By receiving environmental noise data at the surface and performing cross-correlation and tomography inversion, a two-dimensional or three-dimensional model of the target area under the surface is generated, which solves the problems of high cost, time-consuming and insufficient accuracy of intrusive technology, and achieves rapid, low-cost and environmentally friendly land characteristics analysis.

CN120418690APending Publication Date: 2025-08-01FNV IP BV
View PDF 0 Cites 0 Cited by

Patent Information

Application Number
CN202380088155.2
Authority / Receiving Office
CN · China
Patent Type
Applications(China)
Current Assignee / Owner
Priority Date
2022-12-23
Filing Date
2023-12-19
Publication Date
2025-08-01

AI Technical Summary

Technical Problem

Existing intrusive technologies are used to determine the characteristics of subsurface lands that are costly, time-consuming, environmentally unfriendly and inadequately accurate, especially in urban environments.

Method used

By receiving environmental noise data at or near the surface, using cross-correlation and tomography inversion techniques to generate two-dimensional or three-dimensional models of the target area, determining land characteristics such as shear velocity are avoided, drilling and active noise excitation are adopted, and a non-invasive method is adopted.

Benefits of technology

It provides a fast, low-cost, accurate and less environmentally impacted approach, suitable for urban environments, reduces the number of invasive surveys, reduces engineering design uncertainty and material use, and improves operational safety.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN120418690A_ABST
    Figure CN120418690A_ABST
Patent Text Reader

Abstract

A method for determining one or more land characteristics of a target area under the earth surface is disclosed, the method comprising: receiving a data set comprising: a first response signal indicative of ambient noise measured at or near the earth surface by a first receiver arranged at a first location; and a second response signal indicative of ambient noise measured at or near the earth surface by a second receiver arranged at a second location, where the first location and the second location are different; processing the first response signal and the second response signal; cross-correlating the processed first response signal with the processed second response signal; and performing tomography inversion using the cross-correlated response signals to generate a two-dimensional '2D' or three-dimensional '3D' model of the target area from the one or more land characteristics. The invention releases insight from geographic data and also relates to improvements in sustainability and environmental development: we together create a safe and livable world.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] The present disclosure relates to methods and systems for analyzing a target region beneath the Earth's surface. More specifically, the present disclosure relates to a method and system for determining one or more soil properties of a target region beneath the surface based on ambient noise measured at or near the surface. The present invention unlocks insights from geodata and also relates to improvements in sustainability and environmental development: together, we create a safe and livable world. Background Art

[0002] There has long been a general need for systems and methods for determining soil parameters beneath the surface. In particular, there is a need for systems and methods that can be used to model the properties of a target soil mass beneath the surface to provide information useful for infrastructure planning, such as but not limited to foundation design, underground infrastructure, and underground storage facilities. Determining the soil properties beneath the surface during the early planning stages of a construction project reduces uncertainty during the project's site selection, foundation design, and construction phases. This in turn reduces delays, cost overruns, and unnecessary use of material resources (such as concrete) during construction.

[0003] A key parameter for determining the soil properties in a soil mass of interest is the shear modulus and the shear velocity Vs. The shear velocity Vs is the velocity at which a shear wave travels through a material and is controlled by the shear modulus of the material. The relationship between the shear velocity and the shear modulus G is defined by Vs = √G / ρ, where ρ is the density of the material. Thus, measurement of Vs provides valuable insights into the soil properties of a subsurface soil region. The small-strain shear modulus (Gmax) is also important in foundation design, where Gmax = ρ.Vs 2 。

[0004] Spectral analysis of surface waves (SASW) and multichannel analysis of surface waves (MASW) are both examples of techniques for collecting surface wave information that can be used to determine the soil properties in a subsurface soil mass. For completeness, surface waves are waves that occur at or near the surface of the Earth. In SASW and MASW, surface layer vibrations are measured from a passive source (vibrations of the surface caused by ambient noise sources) or from an active source (such as a heavy object falling), and the resulting dispersion of the surface waves is studied. ReMi (refraction microtremor) is another surface layer technique that uses ambient noise and surface waves based on observations of the surface environment noise to infer the soil properties of a subsurface region.

[0005] Downhole and crosshole techniques can also be used to determine the soil properties of a subsurface region. In both methods, a receiver located in a borehole measures waves received from an active source located at another position. In downhole techniques, one of the active source and the receiver is located at a subsurface location within the borehole, while the other of the active source and the receiver is located at the surface. In crosshole techniques, the active source is located in a first borehole and the receiver is located in a second borehole. In both downhole techniques, the propagation and dispersion of the received waves are studied to infer the properties of the material through which the waves from the receiver pass.

[0006] Intrusive techniques for measuring the soil properties of a subsurface region typically pose logistical challenges such as long process durations and thus low cost - efficiency, such as long process setup times, long acquisition times, and / or bulky machines, equipment, and processes. Intrusive techniques are particularly unpopular, especially in urban or inaccessible environments, and are generally very expensive. Intrusive techniques can also be environmentally unfriendly, such as disturbing local wildlife. In contrast, current surface layer techniques may lack the accuracy and reliability of more invasive analysis techniques. Summary of the Invention

[0007] This summary introduces concepts that are more fully described in the detailed description section. This summary should not be used to identify essential features of the claimed subject matter nor to limit the scope of the claimed subject matter.

[0008] According to a first aspect, there is provided a method for determining one or more soil properties of a target region beneath the Earth's surface, the method comprising: receiving a data set comprising: a first response signal indicative of ambient noise measured by a first receiver disposed at a first location at or near the surface; and a second response signal indicative of ambient noise measured by a second receiver disposed at a second location at or near the surface, wherein the first location and the second location are different; processing the first response signal and the second response signal; cross - correlating the processed first response signal with the processed second response signal; and using the cross - correlated response signals to perform tomographic inversion to generate a two - dimensional '2D' or three - dimensional '3D' model of the target region based on the one or more soil properties.

[0009] Advantageously, a non-invasive technique is provided for measuring the soil properties of a target area below the surface (e.g., a possible site for a new underground structure). The ability to generate two-dimensional or three-dimensional models of the soil properties of such a target area in a non-invasive manner and using only ambient noise has many advantages. Compared to existing invasive methods for site analysis, this method is simpler, faster, cheaper, consumes less energy, and has less impact on the local environment, wildlife, and community. This is because less or no drilling is required and less active noise excitation is needed. Additionally, critically, the disclosed method is both accurate and reliable, making it a practical and technically attractive technique. The disclosed method is particularly suitable for use in urban environments because of its low impact and because urban environments provide sufficient ambient noise within a frequency range that is very suitable for the disclosed method.

[0010] The disclosed method is particularly advantageous for providing targeted follow-up surveys. Although some invasive measurements may still be required for calibration and ground-truthing purposes, the method according to the invention allows the detection of soil anomalies and the prediction of target areas where soil anomalies are present, rather than randomly determining soil properties by drilling across the entire site. This significantly reduces the number of invasive surveys required. Thus, a series of benefits according to the invention include providing the desired insight into the soil or soil properties more quickly, and thus: reducing the time required; reducing capital expenditure compared to conventional soil surveys; allowing the use of less heavy machinery and equipment through lighter, more sustainable engineering designs, thereby improving operational safety and minimizing environmental exposure risks. Thus, in summary, compared to conventional soil survey methods, the method according to the invention allows for better risk management and / or risk transfer.

[0011] The target area can be from 0 to 100 meters below the surface, from 0 to 45 meters below the surface, or from 50 to 100 meters below the surface. Deeper depths can also be targeted because there are few technical limitations to the method according to the invention. Advantageously, the illustrated method is thus suitable for analyzing a range of target area depths and profiles. This makes it very versatile.

[0012] The target area may be a potential site for a new construction or a tunnel. Land risk management framework decisions can be made based on two - dimensional or three - dimensional models. In the early stages of a construction project, the method according to the present invention provides very important information and insights into the land characteristics, thus providing the best chance to influence the project outcome with minimal change costs. For example, the plan for construction in or near the target area can be adjusted based on two - dimensional or three - dimensional models. The method can be implemented as part of a feasibility study for construction in or near the target area, such as early site selection. By placing more emphasis on early site selection, the uncertainty of the construction project can be reduced. This means significant savings in materials, time, and cost. For example, reducing uncertainty leads to better risk management, which in turn leads to avoiding excessive engineering design in the project. This helps to reduce the use of materials, such as concrete and steel, which has significant environmental benefits. In addition, early site selection helps to inform further analysis of the target area, enabling any further analysis that may be invasive to focus only on those areas that are allegedly identified as problematic or only on certain areas of particular interest.

[0013] Ambient noise or seismic ambient noise may be generated by one or more sources. These sources may be natural (i.e., naturally occurring vibrations) or anthropogenic (i.e., vibrations generated by human activities). For example, these sources can include ocean noise (such as tidal or wave noise), wind noise, industrial noise, industrial machinery noise, noise from transportation vehicles such as cars or trains, and human noise (such as footsteps). Analyzing one or more land characteristics of the target area based on the ambient noise that has occurred is a practical, convenient, energy - efficient, and low - cost method.

[0014] The receivers used in the present method can be geophones, accelerometers, seismographs, vibration sensors, and / or transducers. The receivers can collect data over a relatively long period of time. For example, the ambient noise can be continuously measured over a five - day period. This longer recording time results in a sufficient recovery of surface - wave information from the ambient noise recorded at or near the surface of the target area. In turn, the (processed) surface - wave information can be used in tomographic inversion, as described in more detail below, to obtain a shear - wave velocity model beneath the surface.

[0015] Each response signal can indicate Rayleigh waves and / or Love waves caused by the ambient noise.

[0016] Each of the first response signal and the second response signal can respectively indicate the vertical component and / or the horizontal component of the ambient noise measured by the first receiver and the second receiver at or near the surface.

[0017] The method may include determining the locations of a first receiver and a second receiver prior to the step of receiving the data set. These locations may be at or near the surface of the ground. These locations may be based on the minimum and / or maximum depth of the target area and / or the desired resolution of waves caused by ambient noise in the target area. By determining the locations of the receivers based on the minimum and / or maximum depth of the target area and / or the desired resolution of the noise waves in the target area, the first and second receivers are optimally placed to measure the noise waves at their most sensitive points. This contributes to the accuracy of the response signal and thus to the resulting two-dimensional or three-dimensional model.

[0018] There may be more than two receivers. The receivers may be arranged in an array (row or grid).

[0019] The method may include selecting a recording frequency for the first receiver and the second receiver prior to the step of receiving the data set. The recording frequency may be selected based on the depth of the target area and / or the expected wavelength of the noise waves in the target area. The recording frequency may be referred to as the target recording frequency and may be a range.

[0020] The one or more soil properties may include one or more elastic properties of the target area, such as the shear velocity Vs. As described above, obtaining the shear velocity of the target area provides valuable insight into the soil properties of the target area, such as the small-strain shear modulus of the target area. This enables engineers to identify weak areas, or lateral geological variations, in the target area below the surface. For the reasons above, it is better to identify these features earlier in the life of the project.

[0021] The processing step may include processing the first response signal and the second response signal to enhance the representation of the ambient noise. In other words, enhancing the broadband characteristics of the ambient noise. For example, this is achieved by removing the instrument response and / or by filtering out large amplitudes. Such large amplitudes may be caused by seismic (undesired) signals. Advantageously, this prevents large amplitude events from overwhelming the ambient noise response of interest.

[0022] Additionally or alternatively, the processing step may include one or more of the following operations: dividing each response signal into segments; tuning each (optionally segmented) response signal to the nearest second; applying low-pass filtering to each response signal; and / or downsampling each response signal. The method may include downsampling each response signal by an integer factor. For example, the method may include downsampling at least two response signals by approximately 10 times (e.g., 10 times). The advantages of segmenting and tuning the data are that this allows the cross-correlation and superposition of response signals measured by different receivers to obtain an estimate of the Green's function between the two receivers. The advantage of applying low-pass filtering is that this reduces the frequency components of each response signal to only those frequencies that can be acquired without aliasing after subsequent cross-correlation steps. Aliasing occurs when the sampling of a signal is insufficient to reconstruct the waveform of a particular frequency. The advantage of downsampling is that it reduces the computational cost of the method as well as the total time required to generate the model.

[0023] The step of cross-correlating the processed response signals may include estimating the Green's function between a first receiver and a second receiver. In other words, the step of cross-correlating at least two processed response signals may include simulating a wave field that would be recorded at one of the first and second receivers if one of the first and second receivers were a virtual source. For completeness, the cross-correlation step may be applied to segmented and / or downsampled and / or tuned response signals. Advantageously, by cross-correlating the response signals, the similarity between the response signals can be analyzed. This is an important step in generating the model.

[0024] The step of performing tomographic inversion may include using only fundamental-mode surface waves in the inversion. Advantageously, using only fundamental-mode surface waves results in a highly accurate model.

[0025] The inversion may be performed in the time domain or the frequency domain. Advantageously, choosing between the time domain and the frequency domain means that the domain that results in a shorter processing time can be selected.

[0026] The step of performing tomographic inversion may include using a combination of tomography and inversion.

[0027] The step of performing tomographic inversion may follow a one-step method or a two-step method.

[0028] The steps for tomographic inversion may include: providing an initial model of one or more subsurface properties of a target area; providing a noise input to the initial model; using the initial model to calculate the travel times of response signals that will be measured at a first receiver and a second receiver; comparing the travel times calculated using the initial model with the cross-correlated first and second response signals; and updating the initial model based on the result of the comparison to generate a two-dimensional or three-dimensional model.

[0029] The step of providing the initial model may include determining a phase velocity dispersion curve between a virtual source and one of the first and second receivers.

[0030] The step of providing the initial model may include determining the (e.g., two-dimensional or three-dimensional) structure of the initial model by computer tomography of two or more phase velocity dispersion curves selected between different virtual sources and different receivers.

[0031] The two-dimensional model may be a two-dimensional representation of the shear velocity of the target area. The three-dimensional model may be a three-dimensional representation of the shear velocity of the target area.

[0032] The method may be a computer-implemented method.

[0033] The method may be used in combination with one or more other techniques for analyzing surface waves, such as in combination with SASW or MASW. Advantageously, this allows for the analysis of the subsurface over a greater depth range.

[0034] According to another aspect, there is provided a computer program product comprising instructions that, when executed by a computer, cause the computer to perform the method of the first aspect.

[0035] According to another aspect, there is provided an apparatus configured to perform the method of the first aspect.

[0036] The apparatus may be a system that includes: one or more processors; and a memory having stored thereon computer-readable instructions configured to cause the one or more processors to perform operations including the steps of the first aspect.

[0037] According to yet another aspect, there is provided a computer-readable medium comprising instructions that, when executed by a computer, cause the computer to perform the method of the first aspect. BRIEF DESCRIPTION OF THE DRAWINGS

[0038] Exemplary embodiments of the present disclosure will now be illustrated by way of example with reference to the accompanying drawings. In the drawings:

[0039] Figure 1 A cross-sectional view of a target area below the ground surface is shown;

[0040] Figure 2 Shows a plurality of geophones arranged in a two-dimensional array at the surface above a target area below the surface;

[0041] Figure 3 Is a flowchart showing an exemplary method for analyzing one or more soil properties of a target area;

[0042] Figure 4 Is a flowchart showing an exemplary method for performing tomographic inversion that can be used in the method of Figure 3 ;

[0043] Figure 5 Is a perspective view of a shear wave velocity model;

[0044] Figure 6 Is an embodiment of a shear wave velocity model obtained from the method described herein; and

[0045] Figure 7 Shows an embodiment of a computing device 700 that can be used to perform the methods of Figure 3 and Figure 4 ;

[0046] Throughout the specification and the drawings, like reference numerals denote like features. DETAILED DESCRIPTION

[0047] A method for analyzing one or more soil properties of a target area below the surface will now be described. The soil property is, for example, shear velocity. Briefly, the method involves receiving a data set indicating ambient noise at the surface of the subsurface area, which is measured by a receiver at or near the surface; performing various operations on the data; and ultimately generating a two-dimensional or three-dimensional model of the target area based on the one or more soil properties. As previously described, the ability to generate a two-dimensional or three-dimensional model of a target area, such as a candidate site for a new structure, based on shear velocity provides valuable insight into the composition of the target area. This, in turn, can be used to influence subsequent site investigation and construction decisions.

[0048] Ambient noise is generated by various sources. These sources are divided into two categories: natural and anthropogenic. Natural sources are causative sources that cause naturally occurring vibrations, such as the ocean or wind. Anthropogenic sources are noise sources resulting from human activities, such as industrial noise, industrial machinery noise, noise from transportation vehicles such as cars or trains, wire noise, and human noise.

[0049] More specifically, the received data set contains response signals. Each response signal is measured by a corresponding receiver at or near the surface of the earth. The response signal indicates the amplitude of ambient noise that is transmitted through the subsurface from various noise sources (which can be anthropogenic and / or natural sources, as explained in more detail below) and measured by the corresponding receiver. In particular, vibrations caused by surface waves are measured. Surface waves are generated by anthropogenic or natural processes occurring at or near the surface of the earth. The elliptical motion and dispersion characteristics of surface waves allow us to obtain information about the land, such as the shear characteristics of the subsurface region.

[0050] The receivers can also be considered as sensors, and they can be geophones, accelerometers, seismographs, vibration sensors, and / or transducers.

[0051] Now, information that helps in understanding the present invention will be provided.

[0052] First, an explanation of the shear modulus and its use in construction and infrastructure projects will be given. The shear modulus is a measure of the elastic shear stiffness of a material, representing the deformation of a solid when it is subjected to a force parallel to one of its surfaces while the opposite surface is subjected to an opposing force. This force and its effect in the soil mass beneath the surface (such as the target area of the method described in conjunction Figure 3 with the illustration) are important parameters for studies before and during the design of construction and infrastructure projects. To determine the shear modulus of the soil mass, the shear velocity Vs is determined. This, in turn, gives an indication of the stiffness of the subsurface material and its ability to support structures that are above and / or extend through the soil mass.

[0053] Next, an explanation of the wave types will be given. In the context of land studies, two types of waves are typically distinguished: P-waves and S-waves. In P-waves, the particles in the soil oscillate along the direction of wave motion, and this oscillation causes compression and restoration of the land as the wave propagates through it. While S-waves are shear waves, in which the particles oscillate along a direction perpendicular to the direction of wave propagation.

[0054] P-waves and S-waves are body waves and propagate through the soil mass in all directions. The interaction of P-waves and S-waves with the earth's surface generates surface waves that propagate along the surface of the earth. Multiple types of surface waves can be distinguished. In the systems and methods described herein, Rayleigh waves are measured and studied because it is convenient to measure the vertical component of surface vibrations. However, it should be understood that in the systems and methods described herein, other surface waves (such as Love waves) can be measured and utilized. The systems and methods described herein can be adapted to measure Rayleigh waves and / or Love waves using single-component receivers or multi-component receivers.

[0055] Since surface waves propagate in two dimensions (at the Earth's surface), they decay more slowly than body waves (which propagate in three dimensions). Surface waves typically occur within a depth range of one wavelength from the Earth's surface, usually propagate more slowly, and have frequencies significantly lower than body waves.

[0056] This lower attenuation, slower travel times, and lower frequencies of surface waves make their study particularly attractive for determining shear velocity. Since surface waves have lower attenuation, the signal strength remains better over longer propagation distances. Thus, the resulting measurements typically have a higher signal quality (signal-to-noise ratio) than body wave studies.

[0057] Some of the subsequent embodiments will be described in the context of geophones and geophone arrays. However, it should be understood that the disclosed systems and methods are applicable to a variety of receiver (i.e., sensor) types, including but not limited to geophones, accelerometers, seismographs, vibration sensors, and / or transducers.

[0058] Now reference will be made to Figure 1 describe an exemplary physical setup for collecting data indicative of ambient noise at the surface of a target area. Figure 1 A cross-sectional view of a subsurface soil mass 100 is shown, which subsurface soil mass is, for example, the target area of the method described in Figure 3 P-waves and S-waves propagate through the soil mass 100 as body waves. The Earth's surface 102 extends above the subsurface soil mass. Surface waves propagate along the Earth's surface 102.

[0059] At point A, a schematic illustration of particle oscillation (due to Rayleigh wave propagation) at the surface above the target subsurface soil mass is shown. As shown, the oscillation of particle P is partially vertical and partially in the direction of propagation. Thus, the resulting particle motion is substantially elliptical.

[0060] At the Earth's surface 102 above the soil mass 100, a plurality of geophones 104a, 104b are arranged. The geophones 104a, 104b located at the Earth's surface 102 may be configured to measure the vertical component of the oscillation schematically shown at point A.

[0061] The geophones 104a, 104b are arranged in a grid array at the surface, which grid extends in two directions. It should be noted that in many cases, the surface above the target area may not be planar. The array of geophones 104a, 104b may thus not be truly "two-dimensional" since each geophone may deviate in the z-direction from its neighbors in the grid. However, this grid arrangement of geophones will be referred to herein as a two-dimensional array.

[0062] As the surface wave propagates across surface 102, the surface wave across surface 102 will cause vertical motion of multiple geophones 104a, 104b.

[0063] Now, how to determine the shear velocity Vs from the observation of surface waves, especially Rayleigh waves, will be generally described. This is achieved by measuring the dispersion behavior of surface waves. Surface waves are dispersive, which means their velocity depends on frequency. Generally, seismic velocity increases with the increase of the depth of the earth. Therefore, conventional surface wave dispersion indicates that the surface wave velocity decreases with the increase of frequency. By studying the behavior of surface waves at the ground surface above the soil mass, the soil properties of the soil mass can be determined.

[0064] There are two ways to measure the velocity of dispersive surface waves, and there are differences in determining the group velocity or the phase velocity.

[0065] The group velocity of a wave is the velocity at which the overall envelope shape of the wave amplitude propagates in space, and this overall envelope shape is called the modulation or envelope of the wave. The group velocity is equivalent to the velocity at which the energy of the wave propagates in the soil mass. The group velocity is measured by determining the wave propagation between a pair of (synthetic) transducers, and the group velocity is a frequency-dependent point property in the soil mass, which depends on the depth. The group velocity is obtained as a measurement of the flight time (i.e., travel time) between a virtual source and a receiver.

[0066] At the same time, the phase velocity is the velocity at which the phase of any one frequency component of the wave propagates. The phase velocity is the velocity at which each frequency component of the wave propagates. Therefore, the phase velocity is expressed as a function of frequency. To measure the phase velocity, at least two measurement nodes are selected to measure the wave propagating through the soil mass to determine the relative flight time between geophones at different frequencies. The result is the phase velocity as a function of frequency, which is the average of the soil mass between the two measurement nodes. The phase velocity is obtained as a point in a two-dimensional phase space (dispersion spectrum), which is obtained by a two-dimensional transformation (such as slant stack, Radon, FK, etc.) of an array of recorded waveforms (time-distance space).

[0067] Please refer again to Figure 1 , each geophone 104a, 104b provides a measurement node for measuring the vertical component of the passing surface wave. The geophones 104a, 104b can be configured to measure vibrations caused by ambient noise. For completeness, as described above, ambient noise is the background wave field generated by natural or man-made noise (rather than impulse points such as explosions or drop hammers used in active methods).

[0068] By cross-correlating the passive noise signals measured at a pair of receivers, the Green's function for that pair of receivers can be obtained, which represents the wave field as if one of the pair of receivers were a virtual (e.g., noise) source and the other of the pair of receivers were the receiver.

[0069] Figure 2 Multiple pairs of virtual source-receivers across the surface above the subsurface region of interest are shown. The ray paths 206 between the source-receiver pairs are indicated. In particular, the ray paths from a single central geophone near the center of the array and each other geophone in the array are shown. The background shading and contour rings represent the travel time field from the central geophone to the other geophones. In fact, there are also corresponding ray paths between each geophone and all other geophones (i.e., between each pair of geophones), which are not shown in Figure 2 for simplicity.

[0070] Each pair of geophones can provide a signal at location A, which will be cross-correlated with the signal at location B to reproduce a pair of virtual source-receivers using the interferometric principle. Specifically, the cross-correlation of the passive noise measured at each pair of geophones on the surface shown in Figure 1 can be used to reproduce the response from the subsurface target soil mass as if it were induced by a pulsed point source, which is equal to the Green's function.

[0071] In other words, the response received by cross-correlating the recordings of two receivers can be interpreted as the response measured at one receiver location as if there were a sound source at the other receiver location. There are various ways to determine the Green's function for a virtual source-receiver pair, which are outlined in the article "Tutorial on Seismic Interferometry: Part 1 - Basic Principles and Applications" in Volume 75, Issue 5 of GEOPHYSICS (September - October 2010; P. 75A195075A209; Wapenaar et al.).

[0072] Figure 3 An exemplary method 300 for determining one or more soil properties of a target region (e.g., soil mass 100) beneath the Earth's surface is shown. Method 300 can be a computer-implemented method. The computer device 700 that can be used to execute this method will be described later with reference to Figure 7 below.

[0073] Method 300 begins at step 310. At step 310, a data set is received. The data set includes a first response signal. The first response signal indicates ambient noise measured at or near the surface 102 by a first receiver disposed at a first location. The data set also includes a second response signal. The second response signal indicates ambient noise measured at or near the surface by a second receiver disposed at a second location. The first location and the second location are different. The first receiver and the second receiver measure vibrations of the surface and output a voltage response.

[0074] The response signals collected in each measurement can indicate the amplitude of ambient noise propagated through the subsurface from various anthropogenic or natural ambient noise sources and measured by the first and second receivers (e.g., vibration sensors / transducers). In particular, vibrations caused by surface waves are measured. As already discussed, surface waves are generated by anthropogenic or natural processes occurring at or near the surface. Their elliptical motion and dispersion characteristics allow us to obtain information about the shear properties in the target area. In this embodiment, the first response signal and the second response signal each separately indicate the vertical component of the ambient noise measured at or near the surface by the first receiver and the second receiver. In other embodiments, the first response signal and the second response signal each separately indicate the horizontal component of the ambient noise measured at or near the surface by the first receiver and the second receiver, and optionally also an additional vertical component.

[0075] The first signal and the second signal can be received directly from the first receiver and the second receiver, or via one or more intermediate devices. For example, a computing device (as described hereinafter with reference to Figure 7 is) can be located at the local location of the receiver for transmitting the signal by any form of wired or wireless communication. Alternatively, the signal can be received at a location remote from the receiver for remotely performing the methods disclosed herein (i.e., remotely), e.g., by Internet communication or by transferring a physical computer-readable medium on which a signal record is stored. In some embodiments, the first signal and the second signal are received in real time. In other embodiments, the first signal and the second signal are received after the recording is completed, or at regular time intervals during the entire recording session.

[0076] Optionally, method 300 may include selecting a recording frequency for the first and second receivers prior to step 310 of receiving the data set. The recording frequency may be selected based on the depth of the target area and / or the expected wavelength of the noise waves in the target area. Additionally or alternatively, the recording frequency for a given receiver may be selected based on its characteristics, such as its optimal recording frequency or frequency range. The recording frequency may be referred to as the target recording frequency and may be a range. Advantageously, this means that the recording frequency can be optimized for the characteristics of the target area and the ambient noise in the target area, which means that the use of overly high recordings can be avoided. This helps to reduce the processing burden associated with an excessive number of data points.

[0077] Optionally, method 300 may include determining locations for the first and second receivers prior to step 310 of receiving the data set. These locations may be locations on or near the surface of the ground. These locations may be based on the minimum and / or maximum depth of the target area and / or the desired resolution of the waves caused by ambient noise in the target area. By determining the locations of the receivers based on the minimum and / or maximum depth of the target area and / or the desired resolution of the noise waves in the target area, the first and second receivers are optimally positioned to measure the noise waves at their most sensitive points. This helps the accuracy of the response signal and thus the resulting two-dimensional or three-dimensional model.

[0078] Next, in step 320, the first response signal and the second response signal are processed. The first response signal and the second response signal are processed to optimize them for subsequent steps of the method. The processing of the first response signal and the second response signal will be described later in this specification.

[0079] Next, in step 330, the first response signal is cross-correlated with the second response signal. The purpose of the cross-correlation is to simulate the following wave field: the wave field that would be recorded at the other of the first and second receivers if one of the first and second receivers were a virtual noise source, thereby enabling easier comparison of the cross-correlated response signals. The simulated wave field is a necessary input for the next step of the method. The cross-correlation step will be described in more detail later in this specification.

[0080] Next, in step 340, tomographic inversion is performed using the cross-correlated response signals. Tomographic inversion is the combination of tomography and inversion. Tomography is the name of a technique that uses penetrating waves to display a representation of the cross-section of an object. In other words, tomography is the name of an imaging technique based on penetrating waves. In this method, tomography and inversion are combined in one or two steps. In this context, inversion is the name of the process of obtaining the cross-correlated response signals and converting them into a predicted two-dimensional or three-dimensional model of the target area based on one or more subsurface properties (and thus, in this embodiment, based on the shear velocity Vs). Similarly, the steps of tomographic inversion will be described in more detail later in this application.

[0081] Finally, in step 350, the model of the shear velocity of the target area is generated. This is the output of the tomographic inversion and is effectively the overall output of the method. It is this model that provides valuable insights into the composition of the target area and thus confirms its suitability as a potential construction site.

[0082] In some embodiments, the model is sent to an output device. For example, the model can be displayed on a monitor or other user interface, or sent via wired or wireless communication to another computing device for further processing or display.

[0083] The steps of method 300 will now be described in more detail.

[0084] In step 320, the first response signal and the second response signal are processed. The overall objective of this processing is to optimize the first response signal and the second response signal for the subsequent steps of method 300, and most particularly to optimize the first response signal and the second response signal for the cross-correlation step 330. This involves obtaining broadband signals. In this embodiment, each response signal is processed separately and a number of different operations can be performed (on each response signal). These operations will now be described.

[0085] The first operation is to remove the instrument response from the response signal. This operation is sometimes referred to as deconvolution and can be performed in a variety of ways that will be apparent to those skilled in the art.

[0086] The second operation is to remove the linear trend and the mean from the response signal. This is called detrending and involves removing aspects of the response signal that cause distortion (such as an overall linear increase in the mean) to reveal sub-trends. It is these sub-trends that generally better represent the ambient noise (and the shear velocity) at the surface.

[0087] The third operation is to reduce spectral leakage. In the context of the present disclosure and this embodiment, spectral leakage refers to the effect that occurs when waves caused by ambient noise do not have a frequency that is periodic with the sampling interval of the receiver measuring them. The effect is that the frequency distribution in the measured signal is not entirely accurate: for example, a particular frequency in the original wave may leak into (i.e., fall into) two adjacent frequencies in the measured data, thereby giving an inaccurate representation of the frequency profile (and its amplitude values) of the wave. To reduce this unwanted effect, the edges of the seismogram, which is a function of the noise measured at a particular receiver over time, are tapered off with a cosine taper. In this embodiment, the cosine taper is 5% of the trace length. Additionally, in this embodiment, the cosine taper is applied before applying the low-pass filter (the latter will now be described). In other embodiments, the taper value may be different.

[0088] The fourth operation is to anonymously filter out large amplitude values, such as values above an amplitude value threshold. The filtering is performed such that noise originating from transient events such as seismic or instrumental events, which typically result in shear velocity waves in the surface having larger amplitude values than the waves caused by ambient noise, is filtered out. If these large amplitude events are not filtered out, then they may overwhelm the ambient noise response.

[0089] More specifically, in certain embodiments, each response signal is plotted on a graph (over time), and any obvious outliers - i.e., data points with significantly abnormal amplitude values - are identified. This can be done either automatically or by visual inspection (say, by a person). If there are one or more obvious outliers, then further investigation is carried out. The further investigation may include inferring whether the outlier belongs to a certain frequency band, inferring whether it is instrumental noise, or inferring whether it is a transient signal (the signal we want to remove). If an obvious abnormal amplitude associated with an event exceeding the ambient noise amplitude is found, then the next step is to experiment with filtering values and techniques known to those skilled in the art to remove the abnormal data points. This is done to leave the rest of the signal as intact as possible.

[0090] The fifth operation is to divide each response signal into shorter segments. The response signals received in step 310 of method 300 may indicate ambient noise measured at or near the surface over a period of several days. In this embodiment, the response signals represent measurements over five days. It is beneficial to divide these potentially large data packets into shorter time segments in order to: allow the response signals to be synchronized (i.e., time-aligned); reduce the computational burden (thereby speeding up the processing); and allow the signals to be more easily manipulated in subsequent analysis to focus on good blocks (where the recording was done well) and avoid bad data blocks (where the recording was not done as well, e.g., due to the receiver being disturbed by wildlife). A consistent segmentation method is applied to the response signals. In this embodiment, the response signals are divided into 24-hour periods. In other embodiments, the duration of the segments can be shorter, longer, or in fact each response signal can be segmented into segments of different durations.

[0091] Optionally, step 320 may also include "chunking" the segmented response signals together. This means combining multiple segmented response signals into larger and longer signal blocks. These longer signal blocks can then be used in the cross-correlation of step 330 of method 300. For completeness, after the cross-correlation, the longer signal blocks resulting from the cross-correlation can be stacked to represent an even longer time period. In some embodiments, multiple segmented response signals are chunked together to represent a one-hour time window. Such a chunking process is sometimes referred to as concatenation. Chunking is useful for optimization, including making full use of computer resources.

[0092] The sixth operation is to adjust each response signal to the nearest second. This allows the response signals collected at different receivers to be synchronized (i.e., for time alignment), which is necessary in the subsequent parts of method 300. In this embodiment, the segmented response signals are adjusted to the nearest whole second. In other embodiments, the original or amplitude-filtered response signals can be adjusted (segmentation optionally occurring afterwards). And different adjustment values can be used: for example, the response signals can be adjusted to the nearest whole minute or the nearest millisecond.

[0093] The seventh operation is to apply low-pass filtering to each response signal (or actually to each segment of each response signal). This is to reduce the high-frequency components in the signal that may cause aliasing. For completeness, aliasing is an unwanted effect that occurs when the sampling frequency is not high enough to accurately sample the high-frequency components of a signal. In short, aliasing can lead to inaccurate data. However, by applying low-pass filtering and removing the high-frequency components, this effect can be minimized. In this embodiment, the cut-off value (i.e., the corner frequency) is set to one-quarter of the desired sampling frequency. In this embodiment, the desired sampling frequency is 1 Hz, so the exemplary cut-off value is 0.25 Hz.

[0094] The eighth operation can be to downsample the response signal. This reduces the amount of data to be processed in subsequent steps of method 300, thereby reducing memory usage, computational burden, and critically, the time taken to execute method 300. This time saving makes method 300 a practical option for performing site analysis even in the most time-pressured construction projects. More specifically, the response signal is downsampled by an integer n, such that only every nth sample is retained. In this embodiment, n is 100. In other embodiments, the value of n may be different. The value of n selected and used depends on the highest non-aliasing frequency that can be obtained in a particular setting.

[0095] Optionally, the downsampling can be done in integer steps, for example by a factor of 4 or 5. For example, if n is 100, then the data can be downsampled from 100 Hz to 20 Hz, from 20 Hz to 5 Hz, and finally to 1 Hz. For instance, this method is useful in embodiments where n is high, as the data may require a large amount of downsampling.

[0096] In some embodiments, the order of the operations is different from that described above. Additionally or alternatively, in some embodiments, the number of operations performed is different. For example, a sub-selection set of the operations that is most suitable for the particular project underway can be executed. Such a sub-selection set is useful for adapting the analysis to the project requirements, which may include one or more of time, cost, the granularity of the final model required, and the size of the target area.

[0097] In step 330, the processed response signals are cross-correlated. Here, the processed response signals are the first and second processed response signals. In this section, the processed response signals can be simply referred to as response signals. The purpose of cross-correlation is to simulate the wave field that would be recorded at the receiver of one response signal if the receiver of the other response signal were a virtual source. In other words, for the case where one of the first and second receivers is a virtual source, the wave field between the first receiver and the second receiver is simulated. In other words, relative to the response signal recorded at one of the first and second receivers, the response signal recorded at the other of the first and second receivers is measured. This gives the surface wave field propagating between the receivers. Advantageously, by cross-correlating the response signals, the travel time difference in the noise wave field between the first receiver and the second receiver can be obtained. Step 330 is an important step towards generating the model in step 340 of method 300.

[0098] In this embodiment, the response signals have been segmented and adjusted, so the first step in cross-correlation is to select the adjusted segments from two (processed) response signals corresponding to the same time period. The selection of a particular segment (and time period) is based on various factors, such as the data quality within that time period. In other embodiments, the response signals may not have been adjusted and / or may not have been processed at all.

[0099] The next step is to cross-correlate the response signals (or segments) to obtain the Green's function representing the wave field between the two receivers, as if one of the receivers were a virtual source.

[0100] In this embodiment, at least two response signals are cross-correlated in the frequency domain; however, cross-correlation can also be performed in the time domain. Conveniently, the user is able to select the domain with a shorter expected processing time. This helps to reduce the total time required for method 300.

[0101] Before proceeding to step 340, the result from step 330 can be normalized.

[0102] In other embodiments, in the case where there are more than two receivers in the dataset input to method 300 and thus more than two response signals, each response signal (or processed response signal) is cross-correlated with each other response signal, and a Green's function is obtained from each cross-correlation result. In other words, each receiver is paired with each other receiver, and one Green's function (and wavefield simulation) is obtained for each pair of receivers or each pair of virtual source-receivers. The output of the cross-correlation between the response signal measured at a particular receiver chosen as the virtual source and the response signals measured at each other receiver among a plurality of receivers is a virtual source gather that shows the Green's function between the virtual source and each other receiver among the plurality of receivers, thereby generating a virtual shot gather. As described above, Figure 2 shows an exemplary wavefield simulation of a plurality of receivers from a single central geophone (a type of receiver).

[0103] Generally speaking, the result of step 330 is a representation of the surface waves of the target area, which can subsequently be used as an indirect measurement input for generating a final model of the land properties of interest (obtained by solving the inverse problem of the known response signals from a ground area with unknown properties). One way to solve the inverse problem is tomographic inversion, which will now be described.

[0104] In step 340, the cross-correlated response signals are used in the tomographic inversion, which here are the result of the cross-correlation of the first and second processed response signals. The result of the tomographic inversion in step 350 is a model of the target area.

[0105] Given the fact that the ambient noise sources measured at the surface of the target area are not likely to be evenly distributed, in this embodiment, only the fundamental mode (i.e., the first resonance) surface waves are included in the tomographic inversion. This is because, despite the lack of uniformity of the ambient noise, these surface waves are assumed to be well reconstructed for a pair of cross-correlated receivers.

[0106] The preparatory step in the exemplary tomographic inversion operation described herein is to determine the average phase velocity between the virtual source-receiver pairs. Then the average group velocity is calculated from the average phase velocity.

[0107] Now will first refer to Figure 4 describe the exemplary tomographic inversion operation. Figure 4 Figure 400 shows a method for performing tomographic inversion, which can be used in step 340 of method 300. There are two different ways to perform method 400: a one-step method or a two-step method. Now the steps of method 400 that are common to both ways will be described.

[0108] In step 410, an initial model of one or more soil properties of the target area is provided. Advantageously, by using the initial model, the model generated in method 300 can be reached more quickly.

[0109] In step 420, a noise input is provided to the initial model. This noise input is a theoretical noise input.

[0110] In step 430, the travel times of the response symbols to be measured at the first receiver and the second receiver are calculated as a result of the noise input.

[0111] In step 440, the travel times calculated in step 430 are compared with the first response signal and the second response signal that are cross-correlated (i.e., the output from step 330).

[0112] In step 450, the initial model is updated based on the result of the comparison.

[0113] In step 460, based on the update of the initial model, a two-dimensional or three-dimensional model of the target area is generated according to the one or more soil properties.

[0114] This inversion can be performed in the time domain or the frequency domain. Advantageously, this means that the user can choose the method that results in a shorter processing time.

[0115] Now, a one-step method and a two-step method for performing method 400 will be described.

[0116] For the one-step method and the two-step method, it is important to understand: (1) the unit grid that makes up the model of the target area, and (2) the ray paths. Therefore, descriptions of (1) and (2) will now be given.

[0117] Please refer to Figure 5 , the exemplary shear wave velocity model 500 has a unit grid mxy that spans the surface above the target soil mass below the surface, including columns mx1, mx2, mx3, etc. that extend along the x direction, and rows m1y, m2y, m3y, etc. that extend along the y direction. Each unit defines an area of the surface. For example, each unit can define a 5m×5m square; other exemplary options include a 1m×1m square or a 10m×10m square. In other words, the shear wave velocity model includes a plurality of units arranged in a two-dimensional grid 500. This two-dimensional grid spans at least the surface area above the target area below the surface. Selecting a smaller unit area can improve the resolution of the model. Each unit also includes a soil mass that extends vertically below the surface area.

[0118] In some embodiments, the shear wave velocity model 500 is not a grid of square cells, but includes cells having other tessellation shapes including rectangular shapes, diamond shapes, or combinations including irregular shapes or different shapes.

[0119] Each cell is associated with a shear wave velocity value, which is an instance of a soil property value that represents the expected value of the shear wave velocity in the actual subsurface target area of interest. The shear wave value carries depth information, either because the shear wave value is constant throughout the target area below the cell area, or because the shear wave value varies with depth in a certain way. For example, the shear wave value defined for each cell can be an explicit function of depth, or a continuous function, or a series of values, each having an associated depth range. In another example, the shear wave value can be a function of frequency, which corresponds to depth information because surface wave propagation is affected by the physical properties of the subsurface soil within approximately one wavelength depth range. In other words, low-frequency surface waves are affected by the physical properties at deeper depths compared to high-frequency surface waves.

[0120] At the same time, the ray path is defined as the line between a first receiver and a second receiver. According to ray theory, the ray path represents the motion of the surface wave as it propagates from the first receiver to the second receiver (or vice versa). In particular, the ray path is defined as the propagation direction of the surface wave, i.e., the direction perpendicular to the wavefront in wave theory or perpendicular to the travel time contour lines. More generally, the ray path is the line between the receivers in a pair of receivers.

[0121] Also important for the one-step and two-step methods is the initial model of step 410 of method 400. This initial model sets an initial physical property value for each cell of the initial model, and method 400 will refine this initial physical property value using the response signal. Therefore, it is not necessary for the initial model and the initial physical property values to be a high-precision or high-resolution model of the subsurface target area, although a more accurate initial model can improve the expected accuracy of the final result of the method or reduce the computational time to reach the final result.

[0122] In some embodiments, an initial model is determined based on user input (e.g., historical data or map information indicating possible land property values for a target area beneath the entire surface). Alternatively, the first response signal and the second response signal (and other response signals, if there are more than two) can be used to determine the initial model. This is done by inverting the group velocity or phase velocity dispersion curve between the first response signal and the second response signal (and other signals, if present), which can be calculated from the cross-correlation of the response signals as described above, to find using a coarser grid or a faster method. As a last resort, arbitrarily selected typical land property values can be used as a starting point for each cell. In these embodiments, an arbitrary model can be selected based on an estimate of the land properties of the target area beneath the surface.

[0123] The two-step method will now be described. This is described in the context of having multiple receivers. In other embodiments, there can be only a first receiver and a second receiver.

[0124] In this embodiment, the two-step method includes selecting a subset of multiple pairs of virtual source-receiver, where the surface wave ray path between each pair of virtual source-receiver (i.e., the ray path between the virtual source and the receiver in each pair) passes through two or more cells. Advantageously, this selection process reduces the computational amount required to generate a model of the target area, thus accelerating the method.

[0125] Next, the method includes performing tomographic inversion using the cross-correlated response signals of each pair of virtual source-receiver in the subset. In this way, the land property value of each cell can be obtained.

[0126] The tomographic inversion of the two-step method involves a method of performing travel time tomography on the group travel time information obtained for each selected pair of source-receiver. This process maps the group travel time information from each pair of source-receiver tomographically into the cells of a physical specific model. Thus, the result of this process is an empirical model of the group or phase velocity (depending on the travel time data used), with each cell of the model having its respective group or phase velocity value. This process is performed for each of multiple frequencies to obtain the group or phase velocity value of each cell for each frequency. This process is the first stage (i.e., the tomographic stage) of the two-step method of tomographic inversion.

[0127] More specifically, this first stage involves first obtaining an initial group velocity model. This can be done in substantially the same manner as described above.

[0128] Next, the first stage involves determining the simulated travel times for each pair of selected source - receivers using an initial group velocity model. This is achieved by identifying the wave path from the source to the receiver (e.g., a direct ray, a curved ray, or a Fresnel zone) and identifying the cells through which the wave path passes. Then, the simulated travel time based on the initial velocity model is determined based on the known distance between the source and the receiver and the velocity values of each cell through which the wave path passes.

[0129] Next, the first stage involves determining an error value indicative of the difference between the simulated travel time and the empirical travel time for each pair of selected source - receivers.

[0130] Next, the first stage involves using this error value to determine an updated initial group velocity model. The updated model is generally the same in all respects as the initial group velocity model except for the new group velocity values associated with at least some of the cells. In other words, the updated model is an updated version of the initial group velocity model that takes into account the determined error value between the empirical group travel time and the simulated group travel time for each pair of source - receivers. This feedback process can involve least - squares methods, Markov chain Monte Carlo methods, or other inversion techniques to iteratively update the initial model based on the updated model.

[0131] The steps illustrated are generally part of a sub - routine of the first - stage tomography process, and this sub - routine is then iterated according to an inversion method such as least - squares inversion. The iteration continues until the error value reaches an end condition, such as the error value being below a threshold.

[0132] Once the tomography of the first - stage process has been performed, this two - step tomography inversion method can proceed to the second stage: inversion. This second stage is to obtain a resulting model of the subsurface properties of the target area. As in the first - stage tomography process, the second - stage inversion process is performed for each of a plurality of frequencies to obtain the group or phase velocity values for each cell at each frequency.

[0133] The second stage includes a first step of obtaining an initial model of the subsurface properties of the underground area. This can be done in substantially the same manner as described above. The initial model sets initial physical property values for each cell of the model, and this method will use a further iterative process to improve these initial physical property values using the group velocity model obtained from the first - stage (iterative) tomography process.

[0134] Next, the second stage includes determining the simulated surface - wave velocity for each cell based on the initial model. For this purpose, a forward - simulation method is used.

[0135] Next, the second stage includes determining an error value indicative of the difference between the simulated velocity (determined based on the initial model) and the empirical velocity (obtained from the first stage tomography process) for each cell of the model. The error value for each cell can be determined in a manner similar to that described for the first stage with respect to determining the error value between the simulated travel time and the empirical travel time for each source-receiver pair. In some embodiments, the initial model can be determined based on empirical phase dispersion data between source-receiver pairs.

[0136] Next, the second stage includes determining an updated physical model based on the determined error values. The updated physical model is generally the same in all respects as the initial physical model except for the new physical property values associated with at least some of the cells. In other words, the updated model is an updated version of the initial physical model that takes into account the determined error values between the empirical surface wave velocity and the simulated surface wave velocity for each cell of the model. This feedback process can involve least squares, Markov chain Monte Carlo, or other inversion techniques to iteratively update the initial physical model based on the updated physical model.

[0137] The steps described are generally part of a subroutine of the second stage tomography process, and then the subroutine is iterated according to an inversion method such as least squares inversion. The iteration continues until the error value reaches an end condition, such as the error value being below a threshold.

[0138] Finally, in step 460, based on the one or more generated land characteristics, the final model is a two-dimensional or three-dimensional model of the target area (and optionally output to a user device).

[0139] The one-step method will now be described in more detail. As for the two-step method, this is described in the context of having multiple receivers. In other embodiments, there can be only a first receiver and a second receiver.

[0140] In short, the two main differences between the one-step method and the two-step method are: (1) in the two-step method, the analysis is performed cell by cell, while in the one-step method, the analysis is based on each selected ray path; and (2) in the two-step method, tomography and inversion are combined because they involve the same step of the method but are performed sequentially as two separate steps, while in the one-step method, tomography and inversion are combined into a single step.

[0141] The method of the one-step method includes selecting a plurality of ray paths from the total number of possible ray paths between receiver pairs. To reduce the processing burden, it is generally beneficial to select a subset of the ray paths between the receivers. The end result of this selection step is the selection of ray paths, and the dispersion function (e.g., group velocity dispersion function and / or phase velocity dispersion function) of the ray paths can be determined using the response signals of the receivers at each end of the corresponding ray paths.

[0142] Next, the one-step method includes determining an empirical dispersion function for each ray path, particularly for each ray path selected in the selection part of the method. The dispersion function can be a group velocity dispersion function or a phase velocity dispersion function (or both can be used). According to any common method in the art, the group velocity dispersion function can be determined as described above with reference to Figure 1 the group velocity dispersion function described.

[0143] Next, the method includes determining a simulated dispersion function for each ray path using an initial model.

[0144] The steps of determining the empirical dispersion function and the simulated dispersion function can be performed in any order (or actually simultaneously).

[0145] Next, the method includes determining an error value indicative of the difference between the simulated dispersion function (determined according to the model) and the empirical dispersion function (determined according to the detected response signal) for each ray path.

[0146] Finally, the method includes updating the initial model using the error value. The updated model is generally the same in all aspects as the initial model except for the new shear wave velocity values (or any physical property values used) associated with at least some of the cells. In other words, the second model is an updated version of the first model that takes into account the determined error value between the empirical dispersion function and the simulated dispersion function. This feedback process is typically built into inversion methods such as the least squares method, the Markov chain Monte Carlo method, etc., and is performed as part of the inversion program.

[0147] The parts of the method for determining the simulated dispersion function, determining the error value, and determining the updated model using the error value are generally part of a subroutine of the one-step method, and then the subroutine is iterated according to an inversion method such as the least squares gradient descent inversion. In other words, each time an updated model is determined using the error value generated from the initial model and the simulated dispersion function for each ray path obtained, and then the updated model obtained is used as the initial model for the next iteration. The iteration continues until the error value reaches an end condition, such as the error value being below a threshold.

[0148] Finally, in step 460, based on the one or more land characteristics generated, the final updated model is a two-dimensional or three-dimensional model of the target area (and optionally output to a user device).

[0149] The advantage of the one-step method is that by focusing on the ray paths and performing single-step tomography and inversion based on the ray paths (say, not for each cell), the computational burden is reduced, thus reducing the time spent.

[0150] For completeness, the tomographic inversion operations used in method 400 may be similar but different from those used in multi-channel analysis of surface waves (MASW). MASW has been described and it is prior art for collecting surface wave information. Key differences between the described tomographic inversion and the MASW method include that MASW does not use tomographic inversion (but rather one-dimensional (1D) inversion); MASW is used in conjunction with an active noise source (such as a sledgehammer or a falling weight), and MASW studies a two-dimensional line of interest along the surface. The target depth in MASW is approximately 5 to 30 m, so MASW cannot penetrate as deep as the method of the present disclosure.

[0151] Please refer to Figure 6 , an exemplary output of the method described herein is a final shear wave velocity model showing the target area beneath the surface. The target area beneath the surface extends to a depth of 100 m (z-direction) along the x and y directions. The values of the shear wave velocity are shown by the shading in the figure, and the transitions between regions of different shear wave velocities are visible, indicating different compositions or structures of the various parts of the target area beneath the surface. Using the method described herein, the final shear wave velocity model can be determined with higher resolution and accuracy without resulting in an infeasible computation time. Such a subsurface model can be used to better understand the suitability of the target soil mass beneath the surface for supporting man-made structures above or within the target area beneath the surface.

[0152] As previously described, there may be more than two receivers. For example, there may be 100 receivers, 100 to 5000 receivers, 5000 receivers, or more than 5000 receivers.

[0153] There are many different types of receivers suitable for use in method 300. For example, the receivers in method 300 can be any one of (or any combination of) geophones, accelerometers, seismographs, vibration sensors, and transducers.

[0154] The target area can have different depths and positions relative to the surface. For example, in terms of subsurface depth, method 300 may be suitable for determining one or more soil properties of a target area with a depth of up to approximately 100 meters. For example, the target area can be 0 to 100 m beneath the surface; it can be 0 to 45 m beneath the surface; it can be 50 to 100 m beneath the surface. Thus, the target area can be entirely underground. This makes the described method well-suited for a variety of different construction projects, including underground projects.

[0155] In some embodiments, the processing of the response signal in step 320 includes applying spectral whitening. This technique helps to enhance the representation of the frequencies of interest in the ambient noise, thereby avoiding signal dominance in the microseismic frequency band in the cross-correlation performed in step 330.

[0156] The model and the one or more soil properties may be soil properties other than the shear wave velocity, or additional soil properties, such as the compressional wave velocity, density, elastic modulus, shear modulus, or, if a viscoelastic model is used, optionally also the viscosity quality factors Qs and Qp. Generally, the model may define multiple soil property values.

[0157] Figure 7 A block diagram of one embodiment of a computing device 700 is shown, in which a set of instructions may be executed to cause the computing device to perform any one or more of the methods discussed herein. In alternative embodiments, the computing device may be connected (e.g., networked) to other machines on a local area network (LAN), intranet, extranet, or the Internet. The computing device may operate as a server or a client in a client - server network environment, or as a peer machine in a peer - to - peer (or distributed) network environment. The computing device may be a personal computer (PC), a tablet computer, a set - top box (STB), a personal digital assistant (PDA), a cellular phone, a network device, a server, a network router, a switch, or a bridge, or any machine capable of executing a set of instructions (sequentially or otherwise) that specify actions to be taken by that machine. Further, although only a single computing device is shown, the term "computing device" should also be understood to include any collection of machines (e.g., computers) that individually or jointly execute a set (or multiple sets) of instructions to perform any one or more of the methods discussed herein.

[0158] Exemplary computing device 700 includes a processor 702, a main memory 704 (e.g., read - only memory (ROM), flash memory, dynamic random access memory (DRAM) (e.g., synchronous DRAM (SDRAM) or Rambus DRAM (RDRAM), etc.), a static memory 706 (e.g., flash memory, static random access memory (SRAM), etc.), and auxiliary memory (e.g., data storage device 718), which communicate with each other via a bus 730.

[0159] The processor 702 represents one or more general-purpose processors, such as a microprocessor, a central processing unit, etc. More specifically, the processor 702 can be a complex instruction set computing (CISC) microprocessor, a reduced instruction set computing (RISC) microprocessor, a very long instruction word (VLIW) microprocessor, a processor implementing other instruction sets, or a processor implementing a combination of instruction sets. The processor 702 can also be one or more special-purpose processors, such as an application-specific integrated circuit (ASIC), a field-programmable gate array (FPGA), a digital signal processor (DSP), a network processor, etc. The processor 702 is configured to execute processing logic (instructions 722) for performing the operations and steps discussed herein.

[0160] The computing device 700 may also include a network interface device 708. The computing device 700 may also include a video display unit 710 (such as a liquid crystal display (LCD) or a cathode ray tube (CRT)), an alphanumeric input device 712 (such as a keyboard or a touch screen), a cursor control device 714 (such as a mouse or a touch screen), and an audio device 716 (such as a speaker).

[0161] Obviously, Figure 7 Some of the features of the computer device 700 shown may not exist. For example, one or more computing devices 700 may not require a display device 710 (or any associated adapter). For example, this may be the case for a specific server-side computer device 700 that only utilizes its processing power and does not need to display information to the user. Similarly, the user input device 712 may not be required. In its simplest form, the computer device 700 includes a processor 702 and a memory 704.

[0162] The data storage device 718 may include one or more machine-readable storage media (or more specifically, one or more non-transitory computer-readable storage media) 728, on which a set or multiple sets of instructions 722 are stored, and the instructions 1122 embody any one or more of the methods or functions described herein. During the execution of the instructions 722 by the computer system 704, the instructions 1122 may also reside entirely or at least partially in the main memory 702 and / or the processor 700, and the main memory 704 and the processor 702 also constitute computer-readable storage media.

[0163] The various methods described above can be implemented by a computer program. The computer program can include computer code that is arranged to direct a computer to perform the functions of one or more of the various methods described above. The computer program and / or code for performing such a method can be provided to a device, such as a computer, on one or more computer-readable media, or more generally on a computer program product. The computer-readable media can be transient or non-transient. The one or more computer-readable media can be, for example, an electronic, magnetic, optical, electromagnetic, infrared, or semiconductor system, or a propagation medium for data transmission, such as for downloading code over the Internet. Alternatively, the one or more computer-readable media can take the form of one or more physical computer-readable media, such as semiconductor or solid state memory, magnetic tape, removable computer disk, random access memory (RAM), read-only memory (ROM), hard disk, and optical disk, such as CD-ROM, CD-R / W, or DVD.

[0164] In one embodiment, the modules, components, and other features described herein can be implemented as discrete components or integrated into the functionality of hardware components such as ASICS, FPGAs, DSPs, or similar devices.

[0165] A “hardware component” is a tangible (e.g., non-transient) physical component (e.g., a collection of one or more processors) capable of performing certain operations and can be configured or arranged in a particular physical manner. A hardware component can include dedicated circuitry or logic permanently configured to perform certain operations. A hardware component can be or include a dedicated processor, such as a field programmable gate array (FPGA) or an ASIC. A hardware component can also include programmable logic or circuitry temporarily configured by software to perform certain operations.

[0166] Thus, the phrase “hardware component” should be understood to encompass a tangible entity capable of being physically constructed, permanently configured (e.g., hardwired), or temporarily configured (e.g., programmed) to operate in a certain manner or to perform certain operations described herein.

[0167] In addition, the modules and components can be implemented as firmware or functional circuitry within a hardware device. Further, the modules and components can be implemented in any combination of hardware devices and software components, or solely in software (e.g., code stored or otherwise embodied in a machine-readable medium or a transmission medium).

[0168] Except as otherwise specifically stated, it will be apparent from the following discussion that, throughout the specification, discussions using terms such as "provide", "calculate", "update", "generate", "receive", "process", "determine", "select", "compare", and "identify" refer to actions and processes of a computer system or similar electronic computing device that manipulate and transform data represented as physical (electronic) quantities in the registers and memories of the computer system into other data similarly represented as physical quantities in the memories or registers of the computer system or other such information storage, transmission, or display devices.

[0169] It should be understood that the above description is merely exemplary and not restrictive. After reading and understanding the above description, many other embodiments will be apparent to those skilled in the art. Although the present disclosure is described with reference to specific exemplary embodiments, it should be recognized that the present disclosure is not limited to the embodiments described, but can be implemented with modifications and variations within the scope of the appended claims. Therefore, the specification and drawings should be regarded as exemplary rather than restrictive. The scope of the present disclosure should be determined with reference to the appended claims and the full scope of equivalents of these claims.

Claims

1. A method for determining one or more land characteristics of a target area beneath the Earth's surface, the method comprising: Receiving a data set, the data set comprising: A first response signal indicative of ambient noise measured at or near the Earth's surface by a first receiver disposed at a first location; And A second response signal indicative of ambient noise measured at or near the Earth's surface by a second receiver disposed at a second location, wherein the first location and the second location are different; Processing the first response signal and the second response signal; Cross-correlating the processed first response signal with the processed second response signal; And Performing tomographic inversion using the cross-correlated response signals to generate a two-dimensional '2D' model or a three-dimensional '3D' model of the target area based on the one or more land characteristics.

2. The method according to claim 1, wherein, The target area is at most 100 meters beneath the Earth's surface.

3. The method according to any one of the preceding claims, wherein, Making a land risk management framework decision based on the two-dimensional model or the three-dimensional model.

4. The method according to any one of the preceding claims, wherein, Each of the first response signal and the second response signal respectively indicates a vertical component and / or a horizontal component of ambient noise measured at or near the Earth's surface by the first receiver and the second receiver.

5. The method according to any one of the preceding claims, comprising determining the locations of the first receiver and the second receiver based on the minimum depth and / or the maximum depth of the target area and / or the desired resolution of waves in the target area caused by ambient noise, prior to the step of receiving the data set.

6. The method according to any one of the preceding claims, including selecting a recording frequency for the first receiver and the second receiver before the step of receiving the data set, wherein, Selecting the recording frequency based on the depth of the target area and / or the expected wavelength of noise waves in the target area.

7. The method according to any one of the preceding claims, wherein, The one or more land characteristics include one or more elastic characteristics of the target area, such as shear velocity Vs.

8. The method according to any one of the preceding claims, wherein, The step of processing includes processing the first response signal and the second response signal to enhance the representation of the ambient noise.

9. The method according to any one of the preceding claims, wherein, The step of processing includes one or more of the following operations: Dividing each response signal into segments; Aligning each segmented response signal to the nearest second; Applying low-pass filtering to each response signal; and / or Downsampling each response signal.

10. The method according to any one of the preceding claims, wherein, The step of cross-correlating at least two processed response signals includes estimating the Green's function between the first receiver and the second receiver.

11. The method according to any one of the preceding claims, wherein, The step of performing tomographic inversion includes: Providing an initial model of the one or more land characteristics of the target area; Providing a noise input for the initial model; Calculating, using the initial model, the travel times of response signals that would be measured at the first receiver and the second receiver; Comparing the travel times calculated using the initial model with the cross-correlated first and second response signals; and Updating the initial model based on the result of the comparison to generate the two-dimensional model or the three-dimensional model.

12. The method according to claim 11, wherein The step of providing the initial model includes determining a phase velocity dispersion curve between a virtual source and a selected one of the first receiver and the second receiver.

13. A computer program product comprising instructions which, when the program is executed by a computer, cause the computer to perform the method according to any of the preceding claims.

14. A system comprising: one or more processors; one or more memories storing computer-readable instructions thereon, the computer-readable instructions being configured to cause the one or more processors to perform operations including the steps according to any one of claims 1 to 12.

15. A computer-readable medium comprising instructions which, when executed by a computer, cause the computer to perform the method according to any one of claims 1 to 12.