diff --git a/src/main/java/gov/usgs/earthquake/nshmp/calc/Transforms.java b/src/main/java/gov/usgs/earthquake/nshmp/calc/Transforms.java index 92120ea7292d98b464d046a95c632e88e245f3d6..7055195cf0e005e7dea32498dd281410845e4b36 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/calc/Transforms.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/calc/Transforms.java @@ -103,8 +103,8 @@ final class Transforms { site.z1p0, site.z2p5); hazardInputs.add(input); + // System.out.println(input); } - return hazardInputs; } } diff --git a/src/main/java/gov/usgs/earthquake/nshmp/fault/surface/ApproxGriddedSurface.java b/src/main/java/gov/usgs/earthquake/nshmp/fault/surface/ApproxGriddedSurface.java index c22ea52128a29c989bd0d925b59a51bc04d0a948..5d6e2e2a7cf2742b125a07c23b6e6cbaec920a25 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/fault/surface/ApproxGriddedSurface.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/fault/surface/ApproxGriddedSurface.java @@ -77,6 +77,7 @@ public class ApproxGriddedSurface extends AbstractGriddedSurface { public ApproxGriddedSurface(LocationList upperTrace, LocationList lowerTrace, double spacing) { + // TODO need to validate like builder methods in DefaultGriddedSurface strikeSpacing = spacing; dipSpacing = spacing; // sameGridSpacing = true; diff --git a/src/main/java/gov/usgs/earthquake/nshmp/geo/LocationList.java b/src/main/java/gov/usgs/earthquake/nshmp/geo/LocationList.java index 87100f6a219201ec6895c5cf1df1634eafe27347..fd0cec85474cd8fcccd830afb79f4757c0d1d933 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/geo/LocationList.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/geo/LocationList.java @@ -71,6 +71,14 @@ public interface LocationList extends List<Location> { prev = loc; } return sum; + + // TODO use this instead; check tests + // if (size() == 1) { + // return 0.0; + // } + // return IntStream.range(0, size() - 1) + // .mapToDouble(i -> Locations.horzDistanceFast(get(i), get(i + 1))) + // .sum(); } /** diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/Deserialize.java b/src/main/java/gov/usgs/earthquake/nshmp/model/Deserialize.java index 84b8b0145bfd73757ed0b44a82085eb347aecc2e..da667ad2167df8658bc39c5eb506c4ab4199a0d3 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/Deserialize.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/Deserialize.java @@ -122,6 +122,11 @@ class Deserialize { builder.aleatoryProperties(GSON.fromJson(aleaProps, AleatoryProperties.class)); } + JsonElement nshmBinModel = config.get(MfdConfig.NSHM_BIN_MODEL_KEY); + if (!nshmBinModel.isJsonNull()) { + builder.nshmBinModel(nshmBinModel.getAsBoolean()); + } + return builder.build(); } diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/FaultRuptureSet.java b/src/main/java/gov/usgs/earthquake/nshmp/model/FaultRuptureSet.java index f2d829c0b8a7f3a5e591c6f1796b447435a7a1f1..e20b78c07d5c760f3f6a0ce123ae581a2385bd7e 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/FaultRuptureSet.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/FaultRuptureSet.java @@ -105,6 +105,7 @@ public class FaultRuptureSet implements RuptureSet { public int size() { // TODO what should this be ?? // a value based on the size of the Mfds?? + // this can be zero if all rate branches are null return 1; } @@ -206,10 +207,10 @@ public class FaultRuptureSet implements RuptureSet { feature = createFeature(); - // TODO test length property rounding + MfdConfig mfdConfig = data.mfdConfig().orElseThrow(); /* Filter of unlikely combination of SINGLEs epistemic uncertainty. */ - checkEpistemic(data.mfdConfig().orElseThrow(), mfdPropsTree); + checkEpistemic(mfdConfig, mfdPropsTree); /* * Fault MFDs are constructed in a variety of ways. Single fault MFDs may @@ -218,7 +219,6 @@ public class FaultRuptureSet implements RuptureSet { * */ - MfdConfig mfdConfig = data.mfdConfig().orElseThrow(); Optional<LogicTree<MoRate>> moRateTree = Optional.empty(); /* Build mfd-tree from fault slip or recurrence rates. */ @@ -238,9 +238,17 @@ public class FaultRuptureSet implements RuptureSet { mfdPropsTree = updateMfdsForRecurrence(); } + // moRateTree.get().stream() + // .map(Branch::value) + // .forEach(System.out::println); + mfdTree = createMfdTree(mfdPropsTree, moRateTree, mfdConfig); mfdTotal = MfdTrees.reduceMfdTree(mfdTree); + // mfdTree.stream() + // .map(Branch::id) + // .forEach(System.out::println); + // for (Branch<Double> b : moRateTree) { // System.out.println(b.id() + " " + b.value()); // double eqMo = Earthquakes.magToMoment(feature.magnitude.getAsDouble()); @@ -366,10 +374,11 @@ public class FaultRuptureSet implements RuptureSet { for (Branch<Mfd.Properties> branch : mfdPropsTree) { switch (branch.value().type()) { case GR: - double mCutoff = data.mfdConfig().orElseThrow().minMagnitude.orElseThrow(); + MfdConfig mfdConfig = data.mfdConfig().orElseThrow(); + double mCutoff = mfdConfig.minMagnitude.orElseThrow(); propsTree.addBranch( branch.id(), - updateGr(branch.value(), Mw, mCutoff), + updateGr(branch.value(), Mw, mCutoff, mfdConfig.nshmBinModel), branch.weight()); break; case SINGLE: @@ -393,7 +402,8 @@ public class FaultRuptureSet implements RuptureSet { private static Mfd.Properties updateGr( Mfd.Properties props, double mMax, - double mCutoff) { + double mCutoff, + boolean nshmBinModel) { Mfd.Properties.GutenbergRichter grProps = props.getAsGr(); @@ -413,10 +423,11 @@ public class FaultRuptureSet implements RuptureSet { /* * USGS NSHM Gutenberg-Richter discretization algorithm. Cutoff for nm = 4 * in fltrate*.f is >0.3, but in reality M6.8 events have 4 GR magnitudes - * in input files. + * in input files. This applies to most WUS faults, however a few use the + * discretization specified in the MFD logic tree (e.g. Seattle) */ int nm = Rm > 0.7 ? 8 : Rm > 0.5 ? 6 : Rm >= 0.3 ? 4 : 2; - double Δm = Maths.round(Rm / nm, 6); + double Δm = nshmBinModel ? Maths.round(Rm / nm, 6) : grProps.Δm(); // System.out.println(Rm + " " + grProps.mMin + " " + mMax + " " + Δm); return new Mfd.Properties.GutenbergRichter(1.0, grProps.b(), Δm, grProps.mMin(), mMax); @@ -485,7 +496,6 @@ public class FaultRuptureSet implements RuptureSet { } } } - // /* * Container class for slip-rate based moRate. @@ -581,10 +591,11 @@ public class FaultRuptureSet implements RuptureSet { RateType rateType = feature.rateType; LogicTree.Builder<MoRate> moRateTree = LogicTree.builder(RATE_TREE); - // TODO clean - if (rateMap.containsValue(null)) { - System.out.println(" " + rateMap); - } + // // TODO clean + // System.out.println(rateMap); + // if (rateMap.containsValue(null)) { + // System.out.println(" " + rateMap); + // } for (Branch<String> branch : data.faultRateTree().orElseThrow()) { @@ -661,6 +672,7 @@ public class FaultRuptureSet implements RuptureSet { * for M7, M + epi 7.2 as upper limit of sequence). */ double mMaxEpi = gr.mMax() + epiBranch.value() + gr.Δm() / 2; + // System.out.println(gr.Δm() + " " + mMaxEpi); double weightEpi = mfdWt * epiBranch.weight(); // TODO check if ever thrown @@ -681,9 +693,6 @@ public class FaultRuptureSet implements RuptureSet { .scaleToMomentRate(epiOk ? moRate : 0) .build(); - // Mfd mfd = newGrBuilder(gr, mMaxEpi) - // .scaleToMomentRate(epiOk ? moRate : 0) - // .build(); mfdTree.addBranch(id, mfd, weightEpi); } } else { diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/FaultSource.java b/src/main/java/gov/usgs/earthquake/nshmp/model/FaultSource.java index acfa83cf6d513f90420cd0e7010f2da72162df33..bac217fe8ab60f1b8732716491cb1762a7c73fcc 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/FaultSource.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/FaultSource.java @@ -3,13 +3,7 @@ package gov.usgs.earthquake.nshmp.model; import static com.google.common.base.Preconditions.checkArgument; import static com.google.common.base.Preconditions.checkNotNull; import static com.google.common.base.Preconditions.checkState; -import static gov.usgs.earthquake.nshmp.Earthquakes.checkCrustalDepth; -import static gov.usgs.earthquake.nshmp.Earthquakes.checkCrustalWidth; -import static gov.usgs.earthquake.nshmp.Faults.checkDip; -import static gov.usgs.earthquake.nshmp.Faults.checkRake; -import static gov.usgs.earthquake.nshmp.Faults.checkTrace; -import static gov.usgs.earthquake.nshmp.Text.checkName; -import static gov.usgs.earthquake.nshmp.data.DoubleData.checkInRange; +import static gov.usgs.earthquake.nshmp.Earthquakes.MAG_RANGE; import static gov.usgs.earthquake.nshmp.model.SourceType.FAULT; import static java.util.stream.Collectors.toUnmodifiableList; @@ -19,15 +13,12 @@ import java.util.List; import java.util.Map; import com.google.common.collect.Iterables; -import com.google.common.collect.Range; import gov.usgs.earthquake.nshmp.data.XyPoint; import gov.usgs.earthquake.nshmp.fault.surface.DefaultGriddedSurface; import gov.usgs.earthquake.nshmp.fault.surface.GriddedSurface; import gov.usgs.earthquake.nshmp.fault.surface.RuptureFloating; -import gov.usgs.earthquake.nshmp.fault.surface.RuptureScaling; import gov.usgs.earthquake.nshmp.geo.Location; -import gov.usgs.earthquake.nshmp.geo.LocationList; import gov.usgs.earthquake.nshmp.geo.Locations; import gov.usgs.earthquake.nshmp.mfd.Mfd; import gov.usgs.earthquake.nshmp.tree.Branch; @@ -48,46 +39,48 @@ import gov.usgs.earthquake.nshmp.tree.LogicTree; */ public class FaultSource implements Source { - final String name; - final int id; - final LocationList trace; - final double length; - final double dip; - final double width; - final double rake; - final LogicTree<Mfd> mfdTree; - final double spacing; - final RuptureScaling rupScaling; - final RuptureFloating rupFloating; + final SourceFeature.NshmFault feature; + final SourceConfig.Fault config; final GriddedSurface surface; - final boolean floatSingleRuptures; + + final LogicTree<Mfd> mfdTree; final double branchWeight; // rupture set branch weight private final List<Mfd> scaledMfds; // includes ruptureSet weight private final List<List<Rupture>> ruptureLists; // 1:1 with Mfds FaultSource(Builder builder) { - this.name = builder.name; - this.id = builder.id; - this.trace = builder.trace; - this.length = builder.length; - this.dip = builder.dip; - this.width = builder.width; - this.surface = builder.surface; - this.rake = builder.rake; + this.feature = builder.feature; + this.config = builder.config; this.mfdTree = builder.mfdTree; - this.spacing = builder.spacing; - this.rupScaling = builder.rupScaling; - this.rupFloating = builder.rupFloating; - this.floatSingleRuptures = builder.floatSingleRuptures; this.branchWeight = builder.branchWeight; + this.surface = DefaultGriddedSurface.builder() + .trace(feature.trace) + .depth(feature.upperDepth) + .dip(feature.dip) + .width(feature.width) + .spacing(config.surfaceSpacing) + .build(); + this.scaledMfds = MfdTrees.scaledMfdList(mfdTree, branchWeight); // System.out.println(mfdTree); this.ruptureLists = scaledMfds.stream() .map(this::createRuptureList) .collect(toUnmodifiableList()); + // scaledMfds.stream() + // .forEach(System.out::println); + + boolean noRuptures = Iterables.size(Iterables.concat(ruptureLists)) == 0; + if (noRuptures) { + System.out.println(" No ruptures: " + feature.name); + } + + // TODO refactor to use size() + // checkState(Iterables.size(Iterables.concat(ruptureLists)) > 0, + // "FaultSource has no ruptures"); + // for (Mfd mfd : scaledMfds) { // System.out.println(mfd.data()); // } @@ -101,8 +94,6 @@ public class FaultSource implements Source { // branch) TODO at least for 2014 the MFD tree and total MFD are identical // across dip branches - checkState(Iterables.size(Iterables.concat(ruptureLists)) > 0, - "FaultSource has no ruptures"); // System.out.println(this); // System.out.println(scaledMfds); @@ -116,17 +107,12 @@ public class FaultSource implements Source { @Override public String name() { - return name; - } - - @Override - public int size() { - return Iterables.size(this); + return feature.name; } @Override public int id() { - return id; + return feature.id; } @Override @@ -134,13 +120,20 @@ public class FaultSource implements Source { return FAULT; } + @Override + public int size() { + return ruptureLists.stream() + .mapToInt(List::size) + .sum(); + } + /** * The closest point on the fault trace, relative to the supplied site * {@code Location}. */ @Override public Location location(Location site) { - return Locations.closestPoint(site, trace); + return Locations.closestPoint(site, feature.trace); } @Override @@ -158,16 +151,16 @@ public class FaultSource implements Source { @Override public String toString() { Map<String, ?> data = Map.of( - "name", name, - "dip", dip, - "width", width, - "rake", rake, + "name", feature.name, + "dip", feature.dip, + "width", feature.width, + "rake", feature.rake, "mfds", mfdTree.size(), - "top", trace.first().depth); + "top", feature.trace.first().depth); return getClass().getSimpleName() + " " + data; } - List<Rupture> createRuptureList(Mfd mfd) { + private List<Rupture> createRuptureList(Mfd mfd) { List<Rupture> ruptureList = new ArrayList<>(); @@ -214,13 +207,40 @@ public class FaultSource implements Source { * 7/30/2020 can specify based on source type whether to float single * ruptures; ruptures will always float if GR */ - // Mfd.Type type = props.type(); - double length = trace.length(); + + // double mMaxGlobal = Earthquakes.MAG_RANGE.upperEndpoint(); + + double magRef = feature.magnitude.orElse(MAG_RANGE.upperEndpoint()); + + double refLength = config.ruptureScaling.dimensions(magRef, feature.width).length; + double traceLength = feature.trace.length(); + // double mRef = feature.magnitude.orElse(Double.MAX_VALUE); + + // System.out.println(" mRef: " + magRef); + // System.out.println(" refLength: " + refLength); + // System.out.println("traceLength: " + traceLength); Mfd.Type mfdType = mfd.properties().type(); - boolean floatSingle = mfdType == Mfd.Type.SINGLE && floatSingleRuptures; - boolean floats = (mfdType == Mfd.Type.GR) || floatSingle; + // if (feature.magnitude.isPresent()) { + // double mLength = config.ruptureScaling.dimensions(feature.magnitude, + // feature.width).length; + // } + // TODO if the original SINGLE magnitude w/o any uncertainty is shorter than + // than the trace length, then we float + + // System.out.println(feature.magnitude.); + + // System.out.println("Mw: " + magRef); + + /* + * First float check; use reference length from source to determine if true + * for SINGLE MFDs, otherwise epistemic branches end up floating. + */ + boolean otherNoFloat = feature.magnitudeNote.isPresent() && + feature.magnitudeNote.orElseThrow().equals(MagnitudeNote.OTHER_NO_FLOAT); + boolean floats = (mfdType == Mfd.Type.GR) || + (mfdType == Mfd.Type.SINGLE && refLength < traceLength && !otherNoFloat); for (XyPoint xy : mfd.data()) { double mag = xy.x(); @@ -238,26 +258,42 @@ public class FaultSource implements Source { continue; // shortcut low rates } - double msrLength = rupScaling.dimensions(mag, surface.width()).length; + // TODO should width here pull directly from feature.width + double mMfdLength = config.ruptureScaling.dimensions(mag, surface.width()).length; // TODO we want to get the 'floats' attribute out of MFDs // the only reason it is there is to allow SINGLE to flip-flop // it should just be a SourceProperty - if (floats && msrLength < length) { - // TODO change to using RuptureFloating.OFF - - // AbstractGriddedSurface surf = (AbstractGriddedSurface) - // surface; - // - List<Rupture> floaters = rupFloating.createFloatingRuptures( - surface, rupScaling, mag, rate, rake); - ruptureList.addAll(floaters); - - } else { - Rupture rup = Rupture.create(mag, rate, rake, surface); - ruptureList.add(rup); - } + RuptureFloating floating = (floats && mMfdLength < traceLength) + ? config.ruptureFloating + : RuptureFloating.OFF; + + List<Rupture> floaters = floating.createFloatingRuptures( + surface, + config.ruptureScaling, + mag, + rate, + feature.rake); + + ruptureList.addAll(floaters); + + // if (floats && msrLength < traceLength) { + // + // // TODO change to using RuptureFloating.OFF + // + // List<Rupture> floaters = config.ruptureFloating.createFloatingRuptures( + // surface, + // config.ruptureScaling, + // mag, + // rate, + // feature.rake); + // ruptureList.addAll(floaters); + // + // } else { + // Rupture rup = Rupture.create(mag, rate, feature.rake, surface); + // ruptureList.add(rup); + // } } // We now have rates of zero for null branches that preserve @@ -273,56 +309,21 @@ public class FaultSource implements Source { private static final String ID = "FaultSource.Builder"; private boolean built = false; - private static final Range<Double> SURFACE_GRID_SPACING_RANGE = Range.closed(0.01, 20.0); - - // required - String name; - Integer id; - LocationList trace; - Double dip; - Double width; - Double depth; - Double rake; + SourceFeature.NshmFault feature; + SourceConfig.Fault config; LogicTree<Mfd> mfdTree; Mfd mfdTotal; - Double spacing; - RuptureScaling rupScaling; - RuptureFloating rupFloating; - Boolean floatSingleRuptures = false; - Double branchWeight; - Builder name(String name) { - this.name = checkName(name, "FaultSource"); - return this; - } - - Builder id(int id) { - this.id = id; - return this; - } - - Builder trace(LocationList trace) { - this.trace = checkTrace(trace); - return this; - } - - Builder dip(double dip) { - this.dip = checkDip(dip); - return this; - } - - Builder width(double width) { - this.width = checkCrustalWidth(width); - return this; - } + Boolean floatSingleRuptures; + Double branchWeight; - Builder depth(double depth) { - this.depth = checkCrustalDepth(depth); + Builder feature(SourceFeature.NshmFault feature) { + this.feature = feature; return this; } - Builder rake(double rake) { - this.rake = checkRake(rake); + Builder config(SourceConfig.Fault config) { + this.config = config; return this; } @@ -338,69 +339,20 @@ public class FaultSource implements Source { return this; } - Builder surfaceSpacing(double spacing) { - this.spacing = checkInRange(SURFACE_GRID_SPACING_RANGE, "Floater Offset", spacing); - return this; - } - - Builder ruptureScaling(RuptureScaling rupScaling) { - this.rupScaling = checkNotNull(rupScaling, "Rup-Scaling Relation is null"); - return this; - } - - Builder ruptureFloating(RuptureFloating rupFloating) { - this.rupFloating = checkNotNull(rupFloating, "Rup-Floating Model is null"); - return this; - } - - Builder floatSingleRuptures(boolean floatSingleRuptures) { - this.floatSingleRuptures = floatSingleRuptures; - return this; - } - Builder branchWeight(double weight) { this.branchWeight = weight; return this; } - void validateState(String buildId) { - checkState(!built, "This %s instance as already been used", buildId); - checkState(name != null, "%s name not set", buildId); - checkState(id != null, "%s id not set", buildId); - checkState(trace != null, "%s trace not set", buildId); - checkState(dip != null, "%s dip not set", buildId); - checkState(width != null, "%s width not set", buildId); - checkState(depth != null, "%s depth not set", buildId); - checkState(rake != null, "%s rake not set", buildId); - checkState(mfdTree != null, "%s MFD tree not set", buildId); - checkState(mfdTotal != null, "%s total MFD not set", buildId); - checkState(spacing != null, "%s surface grid spacing not set", buildId); - checkState(rupScaling != null, "%s rupture-scaling relation not set", buildId); - checkState(rupFloating != null, "%s rupture-floating model not set", buildId); - checkState(floatSingleRuptures != null, "%s float uptures flag not set", buildId); - checkState(branchWeight != null, "%s source branch weight not set", buildId); + FaultSource build() { + checkState(!built, "This %s instance as already been used", ID); + checkState(feature != null, "%s feature not set", ID); + checkState(config != null, "%s config not set", ID); + checkState(mfdTree != null, "%s MFD tree not set", ID); + checkState(mfdTotal != null, "%s total MFD not set", ID); + checkState(branchWeight != null, "%s source branch weight not set", ID); built = true; - } - - GriddedSurface surface; - double length; - - FaultSource buildFaultSource() { - - validateState(ID); - - surface = DefaultGriddedSurface.builder() - .trace(trace) - .depth(depth) - .dip(dip) - .width(width) - .spacing(spacing) - .build(); - - length = trace.length(); - return new FaultSource(this); } } - } diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/FaultSourceSet.java b/src/main/java/gov/usgs/earthquake/nshmp/model/FaultSourceSet.java index fcee12d4775d57973d26b9147b50ec68db7ed4c4..36e500a0dd0716bd2b256e83aba04a36c339147b 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/FaultSourceSet.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/FaultSourceSet.java @@ -69,7 +69,8 @@ public class FaultSourceSet extends AbstractSourceSet<FaultSource> { @Override public boolean test(FaultSource source) { - return filter.test(source.trace.first()) || filter.test(source.trace.last()); + return filter.test(source.feature.trace.first()) || + filter.test(source.feature.trace.last()); } @Override diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/HazardModel.java b/src/main/java/gov/usgs/earthquake/nshmp/model/HazardModel.java index 6cf6ee8c073ec0e12eb95adf0f3a04ef77e426bd..f0b1368df8e00254f45ed731306438b0a5b53cd3 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/HazardModel.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/HazardModel.java @@ -19,7 +19,6 @@ import com.google.common.collect.Streams; import gov.usgs.earthquake.nshmp.Maths; import gov.usgs.earthquake.nshmp.calc.CalcConfig; -import gov.usgs.earthquake.nshmp.geo.LocationList; import gov.usgs.earthquake.nshmp.gmm.Gmm; import gov.usgs.earthquake.nshmp.gmm.GroundMotionModel; import gov.usgs.earthquake.nshmp.model.SourceTree.Leaf; @@ -251,8 +250,6 @@ public final class HazardModel implements Iterable<SourceSet<? extends Source>> if (tree.type() == SourceType.FAULT_SYSTEM) { for (Leaf leaf : tree.branchMap().keySet()) { RuptureSet rs = leaf.ruptureSet; - // System.out.println("Flt to SystemSrcSet: " + rs.type() + " " + - // rs.name()); double branchWeight = tree.branchWeights().get(leaf); SystemRuptureSet srs = (SystemRuptureSet) rs; SystemSourceSet sss = srs.createSourceSet(branchWeight); @@ -282,14 +279,25 @@ public final class HazardModel implements Iterable<SourceSet<? extends Source>> * TODO we're overwriting the id with the id of each rupture set leaf * for dip variants this is ok, but there are many other sources where * this WILL NOT hold; rupture-sets commonly have unique IDs that are - * distinct from the (poissibly multiple) underlying fault sections. + * distinct from the (possibly multiple) underlying fault sections. */ fss.id(frs.id()); // TODO leaf weight is being applied to total MFD - fss.source(faultRuptureSetToSource(frs, leafWeight), leafWeight); + FaultSource faultSource = faultRuptureSetToSource(frs, leafWeight); + // System.out.println("FS size: " + faultSource.size()); + + if (faultSource.size() > 0) { + fss.source(faultSource, leafWeight); + } else { + System.out.println(" No source ruptures: " + faultSource.name()); + } + } + if (fss.sources.size() > 0) { + addSourceSet(fss.buildFaultSet()); + } else { + System.out.println(" Empty Source Set: " + fss.name); } - addSourceSet(fss.buildFaultSet()); return; } @@ -356,24 +364,13 @@ public final class HazardModel implements Iterable<SourceSet<? extends Source>> FaultRuptureSet frs, double leafWeight) { - SourceFeature.NshmFault feature = frs.feature; - SourceConfig.Fault config = frs.data.sourceConfig().asFault(); - return new FaultSource.Builder() - .name(feature.name) - .id(feature.id) - .trace(feature.trace) - .dip(feature.dip) - .width(feature.width) - .depth(feature.upperDepth) - .rake(feature.rake) + .feature(frs.feature) + .config(frs.data.sourceConfig().asFault()) .mfdTree(frs.mfdTree) .mfdTotal(frs.mfdTotal) - .surfaceSpacing(config.surfaceSpacing) - .ruptureScaling(config.ruptureScaling) - .ruptureFloating(config.ruptureFloating) .branchWeight(leafWeight) - .buildFaultSource(); + .build(); } private void interfaceSourceSetsFromTree(SourceTree tree) { @@ -393,26 +390,17 @@ public final class HazardModel implements Iterable<SourceSet<? extends Source>> addSourceSet(builder.buildSubductionSet()); } - private InterfaceSource interfaceRuptureSetToSource(InterfaceRuptureSet irs, double weight) { - SourceFeature.Interface feature = irs.feature; - List<LocationList> traces = feature.traces; - SourceConfig.Interface config = irs.data.sourceConfig().asInterface(); - - InterfaceSource.Builder isb = new InterfaceSource.Builder() - .lowerTrace(traces.get(traces.size() - 1)); + private InterfaceSource interfaceRuptureSetToSource( + InterfaceRuptureSet irs, + double weight) { - isb.name(feature.name) - .id(feature.id) - .trace(traces.get(0)) - .rake(90.0) + return new InterfaceSource.Builder() + .feature(irs.feature) + .config(irs.data.sourceConfig().asInterface()) .mfdTree(irs.mfdTree) .mfdTotal(irs.mfdTotal) - .surfaceSpacing(config.surfaceSpacing) - .ruptureScaling(config.ruptureScaling) - .ruptureFloating(config.ruptureFloating) - .branchWeight(weight); - - return isb.buildSubductionSource(); + .branchWeight(weight) + .build(); } private void zoneSourceSetsFromTree(SourceTree tree) { diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/InterfaceSource.java b/src/main/java/gov/usgs/earthquake/nshmp/model/InterfaceSource.java index 6edec9af507ff93a0ccdd835a6335dadd8425d09..1f2d76aea44fb21d0f9109d161672e6d118d3f6d 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/InterfaceSource.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/InterfaceSource.java @@ -1,22 +1,28 @@ package gov.usgs.earthquake.nshmp.model; +import static com.google.common.base.Preconditions.checkArgument; +import static com.google.common.base.Preconditions.checkNotNull; +import static com.google.common.base.Preconditions.checkState; import static com.google.common.base.StandardSystemProperty.LINE_SEPARATOR; -import static gov.usgs.earthquake.nshmp.Earthquakes.checkInterfaceDepth; -import static gov.usgs.earthquake.nshmp.Earthquakes.checkInterfaceWidth; -import static gov.usgs.earthquake.nshmp.Faults.checkTrace; import static gov.usgs.earthquake.nshmp.model.SourceType.INTERFACE; +import static java.util.stream.Collectors.toUnmodifiableList; import java.util.ArrayList; +import java.util.Iterator; import java.util.List; +import com.google.common.collect.Iterables; + import gov.usgs.earthquake.nshmp.data.XyPoint; import gov.usgs.earthquake.nshmp.fault.surface.ApproxGriddedSurface; -import gov.usgs.earthquake.nshmp.fault.surface.DefaultGriddedSurface; +import gov.usgs.earthquake.nshmp.fault.surface.GriddedSurface; import gov.usgs.earthquake.nshmp.fault.surface.RuptureFloating; import gov.usgs.earthquake.nshmp.geo.Location; import gov.usgs.earthquake.nshmp.geo.LocationList; import gov.usgs.earthquake.nshmp.geo.Locations; import gov.usgs.earthquake.nshmp.mfd.Mfd; +import gov.usgs.earthquake.nshmp.tree.Branch; +import gov.usgs.earthquake.nshmp.tree.LogicTree; /** * Subduction source representation. This class wraps a model of a subduction @@ -31,18 +37,34 @@ import gov.usgs.earthquake.nshmp.mfd.Mfd; * * @author U.S. Geological Survey */ -public class InterfaceSource extends FaultSource { +public class InterfaceSource implements Source { + + final SourceFeature.Interface feature; + final SourceConfig.Interface config; + final GriddedSurface surface; + // temp TODO remove + final LocationList upperTrace; final LocationList lowerTrace; - private InterfaceSource(Builder builder) { + final LogicTree<Mfd> mfdTree; - super(builder); + final double branchWeight; // rupture set branch weight + + private final List<Mfd> scaledMfds; // includes ruptureSet weight + private final List<List<Rupture>> ruptureLists; // 1:1 with Mfds + + private InterfaceSource(Builder builder) { + this.feature = builder.feature; + this.config = builder.config; - this.lowerTrace = (builder.lowerTrace == null) - ? surface.getEvenlyDiscritizedLowerEdge() - : builder.lowerTrace; + this.mfdTree = builder.mfdTree; + this.branchWeight = builder.branchWeight; + // this.lowerTrace = (builder.lowerTrace == null) + // ? surface.getEvenlyDiscritizedLowerEdge() + // : builder.lowerTrace; + // // TODO lowerTrace may be null and this is bad bad; lowerTrace // is referenced in InterfaceSourceSet distanceFilter and // we should populate this even if the original source only @@ -54,6 +76,33 @@ public class InterfaceSource extends FaultSource { // lower trace in distance filter, but given large width of interface // sources TODO clean up Container2D methods + // TODO I don't think the above notes apply, see IntSourceSet distance + // filter + + upperTrace = feature.traces.get(0); + lowerTrace = feature.traces.get(feature.traces.size() - 1); + + surface = new ApproxGriddedSurface(upperTrace, lowerTrace, config.surfaceSpacing); + + this.scaledMfds = MfdTrees.scaledMfdList(mfdTree, branchWeight); + this.ruptureLists = scaledMfds.stream() + .map(this::createRuptureList) + .collect(toUnmodifiableList()); + + // TODO refactor to use size() + checkState(Iterables.size(Iterables.concat(ruptureLists)) > 0, + "FaultSource has no ruptures"); + + } + + @Override + public String name() { + return feature.name; + } + + @Override + public int id() { + return feature.id; } @Override @@ -61,6 +110,13 @@ public class InterfaceSource extends FaultSource { return INTERFACE; } + @Override + public int size() { + return ruptureLists.stream() + .mapToInt(List::size) + .sum(); + } + /** * The closest point on the upper or lower fault trace, relative to the * supplied site {@code Location}. @@ -69,11 +125,28 @@ public class InterfaceSource extends FaultSource { public Location location(Location site) { return Locations.closestPoint(site, LocationList.builder() - .addAll(trace) + .addAll(upperTrace) .addAll(lowerTrace) .build()); } + @Override + public List<Mfd> mfds() { + return mfdTree.stream() + .map(Branch::value) + .collect(toUnmodifiableList()); + } + + @Override + public Iterator<Rupture> iterator() { + return Iterables.concat(ruptureLists).iterator(); + } + + // TODO we do want to be able to return width and everage dip properties of an + // interface source; checkInterfaceDepth checkInterfaceWidth should + // be used on deserialize, perhaps to SourceFeature.NshmInterface with + // additional derived fields not in the interface-fault-section + @Override public String toString() { // TODO use joiner @@ -81,45 +154,61 @@ public class InterfaceSource extends FaultSource { .append("========== Interface Source ==========") .append(LINE_SEPARATOR.value()) .append(" Source name: ").append(name()) - .append(LINE_SEPARATOR.value()) - .append(" dip: ").append(dip) - .append(LINE_SEPARATOR.value()) - .append(" width: ").append(width) - .append(LINE_SEPARATOR.value()) - .append(" rake: ").append(rake) + // .append(LINE_SEPARATOR.value()) + // .append(" dip: ").append(dip) + // .append(LINE_SEPARATOR.value()) + // .append(" width: ").append(width) + // .append(LINE_SEPARATOR.value()) + // .append(" rake: ").append(rake) .append(LINE_SEPARATOR.value()) .append(" mfds: ").append(mfdTree.size()) .append(LINE_SEPARATOR.value()) - .append(" top: ").append(trace.first().depth) + .append(" top: ").append(feature.traces.get(0).first().depth) .append(LINE_SEPARATOR.value()).toString(); } /* Single use builder. */ - static class Builder extends FaultSource.Builder { + static class Builder { private static final String ID = "InterfaceSource.Builder"; + private boolean built = false; + + SourceFeature.Interface feature; + SourceConfig.Interface config; - // required - private LocationList lowerTrace; + LogicTree<Mfd> mfdTree; + Mfd mfdTotal; - @Override - Builder depth(double depth) { - this.depth = checkInterfaceDepth(depth); + Double branchWeight; + + Builder feature(SourceFeature.Interface feature) { + this.feature = feature; + return this; + } + + Builder config(SourceConfig.Interface config) { + this.config = config; + return this; + } + + Builder mfdTree(LogicTree<Mfd> mfdTree) { + checkArgument(checkNotNull(mfdTree).size() > 0, "MFD list is empty"); + this.mfdTree = mfdTree; return this; } - @Override - Builder width(double width) { - this.width = checkInterfaceWidth(width); + Builder mfdTotal(Mfd mfdTotal) { + checkNotNull(mfdTotal, "MFD total is null"); + this.mfdTotal = mfdTotal; return this; } - Builder lowerTrace(LocationList trace) { - this.lowerTrace = checkTrace(trace); + Builder branchWeight(double weight) { + this.branchWeight = weight; return this; } - InterfaceSource buildSubductionSource() { + InterfaceSource build() { /* * Either upper and lower trace must be set, or upper trace depth, dip, @@ -128,41 +217,25 @@ public class InterfaceSource extends FaultSource { * surface implementation; otherwise the three fields are set to NaN to * satisfy parent builder. */ - - if (lowerTrace != null) { - - // if going dual-trace route, satisfy parent builder with NaNs - this.depth = Double.NaN; - this.dip = Double.NaN; - this.width = Double.NaN; - validateState(ID); - surface = new ApproxGriddedSurface(trace, lowerTrace, spacing); - - } else { - - // otherwise build a basic fault source - validateState(ID); - surface = DefaultGriddedSurface.builder() - .trace(trace) - .depth(depth) - .dip(dip) - .width(width) - .spacing(spacing) - .build(); - } + checkState(!built, "This %s instance as already been used", ID); + checkState(feature != null, "%s feature not set", ID); + checkState(config != null, "%s config not set", ID); + checkState(mfdTree != null, "%s MFD tree not set", ID); + checkState(mfdTotal != null, "%s total MFD not set", ID); + checkState(branchWeight != null, "%s source branch weight not set", ID); + built = true; return new InterfaceSource(this); } } - @Override - List<Rupture> createRuptureList(Mfd mfd) { + private List<Rupture> createRuptureList(Mfd mfd) { List<Rupture> ruptureList = new ArrayList<>(); // TODO be careful when working with combined MFDs RuptureFloating floatModel = (mfd.properties().type() == Mfd.Type.SINGLE) ? RuptureFloating.OFF - : rupFloating; + : config.ruptureFloating; for (XyPoint xy : mfd.data()) { double mag = xy.x(); double rate = xy.y(); @@ -183,7 +256,7 @@ public class InterfaceSource extends FaultSource { // SourceConfig.Fault config = (SourceConfig.Fault) data.sourceConfig(); // rupFloating List<Rupture> floaters = floatModel.createFloatingRuptures( - surface, rupScaling, mag, rate, -90.0); + surface, config.ruptureScaling, mag, rate, -90.0); ruptureList.addAll(floaters); } return List.copyOf(ruptureList); diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/InterfaceSourceSet.java b/src/main/java/gov/usgs/earthquake/nshmp/model/InterfaceSourceSet.java index f702e1aeb092813fb7c18f75d29af38929d46a34..de4ffd401f40626878598b6637f7bee774abb4a7 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/InterfaceSourceSet.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/InterfaceSourceSet.java @@ -72,8 +72,8 @@ public class InterfaceSourceSet extends AbstractSourceSet<InterfaceSource> { @Override public boolean test(InterfaceSource source) { - // List<LocationList> traces = source.traces; - return filter.test(source.trace.first()) || filter.test(source.trace.last()) || + return filter.test(source.upperTrace.first()) || + filter.test(source.upperTrace.last()) || filter.test(source.lowerTrace.first()) || filter.test(source.lowerTrace.last()); } diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/MagnitudeNote.java b/src/main/java/gov/usgs/earthquake/nshmp/model/MagnitudeNote.java new file mode 100644 index 0000000000000000000000000000000000000000..fc42efff6d1519e178cc352e598f72430de6aa2c --- /dev/null +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/MagnitudeNote.java @@ -0,0 +1,38 @@ +package gov.usgs.earthquake.nshmp.model; + +/** + * Rupture magnitude description. + * + * NSHMs include rupture magnitudes that may be constrained by ways other than a + * mgnitude-scaling relationship as defined by this class. + * + * @author U.S. Geological Survey + */ +public enum MagnitudeNote { + + /** + * Magnitude override that should be explained in a README in a source model. + */ + OTHER, + + /** + * Magnitude override that should be explained in a README in a source model + * with a hint to not float SINGLE MFD ruptures. + */ + OTHER_NO_FLOAT, + + /** + * Inexplicable inconsistency between geologic and geodetic magnitdues. + */ + GEODETIC_VARIANT, + + /** + * Constrained by a historic earthquake magnitude. + */ + HISTORIC_M_MAX, + + /** + * Constrained by a regional maximum magnitude. + */ + REGIONAL_M_MAX; +} diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/MfdConfig.java b/src/main/java/gov/usgs/earthquake/nshmp/model/MfdConfig.java index b4283fd1814512c78e7016092d99233cbe9e7708..d7d6e5cdd159045c0d35e5e605e7d8fb6c440b6c 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/MfdConfig.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/MfdConfig.java @@ -21,22 +21,20 @@ class MfdConfig { static final String EPISTEMIC_TREE_KEY = "epistemic-tree"; static final String ALEATORY_PROPS_KEY = "aleatory-properties"; static final String MIN_MAGNITUDE_KEY = "minimum-magnitude"; + static final String NSHM_BIN_MODEL_KEY = "nshm-bin-model"; final Optional<LogicTree<Double>> epistemicTree; final OptionalDouble minEpiOffset; final Optional<AleatoryProperties> aleatoryProperties; final OptionalDouble minMagnitude; + final boolean nshmBinModel; private MfdConfig(Builder builder) { this.epistemicTree = builder.epistemicTree; this.minEpiOffset = builder.minEpiOffset; this.aleatoryProperties = builder.aleatoryProperties; this.minMagnitude = builder.minMagnitude; - } - - @Deprecated - static MfdConfig empty() { - return builder().build(); + this.nshmBinModel = builder.nshmBinModel; } static Builder builder() { @@ -50,6 +48,7 @@ class MfdConfig { private OptionalDouble minEpiOffset = OptionalDouble.empty(); private Optional<AleatoryProperties> aleatoryProperties = Optional.empty(); private OptionalDouble minMagnitude = OptionalDouble.empty(); + private boolean nshmBinModel = true; private boolean built = false; Builder epistemicTree(LogicTree<Double> tree) { @@ -65,7 +64,7 @@ class MfdConfig { } Builder aleatoryProperties(AleatoryProperties props) { - this.aleatoryProperties = Optional.ofNullable(checkNotNull(props)); + this.aleatoryProperties = Optional.of(checkNotNull(props)); return this; } @@ -75,6 +74,11 @@ class MfdConfig { return this; } + Builder nshmBinModel(boolean nshmBinModel) { + this.nshmBinModel = nshmBinModel; + return this; + } + MfdConfig build() { checkState(!built, "This MfdConfig.Builder has already been used"); if (epistemicTree.isPresent() || aleatoryProperties.isPresent()) { diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/MfdTrees.java b/src/main/java/gov/usgs/earthquake/nshmp/model/MfdTrees.java index c0c2a76c47dad31334e6385938a7e520a74a2568..c67635db24cda4d3c276460ce785498c05f38ca5 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/MfdTrees.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/MfdTrees.java @@ -236,8 +236,8 @@ class MfdTrees { /* * Convert any LogicTree to a tree of pseudo-paths based on the id() field of - * the supplied tree. This is used to accomodate rate-tree branches in a - * source-tree. + * the supplied tree. This is used to accomodate rate-tree and dip-tree + * branches in a source-tree. */ static <T> LogicTree<Path> toPathTree(LogicTree<T> tree) { LogicTree.Builder<Path> pathTree = LogicTree.builder(tree.name()); 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 316a4986a6bc5454d9b4d21209ad17f76e95b56b..f52317bf14b569b26eb3244a051ceb529b8fc910 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/ModelLoader.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/ModelLoader.java @@ -75,7 +75,7 @@ import gov.usgs.earthquake.nshmp.tree.LogicTree; abstract class ModelLoader { public static void main(String[] args) { - Path testModel = Paths.get("../nshm-conus-2018"); + Path testModel = Paths.get("../nshm-conus-2018-tmp"); HazardModel model = ModelLoader.load(testModel); System.out.println(model); } @@ -268,63 +268,164 @@ abstract class ModelLoader { SourceFeature.NshmFault feature = SourceFeature.newNshmFault(geojson); + SourceTree.Builder treeBuilder = SourceTree.builder() + .name(feature.name) + .type(FAULT) + .gmms(data.gmms()); + + System.out.println(" source: " + root.relativize(geojson)); + + return loadSingleSource(treeBuilder, Optional.empty(), feature, data).build(); + + // + // /* Possibly split on normal fault dip variants. */ + // SourceConfig.Fault config = data.sourceConfig().asFault(); + // boolean isNormal = FocalMech.fromRake(feature.rake) == + // FocalMech.NORMAL; + // + // if (config.dipTree.isPresent() && isNormal) { + // System.out.println("dip tree: " + root.relativize(geojson)); + // // SourceTree.Builder treeBuilder = SourceTree.builder() + // // .name(feature.name) + // // .type(FAULT) + // // .gmms(data.gmms()); + // + // return loadNormalDipTree(treeBuilder, Optional.empty(), feature, data) + // .build(); + // } + // + // System.out.println(" source: " + root.relativize(geojson)); + // + // /* + // * Create new feature with updated properties. Magnitude may be present + // in + // * properties for regional, historic, or other overrides. + // */ + // // double Mw = feature.magnitude.isPresent() + // // ? feature.magnitude.orElseThrow() + // // : FaultRuptureSet.magnitudeWc94(feature.length); + // // double Δz = feature.lowerDepth - feature.upperDepth; + // // double width = computeWidth(Δz, feature.dip); + // // Feature nshmFeature = Feature.copyOf(feature.source) + // // .properties(Properties.fromFeature(feature.source) + // // .put(Key.MAGNITUDE, Mw) + // // .put(Key.WIDTH, width) + // // .build()) + // // .build(); + // // feature = new SourceFeature.NshmFault(nshmFeature); + // + // // feature = updateFeature(feature); + // + // LogicTree<Path> root = LogicTree.singleton( + // Deserialize.SOURCE_TREE, + // feature.name, + // Path.of(feature.name)); + // + // treeBuilder.root(root); + // + // // SourceTree.Builder treeBuilder = SourceTree.builder() + // // .name(feature.name) + // // .type(FAULT) + // // .gmms(data.gmms()) + // // .root(root); + // + // // RuptureSet ruptureSet = createSingleRuptureSet(data, feature); + // // + // // treeBuilder.addLeaf(root.get(0), ruptureSet); + // // + // // loadSingleRuptureSet( + // // root.get(0), + // // treeBuilder, + // // data, + // // feature); + // + // return loadSingle(treeBuilder, root.get(0), data, feature).build(); + // + // // return treeBuilder.build(); + } + + private SourceTree.Builder loadSingleSource( + SourceTree.Builder treeBuilder, + Optional<Branch<Path>> branch, + SourceFeature.NshmFault feature, + ModelData data) { + /* Possibly split on normal fault dip variants. */ SourceConfig.Fault config = data.sourceConfig().asFault(); boolean isNormal = FocalMech.fromRake(feature.rake) == FocalMech.NORMAL; + if (config.dipTree.isPresent() && isNormal) { - System.out.println("dip tree: " + root.relativize(geojson)); - return loadNormalDipTree(feature, data); + return loadNormalDipTree(treeBuilder, branch, feature, data); } - System.out.println(" source: " + root.relativize(geojson)); + return loadSingle(treeBuilder, branch, data, feature); + } + + private SourceTree.Builder loadSingle( + SourceTree.Builder treeBuilder, + Optional<Branch<Path>> branch, + ModelData data, + SourceFeature.NshmFault feature) { + + /* If supplied branch is empty, we need to set a singleton tree root. */ + if (branch.isEmpty()) { + LogicTree<Path> root = LogicTree.singleton( + Deserialize.SOURCE_TREE, + feature.name, + Path.of(feature.name)); + treeBuilder.root(root); + branch = Optional.of(root.get(0)); + } - /* Create new feature with updated properties. */ - double Mw = FaultRuptureSet.magnitudeWc94(feature.length); + double Mw = feature.magnitude.isPresent() + ? feature.magnitude.orElseThrow() + : FaultRuptureSet.magnitudeWc94(feature.length); double Δz = feature.lowerDepth - feature.upperDepth; double width = computeWidth(Δz, feature.dip); - Feature nshmFeature = Feature.copyOf(feature.source) + Feature newFeature = Feature.copyOf(feature.source) .properties(Properties.fromFeature(feature.source) .put(Key.MAGNITUDE, Mw) .put(Key.WIDTH, width) .build()) .build(); - feature = new SourceFeature.NshmFault(nshmFeature); + feature = new SourceFeature.NshmFault(newFeature); - LogicTree<Path> root = LogicTree.singleton( - Deserialize.SOURCE_TREE, - feature.name, - Path.of(feature.name)); - - SourceTree.Builder treeBuilder = SourceTree.builder() - .name(feature.name) - .type(FAULT) - .gmms(data.gmms()) - .root(root); - - loadSingleRuptureSet( - root.get(0), - treeBuilder, - data, - feature); - - return treeBuilder.build(); + RuptureSet ruptureSet = createSingleRuptureSet(data, feature); + treeBuilder.addLeaf(branch.orElseThrow(), ruptureSet); + return treeBuilder; } - private SourceTree loadNormalDipTree( + private SourceTree.Builder loadNormalDipTree( + SourceTree.Builder treeBuilder, + Optional<Branch<Path>> branch, // if absent, set tree root SourceFeature.NshmFault feature, ModelData data) { + System.out.println("dip tree: "); + SourceConfig.Fault config = data.sourceConfig().asFault(); LogicTree<Double> dipTree = config.dipTree.orElseThrow(); LogicTree<Path> dipSourceTree = MfdTrees.toPathTree(dipTree); DipSlipModel dipSlipModel = config.dipSlipModel.orElseThrow(); - SourceTree.Builder treeBuilder = SourceTree.builder() - .name(feature.name) - .type(FAULT) - .gmms(data.gmms()) - .root(dipSourceTree); + // SourceTree.Builder treeBuilder = SourceTree.builder() + // .name(feature.name) + // .type(FAULT) + // .gmms(data.gmms()) + // .root(dipSourceTree); + + if (branch.isPresent()) { + treeBuilder.addBranches(branch.orElseThrow(), dipSourceTree); + } else { + treeBuilder.root(dipSourceTree); + } - double Mw = FaultRuptureSet.magnitudeWc94(feature.length); + /* + * Create new features with updated properties. Magnitude may be present + * in properties for regional, historic, or other overrides. + */ + double Mw = feature.magnitude.isPresent() + ? feature.magnitude.orElseThrow() + : FaultRuptureSet.magnitudeWc94(feature.length); for (int i = 0; i < dipTree.size(); i++) { Branch<Double> dipBranch = dipTree.get(i); int branchDip = (int) (feature.dip + dipBranch.value()); @@ -342,13 +443,19 @@ abstract class ModelLoader { .build()) .build(); - loadSingleRuptureSet( - dipSourceTree.get(i), - treeBuilder, + RuptureSet ruptureSet = createSingleRuptureSet( data, new SourceFeature.NshmFault(dipFeature)); + + treeBuilder.addLeaf(dipSourceTree.get(i), ruptureSet); + + // loadSingleRuptureSet( + // dipSourceTree.get(i), + // treeBuilder, + // data, + // new SourceFeature.NshmFault(dipFeature)); } - return treeBuilder.build(); + return treeBuilder;// .build(); } /* This exists for consistency with legacy codes */ @@ -357,9 +464,12 @@ abstract class ModelLoader { @Deprecated private static double computeWidth(double Δz, double dip) { - return Maths.round(Δz / Math.sin(dip * TO_RAD), 5, RoundingMode.FLOOR); + // Inexplicable rounding up of 45° where other dips go down + RoundingMode rMode = (dip == 45.0) ? RoundingMode.CEILING : RoundingMode.FLOOR; + return Maths.round(Δz / Math.sin(dip * TO_RAD), 5, rMode); } + @Deprecated private void loadSingleRuptureSet( Branch<Path> branch, SourceTree.Builder treeBuilder, @@ -367,8 +477,8 @@ abstract class ModelLoader { SourceFeature.NshmFault feature) { Map<String, LogicTree<Mfd.Properties>> mfdMap = data.mfdMap().orElseThrow(); + // TODO singleton map is a bad assumtion, LogicTree<Mfd.Properties> mfdTree = mfdMap.values().iterator().next(); - data.faultFeatureMap(Map.of(feature.id, feature)); FaultRuptureSet ruptureSet = FaultRuptureSet.builder() @@ -382,6 +492,28 @@ abstract class ModelLoader { treeBuilder.addLeaf(branch, ruptureSet); } + private FaultRuptureSet createSingleRuptureSet( + // Branch<Path> branch, + // SourceTree.Builder treeBuilder, + ModelData data, + SourceFeature.NshmFault feature) { + + Map<String, LogicTree<Mfd.Properties>> mfdMap = data.mfdMap().orElseThrow(); + // TODO singleton map is a bad assumtion, + LogicTree<Mfd.Properties> mfdTree = mfdMap.values().iterator().next(); + data.faultFeatureMap(Map.of(feature.id, feature)); + + return FaultRuptureSet.builder() + .modelData(data) + .name(feature.name) + .id(feature.id) + .sections(List.of(feature.id)) + .mfdTree(mfdTree) + .build(); + + // treeBuilder.addLeaf(branch, ruptureSet); + } + @Override SourceTree loadSourceTree( Path dir, @@ -390,11 +522,13 @@ abstract class ModelLoader { System.out.println(" tree: " + root.relativize(dir)); LogicTree<Path> tree = readSourceTree(dir).orElseThrow(); + /* Set flag to process tree of single geojson sources. */ + boolean singleSourceTree = tree.get(0).value().toString().endsWith(".geojson"); + /* Set flag to process a logic-tree of fault-system solutions. */ var faultSections = readFaultFeatures(dir); - boolean faultSystemTree = faultSections.isEmpty(); - faultSections.ifPresent(data::faultFeatureMap); + boolean faultSystemTree = faultSections.isEmpty() && !singleSourceTree; SourceTree.Builder treeBuilder = SourceTree.builder() .name(checkNotNull(dir.getFileName()).toString()) @@ -403,7 +537,9 @@ abstract class ModelLoader { .root(tree); for (Branch<Path> branch : tree) { - if (faultSystemTree) { + if (singleSourceTree) { + processSingleBranch(dir, branch, treeBuilder, data); + } else if (faultSystemTree) { processSystemBranch(branch, treeBuilder, data); } else { processBranch(branch, treeBuilder, data); @@ -479,6 +615,20 @@ abstract class ModelLoader { readSystemRuptureSet(dir, data).orElseThrow()); } } + + /* + * Process branhces consisting of individual geojson files as would be + * processed using loadSingleSource(). + */ + private void processSingleBranch( + Path dir, + Branch<Path> branch, + SourceTree.Builder treeBuilder, + ModelData data) { + + SourceFeature.NshmFault feature = SourceFeature.newNshmFault(branch.value()); + loadSingleSource(treeBuilder, Optional.of(branch), feature, data); + } } static List<SourceTree> interfaceSources(Path dir) { diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/Source.java b/src/main/java/gov/usgs/earthquake/nshmp/model/Source.java index 3deed542a60fa9109d90341cf1ae086819ade159..19c01d7ef0585f14749d07c7672e27cd529c466f 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/Source.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/Source.java @@ -25,14 +25,14 @@ public interface Source extends Iterable<Rupture> { int id(); /** - * The total number of {@link Rupture}s this source represents. + * The source type. */ - int size(); + SourceType type(); /** - * The source type. + * The total number of {@link Rupture}s this source represents. */ - SourceType type(); + int size(); /** * The {@code Location} of this source relative to the supplied {@code site} diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/SourceFeature.java b/src/main/java/gov/usgs/earthquake/nshmp/model/SourceFeature.java index 51d35515fd138d442edcdc6d4f2a0eddb2a91187..f226d9820d4b5a8904e0f8e000b22f20cd9c5cd8 100644 --- a/src/main/java/gov/usgs/earthquake/nshmp/model/SourceFeature.java +++ b/src/main/java/gov/usgs/earthquake/nshmp/model/SourceFeature.java @@ -8,8 +8,9 @@ import static gov.usgs.earthquake.nshmp.Faults.checkDip; import static gov.usgs.earthquake.nshmp.Faults.checkDipDirection; import static gov.usgs.earthquake.nshmp.Faults.checkRake; import static gov.usgs.earthquake.nshmp.Faults.checkSlipRate; +import static gov.usgs.earthquake.nshmp.Faults.checkTrace; import static gov.usgs.earthquake.nshmp.Text.checkName; -import static java.util.stream.Collectors.toMap; +import static java.util.stream.Collectors.toUnmodifiableMap; import java.awt.geom.Area; import java.lang.reflect.Type; @@ -81,8 +82,7 @@ public abstract class SourceFeature { id = feature.idAsInt().orElseThrow(); checkArgument(id > 0, "Feature ID [%s] < 1", id); Properties props = feature.properties(); - name = props.getString(Key.NAME).orElseThrow(); - checkName(name, "Feature"); + name = checkName(props.getString(Key.NAME).orElseThrow(), "Feature"); country = props.getString(Key.COUNTRY).orElse("US"); state = props.getString(Key.STATE); states = List.of(props.get(Key.STATES, String[].class).orElse( @@ -195,6 +195,7 @@ public abstract class SourceFeature { static final String RATE_DIP = "rate-dip"; static final String M = "m"; static final String MAGNITUDE = "magnitude"; + static final String MAGNITUDE_NOTE = "magnitude-note"; static final String LENGTH = "length"; static final String WIDTH = "width"; static final String DEPTH = "depth"; @@ -244,15 +245,12 @@ public abstract class SourceFeature { private Fault(Feature section, Optional<LocationList> zone) { super(section); trace = section.asLineString(); + checkTrace(trace); Properties props = section.properties(); - dip = props.getDouble(Key.DIP).orElseThrow(); - checkDip(dip); - upperDepth = props.getDouble(Key.UPPER_DEPTH).orElseThrow(); - checkCrustalDepth(upperDepth); - lowerDepth = props.getDouble(Key.LOWER_DEPTH).orElseThrow(); - checkCrustalDepth(lowerDepth); - rake = props.getDouble(Key.RAKE).orElseThrow(); - checkRake(rake); + dip = checkDip(props.getDouble(Key.DIP).orElseThrow()); + upperDepth = checkCrustalDepth(props.getDouble(Key.UPPER_DEPTH).orElseThrow()); + lowerDepth = checkCrustalDepth(props.getDouble(Key.LOWER_DEPTH).orElseThrow()); + rake = checkRake(props.getDouble(Key.RAKE).orElseThrow()); rateType = RateType.valueOf(props.getString(Key.RATE_TYPE).orElseThrow()); rate = (rateType != RateType.RECURRENCE) ? checkSlipRate(props.getDouble(Key.RATE).orElseThrow()) @@ -305,6 +303,7 @@ public abstract class SourceFeature { */ final OptionalDouble rateDip; final OptionalDouble magnitude; + final Optional<MagnitudeNote> magnitudeNote; NshmFault(Feature section) { super(section); @@ -314,12 +313,9 @@ public abstract class SourceFeature { /* Required properties. */ dip = props.getDouble(Key.DIP).orElseThrow(); checkDip(dip); - upperDepth = props.getDouble(Key.UPPER_DEPTH).orElseThrow(); - checkCrustalDepth(upperDepth); - lowerDepth = props.getDouble(Key.LOWER_DEPTH).orElseThrow(); - checkCrustalDepth(lowerDepth); - rake = props.getDouble(Key.RAKE).orElseThrow(); - checkRake(rake); + upperDepth = checkCrustalDepth(props.getDouble(Key.UPPER_DEPTH).orElseThrow()); + lowerDepth = checkCrustalDepth(props.getDouble(Key.LOWER_DEPTH).orElseThrow()); + rake = checkRake(props.getDouble(Key.RAKE).orElseThrow()); rateType = RateType.valueOf(props.getString(Key.RATE_TYPE).orElseThrow()); rateMap = nshmRateMap(props); props.getString(Key.STATE).orElseThrow(); @@ -331,6 +327,7 @@ public abstract class SourceFeature { /* Optional properties. */ rateDip = props.getDouble(Key.RATE_DIP); magnitude = props.getDouble(Key.MAGNITUDE); + magnitudeNote = props.get(Key.MAGNITUDE_NOTE, MagnitudeNote.class); } private double computeWidth() { @@ -398,7 +395,7 @@ public abstract class SourceFeature { } Map<String, Rates> rateMap = rateElements.orElseThrow().entrySet().stream() .filter(e -> !e.getValue().isJsonNull()) - .collect(toMap( + .collect(toUnmodifiableMap( Entry::getKey, e -> new Rates(e.getValue()))); return Optional.of(rateMap); @@ -438,16 +435,13 @@ public abstract class SourceFeature { Properties props = subsection.properties(); index = props.getInt(Key.INDEX).orElseThrow(); parentId = props.getInt(Key.PARENT_ID).orElseThrow(); - dip = props.getDouble(Key.DIP).orElseThrow(); - checkDip(dip); - dipDirection = props.getDouble(Key.DIP_DIRECTION).orElseThrow(); - checkDipDirection(trace, dipDirection); - upperDepth = props.getDouble(Key.UPPER_DEPTH).orElseThrow(); - checkCrustalDepth(upperDepth); - lowerDepth = props.getDouble(Key.LOWER_DEPTH).orElseThrow(); - checkCrustalDepth(lowerDepth); - aseismicity = props.getDouble(Key.ASEISMICITY).orElseThrow(); - checkAseismicityFactor(aseismicity); + dip = checkDip(props.getDouble(Key.DIP).orElseThrow()); + dipDirection = checkDipDirection( + trace, + props.getDouble(Key.DIP_DIRECTION).orElseThrow()); + upperDepth = checkCrustalDepth(props.getDouble(Key.UPPER_DEPTH).orElseThrow()); + lowerDepth = checkCrustalDepth(props.getDouble(Key.LOWER_DEPTH).orElseThrow()); + aseismicity = checkAseismicityFactor(props.getDouble(Key.ASEISMICITY).orElseThrow()); props.getString(Key.STATE).orElseThrow(); } }