diff --git a/src/main/java/gov/usgs/earthquake/nshmp/calc/HazardCalcs.java b/src/main/java/gov/usgs/earthquake/nshmp/calc/HazardCalcs.java
index ee47a2ceff2386fd00604fc705d8fc7f2c2fbe16..63cb9ed7d2560525caa4cb1985097de87a152973 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/calc/HazardCalcs.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/calc/HazardCalcs.java
@@ -164,6 +164,11 @@ public class HazardCalcs {
                   (SystemRuptureSet) ruptures, site, config, ex));
               break;
 
+            case INTERFACE_SYSTEM:
+              curveSets.add(systemToCurves(
+                  (SystemRuptureSet) ruptures, site, config, ex));
+              break;
+
             default:
               curveSets.add(sourcesToCurves(ruptures, site, config, ex));
               break;
@@ -235,6 +240,12 @@ public class HazardCalcs {
               log(log, MSSG_COMPLETED, ruptures.name(), duration(swSource));
               break;
 
+            case INTERFACE_SYSTEM:
+              curveSets.add(systemToCurves(
+                  (SystemRuptureSet) ruptures, site, config));
+              log(log, MSSG_COMPLETED, ruptures.name(), duration(swSource));
+              break;
+
             default:
               curveSets.add(sourcesToCurves(ruptures, site, config));
               log(log, MSSG_COMPLETED, ruptures.name(), duration(swSource));
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 aedbd4cbdea27f36e4fb48d10bfb246de932cf15..2cd5364c8520f310cd06581820030850d9577859 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/model/Deserialize.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/model/Deserialize.java
@@ -4,6 +4,8 @@ import static com.google.common.base.Preconditions.checkNotNull;
 import static com.google.common.base.Preconditions.checkState;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.checkRupturesFile;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readSystemFeatures;
+import static gov.usgs.earthquake.nshmp.model.SourceType.FAULT_SYSTEM;
+import static gov.usgs.earthquake.nshmp.model.SourceType.INTERFACE_SYSTEM;
 import static java.util.stream.Collectors.toUnmodifiableMap;
 
 import java.io.BufferedReader;
@@ -389,13 +391,25 @@ class Deserialize {
         interfaceRuptureSets);
   }
 
