[PATCH] Sediment load: refactored epoch calculation to not average over all but per year

[Wald Commits]

# HG changeset patch
# User Sascha L. Teichmann <teichmann at intevation.de>
# Date 1406819081 -7200
# Node ID 6d24ba2ac96409cc2eb0310dc2141bf6d1abd7b9
# Parent  fe5ef780f8b1934175bc9bb6a7973c770d71c644
Sediment load: refactored epoch calculation to not average over all but per year.

diff -r fe5ef780f8b1 -r 6d24ba2ac964 artifacts/src/main/java/org/dive4elements/river/artifacts/model/minfo/SedimentLoadDataCalculation.java
--- a/artifacts/src/main/java/org/dive4elements/river/artifacts/model/minfo/SedimentLoadDataCalculation.java	Wed Jul 30 19:26:20 2014 +0200
+++ b/artifacts/src/main/java/org/dive4elements/river/artifacts/model/minfo/SedimentLoadDataCalculation.java	Thu Jul 31 17:04:41 2014 +0200
@@ -7,8 +7,11 @@
  */
 package org.dive4elements.river.artifacts.model.minfo;
 
+import java.util.ArrayList;
 import java.util.List;
+import java.util.Map;
 import java.util.Set;
+import java.util.TreeMap;
 import java.util.TreeSet;
 
 import org.dive4elements.river.artifacts.access.SedimentLoadAccess;
@@ -154,17 +157,6 @@
         }
     } // class Sum
 
-    public static final class Avg extends Sum {
-        public Avg() {
-        }
-
-        @Override
-        public double getSum() {
-            return n == 0 ? 0.0 : sum/(double)n;
-        }
-    } // class Sum
-
-
     private String   river;
     private String   yearEpoch;
     private String   unit;
@@ -299,32 +291,46 @@
         boolean isKmUp = isKmUp();
         Set<Integer> missingFractions = new TreeSet<Integer>();
 
+        SedimentDensity sd = getSedimentDensity();
+
         // They are not epochs, they are single years!
         Not notEpochs = new Not(IsEpoch.INSTANCE);
 
         for (int [] epoch: epochs) {
-            Value.Filter filter = new And()
-                .add(notEpochs)
-                .add(new TimeRangeIntersects(epoch[0], epoch[1]));
+            List<double [][]> results = new ArrayList<double [][]>();
 
-            Avg avg = new Avg();
+            int min = Math.min(epoch[0], epoch[1]);
+            int max = Math.max(epoch[0], epoch[1]);
 
-            for (GrainFraction gf: GRAIN_FRACTIONS) {
-                double [][] result = sum(
-                    sld, gf.getGrainFractions(), filter, avg, isKmUp,
-                    missingFractions);
+            for (int year = min; year <= max; ++year) {
+                Value.Filter filter = new And()
+                    .add(notEpochs)
+                    .add(new TimeRangeIntersects(year));
 
-                if (result[0].length == 0 || DoubleUtil.isNaN(result[1])) {
-                    // TODO: resolve i18n
-                    addProblem("minfo.sediment.load.no.fractions",
-                        gf.getDescription());
-                    continue;
+                Sum sum = new Sum();
+
+                for (GrainFraction gf: GRAIN_FRACTIONS) {
+                    double [][] result = sum(
+                        sld, gf.getGrainFractions(), filter, sum, isKmUp,
+                        missingFractions);
+
+                    if (result[0].length == 0 || DoubleUtil.isNaN(result[1])) {
+                        // TODO: resolve i18n
+                        addProblem("minfo.sediment.load.no.fractions",
+                            gf.getDescription());
+                        continue;
+                    }
+
+                    transformT2M3(sd, year, result);
+                    results.add(result);
                 }
-                // TODO: Optionally transform units.
-                SedimentLoadDataResult.Fraction sldrf =
-                    new SedimentLoadDataResult.Fraction(gf.getDescription(), result);
-                sldr.addFraction(sldrf);
             }
+
+            double [][] result = average(results);
+            // TODO: Optionally transform units.
+            SedimentLoadDataResult.Fraction sldrf =
+                new SedimentLoadDataResult.Fraction("TODO: nice description", result);
+            sldr.addFraction(sldrf);
         }
         // TODO: Generate messages for missing fractions.
         return new CalculationResult(sldr, this);
@@ -452,5 +458,52 @@
 
         return result;
     }
+
+    private static final class XSum {
+        private double sum;
+        private int n;
+        public XSum() {
+        }
+        public void add(double v) {
+            sum += v;
+            ++n;
+        }
+        public double avg() {
+            return sum/n;
+        }
+    }
+
+    private static double [][] average(List<double [][]> data) {
+
+        TreeMap<Double, XSum> map = new TreeMap<Double, XSum>();
+
+        for (double [][] pair: data) {
+            double [] kms = pair[0];
+            double [] vs = pair[1];
+            for (int i = 0; i < kms.length; ++i) {
+                double km = kms[i];
+                double v = vs[i];
+                if (Double.isNaN(km) || Double.isNaN(v)) {
+                    continue;
+                }
+                XSum xsum = map.get(km);
+                if (xsum == null) {
+                    map.put(km, xsum = new XSum());
+                }
+                xsum.add(v);
+            }
+        }
+
+        double [][] result = new double[2][map.size()];
+        int i = 0;
+        for (Map.Entry<Double, XSum> entry: map.entrySet()) {
+            result[i][0] = entry.getKey();
+            result[i][1] = entry.getValue().avg();
+            ++i;
+        }
+
+        return null;
+    }
+
 }
 // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf-8 :