From 06aaa101f159ffd95eaa6160e830154e1225b460 Mon Sep 17 00:00:00 2001
From: Peter Powers <pmpowers@usgs.gov>
Date: Tue, 22 Oct 2024 15:27:48 -0600
Subject: [PATCH] adapted system source and rupset to handle interface

---
 .../earthquake/nshmp/model/SourceFeature.java | 28 ++++++++++---
 .../nshmp/model/SystemRuptureSet.java         | 41 ++++++++++++++-----
 2 files changed, 53 insertions(+), 16 deletions(-)

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 f4ed3910..03acc3cf 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/model/SourceFeature.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/model/SourceFeature.java
@@ -3,6 +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 gov.usgs.earthquake.nshmp.Earthquakes.checkCrustalDepth;
+import static gov.usgs.earthquake.nshmp.Earthquakes.checkInterfaceDepth;
 import static gov.usgs.earthquake.nshmp.Faults.checkAseismicityFactor;
 import static gov.usgs.earthquake.nshmp.Faults.checkDepths;
 import static gov.usgs.earthquake.nshmp.Faults.checkDip;
@@ -441,8 +442,11 @@ public abstract class SourceFeature {
     /** ID of the parent fault section. */
     final int parentId;
 
-    /** The fault subsection trace. */
-    final LocationList trace;
+    /**
+     * The fault subsection traces. For crustal fault system sections there will
+     * be only one trace, for interface system sections there will be two.
+     */
+    final List<LocationList> traces;
 
     /** The fault subsection dip. */
     final double dip;
@@ -461,16 +465,28 @@ public abstract class SourceFeature {
 
     SystemSection(Feature subsection) {
       super(subsection);
-      trace = subsection.asLineString();
+      List<LocationList> traces = null;
+      try {
+        traces = List.of(subsection.asLineString());
+      } catch (UnsupportedOperationException uoe) {
+        traces = subsection.asMultiLineString();
+      }
+      this.traces = traces;
       Properties props = subsection.properties();
       index = props.getInt(Key.INDEX).orElseThrow();
       parentId = props.getInt(Key.PARENT_ID).orElseThrow();
       dip = checkDip(props.getDouble(Key.DIP).orElseThrow());
       dipDirection = checkDipDirection(
-          trace,
+          traces.get(0),
           props.getDouble(Key.DIP_DIRECTION).orElseThrow());
-      upperDepth = checkCrustalDepth(props.getDouble(Key.UPPER_DEPTH).orElseThrow());
-      lowerDepth = checkCrustalDepth(props.getDouble(Key.LOWER_DEPTH).orElseThrow());
+      double upperDepth = props.getDouble(Key.UPPER_DEPTH).orElseThrow();
+      double lowerDepth = props.getDouble(Key.LOWER_DEPTH).orElseThrow();
+      this.upperDepth = (traces.size() == 1)
+          ? checkCrustalDepth(upperDepth)
+          : checkInterfaceDepth(upperDepth);
+      this.lowerDepth = (traces.size() == 1)
+          ? checkCrustalDepth(lowerDepth)
+          : checkInterfaceDepth(lowerDepth);
       checkDepths(lowerDepth, upperDepth);
       aseismicity = checkAseismicityFactor(props.getDouble(Key.ASEISMICITY).orElseThrow());
       props.getString(Key.STATE).orElseThrow();
diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/SystemRuptureSet.java b/src/main/java/gov/usgs/earthquake/nshmp/model/SystemRuptureSet.java
index 09167bae..a9943527 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/model/SystemRuptureSet.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/model/SystemRuptureSet.java
@@ -4,6 +4,8 @@ import static com.google.common.base.Preconditions.checkArgument;
 import static com.google.common.base.Preconditions.checkNotNull;
 import static gov.usgs.earthquake.nshmp.Earthquakes.checkCrustalDepth;
 import static gov.usgs.earthquake.nshmp.Earthquakes.checkCrustalWidth;
+import static gov.usgs.earthquake.nshmp.Earthquakes.checkInterfaceDepth;
+import static gov.usgs.earthquake.nshmp.Earthquakes.checkInterfaceWidth;
 import static gov.usgs.earthquake.nshmp.Earthquakes.checkMagnitude;
 import static gov.usgs.earthquake.nshmp.Faults.checkDip;
 import static gov.usgs.earthquake.nshmp.Faults.checkRake;
@@ -46,6 +48,7 @@ import gov.usgs.earthquake.nshmp.data.DelimitedData;
 import gov.usgs.earthquake.nshmp.data.DelimitedData.Record;
 import gov.usgs.earthquake.nshmp.data.Indexing;
 import gov.usgs.earthquake.nshmp.data.IntervalArray;
+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.geo.Location;
@@ -112,10 +115,14 @@ public class SystemRuptureSet extends AbstractRuptureSet<SystemRuptureSet.System
     this.sectionMap = builder.sectionMap;
     this.ruptures = builder.ruptures;
 
+    Function<SystemSection, GriddedSurface> surfaceFunction = (type() == SourceType.FAULT_SYSTEM)
+        ? SystemRuptureSet::featureToCrustalSurface
+        : SystemRuptureSet::featureToInterfaceSurface;
+
     this.sections = sectionMap.keySet().stream()
         .sorted()
         .map(sectionMap::get)
-        .map(SystemRuptureSet::featureToSurface)
+        .map(surfaceFunction)
         .toArray(GriddedSurface[]::new);
 
     this.sectionNames = sectionMap.keySet().stream()
@@ -157,10 +164,16 @@ public class SystemRuptureSet extends AbstractRuptureSet<SystemRuptureSet.System
       checkArgument(Doubles.isFinite(rate), "Non-finite rate value");
       this.rates[i] = rate;
 
-      this.depths[i] = checkCrustalDepth(r.getDouble(DEPTH));
       this.dips[i] = checkDip(r.getDouble(DIP));
-      this.widths[i] = checkCrustalWidth(r.getDouble(WIDTH));
       this.rakes[i] = checkRake(r.getDouble(RAKE));
+      if (type() == SourceType.FAULT_SYSTEM) {
+        this.depths[i] = checkCrustalDepth(r.getDouble(DEPTH));
+        this.widths[i] = checkCrustalWidth(r.getDouble(WIDTH));
+      } else { // interface
+        this.depths[i] = checkInterfaceDepth(r.getDouble(DEPTH));
+        this.widths[i] = checkInterfaceWidth(r.getDouble(WIDTH));
+
+      }
     }
     this.stats = new Statistics(mMin, mMax);
 
@@ -206,17 +219,25 @@ public class SystemRuptureSet extends AbstractRuptureSet<SystemRuptureSet.System
     return IntStream.rangeClosed(min, max);
   }
 
-  private static GriddedSurface featureToSurface(SourceFeature.SystemSection section) {
+  private static GriddedSurface featureToCrustalSurface(SourceFeature.SystemSection section) {
     return DefaultGriddedSurface.builder()
         .depth(section.upperDepth)
         .lowerDepth(section.lowerDepth)
         .aseis(section.aseismicity)
         .dip(section.dip)
         .dipDir(section.dipDirection)
-        .trace(section.trace)
+        .trace(section.traces.get(0))
+        .spacing(1.0)
         .build();
   }
 
+  private static GriddedSurface featureToInterfaceSurface(SourceFeature.SystemSection section) {
+    return new ApproxGriddedSurface(
+        section.traces.get(0),
+        section.traces.get(1),
+        5.0);
+  }
+
   @Override
   public int size() {
     return bitsets.length;
@@ -348,7 +369,7 @@ public class SystemRuptureSet extends AbstractRuptureSet<SystemRuptureSet.System
 
     @Override
     public SourceType type() {
-      return SourceType.FAULT_SYSTEM;
+      return SystemRuptureSet.this.type();
     }
 
     /**
@@ -431,8 +452,8 @@ public class SystemRuptureSet extends AbstractRuptureSet<SystemRuptureSet.System
     // return null; // Locations.closestPoint(site, feature.trace);
   }
 
-  static Builder builder() {
-    return new Builder();
+  static Builder builder(SourceType type) {
+    return new Builder(type);
   }
 
   // Temporary name generating method, called externally
@@ -463,8 +484,8 @@ public class SystemRuptureSet extends AbstractRuptureSet<SystemRuptureSet.System
     private double mMin = Double.POSITIVE_INFINITY;
     private double mMax = Double.NEGATIVE_INFINITY;
 
-    private Builder() {
-      super(SourceType.FAULT_SYSTEM);
+    private Builder(SourceType type) {
+      super(type);
     }
 
     @Override
-- 
GitLab