-  /* Create a fault system rupture set. */
-  static SystemRuptureSet systemRuptureSet(Path json, ModelData data) {
+  /* Create a crustal fault system rupture set. */
+  static SystemRuptureSet crustalSystemRuptureSet(Path json, ModelData data) {
+    return systemRuptureSet(json, data, FAULT_SYSTEM);
+  }
+
+  /* Create an interface system rupture set. */
+  static SystemRuptureSet interfaceSystemRuptureSet(Path json, ModelData data) {
+    return systemRuptureSet(json, data, INTERFACE_SYSTEM);
+  }
+
+  private static SystemRuptureSet systemRuptureSet(
+      Path json,
+      ModelData data,
+      SourceType type) {
 
     Path parent = checkNotNull(json.getParent());
     JsonObject obj = jsonObject(json);
 
-    return SystemRuptureSet.builder()
+    return SystemRuptureSet.builder(type)
         .name(obj.get(NAME).getAsString())
         .id(obj.get(ID).getAsInt())
         .modelData(data)
diff --git a/src/main/java/gov/usgs/earthquake/nshmp/model/ModelFiles.java b/src/main/java/gov/usgs/earthquake/nshmp/model/ModelFiles.java
index 2ae613575d860c04963a5bc0ca7a6bd6d38eced0..3517afcc5c2eff6ba07dfbc01445da311a8d388d 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/model/ModelFiles.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/model/ModelFiles.java
@@ -255,6 +255,11 @@ class ModelFiles {
     return read(dir, CLUSTER_SET, data, Deserialize::interfaceClusterSet);
   }
 
+  /* Interface system rupture-set. */
+  static Optional<SystemRuptureSet> readInterfaceSystemRuptureSet(Path dir, ModelData data) {
+    return read(dir, RUPTURE_SET, data, Deserialize::interfaceSystemRuptureSet);
+  }
+
   /* Fault rupture-set. */
   static Optional<FaultRuptureSet> readFaultRuptureSet(Path dir, ModelData data) {
     return read(dir, RUPTURE_SET, data, Deserialize::faultRuptureSet);
@@ -266,8 +271,8 @@ class ModelFiles {
   }
 
   /* Fault system rupture-set. */
-  static Optional<SystemRuptureSet> readSystemRuptureSet(Path dir, ModelData data) {
-    return read(dir, RUPTURE_SET, data, Deserialize::systemRuptureSet);
+  static Optional<SystemRuptureSet> readCrustalSystemRuptureSet(Path dir, ModelData data) {
+    return read(dir, RUPTURE_SET, data, Deserialize::crustalSystemRuptureSet);
   }
 
   /* Grid rupture-set. */
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 256962ed86a183387e5c13850990edfc6bd524bb..ad1ee10d39410083bcfa97c3fc5839f9f1585aec 100644
--- a/src/main/java/gov/usgs/earthquake/nshmp/model/ModelLoader.java
+++ b/src/main/java/gov/usgs/earthquake/nshmp/model/ModelLoader.java
@@ -9,6 +9,7 @@ import static gov.usgs.earthquake.nshmp.model.ModelFiles.TREE_INFO;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.checkGridDataDirectory;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.checkSourceTree;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readCalcConfig;
+import static gov.usgs.earthquake.nshmp.model.ModelFiles.readCrustalSystemRuptureSet;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readDecollementConfig;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readFaultClusterSet;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readFaultConfig;
@@ -23,6 +24,7 @@ import static gov.usgs.earthquake.nshmp.model.ModelFiles.readInterfaceClusterSet
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readInterfaceConfig;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readInterfaceFeatures;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readInterfaceRuptureSet;
+import static gov.usgs.earthquake.nshmp.model.ModelFiles.readInterfaceSystemRuptureSet;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readMfdConfig;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readMfdMap;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readModelInfo;
@@ -33,7 +35,6 @@ import static gov.usgs.earthquake.nshmp.model.ModelFiles.readSiteData;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readSlabConfig;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readSlabRuptureSets;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readSourceTree;
-import static gov.usgs.earthquake.nshmp.model.ModelFiles.readSystemRuptureSet;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readTreeInfo;
 import static gov.usgs.earthquake.nshmp.model.ModelFiles.readZoneConfig;
 import static gov.usgs.earthquake.nshmp.model.SourceType.DECOLLEMENT;
@@ -103,7 +104,8 @@ abstract class ModelLoader {
 
   public static void main(String[] args) throws IOException {
 
-    Path testModel = Paths.get("../nshm-conus-2018-5.2.x-maint");
+    Path testModel = Paths.get("../nshm-prvi");
+    // Path testModel = Paths.get("../nshm-conus-2018-5.2.x-maint");
     // Path testModel = Paths.get("../nshm-conus");
     // Path testModel = Paths.get("../nshm-alaska");
     // Path testModel = Paths.get("../nshm-alaska-2007-2.2.x-maint");
@@ -728,7 +730,7 @@ abstract class ModelLoader {
           processSystemBranch(child, treeBuilder, data);
         }
       } else {
-        SystemRuptureSet srs = readSystemRuptureSet(dir, data).orElseThrow();
+        SystemRuptureSet srs = readCrustalSystemRuptureSet(dir, data).orElseThrow();
         treeBuilder.addLeaf(branch, srs);
         // add features to 'data'
         // TODO: revisit, do we want parent sections??
@@ -913,10 +915,15 @@ abstract class ModelLoader {
 
       System.out.println("    tree: [" + info.id + "] " + root.relativize(dir));
 
-      data.interfaceFeatureMap(readInterfaceFeatures(dir).orElseThrow());
+      /* Load features if 'classic' subduction model */
+      var interfaceSections = readInterfaceFeatures(dir);
+      interfaceSections.ifPresent(data::interfaceFeatureMap);
 
-      List<Feature> features = ModelData.featureMapToList(
-          data.interfaceFeatureMap().orElseThrow().values());
+      /* Set flag to process a logic-tree of fault-system solutions. */
+      boolean interfaceSystemTree = interfaceSections.isEmpty();
+
+      // List<Feature> features = ModelData.featureMapToList(
+      // data.interfaceFeatureMap().orElseThrow().values());
 
       SourceTree.Builder treeBuilder = SourceTree.builder()
           .path(dir)
@@ -925,11 +932,23 @@ abstract class ModelLoader {
           .setting(data.tectonicSetting())
           .type(INTERFACE)
           .gmms(data.gmms())
-          .features(features)
+          // .features(features)
           .root(tree);
 
       for (Branch<Path> branch : tree) {
-        processBranch(branch, treeBuilder, ModelData.copyOf(data));
+        if (interfaceSystemTree) {
+          processSystemBranch(branch, treeBuilder, data);
+        } else {
+          processBranch(branch, treeBuilder, data);
+        }
+      }
+
+      if (data.systemFeatureMap().isPresent()) {
+        treeBuilder.features(ModelData.featureMapToList(
+            data.systemFeatureMap().orElseThrow().values()));
+      } else {
+        treeBuilder.features(ModelData.featureMapToList(
+            data.interfaceFeatureMap().orElseThrow().values()));
       }
 
       return treeBuilder.build();
@@ -1036,6 +1055,35 @@ abstract class ModelLoader {
         treeBuilder.addLeaf(pathTree.get(i), crsList.get(i));
       }
     }
+
+    /* Process fault system source tree branches. */
+    private void processSystemBranch(
+        Branch<Path> branch,
+        SourceTree.Builder treeBuilder,
+        ModelData data) {
+
+      /*
+       * No configuration overrides are currently expected for
+       * fault-system-solutions so this is simpler than above.
+       */
+      Path dir = branch.value();
+      System.out.println("  branch:        " + root.relativize(dir));
+      Optional<LogicTree<Path>> tree = readSourceTree(dir);
+      if (tree.isPresent()) {
+        LogicTree<Path> children = tree.orElseThrow();
+        treeBuilder.addBranches(branch, children);
+        for (Branch<Path> child : children) {
+          processSystemBranch(child, treeBuilder, data);
+        }
+      } else {
+        SystemRuptureSet srs = readInterfaceSystemRuptureSet(dir, data).orElseThrow();
+        treeBuilder.addLeaf(branch, srs);
+        // add features to 'data'
+        // TODO: revisit, do we want parent sections??
+        data.systemFeatureMap(srs.sectionMap);
+      }
+    }
+
   }
 
   static List<SourceTree> gridSources(
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 f4ed39102c33890b73b52a52467d7c0f44039c45..03acc3cfbb8f7060471a379979dc64731e194fcd 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 09167bae487a504d1d27d97bcda20b65c13195c7..a9943527037e70d8e5a1bccffe03b125ab92a8e5 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