diff --git a/src/main/java/gov/usgs/earthquake/nshmp/gmm/UsgsDampingScaling.java b/src/main/java/gov/usgs/earthquake/nshmp/gmm/UsgsDampingScaling.java
index eb84a2e617548e8cc3ff6059336f9d16277e28d0..66e8b5d87e555ff521b8845f10e466489a3231a0 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/gmm/UsgsDampingScaling.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/gmm/UsgsDampingScaling.java
@@ -2,20 +2,25 @@ package gov.usgs.earthquake.nshmp.gmm;
import static com.google.common.base.Preconditions.checkArgument;
import static gov.usgs.earthquake.nshmp.data.DoubleData.checkInRange;
+import static gov.usgs.earthquake.nshmp.gmm.Gmm.Type.ACTIVE_CRUST;
+import static gov.usgs.earthquake.nshmp.gmm.Gmm.Type.STABLE_CRUST;
+import static gov.usgs.earthquake.nshmp.gmm.Gmm.Type.SUBDUCTION_INTERFACE;
+import static gov.usgs.earthquake.nshmp.gmm.Gmm.Type.SUBDUCTION_SLAB;
import static java.lang.Math.abs;
import static java.lang.Math.log;
import static java.lang.Math.sqrt;
+import static java.util.stream.Collectors.toList;
+import java.util.Arrays;
import java.util.EnumMap;
import java.util.List;
import java.util.Map;
import java.util.Optional;
+import java.util.Set;
import com.google.common.annotations.Beta;
-import com.google.common.collect.FluentIterable;
import com.google.common.collect.Lists;
import com.google.common.collect.Range;
-import com.google.common.primitives.Doubles;
import gov.usgs.earthquake.nshmp.data.Interpolator;
import gov.usgs.earthquake.nshmp.data.XySequence;
@@ -32,8 +37,7 @@ import gov.usgs.earthquake.nshmp.tree.Trees;
*
* <p><b>Implementation notes:</b>
*
- * <ul><li>Active crust model uses coefficients for RotD50, Table 1 in Rezaeian
- * et al. (2014)
+ * <ul><li>Active crust model is used for stable crust GMMs.
*
* </ul>
*
@@ -62,9 +66,8 @@ public class UsgsDampingScaling {
* Current set up is to use instances by tectonic setting.
*
* The calcDsf() model returns a value much closer to zero for a damping ratio
- * of 5.027590868... The assumption is that only coarser damping ratios will
- * be used (e.g. 3.5, 4.0) such that this inconsistency will not make a
- * difference. Are we going to support a range of fixed values?
+ * of 5.027590868. By supporting only fixed, coarser damping ratios (e.g. 3,
+ * 4) this inconsistency does not make a difference.
*
* For subduction, (1) special handling for periods > 5s, (2) correlation
* coefficient table missing 5% column (is 0.0 in crustal), added column for
@@ -79,24 +82,26 @@ public class UsgsDampingScaling {
/** The range of supported damping ratios: {@code [0.5..30]}. */
public static final Range<Double> DAMPING_RATIO_RANGE = Range.closed(0.05, 30.0);
- private static final CoefficientContainer COEFFS_ACTIVE_CRUST;
- private static final CoefficientContainer COEFFS_INTERFACE;
- private static final CoefficientContainer COEFFS_SLAB;
-
- private static final String[] corrCoeffKeys =
+ private static final String[] DSF_KEYS =
{ "0.5", "1", "2", "3", "5", "7", "10", "15", "20", "25", "30" };
- private static final List<Double> CORR_COEFF_DSFS =
- FluentIterable.from(corrCoeffKeys)
- .transform(Doubles.stringConverter())
- .toList();
+ private static final List<Double> DSF_VALUES = Arrays.stream(DSF_KEYS)
+ .map(Double::valueOf)
+ .collect(toList());
private static final Interpolator CC_INTERP = Interpolator.builder().build();
+ private static final Map<Gmm.Type, UsgsDampingScaling> INSTANCE_MAP;
+
static {
- COEFFS_ACTIVE_CRUST = new CoefficientContainer("Rezaeian14-active-crust.csv");
- COEFFS_INTERFACE = new CoefficientContainer("Rezaeian21-interface.csv");
- COEFFS_SLAB = new CoefficientContainer("Rezaeian21-slab.csv");
+ var dsfActive = new UsgsDampingScaling("Rezaeian14-active-crust.csv");
+ var dsfInterface = new UsgsDampingScaling("Rezaeian21-interface.csv");
+ var dsfSlab = new UsgsDampingScaling("Rezaeian21-slab.csv");
+ INSTANCE_MAP = Map.of(
+ ACTIVE_CRUST, dsfActive,
+ STABLE_CRUST, dsfActive,
+ SUBDUCTION_INTERFACE, dsfInterface,
+ SUBDUCTION_SLAB, dsfSlab);
}
private static final class Coefficients {
@@ -122,39 +127,8 @@ public class UsgsDampingScaling {
private final Map<Imt, Coefficients> coeffs = new EnumMap<>(Imt.class);
private final Map<Imt, XySequence> corrCoeffs = new EnumMap<>(Imt.class);
- // private final SourceType type;
-
- /**
- * Create a new damping scaling factor (DSF) model instance for the specified
- * tectonic settign and source type.
- *
- * @param setting tectonic setting for the DSF
- * @param type source type for the DSF; this argument is only relevant if
- * {@code TectonicSetting == SUBDUCTION}
- */
- private UsgsDampingScaling(Gmm.Type type) {
- // Dallin- Modified the function to only use source type to determine the
- // correct coeffs
- // this.type = type;
- CoefficientContainer cc = (type == Gmm.Type.SUBDUCTION_INTERFACE)
- ? COEFFS_INTERFACE
- : (type == Gmm.Type.SUBDUCTION_SLAB)
- ? COEFFS_SLAB
- : COEFFS_ACTIVE_CRUST;
-
- // COEFFS_ACTIVE_CRUST;
- // cc = ()
- // if (type == Gmm.Type.SUBDUCTION_INTERFACE)
- // if (type == INTERFACE || type == INTERFACE_CLUSTER || type ==
- // INTERFACE_SYSTEM ||
- // type == INTERFACE_GRID) {
- // cc = COEFFS_INTERFACE;
- // } else if (type == SLAB || type == SLAB_GRID) {
- // cc = COEFFS_SLAB;
- // } else {
- // cc = COEFFS_ACTIVE_CRUST;
- // }
-
+ private UsgsDampingScaling(String coeffsFile) {
+ CoefficientContainer cc = new CoefficientContainer(coeffsFile);
for (Imt imt : cc.imts()) {
coeffs.put(imt, new Coefficients(imt, cc));
corrCoeffs.put(imt, initCorrCoeff(imt, cc));
@@ -168,10 +142,15 @@ public class UsgsDampingScaling {
private static XySequence initCorrCoeff(Imt imt, CoefficientContainer cc) {
Map<String, Double> coeffs = cc.get(imt);
List<Double> corrCoeffs = Lists.newArrayList();
- for (String key : corrCoeffKeys) {
+ for (String key : DSF_KEYS) {
corrCoeffs.add(coeffs.get(key));
}
- return XySequence.create(CORR_COEFF_DSFS, corrCoeffs);
+ return XySequence.create(DSF_VALUES, corrCoeffs);
+ }
+
+ /** The set of IMTs supported by this class. */
+ public static final Set<Imt> supportedImts() {
+ return INSTANCE_MAP.get(ACTIVE_CRUST).coeffs.keySet();
}
/**
@@ -186,21 +165,9 @@ public class UsgsDampingScaling {
public double factor(Imt imt, double Mw, double rRup, double dampingRatio) {
checkArgument(coeffs.containsKey(imt), "Unsupported IMT: %s", imt.name());
checkInRange(DAMPING_RATIO_RANGE, "Damping ratio", dampingRatio);
- return (dampingRatio != 5.0) ? calcFactor(imt, Mw, rRup, dampingRatio) : 0;
- // double lnDSF = 0.0;
- // if (dampingRatio != 5.0) {
- // if (imt.period() > 5.0 &&
- // (type == INTERFACE || type == INTERFACE_CLUSTER || type ==
- // INTERFACE_SYSTEM ||
- // type == INTERFACE_GRID || type == SLAB || type == SLAB_GRID)) {
- // // System.out.println("LONGPERIOD");
- // lnDSF = calclnDsflongPeriod(imt, Mw, rRup, dampingRatio);
- // } else {
- // // System.out.println("normal");
- // lnDSF = calclnDsf(imt, Mw, rRup, dampingRatio);
- // }
- // }
- // return lnDSF;
+ return (dampingRatio != 5.0)
+ ? calcFactor(imt, Mw, rRup, dampingRatio)
+ : 0;
}
double calcFactor(Imt imt, double Mw, double rRup, double dampingRatio) {
@@ -281,13 +248,6 @@ public class UsgsDampingScaling {
return sqrt(σPsa * σPsa + σDsf * σDsf + 2.0 * σPsa * σDsf * Ï);
}
- private static final UsgsDampingScaling INSTANCE_DEFAULT =
- new UsgsDampingScaling(Gmm.Type.ACTIVE_CRUST);
- private static final UsgsDampingScaling INSTANCE_INTERFACE =
- new UsgsDampingScaling(Gmm.Type.SUBDUCTION_INTERFACE);
- private static final UsgsDampingScaling INSTANCE_SLAB =
- new UsgsDampingScaling(Gmm.Type.SUBDUCTION_SLAB);
-
public static LogicTree<GroundMotion> apply(
Imt imt,
GmmInput in,
@@ -299,20 +259,8 @@ public class UsgsDampingScaling {
return tree;
}
- /*
- * Changed method signature
- *
- * Switch to set instance to use.
- *
- * factor() sigma() and rho() don't depend on mean and sigma values; only be
- * called once
- */
-
- UsgsDampingScaling instance = (type == Gmm.Type.SUBDUCTION_INTERFACE)
- ? INSTANCE_INTERFACE
- : (type == Gmm.Type.SUBDUCTION_SLAB)
- ? INSTANCE_SLAB
- : INSTANCE_DEFAULT;
+ UsgsDampingScaling instance = INSTANCE_MAP.get(type);
+ checkArgument(instance.coeffs.containsKey(imt), "Unsupported IMT: %s", imt.name());
double factor = instance.calcFactor(imt, in.Mw, in.rRup, dampingRatio);
double sigma = instance.calcSigma(imt, dampingRatio);
@@ -324,52 +272,17 @@ public class UsgsDampingScaling {
gm.mean() + factor,
updateSigma(gm.sigma(), sigma, rho)),
Optional.empty());
-
- // if (type == Gmm.Type.SUBDUCTION_INTERFACE) {
- //
- // // if (gmmInput.srcType == INTERFACE || gmmInput.srcType ==
- // // INTERFACE_CLUSTER ||
- // // gmmInput.srcType == INTERFACE_SYSTEM || gmmInput.srcType ==
- // // INTERFACE_GRID) {
- //
- // return Trees.transform(tree,
- // gm -> GroundMotion.create(
- // gm.mean() +
- // INSTANCE_INTERFACE.factor(imt, gmmInput.Mw, gmmInput.rRup, dampingRatio),
- // updateSigma(gm.sigma(), INSTANCE_INTERFACE.sigma(imt, dampingRatio),
- // INSTANCE_INTERFACE.rho(imt, dampingRatio))),
- // Optional.empty());
- //
- // } else if (type == Gmm.Type.SUBDUCTION_SLAB) {
- // // if (gmmInput.srcType == SLAB || gmmInput.srcType == SLAB_GRID) {
- // return Trees.transform(tree,
- // gm -> GroundMotion.create(
- // gm.mean() +
- // INSTANCE_SLAB.factor(imt, gmmInput.Mw, gmmInput.rRup, dampingRatio),
- // updateSigma(gm.sigma(), INSTANCE_SLAB.sigma(imt, dampingRatio),
- // INSTANCE_SLAB.rho(imt, dampingRatio))),
- // Optional.empty());
- //
- // } else {
- // // active and stable crust
- // return Trees.transform(tree,
- // gm -> GroundMotion.create(
- // gm.mean() +
- // INSTANCE_DEFAULT.factor(imt, gmmInput.Mw, gmmInput.rRup, dampingRatio),
- // updateSigma(gm.sigma(), INSTANCE_DEFAULT.sigma(imt, dampingRatio),
- // INSTANCE_DEFAULT.rho(imt, dampingRatio))),
- // Optional.empty());
- // }
}
+ // TODO clean
public static void main(String[] args) {
- UsgsDampingScaling model = new UsgsDampingScaling(Gmm.Type.SUBDUCTION_INTERFACE);
+ UsgsDampingScaling model = INSTANCE_MAP.get(SUBDUCTION_INTERFACE);
- Imt imt = Imt.SA5P0;
+ Imt imt = Imt.SA1P0;
double Mw = 7.0;
double rRup = 120;
- double ratio = 5.0;
+ double ratio = 3.0;
double scaleFactor = Math.exp(model.factor(imt, Mw, rRup, ratio));
double dsfSigma = model.sigma(imt, ratio);