Skip to content

Commit

Permalink
Merge branch 'master' into 2024_09-point_source_dist_corr_refactor
Browse files Browse the repository at this point in the history
  • Loading branch information
Kevin Milner committed Dec 9, 2024
2 parents f7382dc + f5222a1 commit a4aab64
Show file tree
Hide file tree
Showing 16 changed files with 2,109 additions and 216 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -625,35 +625,35 @@ public static void main(String[] args) throws IOException {
* PRVI25 logic tree
* TODO (this is a just a marker to find this part quickly, not an actual todo)
*/
// List<LogicTreeLevel<? extends LogicTreeNode>> levels = PRVI25_LogicTreeBranch.levelsOnFault;
// dirName += "-prvi25_crustal_branches";
// double avgNumRups = 50000;
// gmpes = new AttenRelRef[] { AttenRelRef.USGS_PRVI_ACTIVE };
//
// // random DM sampling
// levels = new ArrayList<>(levels);
// int origNumLevels = levels.size();
// for (int i=levels.size(); --i>=0;)
// if (levels.get(i).getNodes().get(0) instanceof PRVI25_CrustalDeformationModels)
// levels.remove(i);
// Preconditions.checkState(levels.size() == origNumLevels -1);
// individualRandomLevels.add(new PRVI25_CrustalRandomlySampledDeformationModelLevel());
// samplingBranchCountMultiplier = 5; // 5 for each branch
// dirName += "-dmSample";
// if (samplingBranchCountMultiplier > 1)
// dirName += samplingBranchCountMultiplier+"x";
List<LogicTreeLevel<? extends LogicTreeNode>> levels = PRVI25_LogicTreeBranch.levelsOnFault;
dirName += "-prvi25_crustal_branches";
double avgNumRups = 50000;
gmpes = new AttenRelRef[] { AttenRelRef.USGS_PRVI_ACTIVE };

// random DM sampling
levels = new ArrayList<>(levels);
int origNumLevels = levels.size();
for (int i=levels.size(); --i>=0;)
if (levels.get(i).getNodes().get(0) instanceof PRVI25_CrustalDeformationModels)
levels.remove(i);
Preconditions.checkState(levels.size() == origNumLevels -1);
individualRandomLevels.add(new PRVI25_CrustalRandomlySampledDeformationModelLevel());
samplingBranchCountMultiplier = 5; // 5 for each branch
dirName += "-dmSample";
if (samplingBranchCountMultiplier > 1)
dirName += samplingBranchCountMultiplier+"x";

List<LogicTreeLevel<? extends LogicTreeNode>> levels = PRVI25_LogicTreeBranch.levelsSubduction;
dirName += "-prvi25_subduction_branches";
double avgNumRups = 10000;
gmpes = new AttenRelRef[] { AttenRelRef.USGS_PRVI_INTERFACE, AttenRelRef.USGS_PRVI_SLAB };
// List<LogicTreeLevel<? extends LogicTreeNode>> levels = PRVI25_LogicTreeBranch.levelsSubduction;
// dirName += "-prvi25_subduction_branches";
// double avgNumRups = 10000;
// gmpes = new AttenRelRef[] { AttenRelRef.USGS_PRVI_INTERFACE, AttenRelRef.USGS_PRVI_SLAB };

// levels = new ArrayList<>(levels);
// levels.add(NSHM23_LogicTreeBranch.SUB_SECT_CONSTR);

// dirName += "-proxyGriddedTests";

// Class<? extends InversionConfigurationFactory> factoryClass = PRVI25_InvConfigFactory.class;
Class<? extends InversionConfigurationFactory> factoryClass = PRVI25_InvConfigFactory.class;

// Class<? extends InversionConfigurationFactory> factoryClass = PRVI25_InvConfigFactory.GriddedUseM1Bounds.class;
// dirName += "-grid_bounds_m1";
Expand All @@ -670,8 +670,8 @@ public static void main(String[] args) throws IOException {
// Class<? extends InversionConfigurationFactory> factoryClass = PRVI25_InvConfigFactory.RateBalanceAndLimitCrustalBelowObserved_0p9.class;
// dirName += "-limit_below_obs_constraint-grided_rate_balancing";

Class<? extends InversionConfigurationFactory> factoryClass = PRVI25_InvConfigFactory.GriddedForceSlab2Depths.class;
dirName += "-gridded_use_slab2";
// Class<? extends InversionConfigurationFactory> factoryClass = PRVI25_InvConfigFactory.GriddedForceSlab2Depths.class;
// dirName += "-gridded_use_slab2";

if (!factoryClass.equals(PRVI25_InvConfigFactory.class)) {
// try instantiate it to make sure we get any static modifiers that might change branch weights
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -483,10 +483,12 @@ public GeographicMapMaker buildMapPlot(boolean labelMagnitudes, boolean plotVect
if (result == integrals.get(0)) {
double integralMiddleLat = 0.5*(result.startLoc.lat + result.endLoc.lat);
double integralMiddleLon = 0.5*(result.startLoc.lon + result.endLoc.lon);
double furthestLeftLon = Math.max(lowerLeft.lon+0.25, lowerLeft.lon + 0.05*(upperRight.lon - lowerLeft.lon));
System.out.println("Furthest left="+furthestLeftLon+", integral left="+integralMiddleLon);
XYTextAnnotation ann = new XYTextAnnotation("Summation Direction", integralMiddleLon, integralMiddleLat);
ann.setFont(textAnnFont);
ann.setTextAnchor(TextAnchor.TOP_CENTER);
ann.setRotationAnchor(TextAnchor.TOP_CENTER);
ann.setTextAnchor(TextAnchor.BOTTOM_CENTER);
ann.setRotationAnchor(TextAnchor.BOTTOM_CENTER);
double vectorAngle = result.azimuth;
while (vectorAngle > 180d)
vectorAngle -= 360d;
Expand All @@ -496,6 +498,13 @@ public GeographicMapMaker buildMapPlot(boolean labelMagnitudes, boolean plotVect
ann.setTextAnchor(TextAnchor.BOTTOM_CENTER);
ann.setRotationAnchor(TextAnchor.BOTTOM_CENTER);
}
// } else if (integralMiddleLon < furthestLeftLon) {
//// } else {
// // too far to the left, put it on the other side
//// vectorAngle -= 180;
// ann.setTextAnchor(TextAnchor.BOTTOM_CENTER);
// ann.setRotationAnchor(TextAnchor.BOTTOM_CENTER);
// }
double angleRad = Math.toRadians(vectorAngle + 90);
ann.setRotationAngle(mapMaker.getRotationAngleCorrectedForAspectRatio(angleRad));
mapMaker.addAnnotation(ann);
Expand Down Expand Up @@ -587,7 +596,7 @@ public GeographicMapMaker buildMapPlot(boolean labelMagnitudes, boolean plotVect
arrowLen = 0.2*maxVectDist;
arrowThickness = 4f;
} else {
arrowLen = 0.2*maxVectDist;
arrowLen = 0.15*maxVectDist;
arrowThickness = 2f;
}
mapMaker.setFillArrowheads(true);
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -33,7 +33,8 @@ public static void main(String[] args) throws IOException {
// PRVI25_SubductionFaultModels fm = PRVI25_SubductionFaultModels.PRVI_SUB_FM_SMALL;
PRVI25_SubductionDeformationModels dm = PRVI25_SubductionDeformationModels.FULL;

File outputDir = new File("/home/kevin/Documents/papers/2024_PRVI_Subduction/figures/fault_model");
// File outputDir = new File("/home/kevin/Documents/papers/2024_PRVI_Subduction/figures/fault_model");
File outputDir = new File("/home/kevin/Documents/papers/2024_PRVI_ERF/prvi25-erf-paper/Figures/sub_fm");
Preconditions.checkState(outputDir.exists() || outputDir.mkdir());

// first just write out the fault model
Expand Down
20 changes: 13 additions & 7 deletions src/main/java/scratch/kevin/prvi25/SubductionDefModConvert.java
Original file line number Diff line number Diff line change
Expand Up @@ -95,8 +95,9 @@ public static void main(String[] args) throws IOException {
boolean interpolate = true;
boolean interpSymmetry = true;

File debugDir = new File("/home/kevin/Documents/papers/2024_PRVI_Subduction/figures/def_model");
Preconditions.checkState(debugDir.exists() || debugDir.mkdir());
// File plotOutputDir = new File("/home/kevin/Documents/papers/2024_PRVI_Subduction/figures/def_model");
File plotOutputDir = new File("/home/kevin/Documents/papers/2024_PRVI_ERF/prvi25-erf-paper/Figures/sub_dm");
Preconditions.checkState(plotOutputDir.exists() || plotOutputDir.mkdir());

Table<Boolean, Boolean, String> debugPrefixes = HashBasedTable.create();

Expand Down Expand Up @@ -138,7 +139,7 @@ public static void main(String[] args) throws IOException {
String fmTitle = geomTitle;
String fmDmTitle = geomTitle+", "+rateTitle;
debugPrefixes.put(fullRate, largePolys, prefix);
debugWriteSectOrders(inFeatures, debugDir, prefix, fmTitle, fmDmTitle,
debugWriteSectOrders(inFeatures, plotOutputDir, prefix, fmTitle, fmDmTitle,
ratePropName, rateUncertPropName, rakePropName);

Map<String, int[]> stitchInIDs = new HashMap<>();
Expand Down Expand Up @@ -447,7 +448,7 @@ public static void main(String[] args) throws IOException {

// plot them
String interpPrefix = prefix+"_interp_"+name.replaceAll("\\W+", "_");
debugWriteInterpFuncsOrders(debugDir, interpPrefix, name,
debugWriteInterpFuncsOrders(plotOutputDir, interpPrefix, name,
dasSlipFunc, dasSlipUncertFunc, dasRakeFunc, sectStartDASs,
sectMidDASs, interpRanges);
}
Expand Down Expand Up @@ -531,7 +532,7 @@ public static void main(String[] args) throws IOException {
lines.add(tocIndex, "## Table Of Contents");

// write markdown
MarkdownUtils.writeReadmeAndHTML(lines, debugDir);
MarkdownUtils.writeReadmeAndHTML(lines, plotOutputDir);
}

private static final DecimalFormat moDF = new DecimalFormat("0.00E0");
Expand Down Expand Up @@ -652,8 +653,13 @@ private static void debugWriteSectOrders(List<Feature> features, File outputDir,
for (int s=0; s<sects.size(); s++) {
arrows.addAll(PRVI_SubductionSubSectPlots.buildRakeArrows(sects.get(s), slipCPT.getMaxValue()));
}
System.out.println("Plotting "+arrows.size()+" arrow lines");
mapMaker.plotLines(arrows, Color.BLACK, 2f);
System.out.println("Plotting "+arrows.size()+" arrow lines for "+prefix);
// mapMaker.plotLines(arrows, Color.BLACK, 2f);
mapMaker.plotArrows(arrows, 20d, Color.BLACK, 2f);
mapMaker.setFillArrowheads(true);

mapMaker.plotSectScalars(slips, slipCPT, null);
mapMaker.plot(outputDir, prefix+"_slips_no_cpt", fmDmTitle);
}

mapMaker.plotSectScalars(slips, slipCPT, "Slip Rate (mm/yr)");
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -132,11 +132,11 @@ public static void main(String[] args) throws IOException {
PlotUtils.writePlots(outputDir, "combined_mfds_cml", gp, 800, 750, true, false, false);
}

private static IncrementalMagFreqDist calcFaultMFD(Region region, FaultSystemSolution sol, EvenlyDiscretizedFunc refMFD) {
static IncrementalMagFreqDist calcFaultMFD(Region region, FaultSystemSolution sol, EvenlyDiscretizedFunc refMFD) {
return sol.calcNucleationMFD_forRegion(region, refMFD.getMinX(), refMFD.getMaxX(), refMFD.getDelta(), false);
}

private static IncrementalMagFreqDist calcGriddedMFD(Region region, TectonicRegionType trt,
static IncrementalMagFreqDist calcGriddedMFD(Region region, TectonicRegionType trt,
FaultSystemSolution sol, EvenlyDiscretizedFunc refMFD) {
GridSourceProvider gridProv = sol.getGridSourceProvider();
SummedMagFreqDist ret = new SummedMagFreqDist(refMFD.getMinX(), refMFD.getMaxX(), refMFD.size());
Expand All @@ -149,7 +149,7 @@ private static IncrementalMagFreqDist calcGriddedMFD(Region region, TectonicRegi
return ret;
}

private static IncrementalMagFreqDist average(IncrementalMagFreqDist mfd1, IncrementalMagFreqDist mfd2) {
static IncrementalMagFreqDist average(IncrementalMagFreqDist mfd1, IncrementalMagFreqDist mfd2) {
IncrementalMagFreqDist ret = new IncrementalMagFreqDist(mfd1.getMinX(), mfd1.size(), mfd1.getDelta());
for (int i=0; i<ret.size(); i++)
ret.set(i, 0.5*mfd1.getY(i) + 0.5*mfd2.getY(i));
Expand Down
Original file line number Diff line number Diff line change
@@ -0,0 +1,80 @@
package scratch.kevin.prvi25.figures;

import java.io.File;
import java.io.IOException;
import java.util.HashSet;
import java.util.List;
import java.util.Set;
import java.util.concurrent.TimeUnit;

import org.opensha.commons.geo.Location;
import org.opensha.commons.geo.Region;
import org.opensha.commons.gui.plot.GeographicMapMaker;
import org.opensha.sha.earthquake.faultSysSolution.FaultSystemRupSet;
import org.opensha.sha.earthquake.faultSysSolution.FaultSystemSolution;
import org.opensha.sha.earthquake.faultSysSolution.modules.ConnectivityClusters;
import org.opensha.sha.earthquake.faultSysSolution.reports.plots.FaultSectionConnectionsPlot;
import org.opensha.sha.earthquake.faultSysSolution.ruptures.util.ConnectivityCluster;
import org.opensha.sha.faultSurface.FaultSection;

import com.google.common.base.Stopwatch;

class ConnectivityClusterPlot {

public static void main(String[] args) throws IOException {
FaultSystemSolution sol = FaultSystemSolution.load(new File("/data/kevin/nshm23/batch_inversions/"
+ "2024_11_19-prvi25_crustal_branches-dmSample5x/results_PRVI_CRUSTAL_FM_V1p1_branch_averaged.zip"));

File outputDir = new File("/home/kevin/Documents/papers/2024_PRVI_ERF/prvi25-erf-paper/Figures/crustal_fm");

FaultSystemRupSet rupSet = sol.getRupSet();
if (!rupSet.hasModule(ConnectivityClusters.class)) {
System.out.println("Calculating connection clusters");
Stopwatch watch = Stopwatch.createStarted();
ConnectivityClusters clusters = ConnectivityClusters.build(rupSet);
watch.stop();
System.out.println("Found "+clusters.size()+" connectivity clusters in "+(float)(watch.elapsed(TimeUnit.MILLISECONDS)/1000)+" s");
rupSet.addModule(clusters);
}
List<ConnectivityCluster> clustersUnsorted = rupSet.requireModule(ConnectivityClusters.class).get();

// see if any clusters rupture completely
System.out.println("Looking for clusters that rupture fully...");
boolean anyFound = false;
for (ConnectivityCluster cluster : clustersUnsorted) {
HashSet<Integer> parents = new HashSet<>(cluster.getParentSectIDs());
if (parents.size() == 1)
continue;
Set<Integer> sects = cluster.getSectIDs();
for (List<Integer> rupSects : rupSet.getSectionIndicesForAllRups()) {
if (rupSects.size() == sects.size() && rupSects.containsAll(sects)) {
System.out.println("Cluster ruptures fully: "+sects);
anyFound = true;
System.out.println("Parents:");
for (int sectID: rupSects) {
FaultSection sect = rupSet.getFaultSectionData(sectID);
if (parents.contains(sect.getParentSectionId())) {
System.out.println("\t"+sect.getParentSectionId()+". "+sect.getParentSectionName());
parents.remove(sect.getParentSectionId());
}
}
}
}
}
if (!anyFound)
System.out.println("No clusters rupture fully");

GeographicMapMaker plotter = FaultSectionConnectionsPlot.buildConnectedClustersPlot(rupSet, sol,
SlipRateFigures.CRUSTAL_FAULT_MAP_REG, null, clustersUnsorted);

plotter.setWritePDFs(true);
plotter.setWriteGeoJSON(false);

plotter.plot(outputDir, "connectivity_clusters", " ");

CrustalFaultNamesFigure.addStandardFaultLabels(plotter, rupSet.getFaultSectionDataList());

plotter.plot(outputDir, "connectivity_clusters_with_fault_names", " ");
}

}
Loading

0 comments on commit a4aab64

Please sign in to comment.