Skip to content

Commit

Permalink
simulator experiments separating crustal and subduction rups
Browse files Browse the repository at this point in the history
  • Loading branch information
@kevinmilner
kevinmilner committed Nov 30, 2023
1 parent 7d5be01 commit 23ad85a
Show file tree
Hide file tree
Showing 13 changed files with 1,485 additions and 507 deletions.
36 changes: 36 additions & 0 deletions src/main/java/scratch/kevin/MemTest.java
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,36 @@
package scratch.kevin;

import java.text.DecimalFormat;
import java.util.Random;

import com.google.common.base.Preconditions;

public class MemTest {

public static void main(String[] args) {
long maxMemBytes = Runtime.getRuntime().maxMemory();
long gbToBytes = 1024l*1024l*1024l;
double maxMemGB = (double)maxMemBytes/(double)gbToBytes;
long curMemBytes = Runtime.getRuntime().totalMemory();
double curMemGB = (double)curMemBytes/(double)gbToBytes;
DecimalFormat gbDF = new DecimalFormat("0.00");
System.gc();
System.out.println("Max VM memory: "+gbDF.format(maxMemGB)+" GB");
System.out.println("Starting VM memory: "+gbDF.format(curMemGB)+" GB");
int targetGB = Integer.parseInt(args[0]);
long targetBytes = gbToBytes*targetGB;
long numLongs = targetBytes / 8l;
System.out.println("Will try to allocate "+targetGB+" GB = "+targetBytes+" bytes = "+numLongs+" longs");
Preconditions.checkState(numLongs < (long)Integer.MAX_VALUE);
long[] hugeArray = new long[(int)numLongs];
System.out.println("Done initializing array");
curMemBytes = Runtime.getRuntime().totalMemory();
curMemGB = (double)curMemBytes/(double)gbToBytes;
System.out.println("Current VM memory: "+gbDF.format(curMemGB));
System.gc();
Random r = new Random();
while (true)
hugeArray[r.nextInt(hugeArray.length)] = r.nextLong();
}

}
4 changes: 3 additions & 1 deletion src/main/java/scratch/kevin/bbp/BBP_Module.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,9 @@ public static enum VelocityModel {
LA_BASIN_863("LA Basin 863 (m/s)", "LABasin863",
"$BBP_INSTALL_GF/LABasin863/gp/genslip_nr_generic1d-gp01.vmod", 863, 0.02, 3),
LA_BASIN_500("LA Basin 500 (m/s)", "LABasin500",
"$BBP_INSTALL_GF/LABasin500/gp/nr02-vs500.fk1d", 500, 0.2, 2.5);
"$BBP_INSTALL_GF/LABasin500/gp/nr02-vs500.fk1d", 500, 0.2, 2.5),
CENTRAL_JAPAN("Central Japan", "CentralJapan500",
"$BBP_INSTALL_GF/CentralJapan500/gp/niigata04-vs500.fk1d", 500, 0.1, 0.9);

private String name;
private String xmlName;
Expand Down
22 changes: 14 additions & 8 deletions src/main/java/scratch/kevin/simCompare/MultiRupGMPE_ComparePageGen.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@
import java.util.concurrent.ExecutorService;
import java.util.concurrent.Executors;
import java.util.concurrent.Future;
import java.util.function.Supplier;

import org.apache.commons.math3.stat.descriptive.moment.StandardDeviation;
import org.jfree.chart.annotations.XYPolygonAnnotation;
Expand Down Expand Up @@ -55,6 +56,7 @@
import org.opensha.sha.earthquake.faultSysSolution.ruptures.util.RupSetMapMaker;
import org.opensha.sha.faultSurface.FaultSection;
import org.opensha.sha.imr.AttenRelRef;
import org.opensha.sha.imr.AttenRelSupplier;
import org.opensha.sha.imr.ScalarIMR;
import org.opensha.sha.imr.attenRelImpl.MultiIMR_Averaged_AttenRel;
import org.opensha.sha.imr.param.EqkRuptureParams.MagParam;
Expand Down Expand Up @@ -108,7 +110,7 @@ public abstract class MultiRupGMPE_ComparePageGen<E> {

protected ExecutorService exec;

private Map<AttenRelRef, LinkedList<ScalarIMR>> gmpesInstancesCache;
private Map<Supplier<ScalarIMR>, LinkedList<ScalarIMR>> gmpesInstancesCache;

private boolean replotScatters = true;
private boolean replotZScores = true;
Expand Down Expand Up @@ -335,7 +337,7 @@ public boolean matches(RuptureComparison<E> comp, Site site) {
private static final int MAX_SCATTER_NUM_POINTS = 100000;

public boolean plotScatter(Collection<? extends RuptureComparison<E>> eventComps, Collection<Site> sites, IMT imt,
AttenRelRef gmpe, RuptureComparisonFilter<E> filter, List<String> binDescriptions, File outputDir, String prefix)
AttenRelSupplier gmpe, RuptureComparisonFilter<E> filter, List<String> binDescriptions, File outputDir, String prefix)
throws IOException {
DefaultXY_DataSet xy = new DefaultXY_DataSet();

Expand Down Expand Up @@ -451,7 +453,7 @@ public File call() throws Exception {

}

protected ScalarIMR checkOutGMPE(AttenRelRef gmpeRef) {
protected ScalarIMR checkOutGMPE(Supplier<ScalarIMR> gmpeRef) {
synchronized (gmpesInstancesCache) {
LinkedList<ScalarIMR> gmpes = gmpesInstancesCache.get(gmpeRef);
if (gmpes == null) {
Expand All @@ -461,12 +463,12 @@ protected ScalarIMR checkOutGMPE(AttenRelRef gmpeRef) {
if (!gmpes.isEmpty())
return gmpes.pop();
}
ScalarIMR gmpe = gmpeRef.instance(null);
ScalarIMR gmpe = gmpeRef.get();
gmpe.setParamDefaults();
return gmpe;
}

protected void checkInGMPE(AttenRelRef gmpeRef, ScalarIMR gmpe) {
protected void checkInGMPE(Supplier<ScalarIMR> gmpeRef, ScalarIMR gmpe) {
synchronized (gmpesInstancesCache) {
gmpesInstancesCache.get(gmpeRef).push(gmpe);
}
Expand Down Expand Up @@ -507,7 +509,11 @@ private static List<String> getZscoreTableLine(String siteName, ZScoreResult[] s
return line;
}

protected void checkTransSiteParams(AttenRelRef gmpeRef) {
protected void checkTransSiteParams(Supplier<ScalarIMR> gmpeRef) {
checkTransSiteParams(gmpeRef, sites);
}

public static void checkTransSiteParams(Supplier<ScalarIMR> gmpeRef, List<? extends Site> sites) {
SiteTranslator trans = null;
ScalarIMR gmmInstance = gmpeRef.get();
Site site0 = sites.get(0);
Expand Down Expand Up @@ -543,7 +549,7 @@ protected void checkTransSiteParams(AttenRelRef gmpeRef) {
}
}

public void generateGMPE_Page(File outputDir, List<String> headerLines, AttenRelRef gmpeRef, IMT[] imts,
public void generateGMPE_Page(File outputDir, List<String> headerLines, AttenRelSupplier gmpeRef, IMT[] imts,
List<? extends RuptureComparison<E>> comps, List<Site> highlightSites) throws IOException {
File resourcesDir = new File(outputDir, "resources");
Preconditions.checkState(resourcesDir.exists() || resourcesDir.mkdir());
Expand Down Expand Up @@ -1806,7 +1812,7 @@ private static enum HazardRatioType {
}

public void generateNonErgodicMapPage(File outputDir, List<String> headerLines, List<? extends FaultSection> subSects,
Region siteRegion, Region sourceRegion, AttenRelRef gmpeRef, List<? extends RuptureComparison<E>> comps,
Region siteRegion, Region sourceRegion, AttenRelSupplier gmpeRef, List<? extends RuptureComparison<E>> comps,
Map<E, List<FaultSection>> rupSectMappings, Map<E, FaultSection> rupNuclSects, IMT[] imts,
List<Site> highlightSites) throws IOException {
File resourcesDir = new File(outputDir, "resources");
Expand Down
9 changes: 5 additions & 4 deletions src/main/java/scratch/kevin/simCompare/SimulationHazardPlotter.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -42,6 +42,7 @@
import org.opensha.commons.util.FileUtils;
import org.opensha.commons.util.cpt.CPT;
import org.opensha.sha.imr.AttenRelRef;
import org.opensha.sha.imr.AttenRelSupplier;
import org.opensha.sha.imr.param.IntensityMeasureParams.PGA_Param;
import org.opensha.sha.imr.param.IntensityMeasureParams.PGV_Param;
import org.opensha.sha.imr.param.IntensityMeasureParams.SA_Param;
Expand Down Expand Up @@ -100,7 +101,7 @@ public class SimulationHazardPlotter<E> {
private List<? extends RuptureComparison<E>> comps;
private Site site;
private double curveDuration;
private AttenRelRef gmpeRef;
private AttenRelSupplier gmpeRef;
private Map<String, DiscretizedFunc> xValsMap;

private Table<String, E, Double> sourceRupContribFracts;
Expand All @@ -121,12 +122,12 @@ public class SimulationHazardPlotter<E> {
private Table<IMT, String, DiscretizedFunc> gmpeSourceCurves = HashBasedTable.create();

public SimulationHazardPlotter(SimulationHazardCurveCalc<E> simCalc, List<? extends RuptureComparison<E>> comps,
Site site, double curveDuration, AttenRelRef gmpeRef) {
Site site, double curveDuration, AttenRelSupplier gmpeRef) {
this(simCalc, null, comps, site, curveDuration, gmpeRef);
}

public SimulationHazardPlotter(SimulationHazardCurveCalc<E> simCalc, List<SimulationHazardCurveCalc<?>> compCalcs,
List<? extends RuptureComparison<E>> comps, Site site, double curveDuration, AttenRelRef gmpeRef) {
List<? extends RuptureComparison<E>> comps, Site site, double curveDuration, AttenRelSupplier gmpeRef) {
this.simCalc = simCalc;
if (compCalcs == null)
compCalcs = new ArrayList<>();
Expand All @@ -148,7 +149,7 @@ public Site getSite() {
return site;
}

public AttenRelRef getGMPE() {
public AttenRelSupplier getGMPE() {
return gmpeRef;
}

Expand Down
7 changes: 4 additions & 3 deletions src/main/java/scratch/kevin/simCompare/ZScoreHistPlot.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,7 @@
import org.opensha.commons.mapping.gmt.elements.GMT_CPT_Files;
import org.opensha.commons.util.cpt.CPT;
import org.opensha.sha.imr.AttenRelRef;
import org.opensha.sha.imr.AttenRelSupplier;

import com.google.common.base.Preconditions;
import com.google.common.collect.Table;
Expand Down Expand Up @@ -188,7 +189,7 @@ else if (numMatches < 500)
}

public static <E> boolean plotStandardNormal(SimulationRotDProvider<E> simProv, Collection<? extends RuptureComparison<E>> eventComps,
List<Site> sites, IMT[] imts, AttenRelRef gmpe, RuptureComparisonFilter<E> filter, List<String> binDescriptions,
List<Site> sites, IMT[] imts, AttenRelSupplier gmpe, RuptureComparisonFilter<E> filter, List<String> binDescriptions,
File outputDir, String prefix) throws IOException {
return plotStandardNormal(simProv, eventComps, sites, imts, gmpe, filter, binDescriptions, outputDir, prefix, null, 0, false);
}
Expand All @@ -197,7 +198,7 @@ public static <E> boolean plotStandardNormal(SimulationRotDProvider<E> simProv,
public static final double numStdDev = 3.75;

public static <E> boolean plotStandardNormal(SimulationRotDProvider<E> simProv, Collection<? extends RuptureComparison<E>> eventComps,
List<Site> sites, IMT[] imts, AttenRelRef gmpe, RuptureComparisonFilter<E> filter, List<String> binDescriptions,
List<Site> sites, IMT[] imts, AttenRelSupplier gmpe, RuptureComparisonFilter<E> filter, List<String> binDescriptions,
File outputDir, String prefix, Table<String, E, Double> sourceRupContribFracts, int maxNumSourceContribs, boolean pub)
throws IOException {

Expand All @@ -208,7 +209,7 @@ public static <E> boolean plotStandardNormal(SimulationRotDProvider<E> simProv,
return plotStandardNormal(scores, simProv.getName(), imts, gmpe, binDescriptions, outputDir, prefix, maxNumSourceContribs, pub, maxY);
}

public static <E> boolean plotStandardNormal(ZScoreResult[] scores, String name, IMT[] imts, AttenRelRef gmpe,
public static <E> boolean plotStandardNormal(ZScoreResult[] scores, String name, IMT[] imts, AttenRelSupplier gmpe,
List<String> binDescriptions, File outputDir, String prefix, int maxNumSourceContribs, boolean pub, double maxY)
throws IOException {

Expand Down
74 changes: 74 additions & 0 deletions src/main/java/scratch/kevin/simulators/ruptures/ASK_EventData.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -18,11 +18,16 @@
import org.opensha.commons.data.CSVFile;
import org.opensha.commons.data.function.DiscretizedFunc;
import org.opensha.commons.data.function.HistogramFunction;
import org.opensha.commons.data.xyz.EvenlyDiscrXYZ_DataSet;
import org.opensha.commons.gui.plot.HeadlessGraphPanel;
import org.opensha.commons.gui.plot.PlotCurveCharacterstics;
import org.opensha.commons.gui.plot.PlotLineType;
import org.opensha.commons.gui.plot.PlotSpec;
import org.opensha.commons.gui.plot.PlotUtils;
import org.opensha.commons.gui.plot.jfreechart.xyzPlot.XYZPlotSpec;
import org.opensha.commons.mapping.gmt.elements.GMT_CPT_Files;
import org.opensha.commons.util.DataUtils.MinMaxAveTracker;
import org.opensha.commons.util.cpt.CPT;
import org.opensha.sha.imr.attenRelImpl.ngaw2.FaultStyle;

import com.google.common.base.Preconditions;
Expand Down Expand Up @@ -315,6 +320,64 @@ public static void plotMultiCountHist(File outputDir, String prefix, String titl
gp.saveAsPDF(file.getAbsolutePath()+".pdf");
}

public static void plotMagDistXYZ(File outputDir, String prefix, String title,
Collection<List<ASK_EventData>> data, double minMag, double maxMag, double minDist, double maxDist,
boolean rJB, boolean eventCount, boolean logZ) throws IOException {
int numDist = 20;
int numMag = 15;
double distSpacing = (Math.log10(maxDist) - Math.log10(minDist)) / (double)numDist;
double magSpacing = (maxMag - minMag) / (double)numMag;
EvenlyDiscrXYZ_DataSet xyz = new EvenlyDiscrXYZ_DataSet(numDist, numMag, Math.log10(minDist)+0.5*distSpacing,
minMag+0.5*magSpacing, distSpacing, magSpacing);

for (List<ASK_EventData> eventData : data) {
double mag = eventData.get(0).mag;
if (mag < minMag || mag > maxMag)
continue;
int magBin = xyz.getYIndex(mag);
int[] distCounts = new int[xyz.getNumX()];
for (ASK_EventData event : eventData) {
double dist = rJB ? event.rJB : event.rRup;
if (dist < minDist || dist > maxDist)
continue;
double logDist = Math.log10(dist);
int distBin = xyz.getXIndex(logDist);
if (eventCount)
distCounts[distBin] = 1;
else
distCounts[distBin]++;
}
for (int x=0; x<distCounts.length; x++)
if (distCounts[x] > 0)
xyz.set(x, magBin, xyz.get(x, magBin)+distCounts[x]);
}

CPT cpt = GMT_CPT_Files.BLACK_RED_YELLOW_UNIFORM.instance().reverse();
if (logZ) {
xyz.log10();
cpt = cpt.rescale(0d, Math.ceil(xyz.getMaxZ()));
cpt.setNanColor(Color.WHITE);
} else {
cpt = cpt.rescale(1d, xyz.getMaxZ());
}
cpt.setBelowMinColor(Color.WHITE);

String xLabel = rJB ? "Log10 Distance JB" : "Log10 Distance Rup";
String yLabel = "Magnitude";
String zLabel = eventCount ? "Num Events" : "Num Recordings";
if (logZ)
zLabel = "Lot10 "+zLabel;

XYZPlotSpec xyzSpec = new XYZPlotSpec(xyz, cpt, title, xLabel, yLabel, zLabel);
xyzSpec.setCPTTickUnit(0.5d);
HeadlessGraphPanel xyzGP = PlotUtils.initHeadless();
xyzGP.drawGraphPanel(xyzSpec, false, false, new org.jfree.data.Range(Math.log10(minDist), Math.log10(maxDist)),
new org.jfree.data.Range(minMag, maxMag));
// xyzGP.getChartPanel().getChart().setBackgroundPaint(Color.WHITE);
xyzGP.getChartPanel().setSize(700, 550);
PlotUtils.writePlots(outputDir, prefix, xyzGP, 700, 550, true, true, false);
}

private static final DecimalFormat twoDigit = new DecimalFormat("0.00");

public static void main(String[] args) throws IOException {
Expand Down Expand Up @@ -361,6 +424,17 @@ public static void main(String[] args) throws IOException {
System.out.println("Loaded "+totNum+" total records");
System.out.println();
periodData.put(period, data);

double minDist = 1;
double maxDist = 200;
double minMag = 6.5d;
double maxMag = 8.2d;

String perStr = new DecimalFormat("0.#").format(period);
String title = "ASK (2014) Data, "+perStr+"s";

plotMagDistXYZ(new File("/tmp"), "ask_mag_dist_"+perStr+"s", title, data.values(),
minMag, maxMag, minDist, maxDist, false, false, true);
}

for (Range<Double> magRange : magRanges) {
Expand Down
Loading

0 comments on commit 23ad85a

Please sign in to comment.