cpt_layering

cpt_layering(*args)

A function that uses the agglomerative clustering algorithm by Hudson et al. (2023) to automatically identify spatially contiguous soil layers with similar q_c1Ncs and I_c 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 q_c1Ncs and I_c values. Details are contained in the paper cited below.

Citation

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

Parameters

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

Note

arguments in <brackets> are optional. If Nmax = None, 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 the number of data points, whichever occurs first. If you wish to set the number of layers to a specific value, you can do so by setting Nmin = Nmax.

Returns

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.

Example

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)
dGWT = 1.0
ztop, zbot, qc1Ncs_lay, Ic_lay = smt.cpt_layering(qc1Ncs, Ic, depth, dGWT, 1, 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()