diff --git a/src/main/java/gov/usgs/earthquake/nshmp/gmm/Gmm.java b/src/main/java/gov/usgs/earthquake/nshmp/gmm/Gmm.java
index 74fc8ae09adb1f970f33f3096229e55bbfde6cbc..21abfa7b6d3912e650ddec51efd8e9fb7c193a5c 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/gmm/Gmm.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/gmm/Gmm.java
@@ -288,6 +288,90 @@ public enum Gmm {
NgaSubductionUsgs_2018.COEFFS,
NgaSubductionUsgs_2018.CONSTRAINTS),
+ /** @see ParkerEtAl_2020 NGA-Subduction */
+ PSHAB20_INTERFACE(
+ ParkerEtAl_2020.Interface.class,
+ ParkerEtAl_2020.Interface.NAME,
+ ParkerEtAl_2020.COEFFS_INTERFACE,
+ ParkerEtAl_2020.CONSTRAINTS_INTERFACE),
+
+ /** @see ParkerEtAl_2020 NGA-Subduction */
+ PSHAB20_SLAB(
+ ParkerEtAl_2020.Slab.class,
+ ParkerEtAl_2020.Slab.NAME,
+ ParkerEtAl_2020.COEFFS_SLAB,
+ ParkerEtAl_2020.CONSTRAINTS_SLAB),
+
+ /** @see ParkerEtAl_2020 NGA-Subduction */
+ PSHAB20_CASCADIA_INTERFACE(
+ ParkerEtAl_2020.CascadiaInterface.class,
+ ParkerEtAl_2020.CascadiaInterface.NAME,
+ ParkerEtAl_2020.COEFFS_INTERFACE,
+ ParkerEtAl_2020.CONSTRAINTS_INTERFACE),
+
+ /** @see ParkerEtAl_2020 NGA-Subduction */
+ PSHAB20_CASCADIA_SLAB(
+ ParkerEtAl_2020.CascadiaSlab.class,
+ ParkerEtAl_2020.CascadiaSlab.NAME,
+ ParkerEtAl_2020.COEFFS_SLAB,
+ ParkerEtAl_2020.CONSTRAINTS_SLAB),
+
+ /** @see ParkerEtAl_2020 NGA-Subduction */
+ PSHAB20_CASCADIA_SEATTLE_INTERFACE(
+ ParkerEtAl_2020.CascadiaSeattleInterface.class,
+ ParkerEtAl_2020.CascadiaSeattleInterface.NAME,
+ ParkerEtAl_2020.COEFFS_INTERFACE,
+ ParkerEtAl_2020.CONSTRAINTS_INTERFACE),
+
+ /** @see ParkerEtAl_2020 NGA-Subduction */
+ PSHAB20_CASCADIA_SEATTLE_SLAB(
+ ParkerEtAl_2020.CascadiaSeattleSlab.class,
+ ParkerEtAl_2020.CascadiaSeattleSlab.NAME,
+ ParkerEtAl_2020.COEFFS_SLAB,
+ ParkerEtAl_2020.CONSTRAINTS_SLAB),
+
+ /** @see ParkerEtAl_2020 NGA-Subduction */
+ PSHAB20_CASCADIA_PACNW_INTERFACE(
+ ParkerEtAl_2020.CascadiaPacNWInterface.class,
+ ParkerEtAl_2020.CascadiaPacNWInterface.NAME,
+ ParkerEtAl_2020.COEFFS_INTERFACE,
+ ParkerEtAl_2020.CONSTRAINTS_INTERFACE),
+
+ /** @see ParkerEtAl_2020 NGA-Subduction */
+ PSHAB20_CASCADIA_PACNW_SLAB(
+ ParkerEtAl_2020.CascadiaPacNWSlab.class,
+ ParkerEtAl_2020.CascadiaPacNWSlab.NAME,
+ ParkerEtAl_2020.COEFFS_SLAB,
+ ParkerEtAl_2020.CONSTRAINTS_SLAB),
+
+ /** @see ParkerEtAl_2020 NGA-Subduction */
+ PSHAB20_ALASKA_INTERFACE(
+ ParkerEtAl_2020.AlaskaInterface.class,
+ ParkerEtAl_2020.AlaskaInterface.NAME,
+ ParkerEtAl_2020.COEFFS_INTERFACE,
+ ParkerEtAl_2020.CONSTRAINTS_INTERFACE),
+
+ /** @see ParkerEtAl_2020 NGA-Subduction */
+ PSHAB20_ALASKA_SLAB(
+ ParkerEtAl_2020.AlaskaSlab.class,
+ ParkerEtAl_2020.AlaskaSlab.NAME,
+ ParkerEtAl_2020.COEFFS_SLAB,
+ ParkerEtAl_2020.CONSTRAINTS_SLAB),
+
+ /** @see ParkerEtAl_2020 NGA-Subduction */
+ PSHAB20_ALEUTIANS_INTERFACE(
+ ParkerEtAl_2020.AleutiansInterface.class,
+ ParkerEtAl_2020.AleutiansInterface.NAME,
+ ParkerEtAl_2020.COEFFS_INTERFACE,
+ ParkerEtAl_2020.CONSTRAINTS_INTERFACE),
+
+ /** @see ParkerEtAl_2020 NGA-Subduction */
+ PSHAB20_ALEUTIANS_SLAB(
+ ParkerEtAl_2020.AleutiansSlab.class,
+ ParkerEtAl_2020.AleutiansSlab.NAME,
+ ParkerEtAl_2020.COEFFS_SLAB,
+ ParkerEtAl_2020.CONSTRAINTS_SLAB),
+
/** @see BcHydro_2012 */
BCHYDRO_12_INTERFACE(
BcHydro_2012.Interface.class,
diff --git a/src/main/java/gov/usgs/earthquake/nshmp/gmm/ParkerEtAl_2020.java b/src/main/java/gov/usgs/earthquake/nshmp/gmm/ParkerEtAl_2020.java
new file mode 100644
index 0000000000000000000000000000000000000000..d715e776caf772c13eea392e8b7f214d0ca44133
--- /dev/null
+++ b/src/main/java/gov/usgs/earthquake/nshmp/gmm/ParkerEtAl_2020.java
@@ -0,0 +1,944 @@
+package gov.usgs.earthquake.nshmp.gmm;
+
+import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.MW;
+import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.RRUP;
+import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.VS30;
+import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.ZHYP;
+import static gov.usgs.earthquake.nshmp.gmm.Imt.PGA;
+import static java.lang.Math.exp;
+import static java.lang.Math.log;
+import static java.lang.Math.log10;
+import static java.lang.Math.max;
+import static java.lang.Math.min;
+import static java.lang.Math.pow;
+import static java.lang.Math.sqrt;
+
+import java.util.Map;
+
+import com.google.common.collect.Range;
+
+import gov.usgs.earthquake.nshmp.Maths;
+import gov.usgs.earthquake.nshmp.gmm.GmmInput.Constraints;
+
+/**
+ * Implementation of the PEER NGA-Subduction GMM, Parker et al., 2020 GMM
+ * (PSHAB20)
+ *
+ * Notes:
+ *
+ * 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.
+ *
+ * PSHAB20 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.
+ *
+ * Path term (Eqs. 4.2 and 4.3 in PEER report) appears to differ in the official
+ * Python implementation. Compare implementations of getPathTerm() below
+ *
+ * <p><b>Reference:</b> Parker, G.A, Stewart, J.P., Boore, D.M., Atkinson, G.M.,
+ * Hassani, B. (2020) NGA-Subduction global ground-motion models with regional
+ * adjustment factors. Report no. 2020/03. Berkeley, CA: PEER, 131 p.
+ *
+ * <p><b>Component:</b> TODO
+ *
+ * @author U.S. Geological Survey
+ */
+
+public abstract class ParkerEtAl_2020 implements GroundMotionModel {
+
+ static final String NAME = "Parker et al. (2020) NGA-Subduction";
+
+ static final Constraints CONSTRAINTS_INTERFACE = Constraints.builder()
+ .set(MW, Range.closed(4.5, 9.5))
+ .set(RRUP, Range.closed(20.0, 1000.0))
+ .set(ZHYP, Range.closed(0.0, 40.0))
+ .set(VS30, Range.closed(150.0, 2000.0))
+ .build();
+ static final Constraints CONSTRAINTS_SLAB = Constraints.builder()
+ .set(MW, Range.closed(4.5, 8.5))
+ .set(RRUP, Range.closed(35.0, 1000.0))
+ .set(ZHYP, Range.closed(0.0, 200.0))
+ .set(VS30, Range.closed(150.0, 2000.0))
+ .build();
+
+ /*
+ * Regions: global, Alaska, Cascadia, [CAM (Central America), Japan, SA (South
+ * America), Taiwan]
+ *
+ * Saturation Regions: global, Aleutian, Alaska, Cascadia, [Central America S,
+ * Central America N, Japan Pacific, Japan Phillipines, South America N, South
+ * America S, Taiwan W, Taiwan E]
+ */
+ static final class Region {
+ static final String GLOBAL = "Global";
+ static final String ALASKA = "Alaska";
+ static final String ALEUTIANS = "Aleutians";
+ static final String CASCADIA = "Cascadia";
+ }
+
+ static final class Basin {
+ static final String NONE = "None";
+ static final String PACNW_NONE = "PacNW-no-basin";
+ static final String PACNW = "PacNW";
+ static final String SEATTLE = "Seattle";
+ }
+
+ // aleatory uncertainty coefficients are duplicated in interface and interface
+ // CSV files (as are some other coefficients that are also independent of
+ // event type)
+ private static final String INTERFACE_CSV = "PSHAB20_interface.csv";
+ private static final String SLAB_CSV = "PSHAB20_slab.csv";
+ // 4 epistemic coefficients tied to region and interface/slab
+ private static final String EPISTEMIC_INTERFACE_CSV = "PSHAB20_epistemic_interface.csv";
+ private static final String EPISTEMIC_SLAB_CSV = "PSHAB20_epistemic_slab.csv";
+
+ static final CoefficientContainer COEFFS_INTERFACE;
+ static final CoefficientContainer COEFFS_SLAB;
+
+ static {
+ COEFFS_INTERFACE = new CoefficientContainer(INTERFACE_CSV);
+ COEFFS_SLAB = new CoefficientContainer(SLAB_CSV);
+ }
+
+ // private static final double B4 = 0.1;
+ // private static final double V1 = 270.0; // m/s
+ // private static final double VREF = 760.0; // m/s
+ private static final double DB1 = 20.0; // km
+ private static final double DB2 = 67.0; // km
+
+ /* Ï• parameters */
+ private static final double R1 = 200; // km
+ private static final double R2 = 500; // km
+ private static final double V1 = 200; // m/s
+ private static final double V2 = 500; // m/s
+ /* Ï•S2S and Ï•SS parameters */
+ private static final double R3 = 200; // km
+ private static final double R4 = 500; // km
+ private static final double R5 = 500; // km
+ private static final double R6 = 800; // km
+ private static final double V3 = 200; // m/s
+ private static final double V4 = 800; // m/s
+
+ private static final class Coefficients {
+ final Imt imt;
+ final double c0; // model constant
+ final double c1, b4, a0; // path scaling
+ final double mc; // saturation corner magnitude
+ final double c4, c5, c6; // magnitude scaling
+ final double d, m, db; // source depth scaling (slab only)
+ final double v1, v2, vRef, s1, s2; // linear site amplification
+ final double f4, f5; // non-linear site amplification
+ final double c_e1, c_e2, c_e3, del_None, del_Seattle; // Cascadia basin
+ final double Ï„, Ï•21, Ï•22, Ï•2v, vm, Ï•2s2s0, a1, Ï•2ss1, Ï•2ss2, a2; // alleatory
+ // uncertainty
+ final double σε1, σε2, t1, t2;
+
+ Coefficients(Imt imt, CoefficientContainer cc, boolean isSlab, String region) {
+ this.imt = imt;
+ Map<String, Double> coeffs = cc.get(imt);
+ // Model constant
+ c0 = coeffs.get(region + "_c0");
+ // Path
+ c1 = coeffs.get("c1");
+ b4 = coeffs.get("b4"); // period-independent
+ a0 = coeffs.get(getA0Col(isSlab, region));
+ // a0 = region.equals(Region.ALASKA) || (region.equals(Region.CASCADIA) &&
+ // isSlab)
+ // ? coeffs.get(region + "_a0") : coeffs.get(Region.GLOBAL + "_a0");
+ // Magnitude saturation corner magnitude
+ mc = getSaturationMagnitude(region, isSlab);
+ // Magnitude-scaling
+ c4 = coeffs.get("c4");
+ c5 = coeffs.get("c5");
+ c6 = coeffs.get("c6");
+ // Source depth scaling (slab only)
+ d = isSlab ? coeffs.get("d") : 0.0;
+ m = isSlab ? coeffs.get("m") : 0.0;
+ db = isSlab ? coeffs.get("db") : 0.0;
+ // Linear site amplification
+ v1 = coeffs.get("V1");
+ v2 = coeffs.get("V2");;
+ vRef = coeffs.get("Vref");
+ s2 = (region.equals(Region.ALEUTIANS)) ? coeffs.get(Region.GLOBAL + "_s2")
+ : coeffs.get(region + "_s2");
+ s1 = s2; // s1 = s2 except for Japan and Taiwan
+ // Non-linear site amplification
+ f4 = coeffs.get("f4");
+ f5 = coeffs.get("f5");
+ // Cascadia basin
+ c_e1 = coeffs.get("C_e1");
+ c_e2 = coeffs.get("C_e2");
+ c_e3 = coeffs.get("C_e3");
+ del_None = coeffs.get("del_None");
+ del_Seattle = coeffs.get("del_Seattle");
+ // aleatory model
+ Ï„ = coeffs.get("Tau");
+ Ï•21 = coeffs.get("phi21");
+ Ï•22 = coeffs.get("phi22");
+ Ï•2v = coeffs.get("phi2V");
+ vm = coeffs.get("VM");
+ Ï•2s2s0 = coeffs.get("phi2S2S0");
+ a1 = coeffs.get("a1");
+ Ï•2ss1 = coeffs.get("phi2SS1");
+ Ï•2ss2 = coeffs.get("phi2SS2");
+ a2 = coeffs.get("a2");
+ // epistemic model (period-independent)
+ σε1 = getσε1(region, isSlab);
+ σε2 = getσε2(region, isSlab);
+ t1 = getT1(region, isSlab);
+ t2 = getT2(region, isSlab);
+ }
+
+ private static final String getA0Col(boolean isSlab, String region) {
+ String a0Col;
+ if (region.equals(Region.ALASKA)) {
+ a0Col = Region.ALASKA + "_a0";
+ } else if (region.equals(Region.CASCADIA) && isSlab) {
+ a0Col = Region.CASCADIA + "_a0";
+ } else {
+ a0Col = Region.GLOBAL + "_a0";
+ }
+ return a0Col;
+ }
+
+ /* Table E4 */
+ private static final double getσε1(String region, boolean isSlab) {
+ // return region.equals(Region.ALASKA) ? 0.15 :
+ // region.equals(Region.ALEUTIANS) ? 0.15 :
+ // region.equals(Region.CASCADIA) ? 0.43 : 0.4; // else GLOBAL
+ switch (region) {
+ case Region.GLOBAL:
+ return isSlab ? 0.35 : 0.4;
+ case Region.ALASKA:
+ return isSlab ? 0.15 : .15;
+ case Region.ALEUTIANS:
+ return isSlab ? 0.25 : .15;
+ case Region.CASCADIA:
+ return isSlab ? 0.35 : 0.43;
+ default:
+ throw new IllegalArgumentException("Region [" + region + "] not supported");
+ }
+ }
+
+ /* Table E4 */
+ private static final double getσε2(String region, boolean isSlab) {
+ // return region.equals(Region.ALASKA) ? 0.1 :
+ // region.equals(Region.ALEUTIANS) ? 0.11 :
+ // region.equals(Region.CASCADIA) ? 0.33 : 0.4; // else GLOBAL
+ switch (region) {
+ case Region.GLOBAL:
+ return isSlab ? 0.22 : 0.4;
+ case Region.ALASKA:
+ return isSlab ? 0.12 : .1;
+ case Region.ALEUTIANS:
+ return isSlab ? 0.18 : .1;
+ case Region.CASCADIA:
+ return isSlab ? 0.16 : 0.33;
+ default:
+ throw new IllegalArgumentException("Region [" + region + "] not supported");
+ }
+ }
+
+ /* Table E4 */
+ private static final double getT1(String region, boolean isSlab) {
+ // return region.equals(Region.ALASKA) ? 1 :
+ // region.equals(Region.ALEUTIANS) ? 1 :
+ // region.equals(Region.CASCADIA) ? 0.2 : 0.2; // else GLOBAL
+ switch (region) {
+ case Region.GLOBAL:
+ return isSlab ? 0.15 : 0.2;
+ case Region.ALASKA:
+ return isSlab ? 0.5 : 1;
+ case Region.ALEUTIANS:
+ return isSlab ? 0.3 : 1;
+ case Region.CASCADIA:
+ return isSlab ? 0.2 : 0.2;
+ default:
+ throw new IllegalArgumentException("Region [" + region + "] not supported");
+ }
+ }
+
+ /* Table E4 */
+ private static final double getT2(String region, boolean isSlab) {
+ switch (region) {
+ case Region.GLOBAL:
+ return isSlab ? 2 : 0.4;
+ case Region.ALASKA:
+ return isSlab ? 1 : 4;
+ case Region.ALEUTIANS:
+ return isSlab ? 0.8 : 4;
+ case Region.CASCADIA:
+ return isSlab ? 3 : 0.5;
+ default:
+ throw new IllegalArgumentException("Region [" + region + "] not supported");
+ }
+ }
+
+ private static final double getSaturationMagnitude(String region, boolean isSlab) {
+ // TODO? Read from coeff table?
+ // Parket et al. (2020) Table 4.1
+ // Regional saturation magnitudes for interface events computed using
+ // seismogenic fault width [Campbell 2020] and for intraslab events
+ // computed
+ // using slab thickness [Ji and Archuleta 2018].
+ switch (region) {
+ case Region.GLOBAL:
+ return isSlab ? 7.6 : 7.9;
+ case Region.ALASKA:
+ return isSlab ? 7.20 : 8.6;
+ case Region.ALEUTIANS:
+ return isSlab ? 7.98 : 8.0;
+ case Region.CASCADIA:
+ return isSlab ? 7.20 : 7.7;
+ default:
+ throw new IllegalArgumentException("Region [" + region + "] not supported");
+ }
+ }
+ }
+
+ private final Coefficients coeffs;
+ private final Coefficients coeffsPGA;
+
+ ParkerEtAl_2020(final Imt imt) {
+ coeffs = initCoeffs(imt, isSlab(), getRegion());
+ coeffsPGA = initCoeffs(PGA, isSlab(), getRegion());
+ }
+
+ private static Coefficients initCoeffs(final Imt imt, final boolean isSlab, final String region) {
+ CoefficientContainer coeffs = isSlab ? COEFFS_SLAB : COEFFS_INTERFACE;
+ return new Coefficients(imt, coeffs, isSlab, region);
+ }
+
+ @Override
+ public final ScalarGroundMotion calc(final GmmInput in) {
+ /*
+ * Equation 4.1 μ_lny = c0 + Fp + Fm + Fd + Fs
+ */
+
+ // double Mc = getSaturationMagnitude(getRegion(), isSlab());
+ double pgaRef = exp(calcMean(coeffsPGA, isSlab(), in.Mw, in.rRup, in.zHyp));
+
+ double μ = calcMean(coeffs, isSlab(), getBasin(), in.Mw, in.rRup, in.zHyp, pgaRef,
+ in.vs30, in.z2p5);
+
+ double ϕTotal = getϕTotal(coeffs, in.rRup, in.vs30);
+ double σ = Maths.hypot(coeffs.τ, ϕTotal);
+
+ // System.out.println("calc for T = " + coeffs.imt.toString());
+ // System.out.printf(" Ï•Total = %.6f\n", Ï•Total);
+ // System.out.printf(" σ = %.6f\n", σ);
+
+ double ϕS2S = getϕS2S(coeffs, in.rRup, in.vs30);
+ double ϕSS = getϕSS(coeffs, in.rRup, in.vs30);
+
+ double σε = getσε(coeffs);
+
+ return DefaultScalarGroundMotion.create(μ, σ);
+ }
+
+ abstract boolean isSlab();
+
+ abstract String getRegion();
+
+ abstract String getBasin();
+
+ // abstract double σε1();
+ //
+ // abstract double σε2();
+ //
+ // abstract double t1();
+ //
+ // abstract double t2();
+
+ /* calc median at 760 reference conditions */
+ private static final double calcMean(
+ final Coefficients c,
+ final boolean isSlab,
+ final double Mw,
+ final double rRup,
+ final double zHyp) {
+ double h = nearSourceSaturation(Mw, c.mc, isSlab);
+ double Fp = getPathTermPy(c, rRup, Mw, h);
+ double Fm = getMagnitudeTerm(c, Mw);
+ double Fd = getDepthTerm(c, isSlab, zHyp);
+
+ // System.out.println("calcMean Reference for T = " + c.imt.toString());
+ // System.out.printf(" mc = %.2f\n", c.mc);
+ // System.out.printf(" h = %.6f\n", h);
+ // System.out.printf(" Fp = %.6f\n", Fp);
+ // System.out.printf(" Fm(ref) = %.6f\n", Fm);
+ // System.out.printf(" Fd = %.6f\n", Fd);
+ // System.out.printf(" Fb = %.6f (reference calc)\n", 0.0);
+ // System.out.printf(" mu = %.6f\n", c.c0 + Fp + Fm + Fd);
+
+ return c.c0 + Fp + Fm + Fd;
+ }
+
+ /* Equation 4.1 */
+ private static final double calcMean(
+ final Coefficients c,
+ final boolean isSlab,
+ final String basin,
+ final double Mw,
+ final double rRup,
+ final double zHyp,
+ final double pgaRef,
+ final double vs30,
+ final double z2p5) {
+
+ double h = nearSourceSaturation(Mw, c.mc, isSlab);
+ double Fp = getPathTermPy(c, rRup, Mw, h);
+ double Fm = getMagnitudeTerm(c, Mw);
+ double Fd = getDepthTerm(c, isSlab, zHyp);
+
+ /* Equation 4.7: Fs = Flin + Fnl + Fb */
+ double Fslin = getSiteTermLinear(c, vs30);
+ double Fsnl = getSiteTermNonLinear(c, pgaRef, vs30);
+ double Fsb = getSiteBasinTerm(c, basin, vs30, z2p5);
+ double Fs = Fslin + Fsnl + Fsb;
+
+ // System.out.println("calcMean for T = " + c.imt.toString());
+ // System.out.printf(" Vs30 = %.2f\n", vs30);
+ // System.out.printf(" h = %.6f\n", h);
+ // System.out.printf(" Fp = %.6f\n", Fp);
+ // System.out.printf(" Fm = %.6f\n", Fm);
+ // System.out.printf(" Fd = %.6f\n", Fd);
+ // System.out.printf(" PGAref = %.6f\n", pgaRef);
+ // System.out.printf(" Fslin = %.6e\n", Fslin);
+ // System.out.printf(" Fsnl = %.6e\n", Fsnl);
+ // System.out.printf(" Fssb = %.6e\n", Fsb);
+ // System.out.printf(" Fs = %.6e\n", Fs);
+ // System.out.printf(" mu = %.6f\n", c.c0 + Fp + Fm + Fd + Fs);
+ // System.out.printf(" y = %.6f\n", exp(c.c0 + Fp + Fm + Fd + Fs));
+
+ return c.c0 + Fp + Fm + Fd + Fs;
+ }
+
+ /* Equations 4.2 and 4.3 */
+ private static final double getPathTerm(final Coefficients c, double rRup, double Mw, double h) {
+ double R = Maths.hypot(rRup, h);
+ return c.c1 * log(R) + c.b4 * Mw * log(rRup / R) + c.a0 * R;
+ }
+
+ /* Equations 4.2 and 4.3, following official Python implementation */
+ private static final double getPathTermPy(final Coefficients c, double rRup, double Mw,
+ double h) {
+ double Rref = Maths.hypot(1, h);
+ double R = Maths.hypot(rRup, h);
+ double R_Rref = R / Rref;
+ // System.out.printf(" c1 = %.6e\n", c.c1);
+ // System.out.printf(" b4 = %.6e\n", c.b4);
+ // System.out.printf(" a0 = %.6e\n", c.a0);
+ // System.out.printf(" Rref = %.6e\n", R_Rref);
+ // System.out.printf(" R = %.6e\n", R);
+ // System.out.printf(" Mw = %.6e\n", Mw);
+
+ return c.c1 * log(R) + c.b4 * Mw * log(R_Rref) + c.a0 * R;
+ }
+
+ private static final double nearSourceSaturation(double Mw, double mc, boolean slab) {
+ if (slab) {
+ /* Equation 4.4b */
+ return (Mw <= mc) ? pow(10, (1.050 / (mc - 4.0)) * (Mw - mc) + 1.544) : 35;
+ } else {
+ /* Equation 4.4a */
+ return pow(10.0, -0.82 + 0.252 * Mw);
+ }
+ }
+
+ /* Equation 4.5 */
+ private static final double getMagnitudeTerm(final Coefficients c, double Mw) {
+ double dM = Mw - c.mc;
+ return (Mw <= c.mc) ? c.c4 * dM + c.c5 * dM * dM : c.c6 * dM;
+ }
+
+ /* Equation 4.6 */
+ private static final double getDepthTerm(final Coefficients c, boolean slab, double zHyp) {
+ if (!slab) {
+ return 0.0;
+ }
+ // Zhyp (see section 4.3.3) can be replaced with mean depth (Zhyp_bar)
+ // that depends on Ztor, fault width, and dip angle
+ double dTerm;
+ if (zHyp < DB1) { // "less than"
+ dTerm = c.m * (DB1 - DB2) + c.d;
+ } else if (zHyp > DB2) {
+ dTerm = c.d;
+ } else { // "db1 <= zHyp" possible typo? (written as "db1 < xHyp")
+ dTerm = c.m * (zHyp - DB2) + c.d;
+ }
+ return dTerm;
+ }
+
+ /* Equation 4.8 */
+ private static final double getSiteTermLinear(final Coefficients c, double vs30) {
+ double Flin;
+ if (vs30 <= c.v1) {
+ Flin = c.s1 * log(vs30 / c.v1) + c.s2 * log(c.v1 / c.vRef);
+ } else if (vs30 > c.v2) {
+ Flin = c.s2 * log(c.v2 / c.vRef);
+ } else {
+ Flin = c.s2 * log(vs30 / c.vRef);
+ }
+ return Flin;
+ }
+
+ /* Equation 4.9 */
+ private static final double getSiteTermNonLinear(final Coefficients c, double pgaRef,
+ double vs30) {
+ double f3 = 0.05;
+ double f2 = c.f4 * (exp(c.f5 * (min(vs30, 760) - 200)) - exp(c.f5 * (760 - 200)));
+ return 0.0 + f2 * log((pgaRef + f3) / f3);
+ }
+
+ /* Equations 4.11/5.5 - 4.13/5.7 */
+ private static final double getSiteBasinTerm(final Coefficients c, String basin, double vs30,
+ double z2p5) {
+ if (basin.equals(Basin.NONE) || Double.isNaN(z2p5)) {
+ return 0.0;
+ }
+ // Cascadia coefficients
+ double θ0 = 3.94;
+ double θ1 = -0.42;
+ double νμ = 200.0;
+ double νσ = 0.2;
+
+ double e1 = c.c_e1;
+ double e2;
+ double e3;
+
+ if (basin.equals(Basin.PACNW_NONE)) {
+ e2 = c.c_e2 + c.del_None;
+ e3 = c.c_e3 + c.del_None;
+ } else if (basin.equals(Basin.SEATTLE)) {
+ e2 = c.c_e2 + c.del_Seattle;
+ e3 = c.c_e3 + c.del_Seattle;
+ } else if (basin.equals(Basin.PACNW)) {
+ e2 = c.c_e2;
+ e3 = c.c_e3;
+ } else {
+ throw new IllegalArgumentException("Basin [" + basin + "] not supported");
+ }
+
+ // Eq. 4.13/5.7
+ double μz2p5 = exp(
+ log(10) * θ1 * (1 + Maths.erf((log10(vs30) - log10(νμ)) / νσ / sqrt(2)) + log(10) * θ0));
+ double δz2p5 = log(z2p5) - log(μz2p5);
+
+ double Fb;
+ if (δz2p5 <= (e1 / e3)) {
+ Fb = e1;
+ } else if (δz2p5 >= (e2 / e3)) {
+ Fb = e2;
+ } else {
+ Fb = e3 * δz2p5;
+ }
+ return Fb;
+ }
+
+ /* Eqs. 6.3, 6.4, 6.5 */
+ private static final double getϕTotal(final Coefficients c, final double rRup,
+ final double vs30) {
+
+ /* Eq. 6.5 */
+ double Δvar;
+ double rPrime = max(R1, min(R2, rRup));
+ if (vs30 <= V1) {
+ Δvar = c.ϕ2v * (log(R2 / rPrime) / log(R2 / R1));
+ } else if (vs30 >= V2) {
+ Δvar = 0;
+ } else {
+ Δvar = c.ϕ2v * (log(V2 / vs30) / log(V2 / V1)) * (log(R2 / rPrime) / log(R2 / R1));
+ }
+
+ /* Eq. 6.4 */
+ double Ï•2;
+ if (rRup <= R1) {
+ Ï•2 = c.Ï•21;
+ } else if (rRup >= R2) {
+ Ï•2 = c.Ï•22;
+ } else {
+ Ï•2 = (c.Ï•22 - c.Ï•21) / log(R2 / R1) * log(rRup / R1) + c.Ï•21;
+ }
+
+ /* Eq. 6.3 */
+ double ϕTotal = sqrt(Δvar + ϕ2);
+
+ // System.out.println("getϕTotal for T = " + c.imt.toString());
+ // System.out.printf(" rRup = %.6f\n",rRup);
+ // System.out.printf(" vs30 = %.6f\n",vs30);
+ // System.out.printf(" rPrime = %.6f\n",rPrime);
+ // System.out.printf(" c.Ï•2v = %.6f\n",c.Ï•2v);
+ // System.out.printf(" c.Ï•21 = %.6f\n",c.Ï•21);
+ // System.out.printf(" c.Ï•22 = %.6f\n",c.Ï•22);
+ // System.out.printf(" Δvar = %.6f\n",Δvar);
+ // System.out.printf(" Ï•2 = %.6f\n",Ï•2);
+ // System.out.printf(" Ï•Total = %.6f\n",Ï•Total);
+
+ return Ï•Total;
+ }
+
+ /* Eqs. 6.6, 6.7 */
+ private static final double getϕS2S(final Coefficients c, final double rRup, final double vs30) {
+
+ /* Eq. 6.7 */
+ double rPrime = max(R3, min(R4, rRup));
+ double ΔvarS2S;
+ if (vs30 <= V3) {
+ ΔvarS2S = c.a1 * log(V3 / c.vm) * log(R4 / rPrime) / log(R4 / R3);
+ } else if (vs30 < c.vm) {
+ ΔvarS2S = c.a1 * log(vs30 / c.vm) * log(R4 / rPrime) / log(R4 / R3);
+ } else if (vs30 < V4) {
+ ΔvarS2S = c.a1 * log(vs30 / c.vm);
+ } else { // vs30 >= V4
+ ΔvarS2S = c.a1 * log(V4 / c.vm);
+ }
+
+ double ϕS2S = sqrt(c.ϕ2s2s0 + ΔvarS2S);
+
+ // System.out.println("getϕS2S for T = " + c.imt.toString());
+ // System.out.printf(" rRup = %.6f\n",rRup);
+ // System.out.printf(" vs30 = %.6f\n",vs30);
+ // System.out.printf(" rPrime = %.6f\n",rPrime);
+ // System.out.printf(" c.a1 = %.6f\n",c.a1);
+ // System.out.printf(" c.vm = %.6f\n",c.vm);
+ // System.out.printf(" ΔvarS2S = %.6f\n",ΔvarS2S);
+ // System.out.printf(" c.Ï•2s2s0 = %.6f\n",c.Ï•2s2s0);
+ // System.out.printf(" Ï•S2S = %.6f\n",Ï•S2S);
+
+ return Ï•S2S;
+ }
+
+ /* Eqs. 6.8, 6.9, 6.10 */
+ private double getϕSS(Coefficients c, double rRup, double vs30) {
+ /* Eq. 6.10 */
+ double rPrime = max(R3, min(R4, rRup));
+ double ΔvarSS;
+ if (vs30 <= V3) {
+ ΔvarSS = c.a2 * log(V3 / c.vm) * log(R4 / rPrime) / log(R4 / R3);
+ } else if (vs30 < c.vm) {
+ ΔvarSS = c.a2 * log(vs30 / c.vm) * log(R4 / rPrime) / log(R4 / R3);
+ } else if (vs30 < V4) {
+ ΔvarSS = c.a2 * log(vs30 / c.vm);
+ } else { // vs30 >= V4
+ ΔvarSS = c.a2 * log(V4 / c.vm);
+ }
+
+ /* Eq. 6.9 */
+ double Ï•SS2rRup;
+ if (rRup <= R5) {
+ Ï•SS2rRup = c.Ï•2ss1;
+ } else if (rRup < R6) {
+ Ï•SS2rRup = (c.Ï•2ss2 - c.Ï•2ss1) / log(R6 / R5) * log(rRup / R5) + c.Ï•2ss1;
+ } else { // rRup >= R6
+ Ï•SS2rRup = c.Ï•2ss2;
+ }
+
+ double ϕSS = sqrt(ϕSS2rRup + ΔvarSS);
+
+ // System.out.println("getϕSS for T = " + c.imt.toString());
+ // System.out.printf(" rRup = %.6f\n",rRup);
+ // System.out.printf(" vs30 = %.6f\n",vs30);
+ // System.out.printf(" rPrime = %.6f\n",rPrime);
+ // System.out.printf(" c.a2 = %.6f\n",c.a1);
+ // System.out.printf(" c.vm = %.6f\n",c.vm);
+ // System.out.printf(" c.Ï•2ss1 = %.6f\n",c.Ï•2ss1);
+ // System.out.printf(" c.Ï•2ss2 = %.6f\n",c.Ï•2ss2);
+ // System.out.printf(" ΔvarSS = %.6f\n",ΔvarSS);
+ // System.out.printf(" Ï•SS2rRup = %.6f\n",Ï•SS2rRup);
+ // System.out.printf(" Ï•S2S = %.6f\n",Ï•SS);
+
+ return Ï•SS;
+ }
+
+ /* Eq. 8.2 */
+ private double getσε(Coefficients c) {
+ double σε;
+ if (c.imt == Imt.PGA || c.imt == Imt.PGV || c.imt.period() < c.t1) {
+ σε = c.σε1;
+ } else if (c.imt.period() < c.t2) {
+ σε = c.σε1 - (c.σε1 - c.σε2) * log(c.imt.period() / c.t1) / log(c.t2 / c.t1);
+ } else { // T > T2
+ σε = c.σε2;
+ }
+ // System.out.println("getσε for T = " + c.imt.toString());
+ // System.out.printf(" T1 = %.6f\n",c.t1);
+ // System.out.printf(" T2 = %.6f\n",c.t2);
+ // System.out.printf(" c.σε1 = %.6f\n",c.σε1);
+ // System.out.printf(" c.σε2 = %.6f\n",c.σε2);
+ // System.out.printf(" σε = %.6f\n",σε);
+
+ return σε;
+ }
+
+ /*
+ * Global interface model
+ */
+ static class Interface extends ParkerEtAl_2020 {
+ static final String NAME = ParkerEtAl_2020.NAME + " : Interface";
+
+ Interface(Imt imt) {
+ super(imt);
+ }
+
+ @Override
+ final boolean isSlab() {
+ return false;
+ }
+
+ @Override
+ String getRegion() {
+ return Region.GLOBAL;
+ }
+
+ @Override
+ String getBasin() {
+ return Basin.NONE;
+ }
+
+ // @Override
+ // double σε1() {
+ // return 0;
+ // }
+ //
+ // @Override
+ // double σε2() {
+ // return 0;
+ // }
+ //
+ // @Override
+ // double t1() {
+ // return 0;
+ // }
+ //
+ // @Override
+ // double t2() {
+ // return 0;
+ // }
+
+ }
+
+ /*
+ * Global intraslab model
+ */
+ static class Slab extends ParkerEtAl_2020 {
+ static final String NAME = ParkerEtAl_2020.NAME + " : Slab";
+
+ Slab(Imt imt) {
+ super(imt);
+ }
+
+ @Override
+ final boolean isSlab() {
+ return true;
+ }
+
+ @Override
+ String getRegion() {
+ return Region.GLOBAL;
+ }
+
+ @Override
+ String getBasin() {
+ return Basin.NONE;
+ }
+
+ // @Override
+ // double σε1() {
+ // return 0;
+ // }
+ //
+ // @Override
+ // double σε2() {
+ // return 0;
+ // }
+ //
+ // @Override
+ // double t1() {
+ // return 0;
+ // }
+ //
+ // @Override
+ // double t2() {
+ // return 0;
+ // }
+ }
+
+ /*
+ * Cascadia interface model
+ */
+ static class CascadiaInterface extends Interface {
+ static final String NAME = Interface.NAME + " : Cascadia";
+
+ CascadiaInterface(Imt imt) {
+ super(imt);
+ }
+
+ @Override
+ String getRegion() {
+ return Region.CASCADIA;
+ }
+
+ @Override
+ String getBasin() {
+ return Basin.PACNW_NONE;
+ }
+ }
+
+ /*
+ * Cascadia intraslab model
+ */
+ static class CascadiaSlab extends Slab {
+ static final String NAME = Slab.NAME + " : Cascadia";
+
+ CascadiaSlab(Imt imt) {
+ super(imt);
+ }
+
+ @Override
+ String getRegion() {
+ return Region.CASCADIA;
+ }
+
+ @Override
+ String getBasin() {
+ return Basin.PACNW_NONE;
+ }
+ }
+
+ /*
+ * Cascadia interface model, Seattle basin
+ */
+ static class CascadiaSeattleInterface extends CascadiaInterface {
+ static final String NAME = CascadiaInterface.NAME + " : Seattle";
+
+ CascadiaSeattleInterface(Imt imt) {
+ super(imt);
+ }
+
+ @Override
+ String getBasin() {
+ return Basin.SEATTLE;
+ }
+ }
+
+ /*
+ * Cascadia intraslab model, Seattle basin
+ */
+ static class CascadiaSeattleSlab extends CascadiaSlab {
+ static final String NAME = CascadiaSlab.NAME + " : Seattle";
+
+ CascadiaSeattleSlab(Imt imt) {
+ super(imt);
+ }
+
+ @Override
+ String getBasin() {
+ return Basin.SEATTLE;
+ }
+ }
+
+ /*
+ * Cascadia interface model, Other PacNW basin
+ */
+ static class CascadiaPacNWInterface extends CascadiaInterface {
+ static final String NAME = CascadiaInterface.NAME + " : PacNW";
+
+ CascadiaPacNWInterface(Imt imt) {
+ super(imt);
+ }
+
+ @Override
+ String getBasin() {
+ return Basin.PACNW;
+ }
+ }
+
+ /*
+ * Cascadia intraslab model, Other PacNW basin
+ */
+ static class CascadiaPacNWSlab extends CascadiaSlab {
+ static final String NAME = CascadiaSlab.NAME + " : PacNW";
+
+ CascadiaPacNWSlab(Imt imt) {
+ super(imt);
+ }
+
+ @Override
+ String getBasin() {
+ return Basin.PACNW;
+ }
+ }
+
+ /*
+ * Alaska interface model
+ */
+ static class AlaskaInterface extends Interface {
+ static final String NAME = Interface.NAME + " : Alaska";
+
+ AlaskaInterface(Imt imt) {
+ super(imt);
+ }
+
+ @Override
+ String getRegion() {
+ return Region.ALASKA;
+ }
+ }
+
+ /*
+ * Alaska intraslab model
+ */
+ static class AlaskaSlab extends Slab {
+ static final String NAME = Slab.NAME + " : Alaska";
+
+ AlaskaSlab(Imt imt) {
+ super(imt);
+ }
+
+ @Override
+ String getRegion() {
+ return Region.ALASKA;
+ }
+ }
+
+ /*
+ * Aleutians interface model
+ */
+ static class AleutiansInterface extends Interface {
+ static final String NAME = Interface.NAME + " : Aleutians";
+
+ AleutiansInterface(Imt imt) {
+ super(imt);
+ }
+
+ @Override
+ String getRegion() {
+ return Region.ALEUTIANS;
+ }
+ }
+
+ /*
+ * Aleutians intraslab model
+ */
+ static class AleutiansSlab extends Slab {
+ static final String NAME = Slab.NAME + " : Aleutians";
+
+ AleutiansSlab(Imt imt) {
+ super(imt);
+ }
+
+ @Override
+ String getRegion() {
+ return Region.ALEUTIANS;
+ }
+ }
+
+ public static void main(String[] args) {
+
+ }
+
+}
diff --git a/src/main/resources/gmm/coeffs/PSHAB20_interface.csv b/src/main/resources/gmm/coeffs/PSHAB20_interface.csv
new file mode 100644
index 0000000000000000000000000000000000000000..d85285311fcaa97b1ae4c311b062c1eb8d8b3523
--- /dev/null
+++ b/src/main/resources/gmm/coeffs/PSHAB20_interface.csv
@@ -0,0 +1,27 @@
+Period (s),Global_c0,Alaska_c0,Aleutian_c0,Cascadia_c0,c1,b4,Global_a0,Alaska_a0,c4,c5,c6,V1,V2,Vref,Global_s2,Alaska_s2,Cascadia_s2,f4,f5,C_e1,C_e2,C_e3,del_None,del_Seattle,Tau,phi21,phi22,phi2V,VM,phi2S2S0,a1,phi2SS1,phi2SS2,a2
+PGV,8.097,9.283796298,8.374796298,7.728,-1.661,0.1,-3.95E-03,-4.04E-03,1.336,-0.039,1.336,270,850,760,-0.601,-1.031,-0.671,-0.31763,-0.0052,0,0.115,0.068,-0.115,0,0.477,0.348,0.288,-0.179,423,0.142,0.047,0.153,0.166,0.011
+PGA,4.082,4.458796298,3.652796298,3.856,-1.662,0.1,-6.57E-03,-5.41E-03,1.246,-0.021,1.128,270,1350,760,-0.498,-0.785,-0.572,-0.44169,-0.0052,0,0,1,0,0,0.48,0.396,0.565,-0.18,423,0.221,0.093,0.149,0.327,0.068
+0.01,3.714,4.094796298,3.288796298,3.488,-1.587,0.1,-6.57E-03,-5.41E-03,1.246,-0.021,1.128,270,1300,760,-0.498,-0.803,-0.571,-0.4859,-0.0052,0,0,1,0,0,0.476,0.397,0.56,-0.18,423,0.223,0.098,0.148,0.294,0.071
+0.02,3.762,4.132796298,3.338796298,3.536,-1.593,0.1,-6.57E-03,-5.41E-03,1.227,-0.021,1.128,270,1225,760,-0.478,-0.785,-0.575,-0.4859,-0.00518,0,0,1,0,0,0.482,0.401,0.563,-0.181,423,0.227,0.105,0.149,0.294,0.073
+0.025,3.859,4.246796298,3.392796298,3.633,-1.607,0.1,-6.57E-03,-5.41E-03,1.221,-0.021,1.128,270,1200,760,-0.464,-0.745,-0.573,-0.4859,-0.00515,0,0,1,0,0,0.49,0.405,0.575,-0.183,423,0.231,0.12,0.15,0.31,0.076
+0.03,4.014,4.386796298,3.535796298,3.788,-1.63,0.1,-6.57E-03,-5.41E-03,1.215,-0.021,1.128,270,1200,760,-0.446,-0.69,-0.565,-0.4908,-0.00511,0,0,1,0,0,0.5,0.413,0.589,-0.188,423,0.239,0.145,0.153,0.313,0.077
+0.04,4.223,4.553796298,3.747796298,3.997,-1.657,0.1,-6.57E-03,-5.41E-03,1.207,-0.021,1.128,270,1200,760,-0.431,-0.636,-0.546,-0.49569,-0.00505,0,0,1,0,0,0.515,0.439,0.616,-0.205,423,0.261,0.177,0.159,0.322,0.08
+0.05,4.456,4.745796298,3.959796298,4.23,-1.687,0.1,-6.57E-03,-5.41E-03,1.201,-0.021,1.128,270,1225,760,-0.42,-0.594,-0.519,-0.49823,-0.00497,0.1,-0.1,-0.063,-0.05,0,0.528,0.473,0.653,-0.23,423,0.285,0.2,0.167,0.33,0.077
+0.075,4.742,4.972796298,4.231796298,4.516,-1.715,0.1,-6.57E-03,-5.41E-03,1.19,-0.021,1.128,270,1350,760,-0.442,-0.586,-0.497,-0.49724,-0.00489,0.3,-0.34,-0.2,-0.075,0.078,0.53,0.529,0.722,-0.262,423,0.339,0.205,0.184,0.299,0.063
+0.1,4.952,5.160796298,4.471796298,4.726,-1.737,0.1,-6.57E-03,-5.41E-03,1.182,-0.021,1.128,270,1450,760,-0.485,-0.629,-0.486,-0.49471,-0.00478,0.333,-0.377,-0.222,-0.081,0.075,0.524,0.517,0.712,-0.239,423,0.347,0.185,0.176,0.31,0.061
+0.15,5.08,5.285796298,4.665796298,4.848,-1.745,0.1,-6.57E-03,-5.41E-03,1.171,-0.021,1.162,270,1500,760,-0.546,-0.729,-0.499,-0.48583,-0.0046,0.29,-0.29,-0.193,-0.091,0.064,0.51,0.457,0.644,-0.185,423,0.313,0.123,0.164,0.307,0.076
+0.2,5.035,5.277796298,4.661796298,4.798,-1.732,0.1,-6.57E-03,-5.41E-03,1.163,-0.021,1.163,270,1425,760,-0.612,-0.867,-0.533,-0.47383,-0.00434,0.177,-0.192,-0.148,-0.092,0.075,0.501,0.432,0.64,-0.138,423,0.277,0.11,0.163,0.301,0.07
+0.25,4.859,5.154796298,4.503796298,4.618,-1.696,0.1,-6.57E-03,-5.41E-03,1.156,-0.021,1.156,270,1350,760,-0.688,-1.011,-0.592,-0.47696,-0.00402,0.1,-0.035,-0.054,0,0,0.492,0.45,0.633,-0.185,423,0.26,0.119,0.169,0.233,0.077
+0.3,4.583,4.910796298,4.276796298,4.34,-1.643,0.1,-6.57E-03,-5.41E-03,1.151,-0.021,1.151,270,1250,760,-0.748,-1.133,-0.681,-0.4845,-0.0037,0,0,1,0,0,0.492,0.436,0.584,-0.158,423,0.254,0.092,0.159,0.22,0.065
+0.4,4.18,4.548796298,3.919796298,3.935,-1.58,0.1,-6.57E-03,-5.41E-03,1.143,-0.022,1.143,270,1150,760,-0.802,-1.238,-0.772,-0.48105,-0.00342,0,0.05,0.2,0,0,0.492,0.433,0.556,-0.19,423,0.23,0.044,0.158,0.222,0.064
+0.5,3.752,4.168796298,3.486796298,3.505,-1.519,0.1,-6.57E-03,-5.41E-03,1.143,-0.023,1.143,270,1025,760,-0.845,-1.321,-0.838,-0.46492,-0.00322,0,0.1,0.2,0,0,0.492,0.428,0.51,-0.186,423,0.225,0.038,0.16,0.243,0.044
+0.75,3.085,3.510796298,2.710796298,2.837,-1.44,0.1,-6.35E-03,-4.78E-03,1.217,-0.026,1.217,270,900,760,-0.911,-1.383,-0.922,-0.43439,-0.00312,0,0.2,0.125,-0.2,0.012,0.492,0.448,0.471,-0.177,422,0.218,0.04,0.175,0.241,0.04
+1,2.644,3.067796298,2.238796298,2.396,-1.419,0.1,-5.80E-03,-4.15E-03,1.27,-0.028,1.24,270,800,760,-0.926,-1.414,-0.932,-0.38484,-0.0031,0,0.245,0.153,-0.245,0.037,0.492,0.43,0.43,-0.166,422,0.227,0.015,0.195,0.195,0.043
+1.5,2.046,2.513796298,1.451796298,1.799,-1.4,0.1,-5.05E-03,-3.42E-03,1.344,-0.031,1.237,270,760,760,-0.888,-1.43,-0.814,-0.32318,-0.0031,0,0.32,0.2,-0.32,0.064,0.492,0.406,0.406,-0.111,422,0.244,-0.047,0.204,0.204,-0.034
+2,1.556,2.061796298,0.906796298,1.31,-1.391,0.1,-4.29E-03,-2.90E-03,1.396,-0.034,1.232,270,760,760,-0.808,-1.421,-0.725,-0.26577,-0.0031,0,0.37,0.239,-0.28,0.14,0.492,0.393,0.393,0,422,0.231,-0.036,0.196,0.196,-0.036
+2.5,1.167,1.709796298,0.392796298,0.922,-1.394,0.1,-3.69E-03,-2.50E-03,1.437,-0.036,1.227,270,760,760,-0.743,-1.391,-0.632,-0.21236,-0.0031,0,0.4,0.264,-0.313,0.19,0.492,0.381,0.381,0,421,0.222,-0.025,0.169,0.169,-0.029
+3,0.92,1.456796298,0.099796298,0.675,-1.416,0.1,-3.21E-03,-2.17E-03,1.47,-0.038,1.223,270,760,760,-0.669,-1.343,-0.57,0,-0.0031,0,0.43,0.287,-0.355,0.165,0.492,0.367,0.367,0,419,0.199,-0.03,0.177,0.177,-0.011
+4,0.595,1.207796298,-0.356203702,0.352,-1.452,0.1,-2.44E-03,-1.65E-03,1.523,-0.044,1.216,270,760,760,-0.585,-1.297,-0.489,0,-0.0031,0,0.44,0.303,-0.417,0.163,0.492,0.33,0.33,0,416,0.191,-0.042,0.158,0.158,0.033
+5,0.465,1.131796298,-0.601203702,0.223,-1.504,0.1,-1.60E-03,-1.25E-03,1.564,-0.048,1.21,270,760,760,-0.506,-1.233,-0.421,0,-0.0031,0,0.45,0.321,-0.45,0.132,0.492,0.298,0.298,0,415,0.181,0.005,0.132,0.132,0.014
+7.5,0.078,0.758796298,-1.137203702,-0.162,-1.569,0.1,-7.66E-04,-5.19E-04,1.638,-0.059,1.2,270,760,760,-0.418,-1.147,-0.357,0,-0.0031,0,0.406,0.312,-0.35,0.15,0.492,0.254,0.254,0,419,0.181,-0.016,0.113,0.113,0.016
+10,0.046,0.708796298,-1.290203702,-0.193,-1.676,0.1,0.00E+00,0.00E+00,1.69,-0.067,1.194,270,760,760,-0.321,-1.06,-0.302,0,-0.0031,0,0.345,0.265,-0.331,0.117,0.492,0.231,0.231,0,427,0.181,0.04,0.11,0.11,0.017
diff --git a/src/main/resources/gmm/coeffs/PSHAB20_slab.csv b/src/main/resources/gmm/coeffs/PSHAB20_slab.csv
new file mode 100644
index 0000000000000000000000000000000000000000..ab56cdbbdc21fb164050798997eabc9a8a2fb5e4
--- /dev/null
+++ b/src/main/resources/gmm/coeffs/PSHAB20_slab.csv
@@ -0,0 +1,27 @@
+T,Global_c0,Alaska_c0,Aleutian_c0,Cascadia_c0,c1,b4,Global_a0,Alaska_a0,Cascadia_a0,c4,c5,c6,d,m,db,V1,V2,Vref,Global_s2,Alaska_s2,Cascadia_s2,f4,f5,C_e1,C_e2,C_e3,del_None,del_Seattle,Tau,phi21,phi22,phi2V,VM,phi2S2S0,a1,phi2SS1,phi2SS2,a2
+PGV,13.194,12.79,13.6,12.874,-2.422,0.1,-1.90E-03,-2.38E-03,-1.09E-03,1.84,-0.05,0.8,0.2693,0.0252,67,270,850,760,-0.601,-1.031,-0.671,-0.31763,-0.0052,0,0.115,0.068,-0.115,0,0.477,0.348,0.288,-0.179,423,0.142,0.047,0.153,0.166,0.011
+PGA,9.907,9.404,9.912,9.6,-2.543,0.1,-2.55E-03,-2.27E-03,-3.54E-03,1.84,-0.05,0.4,0.3004,0.0314,67,270,1350,760,-0.498,-0.785,-0.572,-0.44169,-0.0052,0,0,1,0,0,0.48,0.396,0.565,-0.18,423,0.221,0.093,0.149,0.327,0.068
+0.01,9.962,9.451,9.954,9.802,-2.554,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.84,-0.05,0.4,0.2839,0.0296,67,270,1300,760,-0.498,-0.803,-0.571,-0.4859,-0.0052,0,0,1,0,0,0.476,0.397,0.56,-0.18,423,0.223,0.098,0.148,0.294,0.071
+0.02,10.099,9.587,10.086,9.933,-2.566,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.05,0.415,0.2854,0.0298,67,270,1225,760,-0.478,-0.785,-0.575,-0.4859,-0.00518,0,0,1,0,0,0.482,0.401,0.563,-0.181,423,0.227,0.105,0.149,0.294,0.073
+0.025,10.181,9.667,10.172,10.009,-2.578,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.05,0.43,0.2891,0.0302,67,270,1200,760,-0.464,-0.745,-0.573,-0.4859,-0.00515,0,0,1,0,0,0.49,0.405,0.575,-0.183,423,0.231,0.12,0.15,0.31,0.076
+0.03,10.311,9.808,10.302,10.133,-2.594,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.05,0.445,0.2932,0.0306,67,270,1200,760,-0.446,-0.69,-0.565,-0.4908,-0.00511,0,0,1,0,0,0.5,0.413,0.589,-0.188,423,0.239,0.145,0.153,0.313,0.077
+0.04,10.588,10.086,10.602,10.404,-2.629,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.05,0.46,0.3004,0.0313,67,270,1200,760,-0.431,-0.636,-0.546,-0.49569,-0.00505,0,0,1,0,0,0.515,0.439,0.616,-0.205,423,0.261,0.177,0.159,0.322,0.08
+0.05,10.824,10.379,10.862,10.634,-2.649,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.05,0.475,0.3048,0.0316,67,270,1225,760,-0.42,-0.594,-0.519,-0.49823,-0.00497,0.1,-0.1,-0.063,-0.05,0,0.528,0.473,0.653,-0.23,423,0.285,0.2,0.167,0.33,0.077
+0.075,11.084,10.65,11.184,10.888,-2.65,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.05,0.49,0.2992,0.0321,67,270,1350,760,-0.442,-0.586,-0.497,-0.49724,-0.00489,0.3,-0.34,-0.2,-0.075,0.078,0.53,0.529,0.722,-0.262,423,0.339,0.205,0.184,0.299,0.063
+0.1,11.232,10.816,11.304,11.03,-2.647,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.05,0.505,0.2854,0.032,67,270,1450,760,-0.485,-0.629,-0.486,-0.49471,-0.00478,0.333,-0.377,-0.222,-0.081,0.075,0.524,0.517,0.712,-0.239,423,0.347,0.185,0.176,0.31,0.061
+0.15,11.311,10.883,11.402,11.103,-2.634,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.06,0.52,0.2814,0.0325,67,270,1500,760,-0.546,-0.729,-0.499,-0.48583,-0.0046,0.29,-0.29,-0.193,-0.091,0.064,0.51,0.457,0.644,-0.185,423,0.313,0.123,0.164,0.307,0.076
+0.2,11.055,10.633,11.183,10.841,-2.583,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.068,0.535,0.291,0.0306,67,270,1425,760,-0.612,-0.867,-0.533,-0.47383,-0.00434,0.177,-0.192,-0.148,-0.092,0.075,0.501,0.432,0.64,-0.138,423,0.277,0.11,0.163,0.301,0.07
+0.25,10.803,10.322,10.965,10.583,-2.539,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.075,0.55,0.2758,0.0306,67,270,1350,760,-0.688,-1.011,-0.592,-0.47696,-0.00402,0.1,-0.035,-0.054,0,0,0.492,0.45,0.633,-0.185,423,0.26,0.119,0.169,0.233,0.077
+0.3,10.669,10.116,10.87,10.443,-2.528,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.082,0.565,0.2719,0.0323,67,270,1250,760,-0.748,-1.133,-0.681,-0.4845,-0.0037,0,0,1,0,0,0.492,0.436,0.584,-0.158,423,0.254,0.092,0.159,0.22,0.065
+0.4,10.116,9.561,10.411,9.884,-2.452,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.091,0.58,0.2539,0.0302,67,270,1150,760,-0.802,-1.238,-0.772,-0.48105,-0.00342,0,0.05,0.2,0,0,0.492,0.433,0.556,-0.19,423,0.23,0.044,0.158,0.222,0.064
+0.5,9.579,8.973,9.901,9.341,-2.384,0.1,-2.55E-03,-2.19E-03,-4.01E-03,1.884,-0.1,0.595,0.2482,0.0295,67,270,1025,760,-0.845,-1.321,-0.838,-0.46492,-0.00322,0,0.1,0.2,0,0,0.492,0.428,0.51,-0.186,423,0.225,0.038,0.16,0.243,0.044
+0.75,8.837,8.246,9.335,8.593,-2.338,0.1,-2.11E-03,-1.89E-03,-3.47E-03,1.884,-0.115,0.61,0.2227,0.0266,67,270,900,760,-0.911,-1.383,-0.922,-0.43439,-0.00312,0,0.2,0.125,-0.2,0.012,0.492,0.448,0.471,-0.177,422,0.218,0.04,0.175,0.241,0.04
+1,8.067,7.507,8.68,7.817,-2.267,0.1,-1.87E-03,-1.68E-03,-3.09E-03,1.884,-0.134,0.625,0.1969,0.0231,67,270,800,760,-0.926,-1.414,-0.932,-0.38484,-0.0031,0,0.245,0.153,-0.245,0.037,0.492,0.43,0.43,-0.166,422,0.227,0.015,0.195,0.195,0.043
+1.5,6.829,6.213,7.581,6.573,-2.166,0.1,-1.54E-03,-1.39E-03,-2.54E-03,1.884,-0.154,0.64,0.1452,0.0118,67,270,760,760,-0.888,-1.43,-0.814,-0.32318,-0.0031,0,0.32,0.2,-0.32,0.064,0.492,0.406,0.406,-0.111,422,0.244,-0.047,0.204,0.204,-0.034
+2,5.871,5.206,6.671,5.609,-2.077,0.1,-1.31E-03,-1.18E-03,-2.16E-03,1.884,-0.154,0.655,0.06,0.007,67,270,760,760,-0.808,-1.421,-0.725,-0.26577,-0.0031,0,0.37,0.239,-0.28,0.14,0.492,0.393,0.393,0,422,0.231,-0.036,0.196,0.196,-0.036
+2.5,5.2,4.594,6.047,4.932,-2.015,0.1,-1.13E-03,-1.01E-03,-1.86E-03,1.884,-0.154,0.67,0,0,0,270,760,760,-0.743,-1.391,-0.632,-0.21236,-0.0031,0,0.4,0.264,-0.313,0.19,0.492,0.381,0.381,0,421,0.222,-0.025,0.169,0.169,-0.029
+3,4.83,4.206,5.667,4.556,-2.012,0.1,-9.79E-04,-8.80E-04,-1.61E-03,1.949,-0.154,0.685,0,0,0,270,760,760,-0.669,-1.343,-0.57,-0.17807,-0.0031,0,0.43,0.287,-0.355,0.165,0.492,0.367,0.367,0,419,0.199,-0.03,0.177,0.177,-0.011
+4,4.173,3.517,4.97,3.893,-1.989,0.1,-7.45E-04,-6.70E-04,-1.23E-03,2.031,-0.154,0.7,0,0,0,270,760,760,-0.585,-1.297,-0.489,-0.13729,-0.0031,0,0.44,0.303,-0.417,0.163,0.492,0.33,0.33,0,416,0.191,-0.042,0.158,0.158,0.033
+5,3.833,3.142,4.592,3.547,-1.998,0.1,-5.64E-04,-5.07E-04,-9.29E-04,2.131,-0.154,0.715,0,0,0,270,760,760,-0.506,-1.233,-0.421,-0.07733,-0.0031,0,0.45,0.321,-0.45,0.132,0.492,0.298,0.298,0,415,0.181,0.005,0.132,0.132,0.014
+7.5,3.132,2.391,3.65,2.84,-2.019,0.1,-2.34E-04,-2.10E-04,-3.85E-04,2.185,-0.154,0.73,0,0,0,270,760,760,-0.418,-1.147,-0.357,-0.05443,-0.0031,0,0.406,0.312,-0.35,0.15,0.492,0.254,0.254,0,419,0.181,-0.016,0.113,0.113,0.016
+10,2.72,2.031,2.95,2.422,-2.047,0.1,0.00E+00,0.00E+00,0.00E+00,2.35,-0.154,0.745,0,0,0,270,760,760,-0.321,-1.06,-0.302,-0.03313,-0.0031,0,0.345,0.265,-0.331,0.117,0.492,0.231,0.231,0,427,0.181,0.04,0.11,0.11,0.017
diff --git a/src/test/java/gov/usgs/earthquake/nshmp/gmm/PSHAB20.java b/src/test/java/gov/usgs/earthquake/nshmp/gmm/PSHAB20.java
new file mode 100644
index 0000000000000000000000000000000000000000..5d2f86bcb3ad22e17eea950d0f024e2564c48321
--- /dev/null
+++ b/src/test/java/gov/usgs/earthquake/nshmp/gmm/PSHAB20.java
@@ -0,0 +1,44 @@
+package gov.usgs.earthquake.nshmp.gmm;
+
+import static gov.usgs.earthquake.nshmp.gmm.Gmm.PSHAB20_INTERFACE;
+import static gov.usgs.earthquake.nshmp.gmm.Gmm.PSHAB20_SLAB;
+import static gov.usgs.earthquake.nshmp.gmm.Imt.PGA;
+
+import java.io.IOException;
+import java.util.EnumSet;
+import java.util.Set;
+import java.util.stream.Stream;
+
+import org.junit.jupiter.api.extension.ExtensionContext;
+import org.junit.jupiter.params.ParameterizedTest;
+import org.junit.jupiter.params.provider.Arguments;
+import org.junit.jupiter.params.provider.ArgumentsProvider;
+import org.junit.jupiter.params.provider.ArgumentsSource;
+
+class PSHAB20 implements ArgumentsProvider {
+
+ private static String GMM_INPUTS = "/gmm/PSHAB20_inputs.csv";
+ private static String GMM_RESULTS = "/gmm/PSHAB20_results.csv";
+ // private static String GMM_INPUTS = "/gmm/nga-west-inputs.csv";
+ // private static String GMM_RESULTS = "/gmm/nga-west2-results.csv";
+
+ @ParameterizedTest(name = "({index}) {0} {2} {1}")
+ @ArgumentsSource(PSHAB20.class)
+ void test(int index, Gmm gmm, Imt imt, double exMedian, double exSigma, String inputs) {
+ GmmTest.test(index, gmm, imt, exMedian, exSigma, inputs);
+ }
+
+ @Override
+ public Stream<? extends Arguments> provideArguments(ExtensionContext context) throws Exception {
+ return GmmTest.loadResults(GMM_RESULTS, GMM_INPUTS);
+ }
+
+ /* Result generation sets */
+ private static Set<Gmm> gmms = EnumSet.of(PSHAB20_INTERFACE, PSHAB20_SLAB);
+ private static Set<Imt> imts = EnumSet.of(PGA);// PSHAB20_INTERFACE.supportedIMTs();
+
+ public static void main(String[] args) throws IOException {
+ GmmTest.generateResults(gmms, imts, GMM_INPUTS, GMM_RESULTS);
+ }
+
+}
diff --git a/src/test/resources/gmm/PSHAB20_inputs.csv b/src/test/resources/gmm/PSHAB20_inputs.csv
new file mode 100644
index 0000000000000000000000000000000000000000..5786a8bc3a3823cd2f244a96a30388a611e6479d
--- /dev/null
+++ b/src/test/resources/gmm/PSHAB20_inputs.csv
@@ -0,0 +1,2 @@
+#Mw,rjB,rRup,rX,dip,width,zTop,zHyp,rake,vs30,vsInf,z1p0,z2p5
+8,0,200,0,90,0,0,55,90,760,0,550,NaN
diff --git a/src/test/resources/gmm/PSHAB20_results.csv b/src/test/resources/gmm/PSHAB20_results.csv
new file mode 100644
index 0000000000000000000000000000000000000000..3829c8b2e7bc72ae5b7d25b46c142657735a16c7
--- /dev/null
+++ b/src/test/resources/gmm/PSHAB20_results.csv
@@ -0,0 +1,46 @@
+0-PSHAB20_INTERFACE-SA0P01,2.08965604936e-02,7.89668284788e-01
+0-PSHAB20_INTERFACE-SA0P02,2.12376645069e-02,7.95816561778e-01
+0-PSHAB20_INTERFACE-SA0P03,2.24575034172e-02,8.14248119433e-01
+0-PSHAB20_INTERFACE-SA0P05,2.58276475982e-02,8.67054784890e-01
+0-PSHAB20_INTERFACE-SA0P075,2.96364966819e-02,8.99944442730e-01
+0-PSHAB20_INTERFACE-SA0P1,3.25369383278e-02,8.89705569276e-01
+0-PSHAB20_INTERFACE-SA0P15,3.55651022874e-02,8.46817571854e-01
+0-PSHAB20_INTERFACE-SA0P2,3.64296484198e-02,8.26438745462e-01
+0-PSHAB20_INTERFACE-SA0P25,3.69487943052e-02,8.31903840597e-01
+0-PSHAB20_INTERFACE-SA0P3,3.71144622014e-02,8.23446415986e-01
+0-PSHAB20_INTERFACE-SA0P4,3.46117876964e-02,8.21622784494e-01
+0-PSHAB20_INTERFACE-SA0P5,3.11738351662e-02,8.18574370476e-01
+0-PSHAB20_INTERFACE-SA0P75,2.56085017514e-02,8.30700908872e-01
+0-PSHAB20_INTERFACE-SA1P0,2.06123859520e-02,8.19795096350e-01
+0-PSHAB20_INTERFACE-SA1P5,1.45673139070e-02,8.05024223238e-01
+0-PSHAB20_INTERFACE-SA2P0,1.08967064936e-02,7.96909028685e-01
+0-PSHAB20_INTERFACE-SA3P0,6.26946254992e-03,7.80425524954e-01
+0-PSHAB20_INTERFACE-SA4P0,4.36498415871e-03,7.56349125735e-01
+0-PSHAB20_INTERFACE-SA5P0,3.44133448941e-03,7.34890468029e-01
+0-PSHAB20_INTERFACE-SA7P5,1.95537566087e-03,7.04318109947e-01
+0-PSHAB20_INTERFACE-SA10P0,1.25157789513e-03,6.87796481526e-01
+0-PSHAB20_INTERFACE-PGA,2.02870001012e-02,7.91454357496e-01
+0-PSHAB20_INTERFACE-PGV,1.95216532021e+00,7.58636276486e-01
+0-PSHAB20_SLAB-SA0P01,7.19805699097e-02,7.89668284788e-01
+0-PSHAB20_SLAB-SA0P02,7.78462260065e-02,7.95816561778e-01
+0-PSHAB20_SLAB-SA0P03,8.37771432985e-02,8.14248119433e-01
+0-PSHAB20_SLAB-SA0P05,1.05691366814e-01,8.67054784890e-01
+0-PSHAB20_SLAB-SA0P075,1.35586500017e-01,8.99944442730e-01
+0-PSHAB20_SLAB-SA0P1,1.58689587375e-01,8.89705569276e-01
+0-PSHAB20_SLAB-SA0P15,1.83281681820e-01,8.46817571854e-01
+0-PSHAB20_SLAB-SA0P2,1.93331748629e-01,8.26438745462e-01
+0-PSHAB20_SLAB-SA0P25,1.88104417586e-01,8.31903840597e-01
+0-PSHAB20_SLAB-SA0P3,1.71253449867e-01,8.23446415986e-01
+0-PSHAB20_SLAB-SA0P4,1.49480063250e-01,8.21622784494e-01
+0-PSHAB20_SLAB-SA0P5,1.26491388754e-01,8.18574370476e-01
+0-PSHAB20_SLAB-SA0P75,8.53899125993e-02,8.30700908872e-01
+0-PSHAB20_SLAB-SA1P0,6.18931082597e-02,8.19795096350e-01
+0-PSHAB20_SLAB-SA1P5,3.59081844251e-02,8.05024223238e-01
+0-PSHAB20_SLAB-SA2P0,2.26680607353e-02,7.96909028685e-01
+0-PSHAB20_SLAB-SA3P0,1.25351157474e-02,7.80425524954e-01
+0-PSHAB20_SLAB-SA4P0,7.74659496901e-03,7.56349125735e-01
+0-PSHAB20_SLAB-SA5P0,5.48590963050e-03,7.34890468029e-01
+0-PSHAB20_SLAB-SA7P5,2.61850755964e-03,7.04318109947e-01
+0-PSHAB20_SLAB-SA10P0,1.57671022525e-03,6.87796481526e-01
+0-PSHAB20_SLAB-PGA,7.18617131166e-02,7.91454357496e-01
+0-PSHAB20_SLAB-PGV,5.11504074039e+00,7.58636276486e-01