-
Notifications
You must be signed in to change notification settings - Fork 1
cpt_layering
A function that uses the agglomerative clustering algorithm by Hudson et al. (2023) to automatically identify spatially contiguous soil layers with similar qc1Ncs and Ic values. The algorithm automatically selects the number of layers, and returns the depths at the top and bottom of each layer, along with the corresponding qc1Ncs and Ic values. Details are contained in the paper cited below.
Hudson, K. S., Ulmer, K., Zimmaro, P., Kramer, S. L., Stewart, J. P., and Brandenberg, S. J. (2023). "Unsupervised Machine Learning for Detecting Soil Layer Boundaries from Cone Penetration Test Data." Earthquake Engineering and Structural Dynamics, 52(11). link
| argument | format | description |
|---|---|---|
| qc1Ncs | Numpy array, dtype=float | overburden- and fines-corrected cone tip resistance |
| Ic | Numpy array, dtype=float | soil behavior type index |
| depth | Numpy array, dtype=float | depth to cpt measurement point |
| <dGWT> | float | depth to groundwater table |
| <Nmin> | integer | minimum number of layers. optional. default = 1. |
| <Nmax> | integer | maximum number of layers. optional. default = None |
| <tref> | float | constant used to define layer thickness parameter in cost function. optional. default = 0.5 |
| <wD> | float | weight factor for distortion score portion of cost function |
| <wT> | float | weight factor for layer thickness portion of cost function |
Note
arguments in <brackets> are optional keyword arguments that should be specified like "kwarg=value" as shown in the example. If Nmax = Nmin, the algorithm will use a single number of layers. If Nmax > Nmin, the optimal number of layers will be obtained by incrementing the number of layers, starting at Nmin, until the minimum value of the cost function is reached, or the number of layers is equal to Nmax, whichever occurs first. If Nmax = None, it is set equal to the number of data points.
| argument | format | description |
|---|---|---|
| ztop | Numpy array, dtype=float | depth to top of each layer |
| zbot | Numpy array, dtype=float | depth to bottom of each layer |
| qc1Ncs_lay | Numpy array, dtype=float | qc1Ncs values for each layer |
| Ic_lay | float | Ic values for each layer |
Note
The ztop and zbot values returned by the algorithm result in a spatially contiguous profile, without gaps, such that the sum of layer thicknesses, defined as zbot - ztop, is the same as the profile thickness. Layer boundaries are set at the midpoint between adjacent clusters, as illustrated by the example below.
import numpy as np
import ngl_tools.smt as smt
import matplotlib.pyplot as plt
qc1Ncs = np.asarray([10, 10, 10, 20, 20, 20, 30, 30, 30], dtype=float)
Ic = np.asarray([3, 3, 3, 2, 2, 2, 1, 1, 1], dtype=float)
depth = np.asarray([0, 1, 2, 3, 4, 5, 6, 7, 8], dtype=float)
ztop, zbot, qc1Ncs_lay, Ic_lay = smt.cpt_layering(qc1Ncs, Ic, depth, dGWT=1.0, Nmin=1, Nmax=4)
def interleave(v1, v2):
vout = np.empty(2 * len(v1), dtype=float)
vout[0::2] = v1
vout[1::2] = v2
return vout
fig, ax = plt.subplots(ncols=2, sharey='row')
ax[0].scatter(qc1Ncs, depth, facecolor='white', edgecolor='C0', s=10, zorder=2, label='data')
ax[1].scatter(Ic, depth, facecolor='white', edgecolor='C0', s=10, zorder=2)
ax[0].plot(interleave(qc1Ncs_lay, qc1Ncs_lay), interleave(ztop, zbot), c='C1', zorder=1, label='clustered')
ax[1].plot(interleave(Ic_lay, Ic_lay), interleave(ztop, zbot), c='C1', zorder=1)
for a in ax:
a.grid(True, alpha=0.5)
ax[0].set_ylabel('depth')
ax[0].set_xlabel(r'$q_{c1Ncs}$')
ax[1].set_xlabel(r'$I_c$')
ax[0].invert_yaxis()
ax[0].legend()