Mercurial > dive4elements > river
changeset 655:913b52064449
Refactored version of "Berechnung 4"
flys-artifacts/trunk@2053 c6561f87-3c4e-4783-a992-168aeb5c3f6f
line wrap: on
line diff
--- a/flys-artifacts/ChangeLog Fri Jun 03 10:21:13 2011 +0000 +++ b/flys-artifacts/ChangeLog Sun Jun 05 18:24:46 2011 +0000 @@ -1,3 +1,51 @@ +2011-06-05 Sascha L. Teichmann <sascha.teichmann@intevation.de> + + Refactored version of "Berechnung 4" + + * src/main/java/de/intevation/flys/artifacts/model/Segment.java: + Added instance fields for a reference point (= location of gauge) + and backup of values (needed for naming). + + * src/main/java/de/intevation/flys/artifacts/model/WQCKms.java: + Added a constructor to be created from a WQKms. This is helpful + if a WQKms is replaced by a back jump correction. + + * src/main/java/de/intevation/flys/artifacts/model/Calculation4.java: + New. Outfactored version of "W bei ungleichmaessigen Abflusslaengsschnitt". + Much cleaner now and it should have a better handling of the corner + cases. + + * src/main/java/de/intevation/flys/artifacts/model/WstValueTable.java: + Removed the linear interpolation stuff. It is now in Linear. Removed + the LinearRemap interpolation method because it is not needed any + longer. Added a method to interpolate a given km with a given + function. + + * src/main/java/de/intevation/flys/artifacts/WINFOArtifact.java: + Removed the old calc 4 and used the new one. + + * src/main/java/de/intevation/flys/artifacts/math/LinearRemap.java: + Deleted. Not needed any longer. + + * src/main/java/de/intevation/flys/artifacts/math/Function.java: + New. Interface for a uni-variate real function. + + * src/main/java/de/intevation/flys/artifacts/math/Identity.java: + New. Implements Function with f(x) = x + + * src/main/java/de/intevation/flys/artifacts/math/Linear.java: + New. Implements Function with f(x) = m*x + b + + * src/main/java/de/intevation/flys/artifacts/FLYSArtifact.java: + Factored some stuff out to DoubleUtil. Removed some dead code. + Does some rounding correct. + + * src/main/java/de/intevation/flys/utils/DoubleUtil.java: New. + Centralized utils surrounding common double operations. + + * src/main/java/de/intevation/flys/exports/ChartInfoGenerator.java: + Removed superfluous imports. + 2011-06-03 Ingo Weinzierl <ingo@intevation.de> flys/issue90(Diagramm: Trennung der Diagrammfläche und Achsen aufheben)
--- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/FLYSArtifact.java Fri Jun 03 10:21:13 2011 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/FLYSArtifact.java Sun Jun 05 18:24:46 2011 +0000 @@ -33,6 +33,8 @@ import de.intevation.artifactdatabase.state.StateEngine; import de.intevation.artifactdatabase.transition.TransitionEngine; +import de.intevation.flys.utils.DoubleUtil; + import de.intevation.flys.model.Gauge; import de.intevation.flys.model.Range; import de.intevation.flys.model.River; @@ -79,8 +81,6 @@ /** The default step width between the start end end kilometer.*/ public static final double DEFAULT_KM_STEPS = 0.1; - public static final double DEFAULT_PRECISION = 1e6; - /** The identifier of the current state. */ protected String currentStateId; @@ -611,6 +611,36 @@ } } + public double [] getFromToStep() { + if (!isRange()) { + return null; + } + double [] fromTo = getDistance(); + + if (fromTo == null) { + return null; + } + + StateData dStep = getData("ld_step"); + if (dStep == null) { + return null; + } + + double [] result = new double[3]; + result[0] = fromTo[0]; + result[1] = fromTo[1]; + + try { + String step = (String)dStep.getValue(); + result[2] = DoubleUtil.round(Double.parseDouble(step) / 1000d); + } + catch (NumberFormatException nfe) { + return null; + } + + return result; + } + /** * Returns the gauge based on the current distance and river. @@ -967,36 +997,9 @@ double to, double step ) { - return getExplodedValues(from, to, step, DEFAULT_PRECISION); + return DoubleUtil.explode(from, to, step); } - public static double[] getExplodedValues( - double from, - double to, - double step, - double precision - ) { - double lower = from; - - double diff = to - from; - double tmp = diff / step; - int num = (int)Math.abs(Math.ceil(tmp)) + 1; - - double [] values = new double[num]; - - if (from > to) { - step = -step; - } - - for (int idx = 0; idx < num; idx++) { - values[idx] = Math.round(lower * precision)/precision; - lower += step; - } - - return values; - } - - /** * Method to dump the artifacts state/data. */
--- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/WINFOArtifact.java Fri Jun 03 10:21:13 2011 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/WINFOArtifact.java Sun Jun 05 18:24:46 2011 +0000 @@ -28,16 +28,15 @@ import de.intevation.flys.artifacts.states.DefaultState; import de.intevation.flys.artifacts.context.FLYSContext; -import de.intevation.flys.artifacts.math.BackJumpCorrector; + import de.intevation.flys.artifacts.model.MainValuesFactory; -import de.intevation.flys.artifacts.model.WQCKms; 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 de.intevation.flys.artifacts.model.Calculation4; +import de.intevation.flys.artifacts.model.Segment; import gnu.trove.TDoubleArrayList; @@ -583,204 +582,31 @@ return new WQKms[0]; } - boolean kmUp = river.getKmUp(); - - WstValueTable wst = WstValueTableFactory.getTable(river); - if (wst == null) { + WstValueTable table = WstValueTableFactory.getTable(river); + if (table == null) { logger.error("No wst found for selected river."); return new WQKms[0]; } - double [][] segments = getRanges(); + List<Segment> segments = getSegments(); - if (logger.isDebugEnabled()) { - logger.debug("segments ----------------- enter"); - for (double [] segment: segments) { - logger.debug(" " + joinDoubles(segment)); - } - logger.debug("segments ----------------- leave"); - } - - if (segments.length < 1) { - logger.warn("no segments given"); + if (segments == null) { + logger.error("Cannot create segments."); return new WQKms[0]; } - 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]); - HashSet<Integer> failed = new HashSet<Integer>(); - WQKms [] results = computeWaterlevelData - (kms, qs, wst, kmUp, failed); - setWaterlevelNames( - results, qs, isQ() ? "Q" : "W", failed); - return results; - } - } - - // more than one segment - - double [] boundKms = extractBoundsKm(river, segments); - - - if (logger.isDebugEnabled()) { - logger.debug("bound kms: " + joinDoubles(boundKms)); - } - - double [][] iqs = null; + double [] range = getFromToStep(); - 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])); - } + if (range == null) { + logger.error("Cannot figure out range."); + return new WQKms[0]; } - double [] boundWs = new double[boundKms.length]; - double [] boundQs = new double[boundKms.length]; - - double [] okms = getKms(new double [] { - boundKms[0], boundKms[boundKms.length-1] }); - - ArrayList<WQKms> results = new ArrayList<WQKms>(); - - int referenceIndex = kmUp ? 0 : boundKms.length-1; - - HashSet<Integer> failed = new HashSet<Integer>(); - - for (int i = 0; i < iqs.length; ++i) { - double [] iqsi = iqs[i]; - - QPosition qPosition = wst.interpolate( - iqsi[0], - boundKms[referenceIndex], - boundKms, boundWs, boundQs); - - if (qPosition == null) { - logger.warn("interpolation failed for " + iqsi[i]); - failed.add(i); - continue; - } - - LinearRemap remap = new LinearRemap(); - - 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]); - } - - double [] oqs = new double[okms.length]; - double [] ows = new double[okms.length]; - - 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)); - } - } - - WQKms [] wqkms = results.toArray(new WQKms[results.size()]); - - setDischargeLongitudinalSectionNames( - wqkms, iqs, isQ() ? "Q" : "W", failed); + Calculation4 calc4 = new Calculation4(segments, river, isQ()); - return wqkms; - } - - protected static String joinDoubles(double [] x) { - if (x == null) { - return ""; - } - 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; - } - + WQKms [] results = calc4.calculate(table, range[0], range[1], range[2]); - /** - * Sets the name for discharge longitudinal section curves where each WQKms - * in <i>r</i> represents a column. - */ - public static void setDischargeLongitudinalSectionNames( - WQKms [] wqkms, - double [][] iqs, - String wq, - Set<Integer> failed - ) { - logger.debug("WINFOArtifact.setDischargeLongitudinalSectionNames"); - - // TODO: I18N - - int pos = 0; - - for (int j = 0; j < iqs.length; ++j) { - if (failed.contains(j)) { - continue; - } - StringBuilder sb = new StringBuilder(wq) - .append(" benutzerdefiniert ("); - - double [] iqsi = iqs[j]; - for (int i = 0; i < iqsi.length; i++) { - if (i > 0) { - sb.append("; "); - } - sb.append(iqsi[i]); - } - sb.append(")"); - - wqkms[pos++].setName(sb.toString()); - } + return results; } } // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/math/Function.java Sun Jun 05 18:24:46 2011 +0000 @@ -0,0 +1,6 @@ +package de.intevation.flys.artifacts.math; + +public interface Function { + double value(double x); +} +// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/math/Identity.java Sun Jun 05 18:24:46 2011 +0000 @@ -0,0 +1,15 @@ +package de.intevation.flys.artifacts.math; + +public final class Identity +implements Function +{ + public static final Identity IDENTITY = new Identity(); + + public Identity() { + } + + public double value(double x) { + return x; + } +} +// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/math/Linear.java Sun Jun 05 18:24:46 2011 +0000 @@ -0,0 +1,68 @@ +package de.intevation.flys.artifacts.math; + +public final class Linear +implements Function +{ + private double m; + private double b; + + public Linear( + 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) { + m = 0; + b = 0.5*(y1 + y2); + } + else { + m = (y2 - y1)/(x2 - x1); + b = y1 - m*x1; + } + } + + 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; + } + + @Override + public double value(double x) { + return m*x + 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; + return a + factor*(b-a); + } +} +// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :
--- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/math/LinearRemap.java Fri Jun 03 10:21:13 2011 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,92 +0,0 @@ -package de.intevation.flys.artifacts.math; - -import org.apache.log4j.Logger; - -public class LinearRemap -{ - private static Logger logger = Logger.getLogger(LinearRemap.class); - - 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; - } - - if (from > to) { double t = from; from = to; to = t; } - - head = new Segment(from, to, m, b, head); - - if (logger.isDebugEnabled()) { - logger.debug("LinearRemap.add --------- enter"); - logger.debug(" range: [" + from + ", " + to + "]"); - logger.debug(" " + x1 + " -> " + y1 + " (" + head.eval(x1) + ")"); - logger.debug(" " + x2 + " -> " + y2 + " (" + head.eval(x2) + ")"); - logger.debug("LinearRemap.add --------- leave"); - } - } - - 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 :
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/Calculation4.java Sun Jun 05 18:24:46 2011 +0000 @@ -0,0 +1,287 @@ +package de.intevation.flys.artifacts.model; + +import java.util.List; +import java.util.Arrays; +import java.util.Comparator; +import java.util.Collections; + +import de.intevation.flys.utils.DoubleUtil; + +import de.intevation.flys.model.River; +import de.intevation.flys.model.Gauge; + +import de.intevation.flys.artifacts.model.WstValueTable.QPosition; + +import de.intevation.flys.artifacts.math.Function; +import de.intevation.flys.artifacts.math.Linear; +import de.intevation.flys.artifacts.math.Identity; +import de.intevation.flys.artifacts.math.BackJumpCorrector; + +import org.apache.log4j.Logger; + +public class Calculation4 +{ + private static Logger logger = Logger.getLogger(Calculation4.class); + + public static final double MINIMAL_STEP_WIDTH = 1e-5; + + public static final Comparator<Segment> REF_CMP = + new Comparator<Segment>() { + public int compare(Segment a, Segment b) { + double d = a.referencePoint - b.referencePoint; + if (d < 0d) return -1; + return d > 0d ? +1 : 0; + } + }; + + protected List<Segment> segments; + + protected boolean isQ; + + public Calculation4() { + } + + public Calculation4(List<Segment> segments, River river, boolean isQ) { + + this.segments = segments; + this.isQ = isQ; + + int numResults = 1; + + // assign reference points + for (Segment segment: segments) { + Gauge gauge = river.determineGauge( + segment.getFrom(), segment.getTo()); + + segment.setReferencePoint(gauge != null + ? gauge.getStation().doubleValue() + : 0.5*(segment.getFrom() + segment.getTo())); + + double [] values = segment.values; + + if (numResults == -1) { + numResults = values.length; + } + else if (numResults != values.length) { + throw new IllegalArgumentException("wrong length of values"); + } + + // convert to Q if needed + if (!isQ && gauge != null) { + double [][] table = new DischargeTables( + river.getName(), gauge.getName()).getFirstTable(); + + // need the original values for naming + segment.backup(); + + for (int i = 0; i < values.length; ++i) { + values[i] = DischargeTables.getQForW(table, values[i]); + } + } + } // for all segments + + Collections.sort(segments, REF_CMP); + } + + public WQKms [] calculate( + WstValueTable table, + double from, double to, double step + ) { + boolean debug = logger.isDebugEnabled(); + + if (debug) { + logger.debug( + "calculate from " + from + " to " + to + " step " + step); + logger.debug("# segments: " + segments.size()); + } + + if (segments.isEmpty()) { + logger.debug("no segments found"); + return new WQKms[0]; + } + + int numResults = segments.get(0).values.length; + + if (numResults < 1) { + logger.debug("no values given"); + return new WQKms[0]; + } + + + WQKms [] results = new WQKms[numResults]; + if (debug) { + logger.debug("after last segment -> gleichwertig"); + } + for (int i = 0; i < results.length; ++i) { + results[i] = new WQKms(); + } + + if (Math.abs(step) < MINIMAL_STEP_WIDTH) { + step = MINIMAL_STEP_WIDTH; + } + + if (from > to) { + step = -step; + } + + QPosition [] qPositions = new QPosition[numResults]; + + Function [] functions = new Function[numResults]; + + double [] out = new double[2]; + + Segment sentinel = new Segment(Double.MAX_VALUE); + Segment s1 = sentinel, s2 = sentinel; + + for (double pos = from; + from < to ? pos <= to : pos >= to; + pos = DoubleUtil.round(pos + step) + ) { + if (pos < s1.referencePoint || pos > s2.referencePoint) { + if (debug) { + logger.debug("need to find new interval for " + pos); + } + // find new interval + if (pos <= segments.get(0).referencePoint) { + // before first segment -> "gleichwertig" + if (debug) { + logger.debug("before first segment -> gleichwertig"); + } + Segment first = segments.get(0); + double [] values = first.values; + double refPos = first.referencePoint; + for (int i = 0; i < qPositions.length; ++i) { + qPositions[i] = table.getQPosition( + refPos, values[i]); + } + sentinel.setReferencePoint(-Double.MAX_VALUE); + s1 = sentinel; + s2 = segments.get(0); + Arrays.fill(functions, Identity.IDENTITY); + } + else if (pos >= segments.get(segments.size()-1).referencePoint) { + // after last segment -> "gleichwertig" + if (debug) { + logger.debug("after last segment -> gleichwertig"); + } + Segment last = segments.get(segments.size()-1); + double [] values = last.values; + double refPos = last.referencePoint; + for (int i = 0; i < qPositions.length; ++i) { + qPositions[i] = table.getQPosition( + refPos, values[i]); + } + sentinel.setReferencePoint(Double.MAX_VALUE); + s1 = last; + s2 = sentinel; + Arrays.fill(functions, Identity.IDENTITY); + } + else { // "ungleichwertig" + // find matching interval + if (debug) { + logger.debug("in segments -> ungleichwertig"); + } + s1 = s2 = null; + for (int i = 1, N = segments.size(); i < N; ++i) { + Segment si1 = segments.get(i-1); + Segment si = segments.get(i); + if (debug) { + logger.debug("check " + pos + " in " + + si1.referencePoint + " - " + si.referencePoint); + } + if (pos >= si1.referencePoint + && pos <= si. referencePoint) { + s1 = si1; + s2 = si; + break; + } + } + + if (s1 == null) { + throw new IllegalStateException("no interval found"); + } + + Segment anchor, free; + + if (from > to) { anchor = s1; free = s2; } + else { anchor = s2; free = s1; } + + // build transforms based on "gleichwertiger" phase + for (int i = 0; i < qPositions.length; ++i) { + QPosition qi = table.getQPosition( + anchor.referencePoint, + anchor.values[i]); + + if ((qPositions[i] = qi) == null) { + // TODO: error report + functions[i] = Identity.IDENTITY; + } + else { + double qA = table.getQ(qi, anchor.referencePoint); + double qF = table.getQ(qi, free .referencePoint); + + functions[i] = Double.isNaN(qA) || Double.isNaN(qF) + ? Identity.IDENTITY + : new Linear( + qA, qF, + anchor.values[i], free.values[i]); + } + } // build transforms + } // "ungleichwertiges" interval + } // find matching interval + + for (int i = 0; i < qPositions.length; ++i) { + QPosition qPosition = qPositions[i]; + + if (qPosition == null) { + continue; + } + + if (table.interpolate(pos, out, qPosition, functions[i])) { + results[i].add(out[0], out[1], pos); + } + else { + // TODO: error report + } + } + } + + // Backjump correction + for (int i = 0; i < results.length; ++i) { + BackJumpCorrector bjc = new BackJumpCorrector(); + + double [] ws = results[i].getWs(); + double [] kms = results[i].getKms(); + + if (bjc.doCorrection(kms, ws)) { + results[i] = new WQCKms(results[i], bjc.getCorrected()); + } + // TODO: error report + } + + // name the curves + for (int i = 0; i < results.length; ++i) { + results[i].setName(createName(i)); + } + + return results; + } + + protected String createName(int index) { + // TODO: i18n + StringBuilder sb = new StringBuilder(isQ ? "Q" : "W"); + sb.append(" benutzerdefiniert ("); + for (int i = 0, N = segments.size(); i < N; ++i) { + if (i > 0) { + sb.append("; "); + } + Segment segment = segments.get(i); + sb.append((segment.backup != null + ? segment.backup + : segment.values)[index]); + } + sb.append(')'); + return sb.toString(); + } +} +// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :
--- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/Segment.java Fri Jun 03 10:21:13 2011 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/Segment.java Sun Jun 05 18:24:46 2011 +0000 @@ -9,6 +9,8 @@ import gnu.trove.TDoubleArrayList; +import de.intevation.flys.utils.DoubleUtil; + public class Segment implements Serializable { @@ -17,10 +19,16 @@ protected double from; protected double to; protected double [] values; + protected double [] backup; + protected double referencePoint; public Segment() { } + public Segment(double referencePoint) { + this.referencePoint = referencePoint; + } + public Segment(double from, double to, double [] values) { this.from = from; this.to = to; @@ -43,6 +51,42 @@ return sb.toString(); } + public void setFrom(double from) { + this.from = from; + } + + public void backup() { + backup = (double [])values.clone(); + } + + public double getFrom() { + return from; + } + + public void setTo(double to) { + this.to = to; + } + + public double getTo() { + return to; + } + + public void setValues(double [] values) { + this.values = values; + } + + public double [] getValues() { + return values; + } + + public void setReferencePoint(double referencePoint) { + this.referencePoint = referencePoint; + } + + public double getReferencePoint() { + return referencePoint; + } + public static List<Segment> parseSegments(String input) { ArrayList<Segment> segments = new ArrayList<Segment>(); @@ -62,7 +106,8 @@ vs.clear(); for (String valueStr: parts[2].split(",")) { - vs.add(Double.parseDouble(valueStr.trim())); + vs.add(DoubleUtil.round( + Double.parseDouble(valueStr.trim()))); } double [] values = vs.toNativeArray();
--- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/WQCKms.java Fri Jun 03 10:21:13 2011 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/WQCKms.java Sun Jun 05 18:24:46 2011 +0000 @@ -13,6 +13,16 @@ protected TDoubleArrayList cw; + public WQCKms() { + } + + public WQCKms(WQKms other, double [] cws) { + this.w = other.w; + this.q = other.q; + this.kms = other.kms; + this.cw = new TDoubleArrayList(cws); + } + public WQCKms(double[] kms, double[] qs, double[] ws, double[] cws) { super(kms, qs, ws);
--- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/WstValueTable.java Fri Jun 03 10:21:13 2011 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/WstValueTable.java Sun Jun 05 18:24:46 2011 +0000 @@ -2,8 +2,8 @@ import java.io.Serializable; - -import de.intevation.flys.artifacts.math.LinearRemap; +import de.intevation.flys.artifacts.math.Linear; +import de.intevation.flys.artifacts.math.Function; import java.util.Arrays; import java.util.ArrayList; @@ -109,7 +109,7 @@ double [] ows, WstValueTable table ) { - double kmWeight = factor(km, this.km, other.km); + double kmWeight = Linear.factor(km, this.km, other.km); QPosition qPosition = new QPosition(); @@ -117,7 +117,7 @@ if (table.getQPosition(km, iqs[i], qPosition) != null) { double wt = getW(qPosition); double wo = other.getW(qPosition); - ows[i] = weight(kmWeight, wt, wo); + ows[i] = Linear.weight(kmWeight, wt, wo); } else { ows[i] = Double.NaN; @@ -137,7 +137,7 @@ return new double[2][0]; } - double factor = factor(km, this.km, other.km); + double factor = Linear.factor(km, this.km, other.km); double [] splineQ = new double[W]; double [] splineW = new double[W]; @@ -146,8 +146,8 @@ double maxQ = -Double.MAX_VALUE; for (int i = 0; i < W; ++i) { - double wws = weight(factor, ws[i], other.ws[i]); - double wqs = weight( + double wws = Linear.weight(factor, ws[i], other.ws[i]); + double wqs = Linear.weight( factor, table.getQIndex(i, km), table.getQIndex(i, other.km)); @@ -231,7 +231,7 @@ return weight == 1.0 ? ws[index] - : weight(weight, ws[index-1], ws[index]); + : Linear.weight(weight, ws[index-1], ws[index]); } public double getW( @@ -239,7 +239,7 @@ double km, QPosition qPosition ) { - double kmWeight = factor(km, this.km, other.km); + double kmWeight = Linear.factor(km, this.km, other.km); int index = qPosition.index; double weight = qPosition.weight; @@ -251,11 +251,11 @@ ow = other.ws[index]; } else { - tw = weight(weight, ws[index-1], ws[index]); - ow = weight(weight, other.ws[index-1], other.ws[index]); + tw = Linear.weight(weight, ws[index-1], ws[index]); + ow = Linear.weight(weight, other.ws[index-1], other.ws[index]); } - return weight(kmWeight, tw, ow); + return Linear.weight(kmWeight, tw, ow); } } // class Row @@ -338,60 +338,45 @@ return r1.interpolateWQ(r2, km, steps, this); } - public void interpolate( - double [] kms, - double [] ws, - double [] qs, - QPosition qPosition, - LinearRemap remap - ) { - interpolate(kms, ws, qs, qPosition, remap, 0, kms.length); - } - - public void interpolate( - double [] kms, - double [] ws, - double [] qs, - QPosition qPosition, - LinearRemap remap, - int startIndex, - int length + public boolean interpolate( + double km, + double [] out, + QPosition qPosition, + Function qFunction ) { int R1 = rows.size()-1; - Row kmKey = new Row(); + out[1] = qFunction.value(getQ(qPosition, km)); + + if (Double.isNaN(out[1])) { + return false; + } QPosition nPosition = new QPosition(); - - for (int i = startIndex, end = startIndex+length; i < end; ++i) { - kmKey.km = kms[i]; - - qs[i] = remap.eval(kms[i], getQ(qPosition, kms[i])); - - if (getQPosition(kms[i], qs[i], nPosition) == null) { - ws[i] = Double.NaN; - continue; - } - - int rowIndex = Collections.binarySearch(rows, kmKey); + if (getQPosition(km, out[1], nPosition) == null) { + return false; + } - if (rowIndex >= 0) { - // direct row match - ws[i] = rows.get(rowIndex).getW(nPosition); - continue; - } - - rowIndex = -rowIndex -1; + int rowIndex = Collections.binarySearch(rows, new Row(km)); - if (rowIndex < 1 || rowIndex >= R1) { - // do not extrapolate - ws[i] = Double.NaN; - continue; - } - Row r1 = rows.get(rowIndex-1); - Row r2 = rows.get(rowIndex); - ws[i] = r1.getW(r2, kms[i], nPosition); + if (rowIndex >= 0) { + // direct row match + out[0] = rows.get(rowIndex).getW(nPosition); + return !Double.isNaN(out[0]); } + + rowIndex = -rowIndex -1; + + if (rowIndex < 1 || rowIndex >= R1) { + // do not extrapolate + return false; + } + + Row r1 = rows.get(rowIndex-1); + Row r2 = rows.get(rowIndex); + out[0] = r1.getW(r2, km, nPosition); + + return !Double.isNaN(out[0]); } public QPosition interpolate( @@ -480,7 +465,7 @@ else { qMin = qCurrent; qMax = qLast; } if (q > qMin && q < qMax) { - double weight = factor(q, qLast, qCurrent); + double weight = Linear.factor(q, qLast, qCurrent); return qPosition.set(i, weight); } qLast = qCurrent; @@ -502,41 +487,8 @@ } double q1 = columns[index-1].getQRangeTree().findQ(km); double q2 = columns[index ].getQRangeTree().findQ(km); - return weight(weight, q1, q2); + return Linear.weight(weight, q1, q2); } - 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; - return a + factor*(b-a); - } } // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :
--- a/flys-artifacts/src/main/java/de/intevation/flys/exports/ChartInfoGenerator.java Fri Jun 03 10:21:13 2011 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/exports/ChartInfoGenerator.java Sun Jun 05 18:24:46 2011 +0000 @@ -1,6 +1,5 @@ package de.intevation.flys.exports; -import java.awt.Color; import java.awt.Transparency; import java.io.IOException; import java.io.OutputStream; @@ -9,11 +8,8 @@ import org.apache.log4j.Logger; -import org.jfree.chart.ChartFactory; import org.jfree.chart.ChartRenderingInfo; import org.jfree.chart.JFreeChart; -import org.jfree.chart.plot.PlotOrientation; -import org.jfree.chart.plot.XYPlot; import de.intevation.artifacts.Artifact; import de.intevation.artifacts.CallContext;
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/utils/DoubleUtil.java Sun Jun 05 18:24:46 2011 +0000 @@ -0,0 +1,48 @@ +package de.intevation.flys.utils; + +public class DoubleUtil +{ + public static final double DEFAULT_STEP_PRECISION = 1e6; + + private DoubleUtil() { + } + + public static final double [] explode(double from, double to, double step) { + return explode(from, to, step, DEFAULT_STEP_PRECISION); + } + + public static final double round(double x, double precision) { + return Math.round(x * precision)/precision; + } + + public static final double round(double x) { + return Math.round(x * DEFAULT_STEP_PRECISION)/DEFAULT_STEP_PRECISION; + } + + public static final double [] explode( + double from, + double to, + double step, + double precision + ) { + double lower = from; + + double diff = to - from; + double tmp = diff / step; + int num = (int)Math.abs(Math.ceil(tmp)) + 1; + + double [] values = new double[num]; + + if (from > to) { + step = -step; + } + + for (int idx = 0; idx < num; idx++) { + values[idx] = round(lower, precision); + lower += step; + } + + return values; + } +} +// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :