/** Copyright (C) 2017 by Bundesanstalt für Gewässerkunde
 * Software engineering by
 *  Björnsen Beratende Ingenieure GmbH
 *  Dr. Schumacher Ingenieurbüro für Wasser und Umwelt
 *
 * This file is Free Software under the GNU AGPL (>=v3)
 * and comes with ABSOLUTELY NO WARRANTY! Check out the
 * documentation coming with Dive4Elements River for details.
 */
package org.dive4elements.river.artifacts.sinfo.common;

import org.apache.commons.lang.math.DoubleRange;
import org.apache.commons.math.FunctionEvaluationException;
import org.apache.commons.math.analysis.UnivariateRealFunction;
import org.apache.commons.math.analysis.interpolation.LinearInterpolator;
import org.apache.log4j.Logger;
import org.dive4elements.river.artifacts.model.Calculation;
import org.dive4elements.river.model.Gauge;
import org.dive4elements.river.model.MainValue;
import org.dive4elements.river.model.MainValueType.MainValueTypeKey;

import gnu.trove.TDoubleArrayList;

/**
 * Search/interpolation of the duration main values of a gauge
 *
 * @author Matthias Schäfer
 *
 */
public final class GaugeDurationValuesFinder {

    /***** FIELDS *****/

    private static Logger log = Logger.getLogger(GaugeDurationValuesFinder.class);

    private final Gauge gauge;

    private Calculation problems;

    private UnivariateRealFunction qInterpolator;

    private DoubleRange qRange;

    private UnivariateRealFunction durInterpolator;

    private DoubleRange durRange;


    /***** CONSTRUCTORS *****/

    private GaugeDurationValuesFinder(final Gauge gauge, final Calculation problems) {
        // Load the duration main values from the database (each duration has a Q)
        this.gauge = gauge;
        this.problems = problems;
        final TDoubleArrayList qs = new TDoubleArrayList();
        final TDoubleArrayList durs = new TDoubleArrayList();
        for (final MainValue v : MainValue.getValuesOfGaugeAndType(gauge, MainValueTypeKey.DURATION)) {
            qs.add(v.getValue().doubleValue());
            durs.add(Integer.valueOf(v.getMainValue().getName()).doubleValue());
        }
        // Merge out-of-order Qs
        final TDoubleArrayList q1s = new TDoubleArrayList();
        final TDoubleArrayList dur1s = new TDoubleArrayList();
        boolean writeWarn = true;
        int cnt = 0;
        int i = -1;
        while ((i + 1) <= qs.size() - 1) {
            i++;
            if ((i == 0) || (qs.getQuick(i - 1) < qs.getQuick(i) - 0.001)) {
                q1s.add(qs.getQuick(i));
                dur1s.add(durs.getQuick(i));
                cnt = 1;
            }
            else {
                q1s.set(q1s.size() - 1, (q1s.getQuick(q1s.size() - 1) * cnt + qs.getQuick(i)) / (cnt + 1));
                dur1s.set(dur1s.size() - 1, (dur1s.getQuick(dur1s.size() - 1) * cnt + durs.getQuick(i)) / (cnt + 1));
                cnt++;
                if (writeWarn) {
                    log.warn("Nicht streng monotone Q in D(Q) des Pegels " + gauge.getName());
                    writeWarn = false;
                }
            }
        }
        // Build the duration-by-Q interpolator
        try {
            this.qInterpolator = new LinearInterpolator().interpolate(q1s.toNativeArray(), dur1s.toNativeArray());
            this.qRange = new DoubleRange(q1s.getQuick(0), q1s.getQuick(q1s.size() - 1));
        }
        catch (final Exception e) {
            this.qInterpolator = null;
            this.qRange = null;
        }
        // Load the Q values by duration from the database
        qs.clear();
        durs.clear();
        for (final MainValue v : MainValue.getDurationDischargesOfGauge(gauge)) {
            durs.add(Integer.valueOf(v.getMainValue().getName()).doubleValue());
            qs.add(v.getValue().doubleValue());
        }
        // Build the Q-by-duration interpolator
        try {
            this.durInterpolator = new LinearInterpolator().interpolate(durs.toNativeArray(), qs.toNativeArray());
            this.durRange = new DoubleRange(durs.getQuick(0), durs.getQuick(durs.size() - 1));
        }
        catch (final Exception e) {
            this.durInterpolator = null;
            this.durRange = null;
        }
        // Report problems
        if (((this.qInterpolator == null) || (this.durInterpolator == null)) && (this.problems != null)) {
            this.problems.addProblem("gauge_duration.missing", gauge.getName());
            // Report only once
            this.problems = null;
        }
    }


    /***** METHODS *****/

    /**
     * Loads the the discharge-duration table of a gauge ({gauge}.sta)
     *
     * @return The main values finder of a a gauge, or null
     */
    public static GaugeDurationValuesFinder loadValues(final Gauge gauge, final Calculation problems) {
        return new GaugeDurationValuesFinder(gauge, problems);
    }

    /**
     * If this provider may return valid data at all.
     */
    public boolean isValid() {
        return (this.qInterpolator != null);
    }

    /**
     * Discharge for a duration
     *
     * @return Q, or NegInf for duration less than all, or PosInf for duration greater then all, or NaN in case of exception
     */
    public double getDischarge(final int duration) {
        try {
            if (this.durInterpolator == null)
                return Double.NaN;
            else if (duration < this.durRange.getMinimumDouble())
                return Double.NEGATIVE_INFINITY;
            else if (duration > this.durRange.getMaximumDouble())
                return Double.POSITIVE_INFINITY;
            else
                return this.durInterpolator.value(duration);
        }
        catch (@SuppressWarnings("unused") final FunctionEvaluationException e) {
            // ignore exception because this can/will happen regularly
            return Double.NaN;
        }
    }

    /**
     * Duration for a discharge
     *
     * @return duration, or 0 for Q less than all, or 365 for duration greater then all, or NaN in case of missing
     *         interpolator or exception
     */
    public double getDuration(final double q) {
        try {
            if (this.qInterpolator == null)
                return Double.NaN;
            else if (q < this.qRange.getMinimumDouble())
                return 0;
            else if (q > this.qRange.getMaximumDouble())
                return 365;
            else
                return this.qInterpolator.value(q);
        }
        catch (@SuppressWarnings("unused") final FunctionEvaluationException e) {
            // ignore exception because this can/will happen regularly
            return Double.NaN;
        }
    }
}