# HG changeset patch
# User Sascha L. Teichmann <sascha.teichmann@intevation.de>
# Date 1261418454 0
# Node ID 21fbd254db71afc4b4ace7422a41e78b9324116e
# Parent  2402173a149086a7048c21668d9d42297f956fce
Added support for converting 2D rasters into polygons.


gnv-artifacts/trunk@472 c6561f87-3c4e-4783-a992-168aeb5c3f6f

diff -r 2402173a1490 -r 21fbd254db71 gnv-artifacts/ChangeLog
--- a/gnv-artifacts/ChangeLog	Mon Dec 21 17:03:10 2009 +0000
+++ b/gnv-artifacts/ChangeLog	Mon Dec 21 18:00:54 2009 +0000
@@ -1,3 +1,18 @@
+2009-12-21	Sascha L. Teichmann	<sascha.teichmann@intevation.de>
+
+	* src/main/java/de/intevation/gnv/raster/Raster.java: New. Models 2D double 
+	  rasters. Has some support for filtering with gauss kernels, building
+	  iso classes, etc.
+
+	* src/main/java/de/intevation/gnv/raster/Palette.java: New. Maps double
+	  values to integer indices and colors.
+
+	* src/main/java/de/intevation/gnv/raster/Vectorizer.java: New. Simple
+	  vectorizer which traces regions in integer rasters.
+
+	* pom.xml: Added dependency to GNU Trove 2.1.1 which is needed by the
+	  vectorizer.
+
 2009-12-21  Ingo Weinzierl <ingo.weinzierl@intevation.de>
 
 	* src/main/java/de/intevation/gnv/utils/WKTUtils.java,
diff -r 2402173a1490 -r 21fbd254db71 gnv-artifacts/pom.xml
--- a/gnv-artifacts/pom.xml	Mon Dec 21 17:03:10 2009 +0000
+++ b/gnv-artifacts/pom.xml	Mon Dec 21 18:00:54 2009 +0000
@@ -118,6 +118,11 @@
       <artifactId>jts</artifactId>
       <version>1.9</version>
     </dependency>
+    <dependency>
+        <groupId>trove</groupId>
+        <artifactId>trove</artifactId>
+        <version>2.1.1</version>
+    </dependency>
   </dependencies>
   <repositories>
     <repository>
diff -r 2402173a1490 -r 21fbd254db71 gnv-artifacts/src/main/java/de/intevation/gnv/raster/Palette.java
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/gnv-artifacts/src/main/java/de/intevation/gnv/raster/Palette.java	Mon Dec 21 18:00:54 2009 +0000
@@ -0,0 +1,188 @@
+package de.intevation.gnv.raster;
+
+import org.w3c.dom.Document;
+import org.w3c.dom.Element;
+import org.w3c.dom.NodeList;
+
+import java.util.Arrays;
+
+import java.awt.Color;
+
+import de.intevation.gnv.raster.Raster.ValueToIndex;
+
+/**
+ * @author Sascha L. Teichmann (sascha.teichmann@intevation.de)
+ */
+public class Palette
+implements   ValueToIndex
+{
+    public static final class Entry 
+    implements                Comparable
+    {
+        private Entry  left;
+        private Entry  right;
+
+        private int    index;
+
+        private double from;
+        private double to;
+
+        private String description;
+
+        private Color  color;
+
+        public Entry() {
+        }
+
+        public Entry(
+            int    index, 
+            double from, 
+            double to, 
+            Color  color,
+            String description
+        ) {
+            this.index       = index;
+            this.from        = from;
+            this.to          = to;
+            this.color       = color;
+            this.description = description;
+        }
+
+        public int compareTo(Object other) {
+            Entry e = (Entry)other;
+            if (from < e.from) return -1;
+            if (from > e.from) return +1;
+            return 0;
+        }
+
+        public Entry locateEntry(double value) {
+            Entry current = this;
+            while (current != null) {
+                boolean b = value >= current.from;
+                if (b && value <= current.to) {
+                    return current;
+                }
+                current = b
+                    ? current.right
+                    : current.left;
+            }
+            return null;
+        }
+
+        public int locate(double value) {
+            Entry entry = locateEntry(value);
+            return entry != null
+                ? entry.index
+                : -1;
+        }
+
+        public Entry findByIndex(int index) {
+            Entry current = this;
+            while (current != null) {
+                if (current.index == index) {
+                    break;
+                }
+                current = index < current.index
+                    ? current.left
+                    : current.right;
+            }
+            return current;
+        }
+
+        public String getDescription() {
+            return description;
+        }
+    } // class Entry
+
+    protected Entry [] entries;
+    protected Entry    lookup;
+    protected Color [] rgbs;
+
+    private Entry buildLookup(Entry [] entries, int lo, int hi) {
+        if (lo > hi) {
+            return null;
+        }
+        int mid = (lo + hi)/2;
+        Entry entry = entries[mid];
+        entry.left  = buildLookup(entries, lo, mid-1);
+        entry.right = buildLookup(entries, mid+1, hi);
+        return entry;
+    }
+
+    public Palette() {
+    }
+
+    public Palette(Document document) {
+
+        NodeList rangeNodes = document.getElementsByTagName("range");
+
+        entries = new Entry[rangeNodes.getLength()];
+        rgbs    = new Color[entries.length];
+
+        for (int i = 0; i < entries.length; ++i) {
+            Element rangeElement = (Element)rangeNodes.item(i);
+            double from  = doubleValue(rangeElement.getAttribute("from"));
+            double to    = doubleValue(rangeElement.getAttribute("to"));
+            Color  color = color(rangeElement.getAttribute("rgb"));
+            String description = rangeElement.getAttribute("description");
+            if (from > to) { double t = from; from = to; to = t; }
+            entries[i] = new Entry(i, from, to, color, description);
+            rgbs[i] = color;
+        }
+
+        buildLookup();
+    }
+
+    private static final double doubleValue(String s) {
+        if (s == null || (s = s.trim()).length() == 0) {
+            return 0d;
+        }
+        if ((s = s.toLowerCase()).startsWith("-inf")) {
+            return -Double.MAX_VALUE; // XXX: Not using Double.NEGATIVE_INFINITY!
+        }
+
+        if (s.startsWith("inf")) {
+            return Double.MAX_VALUE; // XXX: Not using Double.NEGATIVE_INFINITY!
+        }
+
+        return Double.parseDouble(s);
+    }
+
+    private static final Color color(String s) {
+        if (s == null || (s = s.trim()).length() == 0) {
+            return null;
+        }
+        return Color.decode(s);
+    }
+
+
+    protected void buildLookup() {
+        Arrays.sort(entries);
+        lookup = buildLookup(entries, 0, entries.length-1);
+    }
+
+    public int getSize() {
+        return rgbs.length;
+    }
+
+    public Color getColor(int index) {
+        return rgbs[index];
+    }
+
+    public int indexToRGB(int index) {
+        return rgbs[index].getRGB();
+    }
+
+    public int toIndex(double value) {
+        return lookup.locate(value);
+    }
+
+    public Entry getEntry(double value) {
+        return lookup.locateEntry(value);
+    }
+
+    public Entry getEntryByIndex(int index) {
+        return lookup.findByIndex(index);
+    }
+}
+// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8:
diff -r 2402173a1490 -r 21fbd254db71 gnv-artifacts/src/main/java/de/intevation/gnv/raster/Raster.java
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/gnv-artifacts/src/main/java/de/intevation/gnv/raster/Raster.java	Mon Dec 21 18:00:54 2009 +0000
@@ -0,0 +1,233 @@
+package de.intevation.gnv.raster;
+
+/**
+ * @author Sascha L. Teichmann (sascha.teichmann@intevation.de)
+ */
+public class Raster
+{
+    public static class Kernel {
+
+        public static final class Coeff {
+            public int    di;
+            public int    dj;
+            public double c;
+
+            public Coeff() {
+            }
+
+            public Coeff(int di, int dj, double c) {
+                this.di = di;
+                this.dj = dj;
+                this.c  = c;
+            }
+
+            public String toString() {
+                return String.valueOf(c);
+            }
+
+        } // class Coeff
+
+        protected Coeff [] coeffs;
+        protected int      width;
+
+        public Kernel() {
+        }
+
+        public Kernel(Coeff [] coeffs, int width) {
+            this.coeffs = coeffs;
+            this.width  = width;
+        }
+
+        public static double gauss(double x, double y, double sigma) {
+            double s2 = sigma*sigma;
+            return 1.0d/(2.0d*Math.PI*s2)*Math.exp(-(x*x + y*y)/(2.0d*s2));
+        }
+
+        public static Kernel createGauss() {
+            return createGauss(1.0d, 5);
+        }
+
+        public static Kernel createGauss(double sigma, int width) {
+            Coeff [] coeffs = new Coeff[width*width];
+            double sum = 0d;
+            for (int j = 0; j < width; ++j) {
+                int y = -width/2 + j;
+                for (int i = 0; i < width; ++i) {
+                    int x = -width/2 + i;
+                    double c = gauss(x, y, sigma);
+                    coeffs[j*width + i] = new Coeff(x, y, c);
+                    sum += c;
+                }
+            }
+            sum = 1.0/sum;
+            for (int i = 0; i < coeffs.length; ++i) {
+                coeffs[i].c *= sum;
+            }
+            return new Kernel(coeffs, width);
+        }
+
+        public double fold(Access access, int i, int j) {
+
+            double s = 0.0d;
+
+            double unused = 0d;
+
+            for (int k = 0; k < coeffs.length; ++k) {
+                Coeff coeff = coeffs[k];
+                double v = access.get(i + coeff.di, j + coeff.dj);
+                if (Double.isNaN(v)) {
+                    unused += coeff.c;
+                }
+                else {
+                    s += v * coeff.c;
+                }
+            }
+
+            return s*(1.0d - unused);
+        }
+
+        public String toString() {
+            StringBuilder sb = new StringBuilder();
+
+            int height = coeffs.length/width;
+
+            for (int j = 0; j < height; ++j) {
+                if (j > 0) {
+                    sb.append(System.getProperty("line.separator"));
+                }
+                for (int i = 0; i < width; ++i) {
+                    if (i > 0) {
+                        sb.append("\t");
+                    }
+                    sb.append(coeffs[j*width + i]);
+                }
+            }
+
+            return sb.toString();
+        }
+    } // class Kernel
+
+    public interface Access {
+
+        double get(int i, int j);
+
+    } // interface Access
+
+    public interface ValueToIndex {
+
+        int toIndex(double value);
+
+    } // interface ValueToIndex;
+
+    public static class IsoClasses
+    implements          ValueToIndex
+    {
+        protected double m;
+        protected double b;
+        protected int    numClasses;
+
+        public IsoClasses(Raster raster, int numClasses) {
+
+            this.numClasses = numClasses;
+
+            double min =  Double.MAX_VALUE;
+            double max = -Double.MAX_VALUE;
+
+            double [] src = raster.raster;
+
+            for (int i = 0; i < src.length; ++i) {
+                double x = src[i];
+                if (!Double.isNaN(x)) {
+                    if (x < min) min = x;
+                    if (x > max) max = x;
+                }
+            }
+
+            /* f(min) = 0, f(max) = numClasses - 1
+
+               I:  0              = m*min + b
+               II: numClasses - 1 = m*max + b
+
+               II - I:
+               numClasses - 1 = m*(max - min)
+
+               m = (numClasses - 1)/(max - min)
+               b = m*min
+            */
+
+            if (max == min) {
+                m = 0;
+                b = 0;
+            }
+            else {
+                m = (numClasses - 1)/(max - min);
+                b = m*min;
+            }
+        }
+
+        public int toIndex(double value) {
+            return Double.isNaN(value)
+                ? -1
+                : Math.min(Math.max(0, (int)Math.round(m*value + b)), numClasses-1);
+        }
+    } // class IsoClasses
+
+    protected double [] raster;
+    protected int       width;
+
+    public Raster() {
+        this(new double[0], 0);
+    }
+
+    public Raster(double [] raster, int width) {
+        this.raster = raster;
+        this.width  = width;
+    }
+
+    public int getWidth() {
+        return width;
+    }
+
+    public int getHeight() {
+        return raster.length / width;
+    }
+
+    public int [] toIndexed(ValueToIndex valueToIndex) {
+        int [] dst = new int[raster.length];
+        for (int i = 0; i < raster.length; ++i) {
+            dst[i] = valueToIndex.toIndex(raster[i]);
+        }
+        return dst;
+    }
+
+    public Raster create(Kernel kernel) {
+        double [] dst = new double[raster.length];
+        Raster r = new Raster(dst, width);
+        r.apply(kernel, continueBorder());
+        return r;
+    }
+
+    public void apply(Kernel kernel, Access access) {
+        int height = getHeight();
+        for (int j = 0; j < height; ++j) {
+            int row = j*width;
+            for (int i = 0; i < width; ++i) {
+                raster[row + i] = kernel.fold(access, i, j);
+            }
+        }
+    }
+
+    public Access continueBorder() {
+        return new Access() {
+            int height = getHeight()-1;
+            public double get(int i, int j) {
+                     if (i < 0)      i = 0;
+                else if (i >= width) i = width-1;
+                     if (j < 0)      j = 0;
+                else if (j > height) j = height;
+                return raster[j*width + i];
+            }
+        };
+    }
+}
+// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8:
diff -r 2402173a1490 -r 21fbd254db71 gnv-artifacts/src/main/java/de/intevation/gnv/raster/Vectorizer.java
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/gnv-artifacts/src/main/java/de/intevation/gnv/raster/Vectorizer.java	Mon Dec 21 18:00:54 2009 +0000
@@ -0,0 +1,289 @@
+package de.intevation.gnv.raster;
+
+import java.util.BitSet;
+import java.util.List;
+import java.util.ArrayList;
+
+import gnu.trove.TIntStack;
+import gnu.trove.TIntObjectHashMap;
+
+/**
+ * @author Sascha L. Teichmann (sascha.teichmann@intevation.de)
+ */
+public class Vectorizer
+{
+    public interface RingsHandler {
+
+        void handleRings(List rings, int value, int width);
+
+    } // interface RingsHandler
+
+    public static final class Edge {
+
+        public Edge prev;
+        public Edge next;
+
+        public int a;
+        public int b;
+
+        public Edge() {
+        }
+
+        public Edge(int a, int b) {
+            this.a = a;
+            this.b = b;
+        }
+
+        public void chain(Edge other, int found) {
+
+            if (found == a) {
+                other.next = this;
+                prev = other;
+                return;
+            }
+
+            if (found == b) {
+                next = other;
+                other.prev = this;
+                return;
+            }
+
+            throw new IllegalStateException("cannot chain");
+        }
+
+        public boolean isComplete() {
+            Edge current = this;
+            do {
+                if (current.prev == null || current.next == null) {
+                    return false;
+                }
+                current = current.next;
+            }
+            while (current != this);
+            return true;
+        }
+
+        public int length() {
+            int length = 0;
+            Edge current = this;
+            do { ++length; } while ((current = current.next) != this);
+            return length;
+        }
+    } // class Edge
+
+    protected static Edge simplify(Edge edge, int width) {
+
+        Edge e1 = edge;
+        Edge e2 = edge.next;
+        Edge e3 = e2.next;
+
+        int length = edge.length();
+
+        int count = 0;
+
+        do {
+            int e1x = e1.a % width;
+            int e1y = e1.a / width;
+            int e2x = e2.a % width;
+            int e2y = e2.a / width;
+            int e3x = e3.a % width;
+            int e3y = e3.a / width;
+
+            if ((e1x == e2x && e2x == e3x && e1x == e3x)
+            ||  (e1y == e2y && e2y == e3y && e1y == e3y)) {
+                e1.next = e3;
+                e3.prev = e1;
+                e2 = e3;
+                e3 = e2.next;
+                count = 0;
+                --length;
+            }
+            else {
+                e1 = e1.next;
+                e2 = e2.next;
+                e3 = e3.next;
+                ++count;
+            }
+        }
+        while (count < length + 2);
+
+        return e1;
+    }
+
+    protected int []            raster;
+    protected int               width;
+    protected TIntObjectHashMap openEdges;
+    protected ArrayList         rings;
+
+    public Vectorizer() {
+        openEdges = new TIntObjectHashMap();
+        rings     = new ArrayList();
+    }
+
+    public Vectorizer(int [] raster, int width) {
+        this();
+        this.raster = raster;
+        this.width  = width;
+    }
+
+    public static final int tl(int i, int w) {
+        int x = i % w;
+        int y = i / w;
+        return x + (w + 1)*y;
+    }
+
+    public static final int tr(int i, int w) {
+        return tl(i, w) + 1;
+    }
+
+    public static final int bl(int i, int w) {
+        return tl(i, w) + w + 1;
+    }
+
+    public static final int br(int i, int w) {
+        return bl(i, w) + 1;
+    }
+
+    protected void resetRegion() {
+        openEdges.clear();
+        rings.clear();
+    }
+
+    protected void emit(Edge edge) {
+
+        Edge otherA = (Edge)openEdges.remove(edge.a);
+
+        if (otherA != null) {
+            otherA.chain(edge, edge.a);
+        }
+
+        Edge otherB = (Edge)openEdges.remove(edge.b);
+        if (otherB != null) {
+            otherB.chain(edge, edge.b);
+        }
+
+        if (edge.isComplete()) {
+            rings.add(simplify(edge, width + 1));
+        }
+        else {
+            if (otherA == null) {
+                openEdges.put(edge.a, edge);
+            }
+            if (otherB == null) {
+                openEdges.put(edge.b, edge);
+            }
+        }
+    }
+
+    public int process(RingsHandler handler) {
+
+        BitSet visited = new BitSet(raster.length);
+
+        TIntStack stack = new TIntStack();
+
+        int regions = 0;
+
+
+        for (int i = 0; i < raster.length; ++i) {
+            if (visited.get(i)) {
+                continue;
+            }
+
+            ++regions;
+
+            int currentValue = raster[i];
+            visited.set(i);
+
+            int current = i;
+            visited.set(current);
+
+            for (;;) {
+                int tl = tl(current, width);
+                int tr = tr(current, width);
+                int bl = bl(current, width);
+                int br = br(current, width);
+
+                int t = current - width;
+
+                if (t < 0) {
+                    emit(new Edge(tr, tl));
+                }
+                else {
+                    if (raster[t] != currentValue) {
+                        emit(new Edge(tr, tl));
+                    }
+                    else {
+                        if (!visited.get(t)) {
+                            visited.set(t);
+                            stack.push(t);
+                        }
+                    }
+                }
+
+                int b = current + width;
+
+                if (b >= raster.length) {
+                    emit(new Edge(bl, br));
+                }
+                else {
+                    if (raster[b] != currentValue) {
+                        emit(new Edge(bl, br));
+                    }
+                    else {
+                        if (!visited.get(b)) {
+                            visited.set(b);
+                            stack.push(b);
+                        }
+                    }
+                }
+
+                int x = current % width;
+
+                if (x == 0) {
+                    emit(new Edge(tl, bl));
+                }
+                else {
+                    int l = current - 1;
+                    if (raster[l] != currentValue) {
+                        emit(new Edge(tl, bl));
+                    }
+                    else {
+                        if (!visited.get(l)) {
+                            visited.set(l);
+                            stack.push(l);
+                        }
+                    }
+                }
+
+                if (x == width - 1) {
+                    emit(new Edge(br, tr));
+                }
+                else {
+                    int r = current + 1;
+                    if (raster[r] != currentValue) {
+                        emit(new Edge(br, tr));
+                    }
+                    else {
+                        if (!visited.get(r)) {
+                            visited.set(r);
+                            stack.push(r);
+                        }
+                    }
+                }
+
+                if (stack.size() == 0) {
+                    break;
+                }
+
+                current = stack.pop();
+            }
+
+            handler.handleRings(rings, currentValue, width + 1);
+
+            resetRegion();
+        }
+
+        return regions;
+    }
+}
+// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8: