diff artifacts/src/main/java/org/dive4elements/river/artifacts/sinfo/tkhcalculation/FlowVelocityModelKmValueFinder.java @ 8915:d9dbf0b74bc2

Refaktoring of flow depth calculation, extracting tkh part. First implementation of tkh calculation.
author gernotbelger
date Wed, 28 Feb 2018 17:27:15 +0100
parents artifacts/src/main/java/org/dive4elements/river/artifacts/sinfo/flowdepth/FlowVelocityModelKmValueFinder.java@89f3c5462a16
children 45f1ad66560e
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/artifacts/src/main/java/org/dive4elements/river/artifacts/sinfo/tkhcalculation/FlowVelocityModelKmValueFinder.java	Wed Feb 28 17:27:15 2018 +0100
@@ -0,0 +1,337 @@
+/* Copyright (C) 2017 by Bundesanstalt für Gewässerkunde
+ * Software engineering by
+ *  Björnsen Beratende Ingenieure GmbH
+ *  Dr. Schumacher Ingenieurbüro für Wasser und Umwelt
+ *
+ * This file is Free Software under the GNU AGPL (>=v3)
+ * and comes with ABSOLUTELY NO WARRANTY! Check out the
+ * documentation coming with Dive4Elements River for details.
+ */
+
+package org.dive4elements.river.artifacts.sinfo.tkhcalculation;
+
+import java.util.ArrayList;
+import java.util.List;
+
+import org.apache.commons.lang.math.DoubleRange;
+import org.apache.log4j.Logger;
+import org.dive4elements.river.artifacts.math.Linear;
+import org.dive4elements.river.artifacts.math.Utils;
+import org.dive4elements.river.artifacts.sinfo.flowdepth.FlowVelocityKmModelValues;
+import org.dive4elements.river.backend.SessionHolder;
+import org.dive4elements.river.model.River;
+import org.hibernate.SQLQuery;
+import org.hibernate.Session;
+import org.hibernate.type.StandardBasicTypes;
+
+import gnu.trove.TDoubleArrayList;
+
+/**
+ * Searchable sorted km array with parallel FlowVelocityKmModelValues array and linear interpolation for km and the
+ * model values between the array elements.<br />
+ * {@link loadValues} loads all the model values for a given km range of a river.<br />
+ * {@link findKmQValues} then searches a km in the values table or the nearest including km interval, resp.
+ * The v and tau values for a given discharge are either found directly or also interpolated linearly.<br />
+ *
+ * (Created based on a copy of FlowVelocityMeasurementFactory.)
+ *
+ * @author Matthias Schäfer
+ */
+final class FlowVelocityModelKmValueFinder {
+    /***** FIELDS *****/
+
+    /**
+     * Private log to use here.
+     */
+    private static Logger log = Logger.getLogger(FlowVelocityModelKmValueFinder.class);
+
+    /**
+     * Query for a range of stations of a river with all their q, main-v and tau values.<br />
+     * (Might be several 10000 rows if many stations and large q range)
+     */
+    private static final String SQL_SELECT_ALL = //
+            "SELECT fvmv.station AS station, fvmv.q AS q, fvmv.main_channel AS vmain, fvmv.shear_stress AS tau"
+            + "  FROM (discharge_zone dz INNER JOIN flow_velocity_model fvm ON dz.id = fvm.discharge_zone_id)"
+            + "    INNER JOIN flow_velocity_model_values fvmv ON fvm.id = fvmv.flow_velocity_model_id"
+            + "  WHERE (dz.river_id = :river_id) AND (fvmv.station BETWEEN :kmfrom - 0.0001 AND :kmto + 0.0001)"
+            /* + "  WHERE (dz.river_id = :river_id) AND (fvmv.q BETWEEN :qmin AND :qmax)" */
+            + "  ORDER BY fvmv.station ASC, fvmv.q ASC";
+
+    /**
+     * Query for a river's max km below a limit with all its q, main-v and tau values.
+     */
+    private static final String SQL_SELECT_KMLOWER = //
+            "SELECT fvmv.station AS station, fvmv.q AS q, fvmv.main_channel AS vmain, fvmv.shear_stress AS tau" + "  FROM flow_velocity_model_values fvmv"
+            + "    INNER JOIN (SELECT MAX(fvmvi.station) AS kmmax"
+            + "      FROM(discharge_zone dz INNER JOIN flow_velocity_model fvm ON dz.id = fvm.discharge_zone_id)"
+            + "        INNER JOIN flow_velocity_model_values fvmvi ON fvm.id = fvmvi.flow_velocity_model_id"
+            + "        WHERE (dz.river_id = :river_id) AND (fvmvi.station < :kmfrom - 0.0001)) finf ON fvmv.station = finf.kmmax"
+            + "  ORDER BY fvmv.q ASC";
+
+    /**
+     * Query for a river's min km above a limit with all its q, main-v and tau values.
+     */
+    private static final String SQL_SELECT_KMUPPER = //
+            "SELECT fvmv.station AS station, fvmv.q AS q, fvmv.main_channel AS vmain, fvmv.shear_stress AS tau" + "  FROM flow_velocity_model_values fvmv"
+            + "    INNER JOIN (SELECT MIN(fvmvi.station) AS kmmin"
+            + "      FROM(discharge_zone dz INNER JOIN flow_velocity_model fvm ON dz.id = fvm.discharge_zone_id)"
+            + "        INNER JOIN flow_velocity_model_values fvmvi ON fvm.id = fvmvi.flow_velocity_model_id"
+            + "        WHERE (dz.river_id = :river_id) AND (fvmvi.station > :kmto + 0.0001)) fsup ON fvmv.station = fsup.kmmin"
+            + "  ORDER BY fvmv.q ASC";
+
+    // /**
+    // * Query to select all km-q-v-tau of a river that are the q maxima below a q limit
+    // */
+    // private static final String SQL_SELECT_QLOWER =
+    // "SELECT fvmv.station AS station, fvmv.q AS q, fvmv.main_channel AS vmain, fvmv.shear_stress AS tau"
+    // + " FROM flow_velocity_model_values fvmv"
+    // + " INNER JOIN (SELECT fv2.station, MAX(fv2.q) AS q"
+    // + " FROM (discharge_zone dz INNER JOIN flow_velocity_model fvm ON dz.id = fvm.discharge_zone_id)"
+    // + " INNER JOIN flow_velocity_model_values fv2 ON fvm.id = fv2.flow_velocity_model_id"
+    // + " WHERE (dz.river_id = :river_id) AND (fv2.q < :qlim) GROUP BY fv2.station) qx"
+    // + " ON (fvmv.station=qx.station) AND (fvmv.q=qx.q)"
+    // + " ORDER BY fvmv.station ASC";
+    //
+    // /**
+    // * Query to select all km-q-v-tau of a river that are the q minima above a q limit
+    // */
+    // private static final String SQL_SELECT_QUPPER =
+    // "SELECT fvmv.station AS station, fvmv.q AS q, fvmv.main_channel AS vmain, fvmv.shear_stress AS tau"
+    // + " FROM flow_velocity_model_values fvmv"
+    // + " INNER JOIN (SELECT fv2.station, MIN(fv2.q) AS q"
+    // + " FROM (discharge_zone dz INNER JOIN flow_velocity_model fvm ON dz.id = fvm.discharge_zone_id)"
+    // + " INNER JOIN flow_velocity_model_values fv2 ON fvm.id = fv2.flow_velocity_model_id"
+    // + " WHERE (dz.river_id = :river_id) AND (fv2.q > :qlim) GROUP BY fv2.station) qx"
+    // + " ON (fvmv.station=qx.station) AND (fvmv.q=qx.q)"
+    // + " ORDER BY fvmv.station ASC";
+
+    /**
+     * Kms of the loaded river range
+     */
+    private final TDoubleArrayList kms = new TDoubleArrayList();
+
+    /**
+     * For each km in kms a list of q-v-tau-tupels
+     */
+    private final List<FlowVelocityKmModelValues> values = new ArrayList<>();
+
+    /**
+     * Searched km of the last findKmValue
+     */
+    private double findKm;
+
+    /**
+     * kms and values index of the interval start found by the last findKmValue
+     */
+    private int leftIndexFound = -1;
+
+    /**
+     * kms and values index of the interval end found by the last findKmValue
+     */
+    private int rightIndexFound = -1;
+
+    /**
+     * Q of the last findKmQValues
+     */
+    private double findQ;
+
+    /***** METHODS *****/
+
+    /**
+     * Discharge of the last {@link findKmQValue}
+     */
+    public double getFindQ() {
+        return this.findQ;
+    }
+
+    /**
+     * Velocity of the last {@link findKmQValues}
+     */
+    public double getFindVmainFound() {
+        if (this.leftIndexFound < 0)
+            return Double.NaN;
+        else if (this.leftIndexFound == this.rightIndexFound)
+            return getLeftValues().getVmainFound();
+        else
+            return Linear.linear(this.findKm, getLeftValues().getKm(), getRightValues().getKm(), getLeftValues().getVmainFound(),
+                    getRightValues().getVmainFound());
+    }
+
+    /**
+     * Shear stress tau of the last {@link findKmQValues}
+     */
+    public double getFindTauFound() {
+        if (this.leftIndexFound < 0)
+            return Double.NaN;
+        else if (this.leftIndexFound == this.rightIndexFound)
+            return getLeftValues().getTauFound();
+        else
+            return Linear.linear(this.findKm, getLeftValues().getKm(), getRightValues().getKm(), getLeftValues().getTauFound(), getRightValues().getTauFound());
+    }
+
+    /**
+     * Whether the discharge has been interpolated in the last {@link findKmQValues}
+     */
+    public boolean getFindIsQInterpolated() {
+        return (getLeftValues() != null) && (getLeftValues().getIsInterpolated() || getRightValues().getIsInterpolated());
+    }
+
+    /**
+     * Static constructor: queries a range of a river's kms with all their q-v-tau values.
+     *
+     * @return Whether the load has been successful the new instance, <code>null</code> otherwise.
+     */
+    public static FlowVelocityModelKmValueFinder loadValues(final River river, final DoubleRange kmRange, final DoubleRange qRange) {
+        // DB session
+        log.debug(String.format("loadValues km %.3f - %.3f / q %.1f - %.1f", kmRange.getMinimumDouble(), kmRange.getMaximumDouble(), qRange.getMinimumDouble(),
+                qRange.getMaximumDouble()));
+
+        final FlowVelocityModelKmValueFinder instance = new FlowVelocityModelKmValueFinder();
+
+        final TDoubleArrayList kms = instance.kms;
+        final List<FlowVelocityKmModelValues> values = instance.values;
+
+        final boolean isDemoValuesCorrection = river.getName().equalsIgnoreCase("beispielfluss");
+        final Session session = SessionHolder.HOLDER.get();
+
+        // Select km infimum
+        SQLQuery sqlQuery = session.createSQLQuery(SQL_SELECT_KMLOWER).addScalar("station", StandardBasicTypes.DOUBLE).addScalar("q", StandardBasicTypes.DOUBLE)
+                .addScalar("vmain", StandardBasicTypes.DOUBLE).addScalar("tau", StandardBasicTypes.DOUBLE);
+        sqlQuery.setParameter("river_id", river.getId());
+        sqlQuery.setParameter("kmfrom", kmRange.getMinimumDouble());
+        instance.addKms(sqlQuery.list(), isDemoValuesCorrection);
+
+        // Select km range
+        sqlQuery = session.createSQLQuery(SQL_SELECT_ALL).addScalar("station", StandardBasicTypes.DOUBLE).addScalar("q", StandardBasicTypes.DOUBLE)
+                .addScalar("vmain", StandardBasicTypes.DOUBLE).addScalar("tau", StandardBasicTypes.DOUBLE);
+        sqlQuery.setParameter("river_id", river.getId());
+        sqlQuery.setParameter("kmfrom", kmRange.getMinimumDouble());
+        sqlQuery.setParameter("kmto", kmRange.getMaximumDouble());
+        // sqlQuery.setParameter("qmin", qRange.getMinimumDouble());
+        // sqlQuery.setParameter("qmax", qRange.getMaximumDouble());
+
+        int kmcount = kms.size();
+        final int rowcount = instance.addKms(sqlQuery.list(), isDemoValuesCorrection);
+        kmcount = kms.size() - kmcount;
+
+        // Select km supremum
+        sqlQuery = session.createSQLQuery(SQL_SELECT_KMUPPER).addScalar("station", StandardBasicTypes.DOUBLE).addScalar("q", StandardBasicTypes.DOUBLE)
+                .addScalar("vmain", StandardBasicTypes.DOUBLE).addScalar("tau", StandardBasicTypes.DOUBLE);
+        sqlQuery.setParameter("river_id", river.getId());
+        sqlQuery.setParameter("kmto", kmRange.getMaximumDouble());
+        final int supcnt = instance.addKms(sqlQuery.list(), isDemoValuesCorrection);
+
+        // Add copy of last km for search of max km value
+        if ((supcnt == 0) && (values.size() >= 1)) {
+            kms.add(kms.getQuick(kms.size()) + 0.0001);
+            values.add(new FlowVelocityKmModelValues(kms.getQuick(kms.size() - 1), values.get(values.size() - 1)));
+        }
+
+        // log.debug
+        if (values.size() - 1 >= 0) {
+            log.debug(String.format("loadValues %d: km %.3f - %d values", 0, values.get(0).getKm(), values.get(0).size()));
+
+            if (values.size() - 1 >= 1) {
+                log.debug(String.format("loadValues %d: km %.3f - %d values", 1, values.get(1).getKm(), values.get(1).size()));
+
+                if (values.size() - 1 >= 2)
+                    log.debug("loadValues ...");
+
+                if (values.size() - 2 >= 3)
+                    log.debug(String.format("loadValues %d: km %.3f - %d values", values.size() - 2, values.get(values.size() - 2).getKm(),
+                            values.get(values.size() - 2).size()));
+
+                if (values.size() - 1 >= 3)
+                    log.debug(String.format("loadValues %d: km %.3f - %d values", values.size() - 1, values.get(values.size() - 1).getKm(),
+                            values.get(values.size() - 1).size()));
+            }
+        }
+
+        log.debug(String.format("loadValues %d kms, %d values loaded", kmcount, rowcount));
+
+        if (kms.size() == 0)
+            return null;
+
+        return instance;
+    }
+
+    /**
+     * Adds the km-q-v-tau values of a query result row to the last km of the list, or a new one resp.
+     *
+     * @return Number of rows
+     */
+    private int addKms(final List<Object[]> rows, final boolean isDemoValuesCorrection) {
+        for (final Object[] row : rows) {
+            if ((this.kms.size() == 0) || !Utils.epsilonEquals(this.kms.get(this.kms.size() - 1), (double) row[0], 0.0001)) {
+                this.kms.add((double) row[0]);
+                this.values.add(new FlowVelocityKmModelValues(this.kms.get(this.kms.size() - 1)));
+            }
+            if (isDemoValuesCorrection)
+                // "Verfremdung" der v-Werte etwas korrigieren (Originalwerte wurden mit Zufallswert zwischen 10 und 20 multipliziert)
+                this.values.get(this.values.size() - 1).addValues((double) row[1], ((double) row[2]) / 10, (double) row[3]);
+            else
+                this.values.get(this.values.size() - 1).addValues((double) row[1], (double) row[2], (double) row[3]);
+        }
+        return rows.size();
+    }
+
+    /**
+     * Searches a km and finds or interpolates the velocity and shear stress values for a discharge<br />
+     * The values may be got via {@link getVmainFound} etc.
+     *
+     * @return Whether values have been found
+     */
+    public boolean findKmQValues(final double km, final double q) {
+        this.findQ = q;
+        if (!searchKm(km))
+            return false;
+        if (this.leftIndexFound == this.rightIndexFound) {
+            // Exact km match
+            final double qfound = getLeftValues().findQ(q);
+            log.debug(String.format("findKmQValues km %.3f q %.0f = %.0f (%d)", km, q, qfound, this.leftIndexFound));
+            return !Double.isNaN(qfound);
+        } else {
+            final double[] qfound = { getLeftValues().findQ(q), getRightValues().findQ(q) };
+            log.debug(String.format("findKmQValues km %.3f q %.0f = %.0f (%d, %.3f) - %.0f (%d, %.3f)", km, q, qfound[0], this.leftIndexFound,
+                    getLeftValues().getKm(), qfound[1], this.rightIndexFound, getRightValues().getKm()));
+            return !Double.isNaN(qfound[0]) && !Double.isNaN(qfound[1]);
+        }
+    }
+
+    /**
+     * Searches a km
+     *
+     * @return Whether the km was within the supported range
+     */
+    private boolean searchKm(final double km) {
+        this.findKm = km;
+        this.leftIndexFound = -1;
+        this.rightIndexFound = -1;
+        if ((this.kms == null) || (this.kms.size() == 0))
+            return false;
+        int i = this.kms.binarySearch(km);
+        if (i >= 0) {
+            // Exact km match
+            this.leftIndexFound = i;
+            this.rightIndexFound = i;
+            return true;
+        } else {
+            // Out of range or within km interval
+            if (i < 0)
+                i = -i - 1;
+            if ((i <= 0) || (i >= this.kms.size()))
+                return false;
+            this.leftIndexFound = i - 1;
+            this.rightIndexFound = i;
+            return true;
+        }
+    }
+
+    private FlowVelocityKmModelValues getLeftValues() {
+        return this.values.get(this.leftIndexFound);
+    }
+
+    private FlowVelocityKmModelValues getRightValues() {
+        return this.values.get(this.rightIndexFound);
+    }
+}
\ No newline at end of file

http://dive4elements.wald.intevation.org