ingo@1115: /*
ingo@1115:  * Copyright (c) 2010 by Intevation GmbH
ingo@1115:  *
ingo@1115:  * This program is free software under the LGPL (>=v2.1)
ingo@1115:  * Read the file LGPL.txt coming with the software for details
ingo@1115:  * or visit http://www.gnu.org/licenses/ if it does not exist.
ingo@1115:  */
ingo@1115: 
sascha@424: package de.intevation.gnv.raster;
sascha@424: 
sascha@779: import gnu.trove.TIntObjectHashMap;
sascha@779: import gnu.trove.TIntStack;
sascha@779: 
sascha@779: import java.util.ArrayList;
sascha@424: import java.util.BitSet;
sascha@424: import java.util.List;
sascha@424: 
sascha@445: import org.apache.log4j.Logger;
sascha@445: 
sascha@424: /**
sascha@801:  * Instances of this class are able to vectorize 2D integer arrays
sascha@801:  * with a kind of flood filling regions detection mechanism to
sascha@801:  * a set of line strings and polygons.
sascha@801:  *
sascha@780:  * @author <a href="mailto:sascha.teichmann@intevation.de">Sascha L. Teichmann</a>
sascha@424:  */
sascha@424: public class Vectorizer
sascha@424: {
sascha@445:     private static Logger log = Logger.getLogger(Vectorizer.class);
sascha@445: 
sascha@801:     /**
sascha@801:      * Callback to process the found line strings and polygons.
sascha@801:      */
sascha@424:     public interface RingsHandler {
sascha@424: 
sascha@801:         /**
sascha@801:          * Called from {@link #process(de.intevation.gnv.raster.Vectorizer.RingsHandler) }
sascha@801:          * to give the found features to the postprocessing backend.
sascha@801:          * @param rings The found line strings and polygons.
sascha@801:          * @param value The integer value of the raster surrounded by
sascha@801:          * the features.
sascha@801:          * @param width The width of the index space of the vectorized data.
sascha@801:          * @param height The height of the index space of the vectorized data.
sascha@801:          */
sascha@436:         void handleRings(
sascha@436:             List<Edge> rings,
sascha@778:             int        value,
sascha@436:             int        width,
sascha@436:             int        height);
sascha@424: 
sascha@424:     } // interface RingsHandler
sascha@424: 
sascha@801:     /**
sascha@801:      * Doubly link list representing a line string or a polygon.
sascha@801:      */
sascha@424:     public static final class Edge {
sascha@424: 
sascha@801:         /**
sascha@801:          * The predecessor.
sascha@801:          */
sascha@424:         public Edge prev;
sascha@801:         /**
sascha@801:          * The successor.
sascha@801:          */
sascha@424:         public Edge next;
sascha@424: 
sascha@801:         /**
sascha@801:          * Index of first vertex. To separate x and y values
sascha@801:          * you have to divide the value by the width of the
sascha@801:          * index space.
sascha@801:          */
sascha@424:         public int a;
sascha@801:         /**
sascha@801:          * Index of second vertex. To separate x and y values
sascha@801:          * you have to divide the value by the width of the
sascha@801:          * index space.
sascha@801:          */
sascha@424:         public int b;
sascha@424: 
sascha@801:         /**
sascha@801:          * Default constructor.
sascha@801:          */
sascha@424:         public Edge() {
sascha@424:         }
sascha@424: 
sascha@801:         /**
sascha@801:          * Constructor to create Edge with two vertices.
sascha@801:          * @param a The first vertex.
sascha@801:          * @param b The second vertex.
sascha@801:          */
sascha@424:         public Edge(int a, int b) {
sascha@424:             this.a = a;
sascha@424:             this.b = b;
sascha@424:         }
sascha@424: 
sascha@801:         /**
sascha@801:          * Copy constructor
sascha@801:          * @param other The edge to clone.
sascha@801:          */
sascha@437:         public Edge(Edge other) {
sascha@437:             a = other.a;
sascha@437:             b = other.b;
sascha@437:         }
sascha@437: 
sascha@801:         /**
sascha@801:          * Chain a given edge segment to this segment. The
sascha@801:          * side of the chaining is determined by a given
sascha@801:          * parameter.
sascha@801:          * @param other The segment to chain.
sascha@801:          * @param found The side to chain.
sascha@801:          */
sascha@424:         public void chain(Edge other, int found) {
sascha@424: 
sascha@424:             if (found == a) {
sascha@424:                 other.next = this;
sascha@424:                 prev = other;
sascha@424:                 return;
sascha@424:             }
sascha@424: 
sascha@424:             if (found == b) {
sascha@424:                 next = other;
sascha@424:                 other.prev = this;
sascha@424:                 return;
sascha@424:             }
sascha@424: 
sascha@424:             throw new IllegalStateException("cannot chain");
sascha@424:         }
sascha@424: 
sascha@801:         /**
sascha@801:          * Tells if the list is complete which means that this
sascha@801:          * edge list is a closed polygon.
sascha@801:          * @return true if edge list is closed else false.
sascha@801:          */
sascha@424:         public boolean isComplete() {
sascha@424:             Edge current = this;
sascha@424:             do {
sascha@424:                 if (current.prev == null || current.next == null) {
sascha@424:                     return false;
sascha@424:                 }
sascha@424:                 current = current.next;
sascha@424:             }
sascha@424:             while (current != this);
sascha@424:             return true;
sascha@424:         }
sascha@424: 
sascha@801:         /**
sascha@801:          * Returns the length of this edge list in next direction.
sascha@801:          * Segments in prev direction are ignored.
sascha@801:          * @return The length of this edge list.
sascha@801:          */
sascha@424:         public int length() {
sascha@424:             int length = 0;
sascha@424:             Edge current = this;
sascha@778:             do { ++length; }
sascha@437:             while ((current = current.next) != null && current != this);
sascha@424:             return length;
sascha@424:         }
sascha@436: 
sascha@801:         /**
sascha@801:          * Returns the head node of this edge list.
sascha@801:          * @return The head node.
sascha@801:          */
sascha@801:         public Edge head() {
sascha@801:             Edge current = this;
sascha@801:             while (current.prev != null) {
sascha@801:                 current = current.prev;
sascha@801:             }
sascha@801:             return current;
sascha@801:         }
sascha@447: 
sascha@801:         /**
sascha@801:          * Hashes the two vertex indices to a common value.
sascha@801:          * @return The hash value.
sascha@801:          */
sascha@801:         @Override
sascha@436:         public int hashCode() {
sascha@436:             return (a << 16) | b;
sascha@436:         }
sascha@436: 
sascha@801:         /**
sascha@801:          * Two edges are considered equal if they have the same
sascha@801:          * a and b vertices.
sascha@801:          * @param other The other edge.
sascha@801:          * @return true if the edges are equal else false.
sascha@801:          */
sascha@801:         @Override
sascha@436:         public boolean equals(Object other) {
sascha@436:             Edge e = (Edge)other;
sascha@436:             return a == e.a && b == e.b;
sascha@436:         }
sascha@424:     } // class Edge
sascha@424: 
sascha@801:     /**
sascha@801:      * Simplifies a given edge list by removing collinear vertices.
sascha@801:      * Attention: The original list is modified.
sascha@801:      * @param edge The edge list to simplify.
sascha@801:      * @param width The width of the vertex index space.
sascha@801:      * @return The simplified list.
sascha@801:      */
sascha@424:     protected static Edge simplify(Edge edge, int width) {
sascha@424: 
sascha@447:         Edge e1 = edge, start = edge;
sascha@424: 
sascha@424:         int length = edge.length();
sascha@424: 
sascha@437:         if (length < 2) {
sascha@437:             return e1;
sascha@437:         }
sascha@437: 
sascha@447:         Edge e2 = edge.next;
sascha@447: 
sascha@424:         int count = 0;
sascha@424: 
sascha@424:         do {
sascha@424:             int e1x = e1.a % width;
sascha@424:             int e1y = e1.a / width;
sascha@447:             int e2x = e1.b % width;
sascha@447:             int e2y = e1.b / width;
sascha@447:             int e3x = e2.b % width;
sascha@447:             int e3y = e2.b / width;
sascha@424: 
sascha@424:             if ((e1x == e2x && e2x == e3x && e1x == e3x)
sascha@447:                     ||  (e1y == e2y && e2y == e3y && e1y == e3y)) {
sascha@447:                 e1.b = e2.b;
sascha@447:                 Edge removed = e1.next;
sascha@447:                 e1.next = e2.next;
sascha@447:                 if (e1.next != null) {
sascha@447:                     e1.next.prev = e1;
sascha@447:                 }
sascha@447:                 e2 = e1.next;
sascha@424:                 count = 0;
sascha@424:                 --length;
sascha@447:                 if (removed == start) {
sascha@447:                     start = e1;
sascha@447:                 }
sascha@424:             }
sascha@424:             else {
sascha@424:                 e1 = e1.next;
sascha@424:                 e2 = e2.next;
sascha@424:                 ++count;
sascha@424:             }
sascha@424:         }
sascha@447:         while (length > 1 && e2 != null && count < length + 2);
sascha@424: 
sascha@447:         return start;
sascha@424:     }
sascha@424: 
sascha@801:     /**
sascha@801:      * The raster to be traced.
sascha@801:      */
sascha@424:     protected int []            raster;
sascha@801:     /**
sascha@801:      * The width of the raster.
sascha@801:      */
sascha@424:     protected int               width;
sascha@801:     /**
sascha@801:      * Map of the currently open edges.
sascha@801:      */
sascha@424:     protected TIntObjectHashMap openEdges;
sascha@801:     /**
sascha@801:      * List of rings already found.
sascha@801:      */
sascha@436:     protected List<Edge>        rings;
sascha@801:     /**
sascha@801:      * Flag to signal if a simplification should be performed
sascha@801:      * after a ring is completed.
sascha@801:      */
sascha@436:     protected boolean           simplify;
sascha@424: 
sascha@801:     /**
sascha@801:      * Default constructor. Simplification is turned on.
sascha@801:      */
sascha@424:     public Vectorizer() {
sascha@436:         this(true);
sascha@436:     }
sascha@436: 
sascha@801:     /**
sascha@801:      * Constructor to create a vectorized with an explicit
sascha@801:      * simplification policy.
sascha@801:      * @param simplify Indicates if simplification should be
sascha@801:      * used on ring completion.
sascha@801:      */
sascha@436:     public Vectorizer(boolean simplify) {
sascha@436:         openEdges     = new TIntObjectHashMap();
sascha@436:         rings         = new ArrayList<Edge>();
sascha@436:         this.simplify = simplify;
sascha@424:     }
sascha@424: 
sascha@801:     /**
sascha@801:      * Constructor to create a vectorizer with a given raster and width.
sascha@801:      * Simplification is turn on.
sascha@801:      * @param raster The raster to be vectorized.
sascha@801:      * @param width The width of the raster.
sascha@801:      */
sascha@424:     public Vectorizer(int [] raster, int width) {
sascha@445:         this(true, raster, width);
sascha@445:     }
sascha@445: 
sascha@801:     /**
sascha@801:      * Constructor to create a vectorizer with a given raster, width
sascha@801:      * and an explicit simplification policy.
sascha@801:      * @param simplify Indicates if simplification should be
sascha@801:      * used on ring completion.
sascha@801:      * @param raster The raster to be vectorized.
sascha@801:      * @param width The width of the raster.
sascha@801:      */
sascha@445:     public Vectorizer(boolean simplify, int [] raster, int width) {
sascha@445:         this(simplify);
sascha@424:         this.raster = raster;
sascha@424:         this.width  = width;
sascha@424:     }
sascha@424: 
sascha@801:     /**
sascha@801:      * Returns (x, y+1) for given vertex in index space.
sascha@801:      * @param i vertex in index space.
sascha@801:      * @param w width of raster.
sascha@801:      * @return (x, y+1) in index space.
sascha@801:      */
sascha@424:     public static final int tl(int i, int w) {
sascha@424:         int x = i % w;
sascha@424:         int y = i / w;
sascha@424:         return x + (w + 1)*y;
sascha@424:     }
sascha@424: 
sascha@801:     /**
sascha@801:      * Returns tl(i, w) + 1.
sascha@801:      * @param i vertex in index space.
sascha@801:      * @param w width of raster.
sascha@801:      * @return tl(i, w) + 1.
sascha@801:      */
sascha@424:     public static final int tr(int i, int w) {
sascha@424:         return tl(i, w) + 1;
sascha@424:     }
sascha@424: 
sascha@801:     /**
sascha@801:      * Returns tl(i, w) + w + 1.
sascha@801:      * @param i vertex in index space.
sascha@801:      * @param w width of raster.
sascha@801:      * @return tl(i, w) + w + 1.
sascha@801:      */
sascha@424:     public static final int bl(int i, int w) {
sascha@424:         return tl(i, w) + w + 1;
sascha@424:     }
sascha@424: 
sascha@801:     /**
sascha@801:      * Returns bl(i, w) + 1.
sascha@801:      * @param i vertex in index space.
sascha@801:      * @param w width of raster.
sascha@801:      * @return bl(i, w) + 1.
sascha@801:      */
sascha@424:     public static final int br(int i, int w) {
sascha@424:         return bl(i, w) + 1;
sascha@424:     }
sascha@424: 
sascha@801:     /**
sascha@801:      * Resets open resources after a set of features were found.
sascha@801:      */
sascha@424:     protected void resetRegion() {
sascha@424:         openEdges.clear();
sascha@424:         rings.clear();
sascha@424:     }
sascha@424: 
sascha@801:     /**
sascha@801:      * Adds an edge to the map of open edges, joins it
sascha@801:      * with its direct neighbors of if complete add the
sascha@801:      * list to the complete features.
sascha@801:      * @param edge
sascha@801:      */
sascha@424:     protected void emit(Edge edge) {
sascha@424: 
sascha@424:         Edge otherA = (Edge)openEdges.remove(edge.a);
sascha@424:         if (otherA != null) {
sascha@424:             otherA.chain(edge, edge.a);
sascha@424:         }
sascha@424: 
sascha@424:         Edge otherB = (Edge)openEdges.remove(edge.b);
sascha@424:         if (otherB != null) {
sascha@424:             otherB.chain(edge, edge.b);
sascha@424:         }
sascha@424: 
sascha@424:         if (edge.isComplete()) {
sascha@436:             rings.add(simplify ? simplify(edge, width + 1) : edge);
sascha@424:         }
sascha@424:         else {
sascha@424:             if (otherA == null) {
sascha@424:                 openEdges.put(edge.a, edge);
sascha@424:             }
sascha@424:             if (otherB == null) {
sascha@424:                 openEdges.put(edge.b, edge);
sascha@424:             }
sascha@424:         }
sascha@424:     }
sascha@424: 
sascha@801:     /**
sascha@801:      * Vectorize the raster. The found features are fed into
sascha@801:      * the given ring handler.
sascha@801:      * @param handler The RingHandler to postprocess the found features.
sascha@801:      * @return The number of regions found.
sascha@801:      */
sascha@424:     public int process(RingsHandler handler) {
sascha@424: 
sascha@424:         BitSet visited = new BitSet(raster.length);
sascha@424: 
sascha@424:         TIntStack stack = new TIntStack();
sascha@424: 
sascha@424:         int regions = 0;
sascha@424: 
sascha@436:         int height = raster.length / width;
sascha@424: 
sascha@424:         for (int i = 0; i < raster.length; ++i) {
sascha@424:             if (visited.get(i)) {
sascha@424:                 continue;
sascha@424:             }
sascha@424: 
sascha@424:             ++regions;
sascha@424: 
sascha@424:             int currentValue = raster[i];
sascha@424:             visited.set(i);
sascha@424: 
sascha@424:             int current = i;
sascha@424:             visited.set(current);
sascha@424: 
sascha@424:             for (;;) {
sascha@424:                 int tl = tl(current, width);
sascha@424:                 int tr = tr(current, width);
sascha@424:                 int bl = bl(current, width);
sascha@424:                 int br = br(current, width);
sascha@424: 
sascha@424:                 int t = current - width;
sascha@424: 
sascha@424:                 if (t < 0) {
sascha@424:                     emit(new Edge(tr, tl));
sascha@424:                 }
sascha@424:                 else {
sascha@424:                     if (raster[t] != currentValue) {
sascha@424:                         emit(new Edge(tr, tl));
sascha@424:                     }
sascha@424:                     else {
sascha@424:                         if (!visited.get(t)) {
sascha@424:                             visited.set(t);
sascha@424:                             stack.push(t);
sascha@424:                         }
sascha@424:                     }
sascha@424:                 }
sascha@424: 
sascha@424:                 int b = current + width;
sascha@424: 
sascha@424:                 if (b >= raster.length) {
sascha@424:                     emit(new Edge(bl, br));
sascha@424:                 }
sascha@424:                 else {
sascha@424:                     if (raster[b] != currentValue) {
sascha@424:                         emit(new Edge(bl, br));
sascha@424:                     }
sascha@424:                     else {
sascha@424:                         if (!visited.get(b)) {
sascha@424:                             visited.set(b);
sascha@424:                             stack.push(b);
sascha@424:                         }
sascha@424:                     }
sascha@424:                 }
sascha@424: 
sascha@424:                 int x = current % width;
sascha@424: 
sascha@424:                 if (x == 0) {
sascha@424:                     emit(new Edge(tl, bl));
sascha@424:                 }
sascha@424:                 else {
sascha@424:                     int l = current - 1;
sascha@424:                     if (raster[l] != currentValue) {
sascha@424:                         emit(new Edge(tl, bl));
sascha@424:                     }
sascha@424:                     else {
sascha@424:                         if (!visited.get(l)) {
sascha@424:                             visited.set(l);
sascha@424:                             stack.push(l);
sascha@424:                         }
sascha@424:                     }
sascha@424:                 }
sascha@424: 
sascha@424:                 if (x == width - 1) {
sascha@424:                     emit(new Edge(br, tr));
sascha@424:                 }
sascha@424:                 else {
sascha@424:                     int r = current + 1;
sascha@424:                     if (raster[r] != currentValue) {
sascha@424:                         emit(new Edge(br, tr));
sascha@424:                     }
sascha@424:                     else {
sascha@424:                         if (!visited.get(r)) {
sascha@424:                             visited.set(r);
sascha@424:                             stack.push(r);
sascha@424:                         }
sascha@424:                     }
sascha@424:                 }
sascha@424: 
sascha@424:                 if (stack.size() == 0) {
sascha@424:                     break;
sascha@424:                 }
sascha@424: 
sascha@424:                 current = stack.pop();
sascha@424:             }
sascha@424: 
sascha@436:             handler.handleRings(
sascha@778:                 rings,
sascha@778:                 currentValue,
sascha@436:                 width + 1,
sascha@436:                 height + 1);
sascha@424: 
sascha@424:             resetRegion();
sascha@424:         }
sascha@424: 
sascha@424:         return regions;
sascha@424:     }
sascha@424: }
sascha@798: // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :