changeset 451:73bc64c4a7b0

Use new logic to calculate "W für ungleichwertige Abfluesse". Not working yet. flys-artifacts/trunk@1946 c6561f87-3c4e-4783-a992-168aeb5c3f6f
author Sascha L. Teichmann <sascha.teichmann@intevation.de>
date Wed, 18 May 2011 17:37:06 +0000
parents c8bb38115290
children 343f248e4c8c
files flys-artifacts/ChangeLog flys-artifacts/src/main/java/de/intevation/flys/artifacts/FLYSArtifact.java flys-artifacts/src/main/java/de/intevation/flys/artifacts/WINFOArtifact.java
diffstat 3 files changed, 174 insertions(+), 47 deletions(-) [+]
line wrap: on
line diff
--- a/flys-artifacts/ChangeLog	Wed May 18 10:59:38 2011 +0000
+++ b/flys-artifacts/ChangeLog	Wed May 18 17:37:06 2011 +0000
@@ -1,3 +1,14 @@
+2011-05-18	Sascha L. Teichmann	<sascha.teichmann@intevation.de>
+
+	* src/main/java/de/intevation/flys/artifacts/FLYSArtifact.java:
+	  Made getExplodedValues static.
+
+	* src/main/java/de/intevation/flys/artifacts/WINFOArtifact.java:
+	  Use new logic to calculate "W für ungleichwertige Abfluesse".
+	  Not working, yet.
+
+	* ChangeLog: Fixed former entry.
+
 2011-05-18  Ingo Weinzierl <ingo@intevation.de>
 
 	* doc/conf/artifacts/winfo.xml: Registered the WST export for discharge
@@ -16,15 +27,15 @@
 	* src/main/java/de/intevation/flys/exports/DischargeLongitudinalSectionExporter.java:
 	  Adapted the column names for the WST export.
 
-2011-05-19	Sascha L. Teichmann	<sascha.teichmann@intevation.de>
-
-    Work on flys/issue69
-
-    * src/main/java/de/intevation/flys/artifacts/WINFOArtifact.java:
-      Use new logic to calculate "Wasserstand/Wasserspiegellage".
-      Compared to desktop FLYS are the results are structurally right 
-      but a bit off in the positions after the decimal points.
-      Maybe a result of the interpolation? Need to debug this.
+2011-05-18	Sascha L. Teichmann	<sascha.teichmann@intevation.de>
+
+	Work on flys/issue69
+
+	* src/main/java/de/intevation/flys/artifacts/WINFOArtifact.java:
+	  Use new logic to calculate "Wasserstand/Wasserspiegellage".
+	  Compared to desktop FLYS are the results are structurally right 
+	  but a bit off in the positions after the decimal points.
+	  Maybe a result of the interpolation? Need to debug this.
 
 2011-05-18  Ingo Weinzierl <ingo@intevation.de>
 
--- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/FLYSArtifact.java	Wed May 18 10:59:38 2011 +0000
+++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/FLYSArtifact.java	Wed May 18 17:37:06 2011 +0000
@@ -906,7 +906,7 @@
      *
      * @return an array of double values.
      */
