changeset 427:909196be11a0

First step to calculate "W fuer ungleichwertige Abfluesse" correctly. flys/issue55 flys-artifacts/trunk@1925 c6561f87-3c4e-4783-a992-168aeb5c3f6f
author Sascha L. Teichmann <sascha.teichmann@intevation.de>
date Sun, 15 May 2011 21:08:41 +0000
parents c6b7ca0febcb
children 5ee2e2c9b122
files flys-artifacts/ChangeLog flys-artifacts/src/main/java/de/intevation/flys/artifacts/math/LinearRemap.java flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/WstValueTable.java
diffstat 3 files changed, 184 insertions(+), 4 deletions(-) [+]
line wrap: on
line diff
--- a/flys-artifacts/ChangeLog	Sun May 15 10:59:08 2011 +0000
+++ b/flys-artifacts/ChangeLog	Sun May 15 21:08:41 2011 +0000
@@ -1,3 +1,16 @@
+2011-05-10	Sascha L. Teichmann	<sascha.teichmann@intevation.de>
+
+	First step to calculate "W fuer ungleichwertige Abfluesse" correctly.
+	flys/issue55
+
+	* src/main/java/de/intevation/flys/artifacts/math/LinearRemap.java:
+	  New. Remaps "gleichwertige" Q values to the corresponding
+	  "ungleichwertige" Q values depending on km.
+	  
+	* src/main/java/de/intevation/flys/artifacts/model/WstValueTable.java:
+	  Remap the Q values "ungleichwertig" depending on the 
+	  "gleichwertige" ones.
+
 2011-05-10	Sascha L. Teichmann	<sascha.teichmann@intevation.de>
 
 	First step to fix flys/issue69
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/math/LinearRemap.java	Sun May 15 21:08:41 2011 +0000
@@ -0,0 +1,78 @@
+package de.intevation.flys.artifacts.math;
+
+public class LinearRemap
+{
+    public static class Segment {
+
+        protected Segment next;
+
+        protected double from;
+        protected double to;
+
+        protected double m;
+        protected double b;
+
+        public Segment() {
+        }
+
+        public Segment(
+            double from, double to, 
+            double m,    double b,
+            Segment next
+        ) {
+            this.from = from;
+            this.to   = to;
+            this.m    = m;
+            this.b    = b;
+        }
+
+        public double eval(double x) {
+            return m*x + b;
+        }
+    } // class Segment
+
+    protected Segment head;
+
+    public LinearRemap() {
+    }
+
+    public void add(
+        double from, double to,
+        double x1,   double y1,
+        double x2,   double y2
+    ) {
+        // y1 = m*x1 + b <=> b = y1 - m*x1
+        // y2 = m*x2 + b
+        // y2 - y1 = m*(x2 - x1)
+        // m = (y2 - y1)/(x2 - x1)
+
+        double m, b;
+
+        if (x2 == x1) {
+            m = 0.0;
+            b = 0.5*(y2 + y1);
+        }
+        else {
+            m = (y2 - y1)/(x2 - x1);
+            b = y1 - m*x1;
+        }
+
+        head = new Segment(from, to, m, b, head);
+    }
+
+    public double eval(double pos, double x) {
+        Segment current = head;
+
+        while (current != null) {
+
+            if (pos >= current.from && pos <= current.to) {
+                return current.eval(x);
+            }
+
+            current = current.next;
+        }
+
+        return Double.NaN;
+    }
+}
+// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :
--- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/WstValueTable.java	Sun May 15 10:59:08 2011 +0000
+++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/WstValueTable.java	Sun May 15 21:08:41 2011 +0000
@@ -6,6 +6,8 @@
 import de.intevation.flys.model.Wst;
 import de.intevation.flys.model.WstColumn;
 
+import de.intevation.flys.artifacts.math.LinearRemap;
+
 import de.intevation.flys.artifacts.cache.CacheFactory;
 
 import de.intevation.flys.backend.SessionHolder;
@@ -441,6 +443,60 @@
             oqs[index] = weight(kmWeight, oq, tq);
             ows[index] = weight(kmWeight, ow, tw);
         }
+
+        public void storeWQ(
+            QPosition   qPosition,
+            LinearRemap remap,
+            double []   ows,
+            double []   oqs,
+            int         index
+        ) {
+            int    qIdx    = qPosition.index;
+            double qWeight = qPosition.weight;
+
+            oqs[index] = remap.eval(km, qWeight == 1.0
+                ? qs[qIdx]
+                : weight(qWeight, qs[qIdx-1], qs[qIdx])); 
+
+            ows[index] = interpolateW(oqs[index]);
+        }
+
+        public void storeWQ(
+            Row         other,
+            double      km,
+            QPosition   qPosition,
+            LinearRemap remap,
+            double []   ows,
+            double []   oqs,
+            int         index
+        ) {
+            double kmWeight = factor(km, this.km, other.km);
+
+            double tq, tw;
+            double oq, ow;
+
+            int    qIdx    = qPosition.index;
+            double qWeight = qPosition.weight;
+
+            if (qWeight == 1.0) {
+                tq = remap.eval(this.km,  qs[qIdx]);
+                oq = remap.eval(other.km, other.qs[qIdx]);
+            }
+            else {
+                tq = remap.eval(
+                    this.km,
+                    weight(qWeight, qs[qIdx-1], qs[qIdx]));
+                oq = remap.eval(
+                    other.km,
+                    weight(qWeight, other.qs[qIdx-1], other.qs[qIdx])); 
+            }
+
+            tw = interpolateW(tq);
+            ow = other.interpolateW(oq);
+
+            oqs[index] = weight(kmWeight, oq, tq);
+            ows[index] = weight(kmWeight, ow, tw);
+        }
     }
     // class Row
 
@@ -515,7 +571,40 @@
         return r1.interpolateWQ(r2, km, steps);
     }
 
-    public boolean interpolate(
+    public void interpolate(
+        double []   kms,
+        double []   ws,
+        double []   qs,
+        QPosition   qPosition,
+        LinearRemap remap
+    ) {
+        int R1 = rows.size()-1;
+
+        Row kmKey = new Row();
+        for (int i = 0; i < kms.length; ++i) {
+            kmKey.km = kms[i];
+            int rowIndex = Collections.binarySearch(rows, kmKey);
+
+            if (rowIndex >= 0) {
+                // direct row match
+                rows.get(rowIndex).storeWQ(qPosition, remap, ws, qs, i);
+                continue;
+            }
+
+            rowIndex = -rowIndex -1;
+
+            if (rowIndex < 1 || rowIndex >= R1) {
+                // do not extrapolate
+                ws[i] = qs[i] = Double.NaN;
+                continue;
+            }
+            Row r1 = rows.get(rowIndex-1);
+            Row r2 = rows.get(rowIndex);
+            r1.storeWQ(r2, kms[i], qPosition, remap, ws, qs, i);
+        }
+    }
+
+    public QPosition interpolate(
         double    q,
         int       referenceIndex,
         double [] kms,
@@ -538,7 +627,7 @@
 
             if (rowIndex < 1 || rowIndex >= R1) {
                 // do not extrapolate
-                return false;
+                return null;
             }
 
             // interpolation case
@@ -550,7 +639,7 @@
 
         if (qPosition == null) {
             // we cannot locate q in reference row
-            return false;
+            return null;
         }
 
         for (int i = 0; i < kms.length; ++i) {
@@ -575,7 +664,7 @@
             r1.storeWQ(r2, kms[i], qPosition, ws, qs, i);
         }
 
-        return true;
+        return qPosition;
     }
 
     public static WstValueTable getTable(River river) {

http://dive4elements.wald.intevation.org