From cf9a330e97ca20536e781d13eecedb58bd1315a5 Mon Sep 17 00:00:00 2001 From: Kevin Milner Date: Thu, 7 Nov 2024 10:59:41 -0800 Subject: [PATCH] improvements to residual plots; XYZ now with log axis --- .../kevin/quantum/ToyRuptureSetBuilder.java | 34 ++++++++-- .../MultiRupGMPE_ComparePageGen.java | 19 ++++-- .../kevin/simCompare/ResidualScatterPlot.java | 64 ++++++++++--------- .../ruptures/CatalogGMPE_Compare.java | 10 +-- 4 files changed, 82 insertions(+), 45 deletions(-) diff --git a/src/main/java/scratch/kevin/quantum/ToyRuptureSetBuilder.java b/src/main/java/scratch/kevin/quantum/ToyRuptureSetBuilder.java index 82365507..edb100d1 100644 --- a/src/main/java/scratch/kevin/quantum/ToyRuptureSetBuilder.java +++ b/src/main/java/scratch/kevin/quantum/ToyRuptureSetBuilder.java @@ -4,7 +4,9 @@ import java.io.File; import java.io.FileOutputStream; import java.io.IOException; +import java.text.SimpleDateFormat; import java.util.ArrayList; +import java.util.Date; import java.util.List; import org.opensha.commons.data.CSVWriter; @@ -35,10 +37,27 @@ public class ToyRuptureSetBuilder { public static void main(String[] args) throws IOException { boolean rebuildRupSet = false; - File outputDir = new File("/home/kevin/markdown/inversions/2024_10_31-quantum_test_rup_set_problem"); + String dirName = new SimpleDateFormat("yyyy_MM_dd").format(new Date()); +// String dirName = "2024_11_06"; + + dirName += "-quantum_test_rup_set_problem"; + +// String[] nameSearches = { "Likely" }; +// dirName += "-70_sects"; + +// String[] nameSearches = { "Cedar", "Mountain", "Mahogany" }; +// dirName += "-124_sects"; + + String[] nameSearches = { "San", "Felipe", "proxy" }; + dirName += "-213_sects"; + + System.out.println(dirName); + + File outputDir = new File("/home/kevin/markdown/inversions/"+dirName); Preconditions.checkState(outputDir.exists() || outputDir.mkdir()); File rsFile = new File(outputDir, "orig_rup_set.zip"); + File exhaustiveFile = new File(outputDir, "exhaustive_rup_set.zip"); FaultSystemRupSet origRupSet; if (rsFile.exists() && !rebuildRupSet) { origRupSet = FaultSystemRupSet.load(rsFile); @@ -51,7 +70,7 @@ public static void main(String[] args) throws IOException { System.out.println("Building connectivity clusters"); ConnectivityClusters clusters = ConnectivityClusters.build(origRupSet); - int targetParentID = FaultSectionUtils.findParentSectionID(origRupSet.getFaultSectionDataList(), "Likely"); + int targetParentID = FaultSectionUtils.findParentSectionID(origRupSet.getFaultSectionDataList(), nameSearches); FaultSection targetSect = null; for (FaultSection sect : origRupSet.getFaultSectionDataList()) { if (sect.getParentSectionId() == targetParentID) { @@ -116,7 +135,14 @@ public static void main(String[] args) throws IOException { // if true, allow splays (using default settings) rsConfig.setSplays(true); // enabled - FaultSystemRupSet rupSet = rsConfig.build(FaultSysTools.defaultNumThreads()); + FaultSystemRupSet rupSet; + if (exhaustiveFile.exists() && !rebuildRupSet) { + rupSet = FaultSystemRupSet.load(exhaustiveFile); + } else { + rupSet = rsConfig.build(FaultSysTools.defaultNumThreads()); + + rupSet.write(exhaustiveFile); + } System.out.println("Original had "+origRupSet.getNumRuptures()); System.out.println("New has "+rupSet.getNumRuptures()); @@ -171,8 +197,6 @@ public static void main(String[] args) throws IOException { // now write the actual ruptures writeRuptures(rupSet.getSectionIndicesForAllRups(), new File(outputDir, "orig_ruptures.csv"), subSects.size()); writeRuptures(rupSet.getSectionIndicesForAllRups(), new File(outputDir, "exhaustive_ruptures.csv"), subSects.size()); - - rupSet.write(new File(outputDir, "exhaustive_rup_set.zip")); } private static void writeMatrix(double[][] mat, File file) throws IOException { diff --git a/src/main/java/scratch/kevin/simCompare/MultiRupGMPE_ComparePageGen.java b/src/main/java/scratch/kevin/simCompare/MultiRupGMPE_ComparePageGen.java index 5ac2ff05..9a7c6a72 100644 --- a/src/main/java/scratch/kevin/simCompare/MultiRupGMPE_ComparePageGen.java +++ b/src/main/java/scratch/kevin/simCompare/MultiRupGMPE_ComparePageGen.java @@ -634,14 +634,18 @@ public void generateGMPE_Page(File outputDir, List headerLines, AttenRel lines.addAll(table.build()); lines.add(""); - if (!siteBundles.isEmpty() && sites.size() > 1) { - System.out.println("All sites aggregated"); + if (!siteBundles.isEmpty() && this.sites.size() > 1) { + System.out.println("Bundled sites aggregated"); lines.add("### Bundled z-scores"); lines.add(topLink); lines.add(""); for (int s=0; s=", "_ge_").replace(">", "_greater_").replaceAll("\\W+", "_"); + while (sanitizedName.contains("__")) + sanitizedName = sanitizedName.replace("__", "_"); + String prefix = "site_bundle_"+sanitizedName+"_std_norm"; System.out.println("Processing site bundle: "+name); @@ -1054,9 +1058,12 @@ public void generateGMPE_Page(File outputDir, List headerLines, AttenRel List residualTypes = new ArrayList<>(); residualTypes.add(ResidualType.MAG); ScalarIMR tempGMPE = checkOutGMPE(gmpeRef); - if (tempGMPE.getPropagationEffectParams().containsParameter(DistanceJBParameter.NAME)) + ScalarIMR subGMPE = tempGMPE instanceof MultiIMR_Averaged_AttenRel ? ((MultiIMR_Averaged_AttenRel)tempGMPE).getIMRs().get(0) : null; + if (tempGMPE.getPropagationEffectParams().containsParameter(DistanceJBParameter.NAME) + || (subGMPE != null && subGMPE.getPropagationEffectParams().containsParameter(DistanceJBParameter.NAME))) residualTypes.add(ResidualType.DIST_JB); - if (tempGMPE.getPropagationEffectParams().containsParameter(DistanceRupParameter.NAME)) + if (tempGMPE.getPropagationEffectParams().containsParameter(DistanceRupParameter.NAME) + || (subGMPE != null && subGMPE.getPropagationEffectParams().containsParameter(DistanceRupParameter.NAME))) residualTypes.add(ResidualType.DIST_RUP); checkInGMPE(gmpeRef, tempGMPE); sites = new ArrayList<>(this.sites); // we previously added a null element for all sites above @@ -1389,7 +1396,7 @@ public Double getValue(RuptureComparison comp, Site site) { return comp.getDistanceJB(site); } }, - VS30("Vs30", "m/s", Vs30_Param.NAME, true, 200d, Double.NaN) { + VS30("Vs30", "m/s", Vs30_Param.NAME, false, 200d, Double.NaN) { @Override public Double getValue(RuptureComparison comp, Site site) { return site.getParameter(Double.class, Vs30_Param.NAME).getValue(); diff --git a/src/main/java/scratch/kevin/simCompare/ResidualScatterPlot.java b/src/main/java/scratch/kevin/simCompare/ResidualScatterPlot.java index 90be2b31..c44db1e6 100644 --- a/src/main/java/scratch/kevin/simCompare/ResidualScatterPlot.java +++ b/src/main/java/scratch/kevin/simCompare/ResidualScatterPlot.java @@ -150,7 +150,7 @@ else if (minX > myMinX+halfDelta) else xRange = new Range(minX, maxX); - SimpleRegression linearRegression = getRegression(xy, logX, xRange); + SimpleRegression linearRegression = getRegression(xy, logX, xRange, true, true); regressionIntercept = linearRegression.getIntercept(); regressionLinearSlope = linearRegression.getSlope(); System.out.println("Linear regression yInt: "+linearRegression.getIntercept()); @@ -168,7 +168,7 @@ else if (minX > myMinX+halfDelta) Range regressRange = new Range(binCenter - 0.5*deltaX, binCenter + 0.5*deltaX); if (logX) regressRange = new Range(Math.pow(10, regressRange.getLowerBound()), Math.pow(10, regressRange.getUpperBound())); - regressionBinned[i] = getRegressionFit(xy, logX, regressRange); + regressionBinned[i] = getRegressionFit(xy, logX, regressRange, i == 0, i == numRegressionBins-1); if (regressionBinned[i] != null) dataSigmaBinned[i] = calcSigma(regressionBinned[i], xy); } @@ -239,11 +239,11 @@ public void plotScatter(File outputDir, String prefix) throws IOException { funcs.add(zeroLine); chars.add(new PlotCurveCharacterstics(PlotLineType.DASHED, 3f, Color.DARK_GRAY)); - addResidualFuncs(funcs, chars, true); + addResidualFuncs(funcs, chars); PlotSpec spec = new PlotSpec(funcs, chars, title, xAxisLabel, residualLabel); spec.setLegendVisible(true); - List anns = buildAnnotations(false); + List anns = buildAnnotations(); if (isFiltered) { // add filtered note Font font = new Font(Font.SANS_SERIF, Font.BOLD, 20); @@ -355,8 +355,9 @@ else if (!Doubles.isFinite(minZ)) ArrayList chars = new ArrayList<>(); Range xRange = this.xRange; - if (logX) - xRange = new Range(Math.log10(xRange.getLowerBound()), Math.log10(xRange.getUpperBound())); + // now plotting with log axis; XYZ is still in log spacing, but range in linear +// if (logX) +// xRange = new Range(Math.log10(xRange.getLowerBound()), Math.log10(xRange.getUpperBound())); DefaultXY_DataSet zeroLine = new DefaultXY_DataSet(); zeroLine.set(xRange.getLowerBound(), 0d); @@ -365,7 +366,7 @@ else if (!Doubles.isFinite(minZ)) funcs.add(0, zeroLine); chars.add(0, new PlotCurveCharacterstics(PlotLineType.DASHED, 3f, Color.DARK_GRAY)); - addResidualFuncs(funcs, chars, false); + addResidualFuncs(funcs, chars); // add custom y axis guides for (int i=(int)Math.round(yRange.getLowerBound()); i<=(int)Math.round(yRange.getUpperBound()); i++) { @@ -381,7 +382,7 @@ else if (!Doubles.isFinite(minZ)) xyzSpec.setXYElems(funcs); xyzSpec.setXYChars(chars); - xyzSpec.setPlotAnnotations(buildAnnotations(true)); + xyzSpec.setPlotAnnotations(buildAnnotations()); PlotPreferences plotPrefs = PlotPreferences.getDefault(); plotPrefs.setTickLabelFontSize(18); @@ -390,7 +391,7 @@ else if (!Doubles.isFinite(minZ)) plotPrefs.setBackgroundColor(Color.WHITE); HeadlessGraphPanel xyzGP = new HeadlessGraphPanel(plotPrefs); System.out.println("Plotting with xRange: "+xRange); - xyzGP.drawGraphPanel(xyzSpec, false, false, xRange, yRange); + xyzGP.drawGraphPanel(xyzSpec, logX, false, xRange, yRange); // write plot xyzGP.getYAxis().setStandardTickUnits(getYTickUnits()); xyzGP.getChartPanel().setSize(plotWidth, plotWidth-100); @@ -399,7 +400,7 @@ else if (!Doubles.isFinite(minZ)) xyzGP.saveAsPDF(new File(outputDir, prefix+".pdf").getAbsolutePath()); } - private List buildAnnotations(boolean doLog) { + private List buildAnnotations() { List anns = new ArrayList<>(); double yPos = yRange.getUpperBound() * 0.95; @@ -412,7 +413,8 @@ private List buildAnnotations(boolean doLog) { xMin = Math.log10(xMin); } double xPos = xMin + 0.95*(xMax - xMin); - if (logX && !doLog) +// if (logX && !doLog) + if (logX) xPos = Math.pow(10, xPos); Font font = new Font(Font.SANS_SERIF, Font.BOLD, 24); @@ -448,7 +450,8 @@ private List buildAnnotations(boolean doLog) { if (plotLinearResidual) { yPos = yRange.getUpperBound() * 0.95; xPos = xMin + 0.05*(xMax - xMin); - if (logX && !doLog) +// if (logX && !doLog) + if (logX) xPos = Math.pow(10, xPos); ann = new XYTextAnnotation("Linear LSQ Fit", xPos, yPos); @@ -478,7 +481,7 @@ private List buildAnnotations(boolean doLog) { private static final DecimalFormat annDF = new DecimalFormat("0.00"); private static final DecimalFormat regressionDF = new DecimalFormat("0.00E0"); - private void addResidualFuncs(List funcs, List chars, boolean scatter) { + private void addResidualFuncs(List funcs, List chars) { Color linearColor = Color.BLACK; Color binnedColor = Color.BLACK; Color gmpeColor = new Color(35, 72, 132); @@ -489,9 +492,9 @@ private void addResidualFuncs(List funcs, List funcs, List= xRange.getUpperBound(); + if (contains) { + if (logX) + pt = new Point2D.Double(Math.log10(pt.getX()), pt.getY()); + regression.addData(pt.getX(), pt.getY()); + } } if (regression.getN() == 0l) return null; return regression; } - private static ArbitrarilyDiscretizedFunc getRegressionFit(XY_DataSet xy, boolean logX, Range xRange) { - SimpleRegression regression = getRegression(xy, logX, xRange); + private static ArbitrarilyDiscretizedFunc getRegressionFit(XY_DataSet xy, boolean logX, Range xRange, boolean includeBelow, boolean includeAbove) { + SimpleRegression regression = getRegression(xy, logX, xRange, includeBelow, includeAbove); if (regression == null) return null; return getRegressionFit(regression, xy, logX, xRange); diff --git a/src/main/java/scratch/kevin/simulators/ruptures/CatalogGMPE_Compare.java b/src/main/java/scratch/kevin/simulators/ruptures/CatalogGMPE_Compare.java index 3b0ecc1d..6b44993e 100644 --- a/src/main/java/scratch/kevin/simulators/ruptures/CatalogGMPE_Compare.java +++ b/src/main/java/scratch/kevin/simulators/ruptures/CatalogGMPE_Compare.java @@ -636,16 +636,16 @@ public static void main(String[] args) throws ZipException, IOException { // RuptureIdentifier loadIden = FocalMechIden.builder().forRake(-105, -75).forDip(35, 55).build(); // String loadIdenPrefix = "mech_normal"; - boolean replotScatters = true; - boolean replotZScores = true; - boolean replotCurves = true; - boolean replotResiduals = true; - // boolean replotScatters = true; // boolean replotZScores = true; // boolean replotCurves = true; // boolean replotResiduals = true; + boolean replotScatters = false; + boolean replotZScores = false; + boolean replotCurves = false; + boolean replotResiduals = true; + if (timeScale != 1d) doRotD = false;