mschaefer@9295: /** Copyright (C) 2017 by Bundesanstalt für Gewässerkunde mschaefer@9295: * Software engineering by mschaefer@9295: * Björnsen Beratende Ingenieure GmbH mschaefer@9295: * Dr. Schumacher Ingenieurbüro für Wasser und Umwelt mschaefer@9295: * mschaefer@9295: * This file is Free Software under the GNU AGPL (>=v3) mschaefer@9295: * and comes with ABSOLUTELY NO WARRANTY! Check out the mschaefer@9295: * documentation coming with Dive4Elements River for details. mschaefer@9295: */ mschaefer@9295: package org.dive4elements.river.artifacts.uinfo.salix; mschaefer@9295: mschaefer@9375: import java.math.BigDecimal; mschaefer@9295: import java.util.ArrayList; mschaefer@9397: import java.util.Collection; mschaefer@9295: import java.util.List; mschaefer@9316: import java.util.Map.Entry; mschaefer@9309: import java.util.NavigableMap; mschaefer@9295: mschaefer@9295: import org.dive4elements.river.artifacts.WINFOArtifact; mschaefer@9295: import org.dive4elements.river.artifacts.access.ComputationRangeAccess; mschaefer@9397: import org.dive4elements.river.artifacts.common.AbstractResultType; mschaefer@9295: import org.dive4elements.river.artifacts.common.GeneralResultType; mschaefer@9295: import org.dive4elements.river.artifacts.common.ResultRow; mschaefer@9295: import org.dive4elements.river.artifacts.model.Calculation; gernotbelger@9499: import org.dive4elements.river.artifacts.model.river.MainWstValuesCalculator; gernotbelger@9499: import org.dive4elements.river.artifacts.model.river.RiverInfoProvider; mschaefer@9295: import org.dive4elements.river.artifacts.sinfo.tkhstate.WinfoArtifactWrapper; mschaefer@9295: import org.dive4elements.river.artifacts.uinfo.UINFOArtifact; gernotbelger@9328: import org.dive4elements.river.artifacts.uinfo.common.UInfoResultType; mschaefer@9309: import org.dive4elements.river.artifacts.uinfo.salix.SalixLineAccess.ScenarioType; mschaefer@9443: import org.dive4elements.river.artifacts.uinfo.vegetationzones.VegetationZoneServerClientXChange; mschaefer@9361: import org.dive4elements.river.utils.Formatter; mschaefer@9295: mschaefer@9295: /** mschaefer@9295: * Calculation of the result rows of the u-info salix line calc mode mschaefer@9295: * mschaefer@9295: * @author Matthias Schäfer mschaefer@9295: */ gernotbelger@9321: final class SalixLineCalculator { mschaefer@9295: gernotbelger@9499: private static final String MAIN_VALUE_MNQ = "mnq"; gernotbelger@9499: gernotbelger@9499: private static final String MAIN_VALUE_MQ = "mq"; gernotbelger@9499: gernotbelger@9499: private static final String MAIN_VALUE_MHQ = "mhq"; gernotbelger@9499: gernotbelger@9499: private static final String MAIN_VALUE_HQ5 = "hq5"; gernotbelger@9499: mschaefer@9375: private static final BigDecimal SALIX_DISTANCE = new BigDecimal("2.31"); mschaefer@9295: private final List rows = new ArrayList<>(); mschaefer@9295: mschaefer@9295: private final RiverInfoProvider riverInfoProvider; mschaefer@9295: gernotbelger@9321: public SalixLineCalculator(final RiverInfoProvider riverInfoProvider) { mschaefer@9295: this.riverInfoProvider = riverInfoProvider; mschaefer@9295: } mschaefer@9295: mschaefer@9295: /** mschaefer@9295: * Calculate the salix line result rows mschaefer@9295: */ mschaefer@9309: public void execute(final Calculation problems, final UINFOArtifact uinfo, final NavigableMap> rangeScenarios, mschaefer@9361: final ScenarioType scenarioType, final String[] scenarioLabels, final String rangeString, final String additionalString, mschaefer@9361: final SalixLineCalculationResults results) { mschaefer@9295: gernotbelger@9499: final MainWstValuesCalculator mainWstValues = fetchGaugeMainValuePositions2(problems); mschaefer@9295: mschaefer@9295: final WINFOArtifact winfo = new WinfoArtifactWrapper(uinfo); mschaefer@9295: winfo.addStringData("ld_mode", "distance"); mschaefer@9295: winfo.addStringData("ld_step", "100"); mschaefer@9295: for (final double station : new ComputationRangeAccess(winfo).getKms()) { gernotbelger@9499: this.rows.add(createRow(mainWstValues, station, rangeScenarios)); mschaefer@9295: } mschaefer@9309: if (scenarioType == ScenarioType.REGIONAL) mschaefer@9361: results.addResult(new SalixLineCalculationRegionalResult("Salix-regional", scenarioLabels, rangeString, additionalString, this.rows), problems); mschaefer@9309: else if (scenarioType == ScenarioType.SUPRAREGIONAL) mschaefer@9361: results.addResult(new SalixLineCalculationSupraRegionalResult("Salix-supra", scenarioLabels, rangeString, additionalString, this.rows), problems); mschaefer@9309: else if (scenarioType == ScenarioType.HISTORICAL) mschaefer@9361: results.addResult(new SalixLineCalculationHistoricalResult("Salix-hist", scenarioLabels, rangeString, additionalString, this.rows), problems); mschaefer@9309: else mschaefer@9361: results.addResult(new SalixLineCalculationResult("Salix-simple", this.rows), problems); mschaefer@9295: } mschaefer@9295: gernotbelger@9499: private MainWstValuesCalculator fetchGaugeMainValuePositions2(final Calculation problems) { gernotbelger@9499: final MainWstValuesCalculator mainWstValues = MainWstValuesCalculator.forRiverInfo(this.riverInfoProvider, MAIN_VALUE_MQ, MAIN_VALUE_MNQ, gernotbelger@9499: MAIN_VALUE_MHQ, MAIN_VALUE_HQ5); mschaefer@9430: gernotbelger@9499: if (!mainWstValues.hasPosition(MAIN_VALUE_MQ)) gernotbelger@9499: problems.addProblem("uinfo_salix_calc.warning.missing_mq"); gernotbelger@9499: else { gernotbelger@9499: if (!mainWstValues.hasPosition(MAIN_VALUE_MHQ)) gernotbelger@9499: problems.addProblem("uinfo_salix_calc.warning.missing_mhq"); gernotbelger@9499: if (!mainWstValues.hasPosition(MAIN_VALUE_MNQ)) gernotbelger@9499: problems.addProblem("uinfo_salix_calc.warning.missing_mnq"); mschaefer@9368: } gernotbelger@9499: gernotbelger@9499: return mainWstValues; mschaefer@9295: } mschaefer@9295: mschaefer@9295: /** mschaefer@9295: * Create a result row for a station and its gauge, and add w-q-values as selected gernotbelger@9499: * gernotbelger@9499: * @param mainWstValues mschaefer@9295: */ gernotbelger@9499: private ResultRow createRow(final MainWstValuesCalculator mainWstValues, final double station, final NavigableMap> rangeScenarios) { mschaefer@9295: mschaefer@9295: final ResultRow row = ResultRow.create(); mschaefer@9295: row.putValue(GeneralResultType.station, station); mschaefer@9394: // Find station's gauge (obsolete version which calculates gauge-wise) mschaefer@9368: // final Gauge gauge = this.riverInfoProvider.getGauge(station, true); mschaefer@9316: // Interpolate mnw, mw, and mhw mschaefer@9368: // final double mnw = interpolateW(station, this.gaugeMnwPos.get(gauge)); mschaefer@9368: // final double mw = interpolateW(station, this.gaugeMwPos.get(gauge)); mschaefer@9368: // final double mhw = interpolateW(station, this.gaugeMhwPos.get(gauge)); gernotbelger@9499: final double mnw = mainWstValues.interpolateW(station, MAIN_VALUE_MNQ); gernotbelger@9499: final double mw = mainWstValues.interpolateW(station, MAIN_VALUE_MQ); gernotbelger@9499: final double mhw = mainWstValues.interpolateW(station, MAIN_VALUE_MHQ); gernotbelger@9499: final double hw5 = mainWstValues.interpolateW(station, MAIN_VALUE_HQ5); gernotbelger@9429: row.putValue(UInfoResultType.waterlevelMNW, mnw); gernotbelger@9429: row.putValue(UInfoResultType.waterlevelMW, mw); gernotbelger@9429: row.putValue(UInfoResultType.waterlevelMHW, mhw); mschaefer@9430: row.putValue(UInfoResultType.waterlevelMH5, hw5); gernotbelger@9429: mschaefer@9316: // Calc salix-line and mw-mnw mschaefer@9490: row.putValue(UInfoResultType.salixline, calcSalix(mhw, mw, 0.0)); mschaefer@9361: row.putValue(UInfoResultType.salix_mw_mnw, calcMwmnw(mw, mnw)); mschaefer@9490: final double salixw = Formatter.roundW(mhw).subtract(SALIX_DISTANCE).doubleValue(); mschaefer@9490: row.putValue(UInfoResultType.salixw, salixw); mschaefer@9316: // Calc scenario values (always all scenario types set, Result variant extracts the fields needed) mschaefer@9316: final List scenarios = new ArrayList<>(); mschaefer@9394: final List deltaws = getDeltaWs(station, rangeScenarios); mschaefer@9394: for (final Double deltaw : deltaws) { mschaefer@9394: if (deltaw != null) { mschaefer@9490: final double salix = calcSalix(mhw, mw, deltaw.doubleValue()); mschaefer@9490: final double scen = calcSalix(mhw, 0.0, deltaw.doubleValue()); mschaefer@9490: scenarios.add(new SalixScenario((int) (deltaw * 100), salix, scen)); gernotbelger@9499: } else { mschaefer@9361: scenarios.add(null); mschaefer@9361: } mschaefer@9316: } mschaefer@9361: row.putValue(UInfoResultType.customMultiRowColSalixScenarios, scenarios); mschaefer@9397: row.putValue(GeneralResultType.gaugeLabel, this.riverInfoProvider.findGauge(station)); mschaefer@9295: return row; mschaefer@9295: } mschaefer@9295: mschaefer@9295: /** mschaefer@9295: * Calculates the salix value mschaefer@9295: */ mschaefer@9490: private double calcSalix(final double mhw, final double mw, final double deltamw) { mschaefer@9382: if (Double.isNaN(mw) || Double.isInfinite(mw) || Double.isNaN(mhw) || Double.isInfinite(mhw)) mschaefer@9490: return mhw - mw; // preserving NaN or Infinity mschaefer@9490: return Formatter.roundW(mhw).subtract(SALIX_DISTANCE).subtract(Formatter.roundW(mw).add(Formatter.roundW(deltamw))).doubleValue(); mschaefer@9295: } mschaefer@9295: mschaefer@9295: /** mschaefer@9295: * Calculates the inverse MW-MNW difference mschaefer@9295: */ mschaefer@9295: private double calcMwmnw(final double mw, final double mnw) { mschaefer@9382: if (Double.isNaN(mw) || Double.isInfinite(mw) || Double.isNaN(mnw) || Double.isInfinite(mnw)) mschaefer@9382: return mnw - mw; // preserving NaN or Inifinity mschaefer@9375: return Formatter.roundW(mnw).subtract(Formatter.roundW(mw)).doubleValue(); mschaefer@9295: } mschaefer@9316: mschaefer@9316: /** mschaefer@9394: * Gets the station-specific list of delta-ws of the active scenario, at least with one null item in any case mschaefer@9316: */ mschaefer@9394: private List getDeltaWs(final double station, final NavigableMap> rangeScenarios) { gernotbelger@9321: final Entry> stationScenarios = rangeScenarios.floorEntry(station); mschaefer@9394: if (stationScenarios != null) { mschaefer@9394: return stationScenarios.getValue(); mschaefer@9394: } mschaefer@9394: final List noScen = new ArrayList<>(); mschaefer@9394: noScen.add(null); mschaefer@9394: return noScen; mschaefer@9316: } mschaefer@9397: mschaefer@9397: /** mschaefer@9397: * Find and return a height (iota, w main value) of a station in a previously calculated result mschaefer@9397: */ mschaefer@9397: public double fetchStationHeight(final Calculation problems, final double station, final AbstractResultType resultType, mschaefer@9397: final SalixLineCalculationResult result) { mschaefer@9397: mschaefer@9397: // Search the station in the previously calculated result rows mschaefer@9397: final ResultRow stationRow = searchStation(station, result.getRows()); mschaefer@9397: if (stationRow != null) mschaefer@9397: return stationRow.getDoubleValue(resultType); mschaefer@9397: return Double.NaN; mschaefer@9397: } mschaefer@9397: mschaefer@9397: /** mschaefer@9397: * Searches the row of a station in a result rows collection mschaefer@9397: */ mschaefer@9397: private ResultRow searchStation(final double station, final Collection rows) { mschaefer@9397: for (final ResultRow row : rows) { mschaefer@9397: if (row.getDoubleValue(GeneralResultType.station) > station + 0.0001) mschaefer@9397: return row; mschaefer@9397: } mschaefer@9397: return null; mschaefer@9397: } mschaefer@9443: mschaefer@9443: /** mschaefer@9443: * Computes the height of a vegetation zone limit for a station and a previously computed salix line result mschaefer@9443: */ mschaefer@9443: public double computeVegetationZoneHeight(final Calculation problems, final double station, final int vegetationZoneType, mschaefer@9443: final SalixLineCalculationResult result) { mschaefer@9443: mschaefer@9443: // Search the station in the previously calculated result rows mschaefer@9443: final ResultRow stationRow = searchStation(station, result.getRows()); mschaefer@9443: if (stationRow == null) mschaefer@9443: return Double.NaN; gernotbelger@9499: mschaefer@9443: // Compute height from overflow duration days mschaefer@9443: final List vzs = VegetationZoneServerClientXChange.getStandardList(null, null); // TODO river, context mschaefer@9443: if ((vegetationZoneType >= 1) && (vegetationZoneType <= vzs.size())) { mschaefer@9443: final int uefd = vzs.get(vegetationZoneType - 1).getMin_day_overflow(); mschaefer@9443: // Üfd = -70,559 ∗ ln((DGM - MW) + 0,5) + 80,711 mschaefer@9443: final double f1 = -70.559; mschaefer@9443: final double f2 = -88.711; mschaefer@9443: final double mw = stationRow.getDoubleValue(UInfoResultType.waterlevelMW); mschaefer@9443: final double dgm = Math.exp((uefd - f2) / f1) + mw - 0.5; mschaefer@9443: return dgm; mschaefer@9443: } mschaefer@9443: return Double.NaN; mschaefer@9443: } mschaefer@9443: }