--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/gnv-artifacts/src/main/java/de/intevation/gnv/math/XYColumn.java	Fri Sep 28 12:14:00 2012 +0200
@@ -0,0 +1,177 @@
+/*
+ * Copyright (c) 2010 by Intevation GmbH
+ *
+ * This program is free software under the LGPL (>=v2.1)
+ * Read the file LGPL.txt coming with the software for details
+ * or visit http://www.gnu.org/licenses/ if it does not exist.
+ */
+
+package de.intevation.gnv.math;
+
+import java.util.ArrayList;
+import java.util.Collections;
+import java.util.List;
+
+import org.apache.commons.math.FunctionEvaluationException;
+import org.apache.commons.math.MathException;
+
+import org.apache.commons.math.analysis.UnivariateRealFunction;
+
+import org.apache.commons.math.analysis.interpolation.SplineInterpolator;
+import org.apache.commons.math.analysis.interpolation.UnivariateRealInterpolator;
+
+import org.apache.log4j.Logger;
+
+/**
+ * A column of discrete attributed height values located at a point.
+ * Values between the discrete height values are interpolated.
+ *
+ * @author <a href="mailto:ingo.weinzierl@intevation.de">Ingo Weinzierl</a>
+ * @author <a href="mailto:sascha.teichmann@intevation.de">Sascha L. Teichmann</a>
+ */
+public class XYColumn
+extends      Point2d
+implements   UnivariateRealFunction
+{
+    private static Logger log = Logger.getLogger(XYColumn.class);
+
+    /**
+     * The list of discrete height values.
+     */
+    protected List<HeightValue> values;
+
+    /**
+     * The curve used to interpolate the points between the
+     * discrete height values.
+     */
+    protected transient UnivariateRealFunction curve;
+
+    /**
+     * Default constructor.
+     */
+    public XYColumn() {
+        values = new ArrayList<HeightValue>();
+    }
+
+    /**
+     * Constructor to create an XYColumn with a given (x, y) coordinate
+     * and an (i, j) index tuple.
+     * @param x The x coordinate.
+     * @param y The y coordinate.
+     * @param i The i component.
+     * @param j The j component.
+     */
+    public XYColumn(double x, double y, int i, int j) {
+        super(x, y, i, j);
+        values = new ArrayList<HeightValue>();
+    }
+
+    /**
+     * Adds a given height value to the list of height values.
+     * @param value The height value.
+     */
+    public void add(HeightValue value) {
+        values.add(value);
+    }
+
+    /**
+     * Returns the list of height values.
+     * @return The list of height values.
+     */
+    public List<HeightValue> getValues() {
+        return values;
+    }
+
+    public double value(double depth) {
+        try {
+            if (curve != null) {
+                HeightValue h = values.get(0);
+                // extrapolate beyond boundaries by repeating
+                if (depth > h.z) return h.v;
+                h = values.get(values.size()-1);
+                if (depth < h.z) return h.v;
+                return curve.value(depth);
+            }
+        }
+        catch (FunctionEvaluationException fee) {
+            log.error("evaluation failed", fee);
+        }
+
+        return Double.NaN;
+    }
+
+    /**
+     * Prepares this XYColumn to be queried. A given XYDepth
+     * object is used to fill the values to a certain depth.
+     * @param xyDepth To figure out the depth a the cordinate.
+     * @return true if preparation succeeds else false.
+     */
+    public boolean prepare(XYDepth xyDepth) {
+        if (curve == null) {
+            int N = values.size();
+            if (N == 0) {
+                log.error("no points for interpolation");
+                return false;
+            }
+
+            if (N == 1) {
+                // only one value -> constant function
+                curve = new ConstantFunction(values.get(0).v);
+            }
+            else { // more than on value
+                double depth = xyDepth.depth(this);
+                Collections.sort(values, HeightValue.INV_Z_COMPARATOR);
+
+                // if there is no value at 0 repeat first value
+                HeightValue first = values.get(0);
+                if (first.z < 0d) {
+                    values.add(0, new HeightValue(0d, first.v, first.k-1));
+                    ++N;
+                }
+
+                // if there is no value at depth repeat last value
+                HeightValue last = values.get(N-1);
+                if (last.z > depth) {
+                    values.add(new HeightValue(depth, last.v, last.k+1));
+                    ++N;
+                }
+                if (N < 3) { // interpolate linear
+                    first = values.get(0);
+                    last  = values.get(N-1);
+                    curve = new LinearFunction.Univariate(
+                        first.z, first.v,
+                        last.z,  last.v);
+                }
+                else { // higher degree interpolation
+                    double [] z = new double[N];
+                    double [] v = new double[N];
+                    for (int i = 0; i < N; ++i) {
+                        HeightValue h = values.get(N-1-i);
+                        z[i] = h.z;
+                        v[i] = h.v;
+                    }
+                    try {
+                        curve = getInterpolator().interpolate(z, v);
+                    }
+                    catch (MathException me) {
+                        log.error("interpolation failed", me);
+                        return false;
+                    }
+                }
+            }
+        }
+        return true;
+    }
+
+    /**
+     * Returns the interpolator used to interpolate the values between
+     * the discrete height values. This class returns an instance of
+     * {@link org.apache.commons.math.analysis.interpolation.SplineInterpolator}.
+     * Override this if you want to use another kind of interpolation.
+     * @return The interpolator to be used.
+     */
+    protected UnivariateRealInterpolator getInterpolator() {
+        return new SplineInterpolator();
+    }
+}
+// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :