From 30bb0e017322e14772edf2fd1a600649d14cd2ba Mon Sep 17 00:00:00 2001
From: Peter Powers <pmpowers@usgs.gov>
Date: Tue, 14 Jun 2022 20:44:20 -0600
Subject: [PATCH] kuehn cleanup
---
.../earthquake/nshmp/gmm/KuehnEtAl_2020.java | 459 +++++-------------
1 file changed, 121 insertions(+), 338 deletions(-)
diff --git a/src/main/java/gov/usgs/earthquake/nshmp/gmm/KuehnEtAl_2020.java b/src/main/java/gov/usgs/earthquake/nshmp/gmm/KuehnEtAl_2020.java
index c8a578f4..352ec391 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/gmm/KuehnEtAl_2020.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/gmm/KuehnEtAl_2020.java
@@ -26,22 +26,28 @@ import gov.usgs.earthquake.nshmp.gmm.GroundMotionTables.GroundMotionTable.Positi
import gov.usgs.earthquake.nshmp.tree.LogicTree;
/**
- * Implementation of the PEER NGA-Subduction GMM, Kuehn et al., 2020 GMM
- * (KBCG20)
+ * Implementation of the Kuehn et al. (2020) ground motion model for subduction
+ * regions developed as part of the <a
+ * href="https://www.risksciences.ucla.edu/nhr3/nga-subduction">NGA
+ * Subduction</a> project.
*
- * Notes:
+ * <p><b>Implementation notes:</b>
*
- * This model is only applicable to sites in the forearc, future residuals
- * analysis is planned to determine if there are attenuation differences between
- * forearc and backarc zones.
+ * <ul><li>Implementations are provided for the global model as well as Cascadia
+ * and Alaska regionalized models, along with a specialized model for use in the
+ * deepest parts of the Seattle basin.
*
- * PSBAH20 uses zHyp (randomized source location on fault), but this can be
- * replaced this event-specific zHyp with the mean depth expected for the fault
- * plane, zHypBar, in Eq. 4.6 without modification to coefficients or
- * uncertainty.
+ * <li>Simplified epistemic uncertainty tables for the model were obtained from
+ * <a href="https://github.com/nikuehn/KBCG20/tree/master/UNCERTAINTY">
+ * https://github.com/nikuehn/KBCG20/tree/master/UNCERTAINTY</a>.
*
- * Path term (Eqs. 4.2 and 4.3 in PEER report) appears to differ in the official
- * Python implementation. Compare implementations of getPathTerm() below
+ * <li>The model, as published does not include a basin term for Alaska, however
+ * if {@code z}-values are supplied then the general Cascadia basin term is
+ * used.</ul>
+ *
+ * <p><b>Note:</b> Direct instantiation of {@code GroundMotionModel}s is
+ * prohibited. Use {@link Gmm#instance(Imt)} to retrieve an instance for a
+ * desired {@link Imt}.
*
* <p><b>Reference:</b> Kuehn, N., Bozorgnia, Y, Campbell, K.W., Gregor, N.
* (2020) Partially non-ergodic ground-motion model for subduction regions using
@@ -49,51 +55,45 @@ import gov.usgs.earthquake.nshmp.tree.LogicTree;
*
* <p><b>Component:</b> average horizontal (RotD50)
*
- * Questions:
- *
- * * δMbreak, defined in Figure 4.2. The paper does not mention regionalization
- * for this parameter, but in the Python code only Japan and South America have
- * non-zero values.
- *
- * * Source Zref(if|slab) are not specified in the paper. Python code specifies
- * 15 and 50 km, respectively
- *
- * * Source-Depth breakpoints
- *
- * * Eq. 4.2 - median PSA at short periods is not allowed to be less than PGA,
- * so, if PSA(T<0.1) < PGA, do we still use sigma(PSA(T<0.1)) or use
- * sigma(PGA)?
- *
- * * Eq. 4.8 in the paper includes `c ln( vs30/k1)^n`, this ln is omitted in the
- * Python and the other papers using this term (CB14, etc.)
- *
- * * I think fSite term in Python code has a typo, using `vs30 < k1` rather
- * than `vs30 ≤ k1` as specified in the paper (and others using the same site
- * term)
- *
- * * Value of k2 (hard-coded in Python) for PGV is a typo, the value should be
- * -1.955 but -1.995 is used
- *
- * * Basin Zref is only calculated (coefficients only provided) for CASC, JP,
- * NZ, and TW? Then is Zref == 0 for all other regions?
- *
- * * Basin term for Cascadia Other PacNW basins, the paper does not allow the
- * basin term to exceed that for the Seattle basin. The Python code does not
- * seem to impose this constraint but I may not be seeing it.
- *
- * * Uncertainty values obtained from:
- * https://github.com/nikuehn/KBCG20/tree/master/UNCERTAINTY
- *
* @author U.S. Geological Survey
*/
public abstract class KuehnEtAl_2020 implements GroundMotionModel {
+ /*
+ *
+ * Developer notes:
+ *
+ * (Some of the issues noted below will likely be cleared up when an
+ * Earthquake Spectra paper is published, as was the case with Parker et al.)
+ *
+ * δMbreak, defined in Figure 4.2. The paper does not mention regionalization
+ * for this parameter, but in the Python code only Japan and South America
+ * have non-zero values.
+ *
+ * Source Zref(if|slab) are not specified in the paper. Python code specifies
+ * 15 and 50 km, respectively
+ *
+ * Eq. 4.2 - median PSA at short periods is not allowed to be less than PGA,
+ * so, if PSA(T<0.1) < PGA, do we still use sigma(PSA(T<0.1)) or use
+ * sigma(PGA)?
+ *
+ * Eq. 4.8 in the paper includes `c ln( vs30/k1)^n`, this ln is omitted in the
+ * Python and the other papers using this term (CB14, etc.)
+ *
+ * I think fSite term in Python code has a typo, using `vs30 < k1` rather
+ * than `vs30 ≤ k1` as specified in the paper (and others using the same site
+ * term)
+ *
+ * Value of k2 (hard-coded in Python) for PGV is a typo, the value should be
+ * -1.955 but -1.995 is used
+ */
+
static final String NAME = "Kuehn et al. (2020)";
/*
- * Constraints for magnitude breakpoint, source depth (as noted below),
- * forearc/backarc subregions, basin depth are region-specific
+ * Note that author-defined constraints for magnitude breakpoint, source
+ * depth, forearc/backarc subregions, basin depth are region-specific.
*/
static final Constraints CONSTRAINTS_INTERFACE = Constraints.builder()
.set(MW, Range.closed(5.0, 9.5))
@@ -101,6 +101,7 @@ public abstract class KuehnEtAl_2020 implements GroundMotionModel {
.set(ZTOR, Range.closed(0.0, 50.0))
.set(VS30, Range.closed(150.0, 1500.0))
.build();
+
static final Constraints CONSTRAINTS_SLAB = Constraints.builder()
.set(MW, Range.closed(5.0, 8.5))
.set(RRUP, Range.closed(10.0, 1000.0))
@@ -108,31 +109,8 @@ public abstract class KuehnEtAl_2020 implements GroundMotionModel {
.set(VS30, Range.closed(150.0, 1500.0))
.build();
- // static final class Basin {
- // // Other region
- // static final String NONE = "None";
- // // PacNW outside of basins
- // static final String PACNW_NONE = "PacNW-no-basin";
- // // Tacoma/Everett basins
- // static final String OTHER = "Other";
- // // Seattle basin
- // static final String SEATTLE = "Seattle";
- // }
+ static final CoefficientContainer COEFFS = new CoefficientContainer("KBCG20.csv");
- private static final String KBCG20_CSV = "KBCG20.csv";
-
- /*
- * Define coefficient containers
- */
- static final CoefficientContainer COEFFS;
-
- static {
- COEFFS = new CoefficientContainer(KBCG20_CSV);
- }
-
- /*
- * Define constants
- */
private static final double δM = 0.1; // p. 20
private static final double MREF = 6.0; // p. 20
private static final double δZ = 1.0; // p. 22
@@ -165,24 +143,14 @@ public abstract class KuehnEtAl_2020 implements GroundMotionModel {
final double φ;
final double Ï„;
- final double Mbreak;
-
- /**
- * <pre>
- * Mb (different from Mbreak???) see getMbDefault() in python code
- * global Aleutian Alaska Cascadia
- * Interface 7.6 8.0 7.2 7.2
- * Intraslab 7.9 8.0 8.6 8.0
- * </pre>
- */
- final double Zbreak; // dzb_if dzb_slab
- final double Zref;
+ final double mBreak;
+ final double zBreak; // dzb_if dzb_slab
+ final double zRef;
Coefficients(Imt imt, CoefficientContainer cc, SubductionZone zone, boolean slab) {
this.imt = imt;
Map<String, Double> coeffs = cc.get(imt);
- // No regionalization or interface/slab
θ3 = coeffs.get("theta_3");
θ5 = coeffs.get("theta_5");
@@ -191,33 +159,17 @@ public abstract class KuehnEtAl_2020 implements GroundMotionModel {
θ4 = coeffs.get("theta_4_slab");
θ9 = coeffs.get("theta_9_slab");
- Zbreak = 80.0 + coeffs.get("dzb_slab");
- Zref = 50.0; // this value is not in the paper
+ zBreak = 80.0 + coeffs.get("dzb_slab");
+ zRef = 50.0; // this value is not in the paper
} else {
- // set interface coefs
θ2 = coeffs.get("theta_2_if");
θ4 = coeffs.get("theta_4_if");
θ9 = coeffs.get("theta_9_if");
- Zbreak = 30.0 + coeffs.get("dzb_if");
- Zref = 15.0; // this value is not in the paper
+ zBreak = 30.0 + coeffs.get("dzb_if");
+ zRef = 15.0; // this value is not in the paper
}
- /**
- * ***** mu_theta_6[b] (global) vs. theta_6xc vs. regional theta_6_*_reg,
- * especially for Alaska/Cascadia, which don't seem to use these coeffs
- *
- * <pre>
- * Global Alaska Cascadia
- * theta_6_x1 mu_theta_6b theta_6_x1_reg_Al theta_6_x1_reg_Ca
- * theta_6_x2 mu_theta_6b theta_6_x2_reg_Al theta_6_x2_reg_Ca
- * theta_6_x3 mu_theta_6b theta_6_x3_reg_Al theta_6_x3_reg_Ca
- * theta_6_1 mu_theta_6b theta_6_1_reg_Al theta_6_1_reg_Ca
- * theta_6_2 mu_theta_6 theta_6_2_reg_Al theta_6_2_reg_Ca <-- only one used ???
- * theta_6_3 mu_theta_6b theta_6_3_reg_Al theta_6_3_reg_Ca
- * </pre>
- */
-
/* For now we apply general Cascadia basin scaling elsewhere. */
θ11 = coeffs.get("intercept_Ca");
θ11S = coeffs.get("mean_residual_Seattle_basin");
@@ -264,59 +216,16 @@ public abstract class KuehnEtAl_2020 implements GroundMotionModel {
k1 = coeffs.get("k1");
k2 = coeffs.get("k2");
- Mbreak = getBreakpointMag(imt, zone, slab);
- // δZb = (isSlab) ? coeffs.get("dzb_slab") : coeffs.get("dzb_if");
-
- φ = coeffs.get("phi");
- Ï„ = coeffs.get("tau");
- }
-
- /*
- * Breakpoint magnitude, Mb, and δMbreak for the magnitude term (Eq. 4.4).
- *
- * The paper does not mention regionalization for δMbreak, but the Python
- * code only applies δMbreak for Japan and South America
- */
- private double getBreakpointMag(Imt imt, SubductionZone zone, boolean slab) {
- // breakpoint magnitudes are hard coded in python code, see getMbDefault()
-
- // global
- double Mb = slab ? 7.6 : 7.9;
+ double mb = slab ? 7.6 : 7.9; // global
if (zone == CASCADIA)
- Mb = slab ? 7.2 : 8.0;
+ mb = slab ? 7.2 : 8.0;
else if (zone == ALASKA) {
- Mb = slab ? 7.2 : 8.6;
+ mb = slab ? 7.2 : 8.6;
}
+ mBreak = mb;
- return Mb;
-
- /*
- * Figure 4.2: (only for Interface) δMbreak = 0.0 for T<=1.0s, -0.4 for
- * T>=3.0s, and log-linear interpolated in between
- *
- * ??? Python code only calculates this for Japan and South America,
- * δMbreak is 0.0 for all other regions. The paper makes no mention of
- * this term being regionalized
- */
- // double δMbreak = 0.0;
- // if (!slab) {
- // if (!imt.isSA() || imt.period() <= 1.0) {
- // δMbreak = 0.0;
- // } else if (imt.period() >= 3.0) {
- // δMbreak = -0.4;
- // } else {
- // δMbreak = 0.0 + (log(imt.period()) - log(1.0)) * (-0.4 - 0.0) /
- // (log(3.0) - log(1.0));
- // }
- // }
- // // To match Python code:
- // switch (region) {
- // case Region.JAPAN:
- // case Region.SA:
- // return Mb + δMbreak;
- // default:
- // return Mb;
- // }
+ φ = coeffs.get("phi");
+ Ï„ = coeffs.get("tau");
}
}
@@ -327,10 +236,8 @@ public abstract class KuehnEtAl_2020 implements GroundMotionModel {
KuehnEtAl_2020(Imt imt, SubductionZone zone) {
coeffs = new Coefficients(imt, COEFFS, zone, slab());
coeffsPGA = new Coefficients(PGA, COEFFS, zone, slab());
- epistemicTable = GroundMotionTables.getKbcg20(
- zone,
- (slab()) ? SubductionZone.Type.SLAB : SubductionZone.Type.INTERFACE,
- imt);
+ SubductionZone.Type szType = slab() ? SLAB : INTERFACE;
+ epistemicTable = GroundMotionTables.getKbcg20(zone, szType, imt);
}
abstract boolean slab();
@@ -340,48 +247,23 @@ public abstract class KuehnEtAl_2020 implements GroundMotionModel {
@Override
public LogicTree<GroundMotion> calc(GmmInput in) {
- // Calc PGAref
- double pgaRef = calcPgaRef(coeffsPGA, in.Mw, in.rRup, in.zTor); // exp(calcMean());
- pgaRef = exp(pgaRef);
-
+ double pgaRef = exp(calcPgaRef(coeffsPGA, in.Mw, in.rRup, in.zTor));
double μ = calcMean(coeffs, basin(), in.Mw, in.rRup, in.zTor, in.vs30, in.z2p5, pgaRef);
- /*
- * Eq. 4.2 - median PSA at short periods is not allowed to be smaller than
- * PGA
- *
- * TODO: If PSA(T<0.1) < PGA, do we still use sigma(PSA(T<0.1)) or use
- * sigma(PGA)?
- */
+ // short periods can't be lower than PGA
if (coeffs.imt.isSA() && coeffs.imt.period() <= 0.1) {
- // calculate PGA
- double μPga =
- calcMean(coeffsPGA, basin(), in.Mw, in.rRup, in.zTor, in.vs30, in.z2p5, pgaRef);
- // μ can't be smaller than PGA
+ double μPga = calcMean(coeffsPGA, basin(), in.Mw, in.rRup, in.zTor, in.vs30, in.z2p5, pgaRef);
μ = max(μ, μPga);
- // sigma is not period-dependent, so no change to sigma
}
- double sigmaAleatory = Maths.hypot(coeffs.φ, coeffs.τ);
-
- // epistemic uncertainty
- // double sigmaEpistemic = coeffs.ψ; // compute by sampling posterior
- // results
-
- // epistemic uncertainty from pre-computed tables
+ double σ = Maths.hypot(coeffs.φ, coeffs.τ);
Position p = epistemicTable.position(in.rRup, in.Mw);
double ψ = epistemicTable.get(p);
- // double sigmaTotal = Maths.hypot(coeffs.φ, coeffs.τ, sigmaEpistemic);
-
- double σ = sigmaAleatory;
-
- // return GroundMotions.createTree(μ, σ);
-
double[] μs = new double[] { μ + ψ, μ, μ - ψ };
-
- // return GroundMotions.createTree(μ, σ);
- return GroundMotions.createTree(GroundMotions.EPI_IDS, μs, GroundMotions.EPI_WTS, σ);
+ return GroundMotions.createTree(
+ GroundMotions.EPI_IDS, μs,
+ GroundMotions.EPI_WTS, σ);
}
private static double calcPgaRef(
@@ -389,26 +271,13 @@ public abstract class KuehnEtAl_2020 implements GroundMotionModel {
double Mw,
double rRup,
double zTor) {
- double fMag = getMagTerm(c, Mw);
- double fGeom = getGeomAttenTerm(c, Mw, rRup);
-
- double fDepth = getDepthTerm(c, zTor);
-
- double fAtten = getArcCrossingTerm(c, rRup);
-
- double fSite = getSiteAmpTerm(c, VSROCK, 0);
-
- double fBasin = 0;
-
- // System.out.printf("pgaRef fMag = %.8f\n", fMag);
- // System.out.printf("pgaRef fGeom = %.8f\n", fGeom);
- // System.out.printf("pgaRef fDepth = %.8f\n", fDepth);
- // System.out.printf("pgaRef fAtten = %.8f\n", fAtten);
- // System.out.printf("pgaRef fSite = %.8f\n", fSite);
- // System.out.printf("pgaRef fBasin = %.8f\n", fBasin);
-
- return c.θ1 + fMag + fGeom + fDepth + fAtten + fSite + fBasin;
+ double fMag = magnitudeTerm(c, Mw);
+ double fGeom = geometricAttenuationTerm(c, Mw, rRup);
+ double fDepth = depthTerm(c, zTor);
+ double fAtten = arcCrossingTerm(c, rRup);
+ double fSite = siteAmpTerm(c, VSROCK, 0);
+ return c.θ1 + fMag + fGeom + fDepth + fAtten + fSite;
}
private static double calcMean(
@@ -421,88 +290,60 @@ public abstract class KuehnEtAl_2020 implements GroundMotionModel {
double z2p5,
double pgaRef) {
- double fMag = getMagTerm(c, Mw);
-
- double fGeom = getGeomAttenTerm(c, Mw, rRup);
-
- double fDepth = getDepthTerm(c, zTor);
-
- double fAtten = getArcCrossingTerm(c, rRup);
-
- double fSite = getSiteAmpTerm(c, vs30, pgaRef);
-
- double fBasin = basin ? getBasinTerm(c, vs30, z2p5) : 0.0;
-
- // System.out.printf("calcMean fMag = %.8f\n", fMag);
- // System.out.printf("calcMean fGeom = %.8f\n", fGeom);
- // System.out.printf("calcMean fDepth = %.8f\n", fDepth);
- // System.out.printf("calcMean fAtten = %.8f\n", fAtten);
- // System.out.printf("calcMean fSite = %.8f\n", fSite);
- // System.out.printf("calcMean fBasin = %.8f\n", fBasin);
+ double fMag = magnitudeTerm(c, Mw);
+ double fGeom = geometricAttenuationTerm(c, Mw, rRup);
+ double fDepth = depthTerm(c, zTor);
+ double fAtten = arcCrossingTerm(c, rRup);
+ double fSite = siteAmpTerm(c, vs30, pgaRef);
+ double fBasin = basin ? basinTerm(c, vs30, z2p5) : 0.0;
return c.θ1 + fMag + fGeom + fDepth + fAtten + fSite + fBasin;
}
- /*
- * Eq. 4.3 used for mangitude and source-depth terms
- */
- private static double loghinge(double x, double x0, double a, double b0, double b1,
+ /* Eq. 4.3 used for Mw and zTor */
+ private static double loghinge(
+ double x,
+ double x0,
+ double a,
+ double b0,
+ double b1,
double δ) {
+
return a + b0 * (x - x0) + (b1 - b0) * δ * log(1 + exp((x - x0) / δ));
}
- /*
- * Eq. 4.4 - magnitude term
- */
- private static double getMagTerm(
+ /* Eq. 4.4 - magnitude term */
+ private static double magnitudeTerm(
Coefficients c,
double Mw) {
- return loghinge(Mw, c.Mbreak, c.θ4 * (c.Mbreak - MREF), c.θ4, c.θ5, δM);
+
+ return loghinge(Mw, c.mBreak, c.θ4 * (c.mBreak - MREF), c.θ4, c.θ5, δM);
}
/*
- * Eq. 4.5 - geometrical attenuation term Eq. 4.6 - finite fault term
- * ("pow(10,...)")
+ * Eq. 4.5 - geometrical attenuation term
+ *
+ * Eq. 4.6 - finite fault term ("pow(10,...)")
*/
- private static double getGeomAttenTerm(
+ private static double geometricAttenuationTerm(
Coefficients c,
double Mw,
double rRup) {
- return (c.θ2 + c.θ3 * Mw) * log(rRup + pow(10, c.θnft1 + c.θnft2 * (Mw - 6.0)));
+
+ double ffExp = c.θnft1 + c.θnft2 * (Mw - 6.0);
+ return (c.θ2 + c.θ3 * Mw) * log(rRup + pow(10, ffExp));
}
- /**
- * Eq. 4.7 - source depth term
- *
- * <pre>
- * fDepth(interface|slab) = loghinge(zTor, Zbreak + δZbreak, θ9 * (Zbreak + δZbreak - Zref), θ9, θ10, δZ)
- * </pre>
- *
- * "The source-depth breakpoints [Zbreak + δZbreak] were initially set to
- * Zbreak,if = 30 km and Zbreak,slab = 80 km. The period-dependent adjustment
- * coefficients δZbreak were determined by regression because the depth
- * breakpoints have no theoretical basis... Therefore, the effective source
- * depth breakpoints become Zbreak + δZbreak..."
- *
- * "The parameter δZ determines the smoothness of transition in the bilinear
- * form and was set to δZ = 1."
- *
- * "The coefficients for the source-depth scaling rates (DSR) up to the depth
- * breakpoints are θ9,if and θ9,slab. Both of the coefficients for the DSRs
- * above the breakpoints were fixed at θ10 = 0."
- *
- *
- */
- private static double getDepthTerm(
+ /* Eq. 4.7 - source depth term */
+ private static double depthTerm(
Coefficients c,
double zTor) {
- return loghinge(zTor, c.Zbreak, c.θ9 * (c.Zbreak - c.Zref), c.θ9, c.θ10, δZ);
+
+ return loghinge(zTor, c.zBreak, c.θ9 * (c.zBreak - c.zRef), c.θ9, c.θ10, δZ);
}
- /*
- * Eq. 4.8 - site amplification term
- */
- private static double getSiteAmpTerm(Coefficients c, double vs30, double pgaRef) {
+ /* Eq. 4.8 - site amplification term */
+ private static double siteAmpTerm(Coefficients c, double vs30, double pgaRef) {
double siteTerm;
if (vs30 <= c.k1) {
siteTerm = c.θ7 * log(vs30 / c.k1) +
@@ -515,51 +356,30 @@ public abstract class KuehnEtAl_2020 implements GroundMotionModel {
/*
* Eq. 4.9 - anelastic attenuation term (arc-crossing term) -- Region specific
- * (regionalization ignored for now)
- *
- * Alaska and Cascadia do not have different anelastic attenuation in forearc
- * and backarc
- *
- * Global does have some backarc travel path effects, but these are not yet
- * implemented and we will need to track Ri, the travel distance through
- * subregion i (of the arc) for a given geographic region.
+ * (regionalization ignored for now). This is accurate/applicable for ALASKA,
+ * CASCADIA, NZ, and TW.
*/
- private static double getArcCrossingTerm(Coefficients c, double rRup) {
- // This is accurate for ALASKA, CASCADIA, NZ, and TW
+ private static double arcCrossingTerm(Coefficients c, double rRup) {
return c.θ6_2 * rRup;
}
- /*
- * Eq. 4.10 - basin depth term
- *
- * Eq. 4.11 - lnZref(vs30)
- *
- * Uses EITHER Z1.0 or Z2.5: Cascadia --> Z2.5, Alaska???
- *
- * TODO: This function follows the paper, it's not clear if it is consistent
- * with the official python.
- */
- private static double getBasinTerm(
+ /* Eq. 4.10 - z2.5 basin depth term; Cascadia model */
+ private static double basinTerm(
Coefficients c,
double vs30,
double z2p5) {
+ // TODO add Seattle flag
+
if (Double.isNaN(z2p5)) {
return 0.0;
}
double z2p5m = z2p5 * 1000.0;
- double δlnZ = log(z2p5m) - getBasinZref(vs30); // getBasinZref(region,
- // vs30);
- /*
- * The Seattle basin, as defined, is relatively limited in range and
- * encloses sites that are all very deep, hence the fixed Seattle basin
- * scale factor. Here, we opt for the general PNW basin (something seems
- * wrong with using the seattle coefficient as a minimum value). (Eqn. 4.14)
- */
+ double δlnZ = log(z2p5m) - basinZref(vs30);
double fb = c.θ11 + c.θ12 * δlnZ;
- return fb; // Math.min(fb, c.θ11S); disabled seattle min cap
+ return Math.min(fb, c.θ11S); // Seattle cap
// switch (region) {
// case Region.CASCADIA:
@@ -584,51 +404,14 @@ public abstract class KuehnEtAl_2020 implements GroundMotionModel {
// }
}
- /*
- * Eq. 4.11
- */
- // private static double getBasinZref(String region, double vs30) {
- private static double getBasinZref(double vs30) {
- /**
- * Eq. 4.11 - lnZref has regional coefs θz1 - θz4 for CASC, JP, NZ, and TW.
- * But then the text says "For CASC, estimating all parameters of Equation
- * 4.11 is not possible due to the very limited number of stations with
- * associated Vs30 and Z2.5 values at hard rock sites in particular. Hence,
- * in this case the model is constrained to approach a value of Zref = 10 m
- * for large values of Vs30 (Vs30 > 2000 m/s)."
- *
- * - Max vs30 for the model is supposed to be 1500 m/s
- *
- * - What coefs θz1 - θz4 for all other regions?
- *
- * - Is Basin-Depth term only applied to CASC, NZ, TW, JP regions?
- *
- * Hard-coded coeffs in Python code, explicitly for calculating delta_ln_z
- * for Cascadia (all basins), for non-PGA-ref calcs, and the first and last
- * values aren't used:
- *
- * <pre>
- * pars_z_casc = [-999, 8.29404964010203, 2.30258509299405, 6.39692965521615, 0.27081459, 1.7381352625]
- * </pre>
- */
- // if (region == Region.CASCADIA) {
+ /* Eq. 4.11 - basin reference depth for vs30; Cascadia model */
+ private static double basinZref(double vs30) {
double θz1 = 8.29404964010203;
double θz2 = 2.30258509299405;
double θz3 = 6.39692965521615;
double θz4 = 0.27081459;
double vsterm = exp((log(vs30) - θz3) / θz4);
return θz1 + (θz2 - θz1) * vsterm / (1.0 + vsterm);
- // } else {
- // // Python code also has hard-coded coeffs for JA, NZ, and TW
- // return 0.0;
- // }
- }
-
- /*
- * Calculate Aleatory Sigma
- */
- private static double getAleatorySigma(Coefficients c) {
- return Maths.hypot(c.φ, c.τ);
}
private static abstract class Global extends KuehnEtAl_2020 {
--
GitLab