diff --git a/src/main/java/gov/usgs/earthquake/nshmp/calc/CalcConfig.java b/src/main/java/gov/usgs/earthquake/nshmp/calc/CalcConfig.java index b1b08fab85a6223ba65e40a7864a16005c794510..51377f378382efc25b5e2633cf8b5e5ce6ad7f09 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/calc/CalcConfig.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/calc/CalcConfig.java @@ -203,6 +203,14 @@ public final class CalcConfig { */ public final Set<Double> vs30s; + /** + * Whether or not consider site-specific data that may accompany a model. + * Set to false to turn off deep basin scaling effects, for example. + * + * <p><b>Default:</b> true + */ + public final boolean useSiteData; + /** * Viscous damping ratio. The default value is consistent with that used in * {@link GroundMotionModel} most (GMM) development. Values other than 5% @@ -221,6 +229,15 @@ public final class CalcConfig { */ public final double gmmDampingRatio; + /** + * Whether or not {@link GroundMotionModel} (GMM) sigmas are also modified + * according to the Rezaeian et al. (2014) model when applying damping + * scaling factors (DSF). + * + * <p><b>Default:</b> false + */ + public final boolean gmmDampingSigma; + /** * Factor by which the aleatory variability (sigma) for all GMMs should be * scaled. @@ -250,7 +267,9 @@ public final class CalcConfig { this.tectonicSettings = b.tectonicSettings; this.sourceTypes = b.sourceTypes; this.vs30s = b.vs30s; + this.useSiteData = b.useSiteData; this.gmmDampingRatio = b.gmmDampingRatio; + this.gmmDampingSigma = b.gmmDampingSigma; this.gmmSigmaScale = b.gmmSigmaScale; this.valueFormat = b.valueFormat; @@ -292,7 +311,9 @@ public final class CalcConfig { Set<TectonicSetting> tectonicSettings; Set<SourceType> sourceTypes; Set<Double> vs30s; + Boolean useSiteData; Double gmmDampingRatio; + Boolean gmmDampingSigma; Double gmmSigmaScale; ValueFormat valueFormat; Map<Imt, double[]> customImls; @@ -307,11 +328,13 @@ public final class CalcConfig { checkNotNull(tectonicSettings); checkNotNull(sourceTypes); checkNotNull(vs30s); + checkNotNull(useSiteData); checkNotNull(gmmDampingRatio); checkInRange( RezaeianDamping_2014.DAMPING_RATIO_RANGE, "gmmDampingRatio", gmmDampingRatio); + checkNotNull(gmmDampingSigma); checkNotNull(gmmSigmaScale); checkInRange( Range.closed(0.5, 1.0), @@ -329,7 +352,9 @@ public final class CalcConfig { this.tectonicSettings = that.tectonicSettings; this.sourceTypes = that.sourceTypes; this.vs30s = that.vs30s; + this.useSiteData = that.useSiteData; this.gmmDampingRatio = that.gmmDampingRatio; + this.gmmDampingSigma = that.gmmDampingSigma; this.gmmSigmaScale = that.gmmSigmaScale; this.valueFormat = that.valueFormat; this.customImls = that.customImls; @@ -354,9 +379,15 @@ public final class CalcConfig { if (that.vs30s != null) { this.vs30s = that.vs30s; } + if (that.useSiteData != null) { + this.useSiteData = that.useSiteData; + } if (that.gmmDampingRatio != null) { this.gmmDampingRatio = that.gmmDampingRatio; } + if (that.gmmDampingSigma != null) { + this.gmmDampingSigma = that.gmmDampingSigma; + } if (that.gmmSigmaScale != null) { this.gmmSigmaScale = that.gmmSigmaScale; } diff --git a/src/main/java/gov/usgs/earthquake/nshmp/calc/Site.java b/src/main/java/gov/usgs/earthquake/nshmp/calc/Site.java index 5f18d83417558face6c02cd857b30ff12f9b7123..8a90c3a02c56acef816033576804e603e481fe00 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/calc/Site.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/calc/Site.java @@ -29,13 +29,16 @@ import gov.usgs.earthquake.nshmp.gmm.GroundMotionModel; * <li><b>{@code z2.5}:</b> Depth to a shear wave velocity of 2.5 km/sec, in * units of km.</li> * + * <li><b>{@code zSed}:</b> Continental margin sediment thickness, in units of + * km.</li> + * * </ul> * * <p>Both {@code z1.0} and {@code z2.5}, collectively referred to as - * <i>basin-terms</i>, have default values of {@code Double.NaN}. When supplied - * with the default, those GMMs that support basin terms will use an author - * defined model, typically based on {@code Vs30}, to compute basin - * amplification (or demplification). + * <i>basin-terms</i>, and {@code zSed} have default values of + * {@code Double.NaN}. When supplied with the default, those GMMs that support + * basin terms will use an author defined model, typically based on + * {@code Vs30}, to compute basin amplification (or demplification). * * @author U.S. Geological Survey */ @@ -44,39 +47,49 @@ public class Site { /** The name used for a {@code Site} with no supplied name. */ public static final String NO_NAME = "Unnamed"; - /** Default {@link #vs30} value: {@code 760 m/s}. */ + /** Default {@link #vs30()} value: {@code 760 m/s}. */ public static final double VS_30_DEFAULT = 760.0; - /** Supported {@link #vs30} values: {@code [150..2000] m/s}. */ + /** Supported {@link #vs30()} values: {@code [150..2000] m/s}. */ public static final Range<Double> VS30_RANGE = Range.closed(150.0, 3000.0); - /** Default {@link #vsInferred} inferred value: {@code true}. */ + /** Default {@link #vsInferred()} inferred value: {@code true}. */ public static final boolean VS_INF_DEFAULT = true; /** - * Default {@link #z1p0} value: {@code NaN} <br>({@link GroundMotionModel}s + * Default {@link #z1p0()} value: {@code NaN} <br>({@link GroundMotionModel}s * will use a default value or model) */ public static final double Z1P0_DEFAULT = Double.NaN; - /** Supported {@link #z1p0} values: {@code [0..5] km}. */ + /** Supported {@link #z1p0()} values: {@code [0..5] km}. */ public static final Range<Double> Z1P0_RANGE = Range.closed(0.0, 5.0); /** - * Default {@link #z2p5} value: {@code NaN} <br>({@link GroundMotionModel}s + * Default {@link #z2p5()} value: {@code NaN} <br>({@link GroundMotionModel}s * will use a default value or model) */ public static final double Z2P5_DEFAULT = Double.NaN; - /** Supported {@link #z2p5} values: {@code [0..10] km}. */ + /** Supported {@link #z2p5()} values: {@code [0..10] km}. */ public static final Range<Double> Z2P5_RANGE = Range.closed(0.0, 10.0); + /** + * Default {@link #zSed()} value: {@code NaN} <br>({@link GroundMotionModel}s + * will use a default value or model) + */ + public static final double ZSED_DEFAULT = Double.NaN; + + /** Supported {@link #zSed()} values: {@code [0..20] km}. */ + public static final Range<Double> ZSED_RANGE = Range.closed(0.0, 20.0); + private final String name; private final Location location; private final double vs30; private final boolean vsInferred; private final double z1p0; private final double z2p5; + private final double zSed; private Site(Builder builder) { this.name = builder.name; @@ -85,6 +98,7 @@ public class Site { this.vsInferred = builder.vsInferred; this.z1p0 = builder.z1p0; this.z2p5 = builder.z2p5; + this.zSed = builder.zSed; } /** The site name. */ @@ -138,18 +152,31 @@ public class Site { return z2p5; } + /** + * Continental margin sediment thickness, in km. + * + * <p>Default: {@code NaN} <br>({@link GroundMotionModel}s will use a default + * value or model) + */ + public final double zSed() { + return zSed; + } + @Override public String toString() { return new StringBuilder(Strings.padEnd(name, 28, ' ')) .append(String.format("%.3f %.3f ", location.longitude, location.latitude)) - .append(paramString(vs30, vsInferred, z1p0, z2p5)) + .append(paramString(vs30, vsInferred, z1p0, z2p5, zSed)) .toString(); } - private static String paramString(double vs30, boolean vsInf, double z1p0, double z2p5) { + private static String paramString( + double vs30, boolean vsInf, + double z1p0, double z2p5, + double zSed) { return String.format( - "Vs30=%s %s Z1.0=%s Z2.5=%s", - vs30, vsInf ? "inferred " : "measured ", z1p0, z2p5); + "vs30=%s %s z1.0=%s z2.5=%s zSed=%s", + vs30, vsInf ? "inferred " : "measured ", z1p0, z2p5, zSed); } /** @@ -175,6 +202,7 @@ public class Site { private boolean vsInferred = VS_INF_DEFAULT; private double z1p0 = Z1P0_DEFAULT; private double z2p5 = Z2P5_DEFAULT; + private double zSed = ZSED_DEFAULT; boolean built = false; @@ -234,6 +262,16 @@ public class Site { return this; } + /** Sediment thickness, in km. */ + public Builder zSed(double zSed) { + if (Double.isNaN(zSed)) { + this.zSed = zSed; + } else { + this.zSed = checkInRange(ZSED_RANGE, Site.Key.ZSED, zSed); + } + return this; + } + /** * Build the {@code Site}. */ @@ -273,6 +311,8 @@ public class Site { public static final String Z1P0 = "z1p0"; /** The site z2.5 key. */ public static final String Z2P5 = "z2p5"; + /** The site zSed key. */ + public static final String ZSED = "zSed"; } } diff --git a/src/main/java/gov/usgs/earthquake/nshmp/calc/Sites.java b/src/main/java/gov/usgs/earthquake/nshmp/calc/Sites.java index a8ebc50d27b4a7182210440f37d3bb642ad3ab94..33c7c541af59b20bda6318b05365e7936ab93fa4 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/calc/Sites.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/calc/Sites.java @@ -118,10 +118,10 @@ public final class Sites { SiteData siteData) { Site.Builder builder = Site.builder().location(loc); - SiteData.Values sdValues = siteData.get(loc); - sdValues.z1p0.ifPresent(builder::z1p0); - sdValues.z2p5.ifPresent(builder::z2p5); - + SiteData.Values values = siteData.get(loc); + values.z1p0.ifPresent(builder::z1p0); + values.z2p5.ifPresent(builder::z2p5); + values.zSed.ifPresent(builder::zSed); return builder; } diff --git a/src/main/java/gov/usgs/earthquake/nshmp/gmm/GmmInput.java b/src/main/java/gov/usgs/earthquake/nshmp/gmm/GmmInput.java index c92df73fb256dcd5c3fe3384e6e44d627bfaf9ed..9f30ac78a1a0f3c6704e19c6d98df0edb881d7c3 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/gmm/GmmInput.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/gmm/GmmInput.java @@ -15,6 +15,7 @@ import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.WIDTH; import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.Z1P0; import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.Z2P5; import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.ZHYP; +import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.ZSED; import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.ZTOR; import static java.lang.Double.NaN; @@ -76,6 +77,8 @@ public class GmmInput { public final double z1p0; /** Depth to 2.5 km/s (in km). */ public final double z2p5; + /** Sediment thickness (in km). */ + public final double zSed; /** * Create a deterministic rupture and site property container with all @@ -112,6 +115,8 @@ public class GmmInput { * * <li>z2p5 depth to V<sub>s</sub>=2.5 km/sec (in km)</li></ul> * + * <li>zSed sediment thickness (in km)</li></ul> + * * @param builder The GmmInput builder */ protected GmmInput(Builder builder) { @@ -131,6 +136,7 @@ public class GmmInput { vsInf = builder.vsInf; z1p0 = builder.z1p0; z2p5 = builder.z2p5; + zSed = builder.zSed; } /** @@ -169,6 +175,7 @@ public class GmmInput { private boolean vsInf; private double z1p0; private double z2p5; + private double zSed; private Builder() {} @@ -210,9 +217,11 @@ public class GmmInput { * * <li>vsInf: true</li> * + * <li>z1p0: NaN ({@code null} when serialized)</li></ul> + * * <li>z2p5: NaN ({@code null} when serialized)</li> * - * <li>z1p0: NaN ({@code null} when serialized)</li></ul> + * <li>zSed: NaN ({@code null} when serialized)</li></ul> * * @throws IllegalStateException if any other builder method has already * been called without first calling {@link #build()} @@ -236,6 +245,7 @@ public class GmmInput { vsInf = model.vsInf; z1p0 = model.z1p0; z2p5 = model.z2p5; + zSed = model.zSed; flags.set(0, SIZE); return this; } @@ -272,10 +282,12 @@ public class GmmInput { return rake(v); case VS30: return vs30(v); - case Z2P5: - return z2p5(v); case Z1P0: return z1p0(v); + case Z2P5: + return z2p5(v); + case ZSED: + return zSed(v); } } catch (NumberFormatException nfe) { // move along @@ -368,14 +380,6 @@ public class GmmInput { return this; } - /** - * Set both the vs30 at site and whether {@code vs30} is inferred or - * measured. - */ - public Builder vs30(double vs30, boolean vsInf) { - return vs30(vs30).vsInf(vsInf); - } - /** Set the depth to 1.0 km/s (in km). */ public Builder z1p0(double z1p0) { this.z1p0 = validateAndFlag(Z1P0, z1p0); @@ -388,6 +392,12 @@ public class GmmInput { return this; } + /** Set the sediment thickness (in km). */ + public Builder zSed(double zSed) { + this.zSed = validateAndFlag(ZSED, zSed); + return this; + } + public GmmInput build() { checkState(flags.cardinality() == SIZE, "Not all fields set"); reset.clear(); @@ -431,6 +441,7 @@ public class GmmInput { .vsInf(VSINF.defaultValue > 0.0) .z1p0(Z1P0.defaultValue) .z2p5(Z2P5.defaultValue) + .zSed(ZSED.defaultValue) .build(); /** @@ -531,6 +542,13 @@ public class GmmInput { "Depth to Vs=2.5 km/s", "Depth to a shear-wave velocity of 2.5 kilometers per second, in kilometers", Optional.of(DISTANCE_UNIT), + NaN), + + ZSED( + "zSed", + "Sediment thickness", + "Thickness of coastal margin sediments, in kilometers", + Optional.of(DISTANCE_UNIT), NaN); public final String id; @@ -583,6 +601,7 @@ public class GmmInput { .add(VSINF.toString(), vsInf) .add(Z1P0.toString(), z1p0) .add(Z2P5.toString(), z2p5) + .add(ZSED.toString(), zSed) .toString(); } @@ -604,6 +623,9 @@ public class GmmInput { boolean z2p5Check = Double.isNaN(in.z2p5) ? Double.isNaN(this.z2p5) : this.z2p5 == in.z2p5; + boolean zSedCheck = Double.isNaN(in.zSed) + ? Double.isNaN(this.zSed) + : this.zSed == in.zSed; return this.Mw == in.Mw && this.rJB == in.rJB && this.rRup == in.rRup && @@ -616,13 +638,14 @@ public class GmmInput { this.vs30 == in.vs30 && this.vsInf == in.vsInf && z1p0Check && - z2p5Check; + z2p5Check && + zSedCheck; } @Override public int hashCode() { return Objects.hash(Mw, rJB, rRup, rX, dip, width, zTor, zHyp, - rake, vs30, vsInf, z1p0, z2p5); + rake, vs30, vsInf, z1p0, z2p5, zSed); } /** @@ -657,6 +680,7 @@ public class GmmInput { constraintMap.put(VSINF, builder.vsInf); constraintMap.put(Z1P0, builder.z1p0); constraintMap.put(Z2P5, builder.z2p5); + constraintMap.put(ZSED, builder.zSed); } public Optional<?> get(Field field) { @@ -741,6 +765,7 @@ public class GmmInput { private Optional<Range<Boolean>> vsInf = Optional.empty(); private Optional<Range<Double>> z1p0 = Optional.empty(); private Optional<Range<Double>> z2p5 = Optional.empty(); + private Optional<Range<Double>> zSed = Optional.empty(); /** * Set {@code Range<Double>} constraint. @@ -784,6 +809,9 @@ public class GmmInput { case Z2P5: z2p5 = Optional.of(constraint); break; + case ZSED: + zSed = Optional.of(constraint); + break; default: throw new IllegalArgumentException( "GmmInput.Constraints.Builder Unsupported field: " + id.name()); @@ -841,6 +869,7 @@ public class GmmInput { set(VSINF); set(Z1P0, Site.Z1P0_RANGE); set(Z2P5, Site.Z2P5_RANGE); + set(ZSED, Site.ZSED_RANGE); return this; } diff --git a/src/main/java/gov/usgs/earthquake/nshmp/gmm/NgaEastUsgs_2017.java b/src/main/java/gov/usgs/earthquake/nshmp/gmm/NgaEastUsgs_2017.java index b5796f4354399274aae11c3ecb006751c9289588..48005c04777fe5acca504c14c4dc358ea9c4041d 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/gmm/NgaEastUsgs_2017.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/gmm/NgaEastUsgs_2017.java @@ -409,9 +409,9 @@ public abstract class NgaEastUsgs_2017 implements GroundMotionModel { @Override public LogicTree<GroundMotion> calc(GmmInput in) { - double cpa = Double.isNaN(in.z2p5) + double cpa = Double.isNaN(in.zSed) ? 0.0 - : log(ChapmanGuo_2021.cpaPsaRatio(imt, in.z2p5, in.Mw, in.rJB)); + : log(ChapmanGuo_2021.cpaPsaRatio(imt, in.zSed, in.Mw, in.rJB)); Position p = tables[0].position(in.rRup, in.Mw); double[] μs = new double[MODEL_COUNT]; for (int i = 0; i < MODEL_COUNT; i++) { @@ -542,9 +542,9 @@ public abstract class NgaEastUsgs_2017 implements GroundMotionModel { @Override public LogicTree<GroundMotion> calc(GmmInput in) { GmmInput inRock = GmmInput.builder().fromCopy(in).vs30(3000).build(); - double cpa = Double.isNaN(in.z2p5) + double cpa = Double.isNaN(in.zSed) ? 0.0 - : log(ChapmanGuo_2021.cpaPsaRatio(imt, in.z2p5, in.Mw, in.rJB)); + : log(ChapmanGuo_2021.cpaPsaRatio(imt, in.zSed, in.Mw, in.rJB)); 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++) { diff --git a/src/main/java/gov/usgs/earthquake/nshmp/gmm/RezaeianDamping_2014.java b/src/main/java/gov/usgs/earthquake/nshmp/gmm/RezaeianDamping_2014.java index ef95722585f2c2f5299ab2311ad69764dbec9c22..212c4dfe281fa822a9d595d51f57acbd29937feb 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/gmm/RezaeianDamping_2014.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/gmm/RezaeianDamping_2014.java @@ -145,7 +145,7 @@ public class RezaeianDamping_2014 { RezaeianDamping_2014(CalcConfig config) { dampingRatio = config.hazard.gmmDampingRatio * 100.0; - updateSigma = false; // config.hazard.gmmDampingSigma; + updateSigma = config.hazard.gmmDampingSigma; } /** diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/ModelLoader.java b/src/main/java/gov/usgs/earthquake/nshmp/model/ModelLoader.java index 5a26ed125fad8871a9d0c2de61850b9573c516b4..2af32b8354966ce1ae66c10405589fd5d92c727d 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/ModelLoader.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/ModelLoader.java @@ -60,6 +60,7 @@ import java.util.logging.Logger; import com.google.common.collect.Iterables; import gov.usgs.earthquake.nshmp.Maths; +import gov.usgs.earthquake.nshmp.Text; import gov.usgs.earthquake.nshmp.calc.CalcConfig; import gov.usgs.earthquake.nshmp.fault.FocalMech; import gov.usgs.earthquake.nshmp.geo.json.Feature; @@ -98,12 +99,14 @@ abstract class ModelLoader { */ public static void main(String[] args) { + // Path testModel = Paths.get("../nshm-conus-2018-tmp"); Path testModel = Paths.get("../nshm-conus"); // Path testModel = Paths.get("../nshm-hawaii"); HazardModel model = ModelLoader.load(testModel); System.out.println(); System.out.println(model); + } private static Logger log; @@ -115,14 +118,17 @@ abstract class ModelLoader { static HazardModel load(Path dir) { try { - System.out.println(dir); + System.out.println(dir + Text.NEWLINE); CalcConfig config = readCalcConfig(dir); + SiteData siteData = readSiteData(dir); + System.out.println(siteData); + HazardModel.Builder builder = HazardModel.builder() .info(readModelInfo(dir).orElseThrow()) .config(config) - .siteData(readSiteData(dir)); + .siteData(siteData); /* Process tectonic setting source directories. */ loadTectonicSettings(dir, builder); diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/SiteData.java b/src/main/java/gov/usgs/earthquake/nshmp/model/SiteData.java index d0aaadf43789fdcc4f5c2ba4053a4a0bb22af89f..54356f7569b304769397c217fef4fe976474e06c 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/SiteData.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/SiteData.java @@ -1,5 +1,6 @@ package gov.usgs.earthquake.nshmp.model; +import static gov.usgs.earthquake.nshmp.Text.NEWLINE; import static gov.usgs.earthquake.nshmp.model.Deserialize.ID; import static gov.usgs.earthquake.nshmp.model.Deserialize.MODEL; import static gov.usgs.earthquake.nshmp.model.Deserialize.NAME; @@ -8,13 +9,18 @@ import static java.util.stream.Collectors.toUnmodifiableMap; import static java.util.stream.Collectors.toUnmodifiableSet; import java.awt.geom.Area; +import java.io.IOException; import java.math.BigDecimal; +import java.math.RoundingMode; +import java.nio.file.DirectoryStream; import java.nio.file.Files; import java.nio.file.Path; import java.util.Map; import java.util.OptionalDouble; import java.util.Set; +import java.util.function.Function; import java.util.stream.Stream; +import java.util.stream.StreamSupport; import gov.usgs.earthquake.nshmp.Maths; import gov.usgs.earthquake.nshmp.data.DelimitedData; @@ -22,15 +28,13 @@ import gov.usgs.earthquake.nshmp.data.DelimitedData.Record; import gov.usgs.earthquake.nshmp.geo.Location; import gov.usgs.earthquake.nshmp.geo.Locations; import gov.usgs.earthquake.nshmp.geo.json.Feature; -import gov.usgs.earthquake.nshmp.geo.json.FeatureCollection; import gov.usgs.earthquake.nshmp.geo.json.GeoJson; import gov.usgs.earthquake.nshmp.geo.json.Properties; /** * Site data wrapper class. This class provides one method {@code get(Location)} * that returns the superset of all site data parameters required by the current - * suite of supported NSHMs. Currently this only includes the {@code z1.0} and - * {@code z2.5} basin depth values. + * suite of supported NSHMs. * * This class assumes that if a location 'hit test' succeeds on an area for * which there are site data, then non null values will exist for the location. @@ -39,52 +43,89 @@ import gov.usgs.earthquake.nshmp.geo.json.Properties; */ public class SiteData { - /* - * Developer notes: This class will soon be expanded to support a variety of - * site data such as sediment thickness. Web services should map into the - * SiteData.Values class. We're only working with the two basin depth values - * at present, but ultimately we'll probably want - * - * Consider having this class deliver defaults. - */ - /** Empty site data instance. */ - public static final SiteData EMPTY = new SiteData(null) { + public static final SiteData EMPTY = new SiteData(null, null) { @Override public Values get(Location loc) { return Values.EMPTY; } }; - private static Path BASIN_FILE = Path.of("basin/basins.geojson"); + private static Path BASIN_DIR = Path.of("basin"); + private static Path MARGIN_DIR = Path.of("margin"); private Set<Basin> basins; + private Set<Margin> margins; - private SiteData(Set<Basin> basins) { + private SiteData(Set<Basin> basins, Set<Margin> margins) { this.basins = basins; + this.margins = margins; } static SiteData create(Path dir) { - Set<Basin> basins = Basin.init(dir); - return new SiteData(basins); + Set<Basin> basins = init(dir.resolve(BASIN_DIR), Basin::new); + Set<Margin> margins = init(dir.resolve(MARGIN_DIR), Margin::new); + return (basins.isEmpty() && margins.isEmpty()) + ? EMPTY + : new SiteData(basins, margins); } /** - * Return the site data values for the specified location. + * Return the site data values for the specified location. Once a polygon + * 'hit' ocurs, method breaks out of search loops. * * @param location of interest */ public Values get(Location location) { Values.Builder builder = new Values.Builder(); for (Basin basin : basins) { - if (basin.contains(location)) { - builder.basinValues(basin.get(location)); + Location snapped = snapToGrid(location, basin.spacing, basin.scale); + if (basin.contains(snapped)) { + builder.basinValues(basin.data.get(snapped)); + break; + } + } + for (Margin margin : margins) { + Location snapped = snapToGrid(location, margin.spacing, margin.scale); + if (margin.contains(snapped)) { + builder.marginValues(margin.data.get(snapped)); break; } } return builder.build(); } + /* Snap location to closest multiple of spacing. */ + private static Location snapToGrid(Location loc, double spacing, int scale) { + return Location.create( + snapToGrid(loc.longitude, spacing, scale), + snapToGrid(loc.latitude, spacing, scale)); + } + + /* + * Snap value to closest multiple of spacing. + * + * Double precision rounding errors results in grid midpoints rounding both up + * and down in the initial Math.round() call. + */ + private static double snapToGrid(double value, double spacing, int scale) { + return Maths.round( + Math.round(value / spacing) * spacing, + scale, RoundingMode.HALF_UP); + } + + @Override + public String toString() { + StringBuilder sb = new StringBuilder(ModelFiles.SITE_DATA_DIR).append("/"); + if (this == EMPTY) { + return sb.append(" (empty)").toString(); + } + sb.append(NEWLINE); + basins.forEach(basin -> sb.append(" basin: ").append(basin.name).append(NEWLINE)); + margins.forEach(margin -> sb.append(" margin: ").append(margin.name).append(NEWLINE)); + return sb.toString(); + } + /** Optional site data values associated with a location. */ public static class Values { @@ -96,15 +137,20 @@ public class SiteData { /** Depth to the shear-wave velocity horizon of 2.5 km/s, in km. */ public final OptionalDouble z2p5; + /** Cpntinental margin sediment thickness, in km. */ + public final OptionalDouble zSed; + private Values(Builder builder) { this.z1p0 = builder.z1p0; this.z2p5 = builder.z2p5; + this.zSed = builder.zSed; } private static class Builder { OptionalDouble z1p0 = OptionalDouble.empty(); OptionalDouble z2p5 = OptionalDouble.empty(); + OptionalDouble zSed = OptionalDouble.empty(); Builder basinValues(Basin.Values values) { this.z1p0 = OptionalDouble.of(values.z1p0); @@ -112,6 +158,11 @@ public class SiteData { return this; } + Builder marginValues(Margin.Values values) { + this.zSed = OptionalDouble.of(values.zSed); + return this; + } + Values build() { return new Values(this); } @@ -124,28 +175,43 @@ public class SiteData { record.getDouble("lat")); } - /* Container for a basin region. */ + private static <T> Set<T> init(Path dir, Function<Path, T> mapper) { + if (Files.exists(dir)) { + try (DirectoryStream<Path> ds = Files.newDirectoryStream(dir, "*.geojson")) { + return StreamSupport.stream(ds.spliterator(), false) + .map(mapper) + .collect(toUnmodifiableSet()); + } catch (IOException ioe) { + throw new RuntimeException(ioe); + } + } + return Set.of(); + } + + /* Container for a basin region and depth data. */ static final class Basin { final String name; final String id; final String model; + final double spacing; final int scale; final Feature feature; final Area area; - final Map<Location, Values> dataMap; + final Map<Location, Values> data; - private Basin(Path file, Feature feature) { + private Basin(Path geojson) { + Feature feature = GeoJson.from(geojson).toFeature(); Properties properties = feature.properties(); this.name = properties.getString(NAME).orElseThrow(); this.id = properties.getString(ID).orElseThrow(); this.model = properties.getString(MODEL).orElseThrow(); - double spacing = properties.getDouble(SPACING).orElseThrow(); + this.spacing = properties.getDouble(SPACING).orElseThrow(); this.scale = BigDecimal.valueOf(spacing).scale(); this.feature = feature; this.area = new Area(Locations.toPath(feature.asPolygonBorder())); - Path dataFile = file.resolveSibling(id + ".csv"); - this.dataMap = loadData(dataFile); + Path csvpath = geojson.resolveSibling(id + ".csv"); + this.data = loadData(csvpath); } /* Return whether the basin contains the supplied location. */ @@ -153,25 +219,6 @@ public class SiteData { return area.contains(loc.longitude, loc.latitude); } - Values get(Location loc) { - Location roundedLoc = Location.create( - Maths.round(loc.longitude, scale), - Maths.round(loc.latitude, scale)); - return dataMap.get(roundedLoc); - } - - static Set<Basin> init(Path dir) { - Path file = dir.resolve(BASIN_FILE); - if (Files.exists(dir)) { - FeatureCollection fc = GeoJson.from(file).toFeatureCollection(); - Set<Basin> basins = fc.features().stream() - .map(feature -> new Basin(file, feature)) - .collect(toUnmodifiableSet()); - return basins; - } - return Set.of(); - } - private Map<Location, Values> loadData(Path file) { DelimitedData csv = DelimitedData.comma(file); try (Stream<Record> stream = csv.records()) { @@ -199,4 +246,59 @@ public class SiteData { } } + /* Container for continental margin sediment thickness data. */ + static final class Margin { + + final String name; + final String id; + final String model; + final double spacing; + final int scale; + final Feature feature; + final Area area; + final Map<Location, Values> data; + + private Margin(Path geojson) { + Feature feature = GeoJson.from(geojson).toFeature(); + Properties properties = feature.properties(); + this.name = properties.getString(NAME).orElseThrow(); + this.id = properties.getString(ID).orElseThrow(); + this.model = properties.getString(MODEL).orElseThrow(); + this.spacing = properties.getDouble(SPACING).orElseThrow(); + this.scale = BigDecimal.valueOf(spacing).scale(); + this.feature = feature; + this.area = new Area(Locations.toPath(feature.asPolygonBorder())); + Path csvpath = geojson.resolveSibling(id + ".csv"); + this.data = loadData(csvpath); + } + + /* Return whether the basin contains the supplied location. */ + boolean contains(Location loc) { + return area.contains(loc.longitude, loc.latitude); + } + + private Map<Location, Values> loadData(Path file) { + DelimitedData csv = DelimitedData.comma(file); + try (Stream<Record> stream = csv.records()) { + return stream.collect(toUnmodifiableMap( + SiteData::toLocation, + Margin::toValues)); + } + } + + private static Values toValues(Record record) { + double value = Maths.round(record.getDouble("zsed") / 1000.0, 3); + return new Values(value); + } + + static final class Values { + + final double zSed; + + Values(double zSed) { + this.zSed = zSed; + } + } + } + } diff --git a/src/main/resources/calc/calc-config-defaults.json b/src/main/resources/calc/calc-config-defaults.json index 35ab5ca37a20f35b4d907bf3dc7d3cbfc40f52d1..5ce1beba5cfb9ee5531a783a49a12832bd29dfb9 100644 --- a/src/main/resources/calc/calc-config-defaults.json +++ b/src/main/resources/calc/calc-config-defaults.json @@ -12,7 +12,9 @@ "tectonicSettings": [], "sourceTypes": [], "vs30s": [], + "useSiteData": true, "gmmDampingRatio": 0.05, + "gmmDampingSigma": false, "gmmSigmaScale": 1.0, "valueFormat": "ANNUAL_RATE", "customImls": {} diff --git a/src/test/java/gov/usgs/earthquake/nshmp/calc/CalcConfigTests.java b/src/test/java/gov/usgs/earthquake/nshmp/calc/CalcConfigTests.java index bc606ae9a4f0e855dbc5f2470d7f50d5beb45ebb..68d7cf9415181d88f6d8dc5776131b5f26a3c94c 100644 --- a/src/test/java/gov/usgs/earthquake/nshmp/calc/CalcConfigTests.java +++ b/src/test/java/gov/usgs/earthquake/nshmp/calc/CalcConfigTests.java @@ -134,8 +134,10 @@ class CalcConfigTests { assertEquals(3.0, def.truncationLevel); assertEquals(IMTS, def.imts); assertEquals(EnumSet.noneOf(TectonicSetting.class), def.tectonicSettings); - assertEquals(EnumSet.noneOf(SourceType.class), def.tectonicSettings); + assertEquals(EnumSet.noneOf(SourceType.class), def.sourceTypes); + assertEquals(true, def.useSiteData); assertEquals(0.05, def.gmmDampingRatio); + assertEquals(false, def.gmmDampingSigma); assertEquals(ValueFormat.ANNUAL_RATE, def.valueFormat); Map<Imt, XySequence> imlMap = def.modelCurves(); @@ -256,7 +258,9 @@ class CalcConfigTests { assertEquals(EnumSet.of(ACTIVE_CRUST, SUBDUCTION), def.tectonicSettings); assertEquals(EnumSet.of(FAULT, ZONE, SLAB), def.sourceTypes); assertEquals(Set.of(260.0, 760.0), def.vs30s); + assertEquals(false, def.useSiteData); assertEquals(0.03, def.gmmDampingRatio); + assertEquals(true, def.gmmDampingSigma); assertEquals(0.85, def.gmmSigmaScale); assertEquals(ValueFormat.POISSON_PROBABILITY, def.valueFormat); @@ -271,7 +275,9 @@ class CalcConfigTests { assertEquals(EnumSet.noneOf(TectonicSetting.class), def.tectonicSettings); assertEquals(EnumSet.noneOf(SourceType.class), def.tectonicSettings); assertEquals(Set.of(), def.vs30s); + assertEquals(true, def.useSiteData); assertEquals(0.05, def.gmmDampingRatio); + assertEquals(false, def.gmmDampingSigma); assertEquals(1.0, def.gmmSigmaScale); assertEquals(ValueFormat.ANNUAL_RATE, def.valueFormat); } diff --git a/src/test/java/gov/usgs/earthquake/nshmp/calc/SiteTests.java b/src/test/java/gov/usgs/earthquake/nshmp/calc/SiteTests.java index f0283b2575f48fa6e9d7b91be0e28586d7dd4fa9..03bb085691f6463b1c15f4db31f013824839f333 100644 --- a/src/test/java/gov/usgs/earthquake/nshmp/calc/SiteTests.java +++ b/src/test/java/gov/usgs/earthquake/nshmp/calc/SiteTests.java @@ -16,11 +16,12 @@ class SiteTests { static final boolean VS_INF = false; static final double Z1P0 = 0.5; static final double Z2P5 = 3.0; + static final double ZSED = 7.5; static final String TO_STRING_VSINF_TRUE = - "Unnamed -118.000 34.000 Vs30=760.0 inferred Z1.0=NaN Z2.5=NaN"; + "Unnamed -118.000 34.000 vs30=760.0 inferred z1.0=NaN z2.5=NaN zSed=NaN"; static final String TO_STRING_VSINF_FALSE = - "TestSite -118.000 34.000 Vs30=530.0 measured Z1.0=0.5 Z2.5=3.0"; + "TestSite -118.000 34.000 vs30=530.0 measured z1.0=0.5 z2.5=3.0 zSed=7.5"; static Site.Builder siteBuilder() { return Site.builder() @@ -29,7 +30,8 @@ class SiteTests { .vs30(VS_30) .vsInferred(VS_INF) .z1p0(Z1P0) - .z2p5(Z2P5); + .z2p5(Z2P5) + .zSed(ZSED); } @Test @@ -43,6 +45,7 @@ class SiteTests { assertEquals(VS_INF, site1.vsInferred()); assertEquals(Z1P0, site1.z1p0()); assertEquals(Z2P5, site1.z2p5()); + assertEquals(ZSED, site1.zSed()); // check single use assertThrows(IllegalStateException.class, () -> { @@ -55,10 +58,12 @@ class SiteTests { Site site3 = siteBuilder() .z1p0(Double.NaN) .z2p5(Double.NaN) + .zSed(Double.NaN) .build(); assertEquals(Double.NaN, site3.z1p0()); assertEquals(Double.NaN, site3.z2p5()); + assertEquals(Double.NaN, site3.zSed()); // named location NamedLocation loc = new NamedLocation() { diff --git a/src/test/java/gov/usgs/earthquake/nshmp/calc/SitesTests.java b/src/test/java/gov/usgs/earthquake/nshmp/calc/SitesTests.java index e30b36cdf7f4222140d3468fc0955c5cf598c6ca..74e00ffdecebfc9c4bc2996f55ed760e378f4653 100644 --- a/src/test/java/gov/usgs/earthquake/nshmp/calc/SitesTests.java +++ b/src/test/java/gov/usgs/earthquake/nshmp/calc/SitesTests.java @@ -49,11 +49,11 @@ class SitesTests { static final String TO_STRING = "Sites [size=6]\n" + - " Hilo HI -155.060 19.700 Vs30=760.0 inferred Z1.0=NaN Z2.5=NaN\n" + - " Salt Lake City UT -111.900 40.750 Vs30=530.0 inferred Z1.0=0.3 Z2.5=3.0\n" + - " San Jose CA -121.900 37.350 Vs30=350.0 measured Z1.0=NaN Z2.5=NaN\n" + - " Seattle WA -122.300 47.600 Vs30=760.0 inferred Z1.0=1.2 Z2.5=6.7\n" + - " Memphis TN -90.050 35.150 Vs30=3000.0 inferred Z1.0=NaN Z2.5=NaN\n" + + " Hilo HI -155.060 19.700 vs30=760.0 inferred z1.0=NaN z2.5=NaN zSed=NaN\n" + + " Salt Lake City UT -111.900 40.750 vs30=530.0 inferred z1.0=0.3 z2.5=3.0 zSed=NaN\n" + + " San Jose CA -121.900 37.350 vs30=350.0 measured z1.0=NaN z2.5=NaN zSed=NaN\n" + + " Seattle WA -122.300 47.600 vs30=760.0 inferred z1.0=1.2 z2.5=6.7 zSed=NaN\n" + + " Memphis TN -90.050 35.150 vs30=3000.0 inferred z1.0=NaN z2.5=NaN zSed=NaN\n" + " ... and 1 more ..."; } diff --git a/src/test/java/gov/usgs/earthquake/nshmp/gmm/GmmInputTest.java b/src/test/java/gov/usgs/earthquake/nshmp/gmm/GmmInputTest.java index e425d6f8c630dc8da6de6542c56e3aec27cf2e53..74a4d0ce682a3cf6f7c692aa9b8ba783e96de41d 100644 --- a/src/test/java/gov/usgs/earthquake/nshmp/gmm/GmmInputTest.java +++ b/src/test/java/gov/usgs/earthquake/nshmp/gmm/GmmInputTest.java @@ -12,6 +12,7 @@ import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.WIDTH; import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.Z1P0; import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.Z2P5; import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.ZHYP; +import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.ZSED; import static gov.usgs.earthquake.nshmp.gmm.GmmInput.Field.ZTOR; import static org.junit.jupiter.api.Assertions.assertEquals; import static org.junit.jupiter.api.Assertions.assertFalse; @@ -46,6 +47,7 @@ class GmmInputTest { private static final double width = WIDTH.defaultValue; private static final double z1p0 = Z1P0.defaultValue; private static final double z2p5 = Z2P5.defaultValue; + private static final double zSed = ZSED.defaultValue; private static final double zHyp = ZHYP.defaultValue; private static final double zTor = ZTOR.defaultValue; @@ -64,6 +66,7 @@ class GmmInputTest { VALUES.put(WIDTH, Double.toString(width)); VALUES.put(Z1P0, Double.toString(z1p0)); VALUES.put(Z2P5, Double.toString(z2p5)); + VALUES.put(ZSED, Double.toString(zSed)); VALUES.put(ZHYP, Double.toString(zHyp)); VALUES.put(ZTOR, Double.toString(zTor)); } @@ -103,6 +106,7 @@ class GmmInputTest { assertEquals(gmmCopy.width, model.width, 0); assertEquals(gmmCopy.z1p0, model.z1p0, 0); assertEquals(gmmCopy.z2p5, model.z2p5, 0); + assertEquals(gmmCopy.zSed, model.zSed, 0); assertEquals(gmmCopy.zHyp, model.zHyp, 0); assertEquals(gmmCopy.zTor, model.zTor, 0); } @@ -124,6 +128,7 @@ class GmmInputTest { gmmBuilder.set(WIDTH, Double.toString(val)); gmmBuilder.set(Z1P0, Double.toString(val)); gmmBuilder.set(Z2P5, Double.toString(val)); + gmmBuilder.set(ZSED, Double.toString(val)); gmmBuilder.set(ZHYP, Double.toString(val)); gmmBuilder.set(ZTOR, Double.toString(val)); GmmInput gmm = gmmBuilder.build(); @@ -139,6 +144,7 @@ class GmmInputTest { assertEquals(gmm.width, val, 0); assertEquals(gmm.z1p0, val, 0); assertEquals(gmm.z2p5, val, 0); + assertEquals(gmm.zSed, val, 0); assertEquals(gmm.zHyp, val, 0); assertEquals(gmm.zTor, val, 0); @@ -163,6 +169,7 @@ class GmmInputTest { assertEquals(gmm.width, width, 0); assertEquals(gmm.z1p0, z1p0, 0); assertEquals(gmm.z2p5, z2p5, 0); + assertEquals(gmm.zSed, zSed, 0); assertEquals(gmm.zHyp, zHyp, 0); assertEquals(gmm.zTor, zTor, 0); } @@ -199,6 +206,7 @@ class GmmInputTest { .vsInf(vsInf) .z1p0(z1p0) .z2p5(z2p5) + .zSed(zSed) .build(); assertEquals(gmm.dip, dip, 0); @@ -212,6 +220,7 @@ class GmmInputTest { assertEquals(gmm.width, width, 0); assertEquals(gmm.z1p0, z1p0, 0); assertEquals(gmm.z2p5, z2p5, 0); + assertEquals(gmm.zSed, zSed, 0); assertEquals(gmm.zHyp, zHyp, 0); assertEquals(gmm.zTor, zTor, 0); } @@ -238,6 +247,7 @@ class GmmInputTest { assertFalse(gmm.equals(gmmBuilder.withDefaults().width(v).build())); assertFalse(gmm.equals(gmmBuilder.withDefaults().z1p0(v).build())); assertFalse(gmm.equals(gmmBuilder.withDefaults().z2p5(v).build())); + assertFalse(gmm.equals(gmmBuilder.withDefaults().zSed(v).build())); assertFalse(gmm.equals(gmmBuilder.withDefaults().zHyp(v).build())); assertFalse(gmm.equals(gmmBuilder.withDefaults().zTor(v).build())); } @@ -271,6 +281,7 @@ class GmmInputTest { assertTrue(toString.contains("width=")); assertTrue(toString.contains("z1p0=")); assertTrue(toString.contains("z2p5=")); + assertTrue(toString.contains("zSed=")); assertTrue(toString.contains("zHyp=")); assertTrue(toString.contains("zTor=")); assertEquals(toString, toStringCopy); diff --git a/src/test/java/gov/usgs/earthquake/nshmp/gmm/GmmTest.java b/src/test/java/gov/usgs/earthquake/nshmp/gmm/GmmTest.java index 13a1f3310a9e7cbb08381c20c4b91d83de3591a1..4c01320d9fde3b22ca87ea74d7b9eb0cc9564278 100644 --- a/src/test/java/gov/usgs/earthquake/nshmp/gmm/GmmTest.java +++ b/src/test/java/gov/usgs/earthquake/nshmp/gmm/GmmTest.java @@ -144,9 +144,11 @@ class GmmTest { .zTor(it.next()) .zHyp(it.next()) .rake(it.next()) - .vs30(it.next(), it.next() > 0.0) + .vs30(it.next()) + .vsInf(it.next() > 0.0) .z1p0(it.next()) .z2p5(it.next()) + .zSed(Double.NaN) .build(); } diff --git a/src/test/resources/calc/calc-config-extends.json b/src/test/resources/calc/calc-config-extends.json index 37f44268de927407d0d93cdcf780aae14a0eaa04..39989c7f9087e2aa3a6a35691199b57bf222056b 100644 --- a/src/test/resources/calc/calc-config-extends.json +++ b/src/test/resources/calc/calc-config-extends.json @@ -8,7 +8,9 @@ "tectonicSettings": [ "ACTIVE_CRUST", "SUBDUCTION" ], "sourceTypes": [ "FAULT", "ZONE", "SLAB" ], "vs30s": [ 260, 760 ], + "useSiteData": false, "gmmDampingRatio": 0.03, + "gmmDampingSigma": true, "gmmSigmaScale": 0.85, "valueFormat": "POISSON_PROBABILITY", "customImls": {