-    public double[] getExplodedValues(double from, double to, double step) {
+    public static double[] getExplodedValues(double from, double to, double step) {
         double lower = from;
 
         double diff = to - from;
--- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/WINFOArtifact.java	Wed May 18 10:59:38 2011 +0000
+++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/WINFOArtifact.java	Wed May 18 17:37:06 2011 +0000
@@ -32,8 +32,12 @@
 import de.intevation.flys.artifacts.model.WQDay;
 import de.intevation.flys.artifacts.model.WQKms;
 import de.intevation.flys.artifacts.model.WstValueTable;
+import de.intevation.flys.artifacts.model.WstValueTable.QPosition;
 import de.intevation.flys.artifacts.model.WstValueTableFactory;
 
+import de.intevation.flys.artifacts.math.LinearRemap;
+
+import gnu.trove.TDoubleArrayList;
 
 /**
  * The default WINFO artifact.
@@ -472,73 +476,185 @@
         return wqkms;
     }
 
-
     /**
      * Returns the data computed by the discharge longitudinal section
      * computation.
      *
      * @return an array of WQKms object - one object for each given Q value.
      */
-    public WQKms[] getDischargeLongitudinalSectionData() {
+    public WQKms [] getDischargeLongitudinalSectionData() {
+
         logger.debug("WINFOArtifact.getDischargeLongitudinalSectionData");
 
         River river = getRiver();
         if (river == null) {
             logger.error("No river selected.");
+            return new WQKms[0];
         }
 
         WstValueTable wst = WstValueTableFactory.getTable(river);
         if (wst == null) {
-            logger.error("No Wst found for selected river.");
+            logger.error("No wst found for selected river.");
+            return new WQKms[0];
         }
 
-        double[][] dist = getSplittedDistance();
-        int        num  = dist != null ? dist.length : 0;
+        double [][] segments = getSplittedDistance();
 
-        double[][] allQs = new double[num][];
-        double[][] allWs = new double[num][];
-        WQKms[][]  wqkms = new WQKms[num][];
+        if (segments.length < 1) {
+            logger.warn("no segments given");
+            return new WQKms[0];
+        }
 
-        boolean qSel = true;
+        if (segments.length == 1) {
+            // fall back to normal "Wasserstand/Wasserspiegellage" calculation
+            double [] qs = toQs(segments[0]);
+            if (qs == null) {
+                logger.warn("no qs given");
+                return new WQKms[0];
+            }
+            if (qs.length == 1) {
+                double [] kms = getKms(segments[0]);
+                return computeWaterlevelData(kms, qs, wst);
+            }
+        }
 
-        for (int i = 0; i < num; i++) {
-            double[] kms = getKms(dist[i]);
-            if (kms == null) {
-                // XXX maybe we should cancel this operation here.
+        // more than one segment
+
+        double [] boundKms;
+
+        if (segments.length == 2) {
+            boundKms = new double [] { segments[0][0], segments[1][1] };
+        }
+        else {
+            TDoubleArrayList bounds = new TDoubleArrayList();
+
+            bounds.add(segments[0][0]);
+
+            for (int i = 1; i < segments.length-1; ++i) {
+                double [] segment = segments[i];
+
+                Gauge gauge = river.determineGauge(segment[0], segment[1]);
+
+                if (gauge == null) {
+                    logger.warn("no gauge found between " + 
+                        segment[0] + " and " + segment[1]);
+                    bounds.add(0.5*(segment[0] + segment[1]));
+                }
+                else {
+                    bounds.add(gauge.getStation().doubleValue());
+                }
+            }
+
+            bounds.add(segments[segments.length-1][1]);
+            boundKms = bounds.toNativeArray();
+        }
+
+        if (logger.isDebugEnabled()) {
+            logger.debug("bound kms: " + joinDoubles(boundKms));
+        }
+
+        double [][] iqs = null;
+
+        for (int i = 0; i < segments.length; ++i) {
+            double [] iqsi = toQs(segments[i]);
+            if (iqsi == null) {
+                logger.warn("iqsi == null");
+                return new WQKms[0];
+            }
+
+            if (iqs == null) {
+                iqs = new double[iqsi.length][boundKms.length];
+            }
+            else if (iqs.length != iqsi.length) {
+                logger.warn("iqsi.logger != iqs.length: "
+                    + iqsi.length + " " + iqsi.length);
+                return new WQKms[0];
+            }
+
+            if (logger.isDebugEnabled()) {
+                logger.debug("segments qs[ " + i + "]: " + joinDoubles(iqsi));
+            }
+
+            for (int j = 0; j < iqs.length; ++j) {
+                iqs[j][i] = iqsi[j];
+            }
+        }
+
+        if (logger.isDebugEnabled()) {
+            for (int i = 0; i < iqs.length; ++i) {
+                logger.debug("iqs[" + i + "]: " + joinDoubles(iqs[i]));
+            }
+        }
+
+        double [] boundWs = new double[boundKms.length];
+        double [] boundQs = new double[boundKms.length];
+
+        // XXX: Is there some state missing?
+
+        double [] okms = getExplodedValues(
+            boundKms[0], boundKms[boundKms.length-1], DEFAULT_KM_STEPS);
+
+        ArrayList<WQKms> results = new ArrayList<WQKms>();
+
+        for (int i = 0; i < iqs.length; ++i) {
+            double [] iqsi = iqs[i];
+
+            QPosition qPosition = wst.interpolate(
+                iqsi[0], 0, boundKms, boundWs, boundQs);
+
+            if (qPosition == null) {
+                logger.warn("interpolation failed for " + iqsi[i]);
                 continue;
             }
 
-            double[] qs = getQs(dist[i]);
-            allQs[i]    = qs;
+            LinearRemap remap = new LinearRemap();
 
-            if (qs == null) {
-                logger.debug("Determine Q values based on a set of W values.");
-                qSel = false;
-
-                allWs[i] = getWs(dist[i]);
-                qs       = getQsForWs(allWs[i]);
+            for (int j = 1; j < boundKms.length; ++j) {
+                remap.add(
+                    boundKms[j-1], boundKms[j],
+                    boundQs[j-1],  iqsi[j-1],
+                    boundQs[j],    iqsi[j]);
             }
 
-            wqkms[i] = computeWaterlevelData(kms, qs, wst);
-        }
-
-        WQKms[] merged    = WQKms.merge(wqkms);
-        int     numMerged = merged.length;
-        WQKms[] computed  = new WQKms[numMerged];
-
-        for (int i = 0; i < numMerged; i++) {
-            computed[i] = computeDischargeLongitudinalSectionData(merged[i]);
+            double [] oqs = new double[okms.length];
+            double [] ows = new double[okms.length];
 
-            setDischargeLongitudinalSectionNames(
-                computed[i],
-                qSel ? allQs : allWs,
-                i,
-                qSel ? "Q" : "W");
+            wst.interpolate(okms, ows, oqs, qPosition, remap);
+
+            BackJumpCorrector bjc = new BackJumpCorrector();
+            if (bjc.doCorrection(okms, ows)) {
+                logger.debug("Discharge longitudinal section has backjumps.");
+                results.add(new WQCKms(okms, oqs, ows, bjc.getCorrected()));
+            }
+            else {
+                results.add(new WQKms(okms, oqs, ows));
+            }
         }
+        // TODO: set names
 
-        // TODO Introduce a caching mechanism here!
+        return results.toArray(new WQKms[results.size()]);
+    }
 
-        return computed;
+    protected static String joinDoubles(double [] x) {
+        if (x == null) {
+            return "null";
+        }
+        StringBuilder sb = new StringBuilder();
+        for (int i = 0; i < x.length; ++i) {
+            if (i > 0) sb.append(", ");
+            sb.append(x[i]);
+        }
+        return sb.toString();
+    }
+
+    protected double [] toQs(double [] range) {
+        double [] qs = getQs(range);
+        if (qs == null) {
+            logger.debug("Determine Q values based on a set of W values.");
+            double [] ws = getWs(range);
+            qs = getQsForWs(ws);
+        }
+        return qs;
     }
 
 

http://dive4elements.wald.intevation.org