Mercurial > dive4elements > river
changeset 335:64cfbd631f29
Implemented the "Wasserstand/Wasserspiegellagen" calculation.
flys-artifacts/trunk@1733 c6561f87-3c4e-4783-a992-168aeb5c3f6f
author | Sascha L. Teichmann <sascha.teichmann@intevation.de> |
---|---|
date | Tue, 19 Apr 2011 17:43:39 +0000 |
parents | b7c8df643dc4 |
children | 7f13ed751277 |
files | flys-artifacts/ChangeLog flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/WstValueTable.java |
diffstat | 2 files changed, 187 insertions(+), 8 deletions(-) [+] |
line wrap: on
line diff
--- a/flys-artifacts/ChangeLog Tue Apr 19 17:37:21 2011 +0000 +++ b/flys-artifacts/ChangeLog Tue Apr 19 17:43:39 2011 +0000 @@ -1,3 +1,16 @@ +2011-04-19 Sascha L. Teichmann <sascha.teichmann@intevation.de> + + * src/main/java/de/intevation/flys/artifacts/model/WstValueTable.java: + add a method interpolateW() method which takes an array of + q values and returns an equal sized array of w values. + This is essentially the "Wasserstand/Wasserspiegellagen" calculation + of desktop FLYS. + + If you want to do a calculation with given w values you have + to convert the w values with DischargeTables.getQForW() first. + + !!! This code needs heavy testing !!! + 2011-04-19 Sascha L. Teichmann <sascha.teichmann@intevation.de> * src/main/java/de/intevation/flys/artifacts/model/DischargeTables.java:
--- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/WstValueTable.java Tue Apr 19 17:37:21 2011 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/WstValueTable.java Tue Apr 19 17:43:39 2011 +0000 @@ -8,6 +8,7 @@ import de.intevation.flys.backend.SessionHolder; +import java.util.Arrays; import java.util.ArrayList; import java.util.Comparator; import java.util.List; @@ -53,17 +54,139 @@ // class Column public static class Row - implements Serializable + implements Serializable, Comparable<Row> { double km; - double [] wq; + double [] ws; + double [] qs; + boolean qSorted; public Row() { } - public Row(double km, double [] wq) { + public Row(double km) { this.km = km; - this.wq = wq; + } + + public Row(double km, double [] ws, double [] qs) { + this(km); + this.ws = ws; + this.qs = qs; + } + + public static final double linear( + double x, + double x1, double x2, + double y1, double y2 + ) { + // y1 = m*x1 + b + // y2 = m*x2 + b + // y2 - y1 = m*(x2 - x1) + // m = (y2 - y1)/(x2 - x1) # x2 != x1 + // b = y1 - m*x1 + + if (x1 == x2) { + return 0.5*(y1 + y2); + } + double m = (y2 - y1)/(x2 - x1); + double b = y1 - m*x1; + return x*m + b; + } + + public static final double factor(double x, double p1, double p2) { + // 0 = m*p1 + b <=> b = -m*p1 + // 1 = m*p2 + b + // 1 = m*(p2 - p1) + // m = 1/(p2 - p1) # p1 != p2 + // f(x) = x/(p2-p1) - p1/(p2-p1) <=> (x-p1)/(p2-p1) + + return p1 == p2 ? 0.0 : (x-p1)/(p2-p1); + } + + public static final double weight(double factor, double a, double b) { + return (1.0-factor)*a + factor*b; + } + + public double getWForKM(Row other, int index, double km) { + double w1 = ws[index]; + double w2 = other.ws[index]; + double km1 = this.km; + double km2 = other.km; + return linear(km, km1, km2, w1, w2); + } + + public void interpolateW( + Row other, + double km, + double [] input, + double [] output + ) { + double factor = factor(km, this.km, other.km); + + for (int i = 0; i < input.length; ++i) { + double q = input[i]; + int idx1 = getQIndex(q); + int idx2 = other.getQIndex(q); + + double w1 = idx1 >= 0 + ? ws[idx1] + : interpolateW(-idx1-1, q); + + double w2 = idx2 >= 0 + ? other.ws[idx2] + : other.interpolateW(-idx2-1, q); + + output[i] = weight(factor, w1, w2); + } + } + + public double interpolateW(int idx, double q) { + return idx < 1 || idx >= qs.length + ? Double.NaN // do not extrapolate + : linear(q, qs[idx-1], qs[idx], ws[idx-1], ws[idx]); + } + + public double interpolateW(double q) { + if (Double.isNaN(q)) return Double.NaN; + int index = getQIndex(q); + return index >= 0 ? ws[index] : interpolateW(-index -1); + } + + public int getQIndex(double q) { + return qSorted ? binaryQIndex(q) : linearQIndex(q); + } + + public int binaryQIndex(double q) { + return Arrays.binarySearch(qs, q); + } + + public int linearQIndex(double q) { + switch (qs.length) { + case 0: break; + case 1: + if (qs[0] == q) return 0; + if (qs[0] < q) return -(1+1); + return -(0+1); + default: + for (int i = 1; i < qs.length; ++i) { + double qa = qs[i-1]; + double qb = qs[i]; + if (qa == q) return i-1; + if (qb == q) return i; + if (qa > qb) { double t = qa; qa = qb; qb = t; } + if (q > qa && q < qb) return -(i+1); + } + return -qs.length; + } + + return -1; + } + + public int compareTo(Row other) { + double d = km - other.km; + if (d < 0.0) return -1; + if (d > 0.0) return +1; + return 0; } } // class Row @@ -81,6 +204,36 @@ this.columns = columns; } + + public double [] interpolateW(double km, double [] qs) { + + double [] ws = new double[qs.length]; + + int rowIndex = Collections.binarySearch(rows, new Row(km)); + + if (rowIndex >= 0) { // direct row match + Row row = rows.get(rowIndex); + for (int i = 0; i < qs.length; ++i) { + ws[i] = row.interpolateW(qs[i]); + } + } + else { // needs bilinear interpolation + rowIndex = -rowIndex -1; + + if (rowIndex < 1 || rowIndex >= rows.size()) { + // do not extrapolate + Arrays.fill(ws, Double.NaN); + } + else { + rows.get(rowIndex-1).interpolateW( + rows.get(rowIndex), + km, qs, ws); + } + } + + return ws; + } + public static WstValueTable getTable(River river) { return getTable(river, 0); } @@ -133,6 +286,9 @@ int lastColumnNo = -1; Row row = null; + Double lastQ = -Double.MAX_VALUE; + boolean qSorted = true; + for (Iterator<Object []> iter = sqlQuery.iterate(); iter.hasNext();) { Object [] result = iter.next(); double km = (Double) result[0]; @@ -142,14 +298,24 @@ if (columnNo > lastColumnNo) { // new row if (row != null) { + row.qSorted = qSorted; valueTable.rows.add(row); } - row = new Row(km, new double[(columnNo+1) << 1]); + row = new Row( + km, + new double[columnNo+1], + new double[columnNo+1]); + lastQ = -Double.MAX_VALUE; + qSorted = true; } - int idx = columnNo << 1; - row.wq[idx ] = w != null ? w : Double.NaN; - row.wq[idx+1] = q != null ? q : Double.NaN; + row.ws[columnNo] = w != null ? w : Double.NaN; + row.qs[columnNo] = q != null ? q : Double.NaN; + + if (qSorted && (q == null || lastQ > q)) { + qSorted = false; + } + lastQ = q; lastColumnNo = columnNo; }