diff --git a/src/main/java/gov/usgs/earthquake/nshmp/mfd/Mfd.java b/src/main/java/gov/usgs/earthquake/nshmp/mfd/Mfd.java
index b0340803b76435a7233f48f5b85d37160370462e..3378ad66b672cd06e7fac9af1d9fb73cb330513c 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/mfd/Mfd.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/mfd/Mfd.java
@@ -27,6 +27,7 @@ import gov.usgs.earthquake.nshmp.data.MutableXySequence;
import gov.usgs.earthquake.nshmp.data.Sequences;
import gov.usgs.earthquake.nshmp.data.XyPoint;
import gov.usgs.earthquake.nshmp.data.XySequence;
+import gov.usgs.earthquake.nshmp.mfd.Mfd.Properties.TaperedGr;
/**
* Entry point for creating a magnitude frequency distribution (MFD). An MFD
@@ -429,6 +430,10 @@ public final class Mfd {
* @return this {@code Builder} object
*/
public Builder scaleToIncrementalRate(double incrementalRate) {
+ // TODO remove (needs sensitivity testing)
+ if (props.type == Type.GR_TAPER) {
+ return TaperedGr.taperCorrection(this, incrementalRate);
+ }
return scale(checkRate(incrementalRate) / mfd.y(0));
}
@@ -453,6 +458,7 @@ public final class Mfd {
mfd.transform(checkNotNull(action));
return this;
}
+
}
/**
@@ -779,6 +785,13 @@ public final class Mfd {
* directly with a given a-value and by creating a generic taperd MFD
* that we then scaleToIncrementalRate() using the given a-value, we
* correct the tapered MFD by the offset.
+ *
+ * Now, this is how we wished the nshmp-haz implementation behaved, but
+ * since it doesn't we need an intermediate solution to match nshmp-haz
+ * until further sensitivity testing through to hazard can be
+ * undertaken. What we want: after calling scaleToIncrementalRate(), the
+ * tapered GR MFD rates should match those that would've resulted from
+ * the direct calculation.
*/
double grRate = Mfds.gutenbergRichterRate(a(), b(), builder.mfd.x(0));
double tgrRate = builder.mfd.y(0);
@@ -798,6 +811,24 @@ public final class Mfd {
type(),
super.propsString() + ", mc=" + mc);
}
+
+ /*
+ * TODO clean after sensitivity testing for tapered rates when calling
+ * scale toIncrementalRate() tapered GR values become exact matches to the
+ * target rate as opposed to the slightly scaled rate resulting from
+ * creating the tapered GR directly with a given a-value. We want to
+ * replicate nshmp-haz behavior for now but get rid of it in the not too
+ * distant future. This is basically reverting the scaling correction
+ * applied in the builder.
+ */
+ private static Builder taperCorrection(Builder b, double rate) {
+ TaperedGr p = b.props.getAsGrTaper();
+ TaperFunction taperFn = new TaperFunction(p.Δm(), p.b(), p.mc);
+ double mMin = b.mfd.x(0);
+ double taperRate = checkRate(rate) * taperFn.scale(mMin);
+ return b.scale(taperRate / b.mfd.y(0));
+ }
+
}
private static final double M_MIN_MOMENT = Earthquakes.magToMoment(4.0);
diff --git a/src/test/java/gov/usgs/earthquake/nshmp/mfd/MfdTests.java b/src/test/java/gov/usgs/earthquake/nshmp/mfd/MfdTests.java
index 6c501842245f6e20e32fca3c0a8a9032c78bac6d..6270622bdf0bdace3ad89574687e77dea941a3fa 100644
--- a/src/test/java/gov/usgs/earthquake/nshmp/mfd/MfdTests.java
+++ b/src/test/java/gov/usgs/earthquake/nshmp/mfd/MfdTests.java
@@ -10,6 +10,7 @@ import java.io.IOException;
import java.nio.file.Files;
import java.nio.file.Path;
import java.nio.file.Paths;
+import java.util.Arrays;
import org.junit.jupiter.api.Test;
@@ -202,6 +203,66 @@ class MfdTests {
grProps.toString());
}
+ /*
+ * TODO clean
+ *
+ * TODO something in the TaperedGR moment based scaling leads to differences
+ * in the 4th to 5th significant figure that should be better understood, and
+ * then this main method cleaned out. The results are still pretty close
+ * overall.
+ */
+ public static void main(String[] args) {
+ double bVal = 0.8;
+ double aVal = Math.log10(4.0);
+ double Δm = 0.1;
+ double mMin = 5.0;
+ double mMax = 8.0;
+ double mc = 7.5;
+
+ // /* GR direct vs scaled builder comparison */
+ // System.out.println("Gutenberg Richter");
+ // Mfd mfd = new GutenbergRichter(aVal, bVal, Δm, mMin, mMax)
+ // .toBuilder()
+ // .build();
+ //
+ // System.out.println(Arrays.toString(mfd.data().xValues().toArray()));
+ // System.out.println(Arrays.toString(mfd.data().yValues().toArray()));
+ //
+ // double incrRate = Mfds.gutenbergRichterRate(aVal, bVal, 5.05);
+ // mfd = Mfd.newGutenbergRichterBuilder(bVal, Δm, mMin, mMax)
+ // .scaleToIncrementalRate(incrRate)
+ // .build();
+ //
+ // System.out.println(Arrays.toString(mfd.data().yValues().toArray()));
+
+ /* GR direct vs scaled builder comparison */
+ System.out.println("Tapered Gutenberg Richter");
+ Mfd mfd = new TaperedGr(aVal, bVal, Δm, mMin, mMax, mc)
+ .toBuilder()
+ .build();
+
+ double[] directRates = mfd.data().yValues().toArray();
+ System.out.println(Arrays.toString(mfd.data().xValues().toArray()));
+ System.out.println(Arrays.toString(directRates));
+
+ double incrRate = Mfds.gutenbergRichterRate(aVal, bVal, 5.05);
+ mfd = Mfd.newTaperedGutenbergRichterBuilder(bVal, Δm, mMin, mMax, mc)
+ .scaleToIncrementalRate(incrRate)
+ .build();
+
+ double[] scaledRates = mfd.data().yValues().toArray();
+ System.out.println(Arrays.toString(scaledRates));
+
+ int size = directRates.length;
+ double[] ratios = new double[size];
+ for (int i = 0; i < size; i++) {
+ double ratio = 100.0 * (directRates[i] / scaledRates[i]);
+ ratios[i] = ratio;
+ }
+ System.out.println(Arrays.toString(ratios));
+
+ }
+
private static final double[] TAPERED_GR_M = {
5.05, 5.15, 5.25, 5.35, 5.45,
5.55, 5.65, 5.75, 5.85, 5.95,
@@ -243,6 +304,39 @@ class MfdTests {
5.437829885686743E-7,
1.5949675112240326E-7 };
+ /* TODO These match nshmp-haz and should be removed */
+ private static final double[] TAPERED_GR_R_TEMP = {
+ 3.648734771264746E-4,
+ 3.0351155247723867E-4,
+ 2.5247694248332143E-4,
+ 2.1003291504182058E-4,
+ 1.7473503715378214E-4,
+ 1.4538201763727166E-4,
+ 1.2097482132130438E-4,
+ 1.0068266229609009E-4,
+ 8.381471718000333E-5,
+ 6.979659287661517E-5,
+ 5.815074326255645E-5,
+ 4.848016071724236E-5,
+ 4.0454775273297684E-5,
+ 3.3800079282566184E-5,
+ 2.8287560844165223E-5,
+ 2.3726576420233305E-5,
+ 1.9957324528955872E-5,
+ 1.6844604598702632E-5,
+ 1.4272075425536471E-5,
+ 1.213680849238659E-5,
+ 1.0344152445466782E-5,
+ 8.803445832444013E-6,
+ 7.426137715686031E-6,
+ 6.1294171417451564E-6,
+ 4.849632254896957E-6,
+ 3.5675342487878286E-6,
+ 2.3362943546550987E-6,
+ 1.2825840184413834E-6,
+ 5.438860517858617E-7,
+ 1.5952698054966957E-7 };
+
@Test
void testTaperedGr() {
@@ -262,7 +356,7 @@ class MfdTests {
assertTrue(xy.size() == 30);
assertArrayEquals(TAPERED_GR_M, xy.xValues().toArray());
/* These end up being very slightly different due to scaling */
- assertArrayEquals(TAPERED_GR_R, xy.yValues().toArray(), 1e-17);
+ assertArrayEquals(TAPERED_GR_R_TEMP, xy.yValues().toArray(), 1e-17);
/* Properties */
Properties props = new TaperedGr(aVal, bVal, Δm, mMin, mMax, mc);