view flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/Calculation4.java @ 2089:0da8874bd378

Added initial state to map artifact to be able to advance and step back. The map artifact overrides describe() to have the complete UI information in the describe response document. flys-artifacts/trunk@3613 c6561f87-3c4e-4783-a992-168aeb5c3f6f
author Raimund Renkert <raimund.renkert@intevation.de>
date Fri, 06 Jan 2012 12:02:10 +0000
parents c09c9e05ecfa
children 2898b1ff6013
line wrap: on
line source
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
extends      Calculation
{
    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.maxOverlap(segment.getFrom(), segment.getTo());

            if (gauge == null) {
                logger.warn("no gauge found. Defaults to mid point.");
                segment.setReferencePoint(
                    0.5*(segment.getFrom()+segment.getTo()));
            }
            else {
                double ref = gauge.getStation().doubleValue();
                logger.debug(
                    "reference gauge: " + gauge.getName() +
                    " (km " + ref + ")");
                segment.setReferencePoint(ref);
            }

            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 CalculationResult 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());
            for (Segment segment: segments) {
                logger.debug("  " + segment);
            }
        }

        if (segments.isEmpty()) {
            logger.debug("no segments found");
            // TODO: I18N
            addProblem("no segments found");
            return new CalculationResult(new WQKms[0], this);
        }

        int numResults = segments.get(0).values.length;

        if (numResults < 1) {
            logger.debug("no values given");
            // TODO: I18N
            addProblem("no values given");
            return new CalculationResult(new WQKms[0], this);
        }


        WQKms [] results = new WQKms[numResults];
        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: I18N
                            addProblem(pos, "cannot find q = " + anchor.values[i]);
                            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]);

                            if (debug) {
                                logger.debug(
                                    anchor.referencePoint + ": " +
                                    qA + " -> " + functions[i].value(qA) +
                                    " / " + free.referencePoint + ": " +
                                    qF + " -> " + functions[i].value(qF));
                            }
                        }
                    } // 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: I18N
                    addProblem(pos, "cannot interpolate w/q");
                }
            }
        }

        // 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, this)) {
                results[i] = new WQCKms(results[i], bjc.getCorrected());
            }
        }

        // name the curves
        for (int i = 0; i < results.length; ++i) {
            results[i].setName(createName(i));
        }

        return new CalculationResult(results, this);
    }

    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 :

http://dive4elements.wald.intevation.org