diff --git a/src/main/java/gov/usgs/earthquake/nshmp/gmm/NgaEast.java b/src/main/java/gov/usgs/earthquake/nshmp/gmm/NgaEast.java
index 5e130b953a54a06d53d9ffe4fc85110c3e7ee322..ae96cbd6cce450cdc7f2c8025096930652ff0ac3 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/gmm/NgaEast.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/gmm/NgaEast.java
@@ -342,7 +342,9 @@ public abstract class NgaEast implements GroundMotionModel {
/* Implementation of sammons map models. */
static class NgaEast_2018 extends NgaEast {
+
static final String NAME = NgaEast.NAME + " (2018)";
+
static final int MODEL_COUNT = 17;
static final String[] MEAN_IDS = IntStream.range(1, MODEL_COUNT + 1)
.mapToObj(i -> "Sammons" + i)
@@ -396,16 +398,30 @@ public abstract class NgaEast implements GroundMotionModel {
super(imt);
}
+ /*
+ * Developer notes: The current implementation assumes a reference site
+ * condition of Vs30 = 590 m/s. This means that we apply the full Chapman &
+ * Guo PSA ratio if Vs30 = 590 and scale it by the ratio of the site term at
+ * the target Vs30 to the reference Vs30.
+ */
@Override
public LogicTree<GroundMotion> calc(GmmInput in) {
- double cpa = Double.isNaN(in.zSed)
- ? 0.0
- : log(ChapmanGuo_2021.cpaPsaRatio(imt, in.zSed, in.Mw, in.rJB));
+ boolean cpa = !Double.isNaN(in.zSed);
+ double fCpa = cpa
+ ? log(ChapmanGuo_2021.cpaPsaRatio(imt, in.zSed, in.Mw, in.rJB))
+ : 0.0;
Position p = tables[0].position(in.rRup, in.Mw);
double[] μs = new double[MODEL_COUNT];
for (int i = 0; i < MODEL_COUNT; i++) {
double μ = tables[i].get(p);
- μs[i] = μ + cpa;
+ double μPga = exp(pgaTables[i].get(p));
+ SiteAmp.Value fSite = siteAmp.calc(μPga, in.vs30);
+ if (cpa) {
+ SiteAmp.Value fCpaRef = siteAmp.calc(μPga, ChapmanGuo_2021.VS_REF);
+ double cpaScale = fSite.siteAmp / fCpaRef.siteAmp;
+ fCpa *= cpaScale;
+ }
+ μs[i] = fSite.apply(μ) + fCpa;
}
double[] σs = new double[] {
sigmaEpri(in.Mw),
@@ -452,6 +468,7 @@ public abstract class NgaEast implements GroundMotionModel {
final Map<Gmm, GroundMotionTable> tables;
final Map<Gmm, GroundMotionTable> pgaTables;
final ShahjoueiPezeshk_2016 sp16;
+ final ShahjoueiPezeshk_2016 sp16pga;
final SiteAmp siteAmp;
static {
@@ -481,6 +498,8 @@ public abstract class NgaEast implements GroundMotionModel {
this.tables = getImtTables(imt);
this.pgaTables = getImtTables(Imt.PGA);
this.sp16 = (ShahjoueiPezeshk_2016) Gmm.NGA_EAST_SEED_SP16.instance(imt);
+ this.sp16pga = (ShahjoueiPezeshk_2016) Gmm.NGA_EAST_SEED_SP16.instance(Imt.PGA);
+
this.siteAmp = new SiteAmp(imt);
}
@@ -499,15 +518,18 @@ public abstract class NgaEast implements GroundMotionModel {
for (int i = 0; i < GMMS.size(); i++) {
Gmm seed = GMMS.get(i);
if (seed == Gmm.NGA_EAST_SEED_SP16) {
- μs[i] = GroundMotions.combine(sp16.calc(in)).mean();
+ double μ = sp16.calcMeanRock(in.Mw, in.rRup);
+ double μPga = exp(sp16pga.calcMeanRock(in.Mw, in.rRup));
+ SiteAmp.Value fSite = siteAmp.calc(μPga, in.vs30);
+ μs[i] = fSite.apply(μ);
} else if (imt == Imt.PGV && noPgvSeeds.contains(seed)) {
// site is included when calling individual seed
μs[i] = UsgsPgvSupport.calcAB20Pgv(seed, in).mean();
} else {
- double μRock = tables.get(seed).get(p);
+ double μ = tables.get(seed).get(p);
double μPga = exp(pgaTables.get(seed).get(p));
SiteAmp.Value fSite = siteAmp.calc(μPga, in.vs30);
- μs[i] = fSite.apply(μRock);
+ μs[i] = fSite.apply(μ);
}
}
double[] σs = new double[] {
@@ -519,6 +541,19 @@ public abstract class NgaEast implements GroundMotionModel {
}
}
+ public static void main(String[] args) {
+ // Gmm gmm = Gmm.NGA_EAST_SEED_SP16;
+ Gmm gmm = Gmm.NGA_EAST_SEEDS_2018;
+ // Gmm gmm = Gmm.NGA_EAST_SEED_SP16;
+ GmmInput in = GmmInput.builder().withDefaults().build();
+ LogicTree<GroundMotion> gms = gmm.instance(Imt.PGA).calc(in);
+ System.out.println(gms);
+ gms.forEach(gm -> System.out.println(gm.id() + ": " + gm.value()));
+ GroundMotion gm = GroundMotions.combine(gms);
+ System.out.println(gm);
+
+ }
+
/*
* Updated NGA-East for use in the 2023 nshm-conus update. This model includes
* (1) the final published nonlinear site amplification model of Hshash et al.
@@ -536,21 +571,34 @@ public abstract class NgaEast implements GroundMotionModel {
@Override
public LogicTree<GroundMotion> calc(GmmInput in) {
GmmInput inRock = GmmInput.builder().fromCopy(in).vs30(3000).build();
- double cpa = Double.isNaN(in.zSed)
- ? 0.0
- : log(ChapmanGuo_2021.cpaPsaRatio(imt, in.zSed, in.Mw, in.rJB));
+ boolean cpa = !Double.isNaN(in.zSed);
+ double fCpa = cpa
+ ? log(ChapmanGuo_2021.cpaPsaRatio(imt, in.zSed, in.Mw, in.rJB))
+ : 0.0;
Position p = tables.values().iterator().next().position(in.rRup, in.Mw);
double[] μs = new double[GMMS.size()];
for (int i = 0; i < GMMS.size(); i++) {
Gmm seed = GMMS.get(i);
+
+ double μ;
+ double μPga;
+
if (seed == Gmm.NGA_EAST_SEED_SP16) {
- μs[i] = GroundMotions.combine(sp16.calc(inRock)).mean() + cpa;
+ μ = sp16.calcMeanRock(in.Mw, in.rRup);
+ μPga = sp16pga.calcMeanRock(in.Mw, in.rRup);
} else if (imt == Imt.PGV && noPgvSeeds.contains(seed)) {
- // use custom input to ensure hard rock from seeds
- μs[i] = UsgsPgvSupport.calcAB20Pgv(seed, inRock).mean() + cpa;
+ // note rock input; TODO this doesn't really work as the PGV
+ // is conditioned on vs30 so site should already be applied
+ μ = UsgsPgvSupport.calcAB20Pgv(seed, inRock).mean();
+ μPga = exp(pgaTables.get(seed).get(p));
} else {
- μs[i] = tables.get(seed).get(p) + cpa;
+ μ = tables.get(seed).get(p);
+ μPga = exp(pgaTables.get(seed).get(p));
}
+ SiteAmp.Value fSite = siteAmp.calc(μPga, in.vs30);
+ double cpaScale = cpa ? cpaScale(μPga, fSite) : 0.0;
+ μs[i] = fSite.apply(μ) + fCpa * cpaScale;
+
}
double[] σs = new double[] {
sigmaEpri(in.Mw),
@@ -559,6 +607,11 @@ public abstract class NgaEast implements GroundMotionModel {
MEAN_IDS, μs, MEAN_WTS,
SIGMA_IDS, σs, SIGMA_WTS);
}
+
+ private double cpaScale(double μPga, SiteAmp.Value fSite) {
+ SiteAmp.Value fCpaRef = siteAmp.calc(μPga, ChapmanGuo_2021.VS_REF);
+ return fSite.siteAmp / fCpaRef.siteAmp;
+ }
}
static abstract class Seed extends NgaEast {