From 020c78f8b0787bc62a10b8b4de04210dedf6af37 Mon Sep 17 00:00:00 2001
From: Peter Powers <pmpowers@usgs.gov>
Date: Mon, 20 Sep 2021 09:24:35 -0600
Subject: [PATCH] zone rupture set MFD tree fix

---
 .../nshmp/model/FaultRuptureSet.java          |  4 +-
 .../usgs/earthquake/nshmp/model/Models.java   |  3 ++
 .../nshmp/model/ZoneRuptureSet.java           | 40 +++++++++++++++----
 3 files changed, 37 insertions(+), 10 deletions(-)

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 c67c2b2a..0c9f6a41 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/model/FaultRuptureSet.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/model/FaultRuptureSet.java
@@ -798,8 +798,8 @@ public class FaultRuptureSet implements RuptureSet {
         Trees.checkTreeTypes(mfdTree, Mfd.Type.SINGLE);
     boolean hasEpistemic = mfdConfig.epistemicTree.isPresent();
     if (multipleSingleBranches && hasEpistemic) {
-      System.err
-          .println("                 WARNING: Multiple SINGLE MFDs with epistemic uncertainty");
+      System.err.println(
+          "                 WARNING: Multiple SINGLE MFDs with epistemic uncertainty");
     }
     // checkState(
     // !(multipleSingleBranches && hasEpistemic),
diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/Models.java b/src/main/java/gov/usgs/earthquake/nshmp/model/Models.java
index bf8dc998..e0b088e4 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/model/Models.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/model/Models.java
@@ -50,6 +50,9 @@ public class Models {
     System.out.println(GSON.toJson(tree(model, 3199)));
     System.out.println();
 
+    System.out.println(GSON.toJson(tree(model, 2539)));
+    System.out.println();
+
     // 2799 is huge wasatch tree
   }
 
diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/ZoneRuptureSet.java b/src/main/java/gov/usgs/earthquake/nshmp/model/ZoneRuptureSet.java
index e820f353..a476764a 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/model/ZoneRuptureSet.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/model/ZoneRuptureSet.java
@@ -1,7 +1,10 @@
 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.mfd.Mfd.Type.GR;
+import static gov.usgs.earthquake.nshmp.mfd.Mfd.Type.SINGLE;
 import static gov.usgs.earthquake.nshmp.model.Deserialize.MFD_TREE;
 import static gov.usgs.earthquake.nshmp.model.Deserialize.RATE_TREE;
 import static java.util.stream.Collectors.toUnmodifiableList;
@@ -43,6 +46,8 @@ class ZoneRuptureSet implements RuptureSet {
    *
    * Rates in rate tree will need to be distributed to grid rates that scale
    * with latitude.
+   *
+   * ZOne sources: all GR or SINGLE, but not both
    */
 
   ZoneRuptureSet(Builder builder) {
@@ -159,6 +164,9 @@ class ZoneRuptureSet implements RuptureSet {
 
     /* Set the MFD tree. */
     Builder mfdTree(LogicTree<Mfd.Properties> mfdTree) {
+      checkArgument(Trees.checkTreeTypes(mfdTree, GR) ||
+          Trees.checkTreeTypes(mfdTree, SINGLE),
+          "Zone MFDs should all be of the same type");
       this.mfdPropertiesTree = mfdTree;
       return this;
     }
@@ -175,7 +183,7 @@ class ZoneRuptureSet implements RuptureSet {
       checkNotNull(mfdPropertiesTree, "%s MFD tree not set", label);
       initLocationAndRateLists();
       mfds = createTotalMfds();
-      mfdsTree = createSingleMfdsTree();
+      mfdsTree = createMfdsTree();
       mfdTree = Trees.reduceMfdListTree(mfdsTree);
     }
 
@@ -200,9 +208,9 @@ class ZoneRuptureSet implements RuptureSet {
     }
 
     private List<Mfd> createTotalMfds() {
-      // TODO hack to get WUS GR zones to work
       double[] rates = reduceRatesTree();
       Mfd.Properties props0 = mfdPropertiesTree.get(0).value();
+      // SINGLE
       if (props0.type().equals(Mfd.Type.SINGLE)) {
         Mfd model = Trees.reduceSingleMfdTree(mfdPropertiesTree);
         return Arrays.stream(rates)
@@ -235,19 +243,34 @@ class ZoneRuptureSet implements RuptureSet {
     }
 
     // SINGLE only
-    private LogicTree<List<Mfd>> createSingleMfdsTree() {
+    private LogicTree<List<Mfd>> createMfdsTree() {
       Mfd.Properties props0 = mfdPropertiesTree.get(0).value();
       // TODO hack short circuit for WUS zones
-      if (props0.type().equals(Mfd.Type.GR)) {
-        return null;
+      if (props0.type().equals(Mfd.Type.SINGLE)) {
+        LogicTree.Builder<List<Mfd>> mfdsTree = LogicTree.builder(MFD_TREE);
+        for (Branch<Mfd.Properties> m : mfdPropertiesTree) {
+          Mfd.Builder builder = m.value().getAsSingle().toBuilder();
+          for (Branch<double[]> r : ratesTree) {
+            List<Mfd> mfds = Arrays.stream(r.value())
+                .mapToObj(builder::scaleToCumulativeRate)
+                .map(Mfd.Builder::build)
+                .collect(toUnmodifiableList());
+            mfdsTree.addBranch(
+                m.id() + "-" + r.id(),
+                mfds,
+                m.weight() * r.weight());
+          }
+        }
+        return mfdsTree.build();
       }
+      // else GR
+      GutenbergRichter grProps = props0.getAsGr();
       LogicTree.Builder<List<Mfd>> mfdsTree = LogicTree.builder(MFD_TREE);
       for (Branch<Mfd.Properties> m : mfdPropertiesTree) {
-        Mfd.Builder builder = m.value().getAsSingle().toBuilder();
         for (Branch<double[]> r : ratesTree) {
           List<Mfd> mfds = Arrays.stream(r.value())
-              .mapToObj(builder::scaleToCumulativeRate)
-              .map(Mfd.Builder::build)
+              .map(Math::log10)
+              .mapToObj(a -> createGrMfd(grProps, a))
               .collect(toUnmodifiableList());
           mfdsTree.addBranch(
               m.id() + "-" + r.id(),
@@ -256,6 +279,7 @@ class ZoneRuptureSet implements RuptureSet {
         }
       }
       return mfdsTree.build();
+
     }
 
   }
-- 
GitLab