/* Copyright (C) 2011, 2012, 2013 by Bundesanstalt für Gewässerkunde
 * Software engineering by Intevation GmbH
 *
 * 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.model.minfo;

import gnu.trove.TDoubleArrayList;

import java.util.ArrayList;
import java.util.List;

import org.apache.log4j.Logger;

import org.dive4elements.river.artifacts.access.SedimentLoadAccess;
import org.dive4elements.river.artifacts.model.Calculation;
import org.dive4elements.river.artifacts.model.CalculationResult;
import org.dive4elements.river.artifacts.model.Range;


/** Calculate sediment load. */
public class SedimentLoadCalculation
extends Calculation
{

    /** Private logger. */
    private static final Logger logger = Logger
        .getLogger(SedimentLoadCalculation.class);

    protected String river;
    protected String yearEpoch;
    protected double kmUp;
    protected double kmLow;
    protected int[] period;
    /** Years of chosen epochs.
     * epoch[0][0] is the start year of first entered epoch, e[0][1] the end-year. */
    protected int[][] epoch;
    protected String unit;

    public SedimentLoadCalculation() {
    }

    /** Returns CalculationResult with array of SedimentLoadResults. */
    public CalculationResult calculate(SedimentLoadAccess access) {
        logger.info("SedimentLoadCalculation.calculate");

        String river = access.getRiverName();
        String yearEpoch = access.getYearEpoch();
        String unit = access.getUnit();
        int[] period = null;
        int[][] epoch = null;
        double kmUp = access.getUpperKM();
        double kmLow = access.getLowerKM();
        if (yearEpoch.equals("year")) {
            period = access.getPeriod();
            epoch = null;
        }
        else if (yearEpoch.equals("epoch") || yearEpoch.equals("off_epoch")) {
            epoch = access.getEpochs();
            period = null;
        }
        else {
            addProblem("minfo.missing.year_epoch");
        }

        if (river == null) {
            // TODO: i18n
            addProblem("minfo.missing.river");
        }

        if (period == null && epoch == null) {
            addProblem("minfo.missing.time");
        }

        if (!hasProblems()) {
            this.river = river;
            this.yearEpoch = yearEpoch;
            this.unit = unit;
            this.period = period;
            this.epoch = epoch;
            this.kmUp = kmUp;
            this.kmLow = kmLow;
            return internalCalculate();
        }

        return new CalculationResult();
    }

    /** Returns CalculationResult with array of SedimentLoadResults. */
    private CalculationResult internalCalculate() {
        logger.debug("internalCalulate; mode:" + yearEpoch);
        if (yearEpoch.equals("year")) {
            List<SedimentLoadResult> results =
                new ArrayList<SedimentLoadResult>();
            for (int i = period[0]; i <= period[1]; i++) {
                SedimentLoadResult res = calculateYear(i);
                results.add(res);
            }
            return new CalculationResult(
                results.toArray(new SedimentLoadResult[results.size()]), this);
        }
        else if (yearEpoch.equals("epoch")) {
            List<SedimentLoadResult> results =
                new ArrayList<SedimentLoadResult>();
            for (int i = 0; i < epoch.length; i++) {
                SedimentLoadResult res = calculateEpoch(i);
                results.add(res);
            }
            return new CalculationResult(
                results.toArray(new SedimentLoadResult[results.size()]), this);
        }
        else if (yearEpoch.equals("off_epoch")) {
            List<SedimentLoadResult> results =
                new ArrayList<SedimentLoadResult>();
            for (int i = 0; i < epoch.length; i++) {
                SedimentLoadResult res = calculateOffEpoch(i);
                results.add(res);
            }
            return new CalculationResult(
                results.toArray(new SedimentLoadResult[results.size()]), this);
        }
        else {
            logger.error("Unknown mode " + yearEpoch);
        }
        return null;
    }


    /** Returns val if not NaN, 0d otherwise. */
    private static double makeNaN0(double val) {
        return Double.isNaN(val)
            ? 0d
            : val;
    }

    /**
     * Take Loads and build average of all fractions at given km.
     * The average fractions value is set in resLoad.
     *
     * @param km km at which to build fractions average.
     * @param epochLoads the loads to build average over.
     * @param[out] resLoad resulting SedimentLoad.
     */
    private void calculateEpochKm(
        List<SedimentLoadLSData> epochLoads,
        SedimentLoadLSData resLoad,
        double km
    ) {
        int cSum = 0;
        int fmSum = 0;
        int sSum = 0;
        int ssSum = 0;
        int ssbSum = 0;
        int sseSum = 0;
        for (SedimentLoadLSData load : epochLoads) {
            SedimentLoadFraction f = load.getFraction(km);
            if (f.getCoarse() > 0d) {
                double c = makeNaN0(resLoad.getFraction(km).getCoarse());
                resLoad.setCoarse(km, c + f.getCoarse(), f.getCoarseRange());
                cSum++;
            }
            if (f.getFineMiddle() > 0d) {
                double fm = makeNaN0(resLoad.getFraction(km).getFineMiddle());
                resLoad.setFineMiddle(km, fm + f.getFineMiddle(), f.getFineMiddleRange());
                fmSum++;
            }
            if (f.getSand() > 0d) {
                double s = makeNaN0(resLoad.getFraction(km).getSand());
                resLoad.setSand(km, s + f.getSand(), f.getSandRange());
                sSum++;
            }
            if (f.getSuspSand() > 0d) {
                double s = makeNaN0(resLoad.getFraction(km).getSuspSand());
                resLoad.setSuspSand(km, s + f.getSuspSand(), f.getSuspSandRange());
                ssSum++;
            }
            if (f.getSuspSandBed() > 0d) {
                double s = makeNaN0(resLoad.getFraction(km).getSuspSandBed());
                resLoad.setSuspSandBed(km, s + f.getSuspSandBed(), f.getSuspSandBedRange());
                ssbSum++;
            }
            if (f.getSuspSediment() > 0d) {
                double s = makeNaN0(resLoad.getFraction(km).getSuspSediment());
                resLoad.setSuspSediment(km, s + f.getSuspSediment(), f.getSuspSedimentRange());
                sseSum++;
            }
        }

        SedimentLoadFraction fr = resLoad.getFraction(km);
        // Prevent divisions by zero.
        if (cSum != 0) {
            resLoad.setCoarse(km, fr.getCoarse()/cSum, fr.getCoarseRange());
        }
        if (fmSum != 0) {
            resLoad.setFineMiddle(km, fr.getFineMiddle()/fmSum,
                fr.getFineMiddleRange());
        }
        if (sSum != 0) {
            resLoad.setSand(km, fr.getSand()/sSum, fr.getSandRange());
        }
        if (ssSum != 0) {
            resLoad.setSuspSand(km, fr.getSuspSand()/ssSum,
                fr.getSuspSandRange());
        }
        if (ssbSum != 0) {
            resLoad.setSuspSandBed(km, fr.getSuspSandBed()/ssbSum,
                fr.getSuspSandBedRange());
        }
        if (sseSum != 0) {
            resLoad.setSuspSediment(km, fr.getSuspSediment()/sseSum, fr.getSuspSedimentRange());
        }
    }

    /**
     * Calculate result for the ith given epoch.
     * @param i index of epoch (if multiple given).
     */
    private SedimentLoadResult calculateEpoch(int i) {
        List<SedimentLoadLSData> epochLoads = new ArrayList<SedimentLoadLSData>();
        for (int j = epoch[i][0]; j <= epoch[i][1]; j++) {
            epochLoads.add(SedimentLoadFactory.getLoadWithData(
                this.river,
                this.yearEpoch,
                this.kmLow,
                this.kmUp,
                j, //syear
                j)); //eyear
        }

        SedimentLoadLSData resLoad = new SedimentLoadLSData();
        TDoubleArrayList kms = new TDoubleArrayList();

        for (SedimentLoadLSData load : epochLoads) {
            for (double km : load.getKms()) {
                if (!kms.contains(km)) {
                    kms.add(km);
                }
            }
        }

        for (int j = 0, J = kms.size(); j < J; j++) {
            calculateEpochKm(epochLoads, resLoad, kms.get(j));
        }
        resLoad.setDescription("");
        resLoad.setEpoch(true);

        SedimentLoadLSData sl = calculateTotalLoad(resLoad, this.epoch[i][0]);

        if (this.unit.equals("m3_per_a")) {
            sl = calculateUnit(sl, this.epoch[i][0]);
        }

        return new SedimentLoadResult(
            this.epoch[i][0],
            this.epoch[i][1],
            sl);
    }

    /**
     * Calculate/Fetch values at off. epochs.
     * @param i index in epochs.
     */
    private SedimentLoadResult calculateOffEpoch(int i) {
        SedimentLoadLSData load = SedimentLoadFactory.getLoadWithData(
            this.river,
            this.yearEpoch,
            kmLow,
            kmUp,
            this.epoch[i][0],
            this.epoch[i][1]);
        SedimentLoadResult result;
        SedimentLoadLSData sl = calculateTotalLoad(load, this.epoch[i][0]);
        if (unit.equals("m3_per_a")) {
            SedimentLoadLSData slu = calculateUnit(sl, epoch[i][0]);
            result = new SedimentLoadResult(
                this.epoch[i][0],
                this.epoch[i][1],
                slu);
        }
        else {
            result = new SedimentLoadResult(
                this.epoch[i][0],
                this.epoch[i][1],
                sl);
        }

        return result;
    }

    /**
     * Fetch loads for a single year, calculate total and
     * return the result containing both.
     * @param y year, e.g. 1980
     */
    private SedimentLoadResult calculateYear(int y) {
        SedimentLoadLSData load = SedimentLoadFactory.getLoadWithData(
            this.river,
            this.yearEpoch,
            this.kmLow,
            this.kmUp,
            y,
            y);

        SedimentLoadResult result;
        SedimentLoadLSData sl = calculateTotalLoad(load, y);
        if (unit.equals("m3_per_a")) {
            SedimentLoadLSData slu = calculateUnit(sl, y);
            result = new SedimentLoadResult(y, 0, slu);
        }
        else {
            result = new SedimentLoadResult(y, 0, sl);
        }
        return result;
    }

    /** Add up the loads of a year. */
    private SedimentLoadLSData calculateTotalLoad(SedimentLoadLSData load, int year) {
        logger.debug("calculateTotalLoad");
        boolean problemThisYear = false;
        if (!load.hasCoarse()) {
            addProblem("missing.fraction.coarse", Integer.toString(year));
            problemThisYear = true;
        }
        if (!load.hasFineMiddle()) {
            addProblem("missing.fraction.fine_middle", Integer.toString(year));
            problemThisYear = true;
        }
        if (!load.hasSand()) {
            addProblem("missing.fraction.sand", Integer.toString(year));
            problemThisYear = true;
        }
        if (!load.hasSuspSand()) {
            addProblem("missing.fraction.susp_sand", Integer.toString(year));
            problemThisYear = true;
        }
        if (!load.hasSuspSediment()) {
            addProblem("missing.fraction.susp_sediment", Integer.toString(year));
            problemThisYear = true;
        }
        if (problemThisYear) {
            logger.warn("Some problem, not calculating total load.");
            return load;
        }
        return partialTotal(load);
    }



    /**
     * Set total values in load.
     *
     * Therefore, run over the sorted kms and find ranges where either all
     * or all Geschiebe or just the Schwebstoff fractions are set.
     * Merge these ranges and add (maybe new) respective fractions to
     * load. In the process, remember any 'unfished' ends from measurements
     * where the km-ranges did not completely match.
     *
     * @param load SedimentLoad to add total values (and ranges) to.
     * @return input param load, with total values set.
     */
    private SedimentLoadLSData partialTotal(SedimentLoadLSData load) {
        // The load with balanced ranges, will be returned.
        SedimentLoadLSData fairLoad = load;

        Range lastOtherRange = null;
        double lastOtherValue = 0d;

        Range lastSuspRange = null;
        double lastSuspValue = 0d;

        for (double km: load.getKms()) { // kms are already sorted!
            logger.debug ("Trying to add at km " + km);
            SedimentLoadFraction fraction = load.getFraction(km);
            if (fraction.isComplete()) {
                double total = fraction.getCoarse() +
                    fraction.getFineMiddle() +
                    fraction.getSand() +
                    fraction.getSuspSand() +
                    fraction.getSuspSediment();
                // Easiest case. Add values up and set'em.
                if (fraction.getCoarseRange().equals(
                    fraction.getSuspSedimentRange())) {
                    lastOtherRange = null;
                    lastSuspRange = null;
                    fairLoad.setTotal(km, total, fraction.getCoarseRange());
                }
                else {
                    // Need to split a range.
                    if (fraction.getCoarseRange().getEnd()
                        < fraction.getSuspSedimentRange().getEnd()) {
                        // Schwebstoff is longer.
                        // Adjust and remember schwebstoffs range and value.
                        lastSuspRange = (Range) fraction.getSuspSedimentRange().clone();
                        lastSuspRange.setStart(fraction.getCoarseRange().getEnd());
                        lastSuspValue = fraction.getSuspSediment();
                        lastOtherRange = null;
                        fairLoad.setTotal(km, total, fraction.getCoarseRange());
                    }
                    else {
                        // Geschiebe is longer.
                        // Adjust and remember other values.
                        lastOtherRange = (Range) fraction.getCoarseRange().clone();
                        lastOtherRange.setStart(fraction.getSuspSedimentRange().getEnd());
                        lastOtherValue = (total - fraction.getSuspSediment());
                        lastSuspRange = null;
                        fairLoad.setTotal(km, total, fraction.getSuspSedimentRange());
                    }
                }
            }
            else if (fraction.hasOnlySuspValues() && lastOtherRange != null) {
                // Split stuff.
                Range suspSedimentRange = fraction.getSuspSedimentRange();
                // if intersects with last other range, cool! merge and add!
                if (lastOtherRange.contains(km)) {
                    double maxStart = 0d;
                    double minEnd = 0d;
                    maxStart = Math.max(suspSedimentRange.getStart(),
                        lastOtherRange.getStart());

                    minEnd = Math.min(suspSedimentRange.getEnd(),
                        lastOtherRange.getEnd());
                    double total = lastOtherValue + fraction.getSuspSediment();
                    Range totalRange = new Range(maxStart, minEnd);
                    if (suspSedimentRange.getEnd() > lastOtherRange.getEnd()) {
                        lastSuspRange = (Range) suspSedimentRange.clone();
                        lastSuspRange.setStart(lastOtherRange.getEnd());
                        lastSuspValue = fraction.getSuspSediment();
                        lastOtherRange = null;
                    }
                    else {
                        // Other is "longer".
                        lastOtherRange.setStart(suspSedimentRange.getEnd());
                        lastSuspRange = null;
                    }
                    if (lastOtherRange != null
                    && Math.abs(suspSedimentRange.getEnd() - lastOtherRange.getEnd()) < 0.1d) {
                        lastOtherRange = null;
                        lastSuspRange = null;
                    }
                    fairLoad.setTotal(km, total, totalRange);
                }
                else {
                    lastSuspRange = suspSedimentRange;
                    lastSuspValue = fraction.getSuspSediment();
                    lastOtherRange = null;
                }
            }
            else if (fraction.hasButSuspValues() && lastSuspRange != null) {
                // If intersects with last suspsed range, merge and add
                double total = fraction.getCoarse() +
                    fraction.getFineMiddle() +
                    fraction.getSand() +
                    fraction.getSuspSand() +
                    lastSuspValue;
                double maxStart = Math.max(fraction.getCoarseRange().getStart(),
                    lastSuspRange.getStart());
                if (lastSuspRange.contains(km)) {
                    double minEnd = Math.min(fraction.getCoarseRange().getEnd(),
                        lastSuspRange.getEnd());
                    Range totalRange = new Range(maxStart, minEnd);
                    if (lastSuspRange.getEnd() > fraction.getCoarseRange().getEnd()) {
                        // SuspSed longer.
                        lastSuspRange.setStart(fraction.getCoarseRange().getEnd());
                        lastOtherRange = null;
                    }
                    else {
                        // Other longer
                        lastOtherRange = (Range) fraction.getCoarseRange().clone();
                        lastOtherRange.setStart(lastSuspRange.getEnd());
                        lastSuspRange = null;
                        lastOtherValue = total - lastSuspValue;
                    }
                    if (lastSuspRange != null
                    &&  lastOtherRange != null
                    && Math.abs(lastSuspRange.getEnd() - lastOtherRange.getEnd()) < 0.1d) {
                        lastOtherRange = null;
                        lastSuspRange = null;
                    }
                    fairLoad.setTotal(km, total, totalRange);
                }
                else {
                    // Ranges are disjoint.
                    lastOtherRange = fraction.getCoarseRange();
                    lastOtherValue = total - fraction.getSuspSediment();
                    lastSuspRange = null;
                }
            }
            else {
                // Some values are missing or no intersection with former values.
                // Stay as we are.
                if (fraction.hasButSuspValues()) {
                    double total = fraction.getCoarse() +
                        fraction.getFineMiddle() +
                        fraction.getSand() +
                        fraction.getSuspSand();
                    lastOtherRange = fraction.getCoarseRange();
                    lastOtherValue = total;
                    lastSuspRange = null;
                }
                else if (fraction.hasOnlySuspValues()) {
                    lastSuspRange = fraction.getSuspSedimentRange();
                    lastSuspValue = fraction.getSuspSediment();
                    lastOtherRange = null;
                }
            }
        }
        return fairLoad;
    }


    /**
     * Transform values in load.
     * Background is to transform values measured in
     * t/a to m^3/a using the specific measured densities.
     *
     * @param load The load of which values should be transformed.
     * @param year The year at which to look at density (e.g. 2003).
     *
     * @return parameter load with transformed values.
     */
    private SedimentLoadLSData calculateUnit(SedimentLoadLSData load, int year) {
        SedimentDensity density =
            SedimentDensityFactory.getSedimentDensity(river, kmLow, kmUp);

        for (double km: load.getKms()) {
            double dens = 1d/density.getDensity(km, year);
            SedimentLoadFraction fraction = load.getFraction(km);
            double coarse = fraction.getCoarse();
            double fineMiddle = fraction.getFineMiddle();
            double sand = fraction.getSand();
            double suspSand = fraction.getSuspSand();
            double bedSand = fraction.getSuspSandBed();
            double sediment = fraction.getSuspSediment();
            double total = fraction.getTotal();
            double loadTotal = fraction.getLoadTotal();
            load.setCoarse(km, (coarse * dens), fraction.getCoarseRange());
            load.setFineMiddle(km, (fineMiddle * dens), fraction.getFineMiddleRange());
            load.setSand(km, (sand * dens), fraction.getSandRange());
            load.setSuspSand(km, (suspSand * dens), fraction.getSuspSandRange());
            load.setSuspSandBed(km, (bedSand * dens), fraction.getSuspSandBedRange());
            load.setSuspSediment(km, (sediment * dens), fraction.getSuspSedimentRange());
            load.setTotal(km, (total * dens), fraction.getTotalRange());
            load.setLoadTotal(km, (loadTotal * dens), fraction.getLoadTotalRange());
        }
        return load;
    }
}
// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf-8 :