view flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/Calculation4.java @ 4187:21f4e4b79121

Refactor GaugeDischargeCurveFacet to be able to set a facet name For adding another output of the GaugeDischargeCurveArtifact it is necessary to provide to facet instances with different names. Therefore the GaugeDischargeCurveFacet is extended to set the facet name in the constructor.
author Björn Ricks <bjoern.ricks@intevation.de>
date Fri, 19 Oct 2012 13:25:49 +0200
parents 22790758b132
children 39885bdfc6fc
line wrap: on
line source
package de.intevation.flys.artifacts.model;

import de.intevation.flys.artifacts.math.BackJumpCorrector;
import de.intevation.flys.artifacts.math.Function;
import de.intevation.flys.artifacts.math.Identity;
import de.intevation.flys.artifacts.math.Linear;

import de.intevation.flys.artifacts.model.WstValueTable.QPosition;

import de.intevation.flys.model.River;

import de.intevation.flys.utils.DoubleUtil;

import java.util.Arrays;
import java.util.List;

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;

    protected List<Segment> segments;

    protected boolean       isQ;

    public Calculation4() {
    }

    public Calculation4(List<Segment> segments, River river, boolean isQ) {

        this.segments = segments;
        this.isQ      = isQ;

        Segment.setReferencePointConvertQ(segments, river, isQ, this);
    }

    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");
            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");
            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) {
                            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 {
                    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