diff --git a/pom.xml b/pom.xml index 519a379..a93ef5e 100644 --- a/pom.xml +++ b/pom.xml @@ -25,6 +25,11 @@ junit 4.1 + + com.vividsolutions + jts + 1.13 + diff --git a/src/main/java/de/rondiplomatico/nds/NDSCoordinate.java b/src/main/java/de/rondiplomatico/nds/NDSCoordinate.java index 343468f..5f4fe06 100644 --- a/src/main/java/de/rondiplomatico/nds/NDSCoordinate.java +++ b/src/main/java/de/rondiplomatico/nds/NDSCoordinate.java @@ -96,6 +96,24 @@ public NDSCoordinate(long ndsMortonCoordinates) { latitude = lat; longitude = lon; } + + /* + * Computes a point along the distance of this and given NDSCoordinate + * + * @param p + * the given coordinate + * @param fraction + * the fraction + */ + public NDSCoordinate fractionBetween(NDSCoordinate p, double fraction) { + assert (0.0 <= fraction) && (fraction <= 1.0); + double lat = p.toWGS84().getLatitude(); + double lon = p.toWGS84().getLongitude(); + WGS84Coordinate lonlat = toWGS84(); + double longitude = lonlat.getLongitude(); + double latitude = lonlat.getLatitude(); + return new NDSCoordinate(longitude + fraction * (lon - longitude), latitude + fraction * (lat - latitude)); + } private void verify(int lon, int lat) { if (lat < MIN_LATITUDE || MAX_LATITUDE < lat) { diff --git a/src/main/java/de/rondiplomatico/nds/NDSEnvelope.java b/src/main/java/de/rondiplomatico/nds/NDSEnvelope.java new file mode 100644 index 0000000..cc4ddce --- /dev/null +++ b/src/main/java/de/rondiplomatico/nds/NDSEnvelope.java @@ -0,0 +1,138 @@ +package de.rondiplomatico.nds; + +import com.vividsolutions.jts.geom.Coordinate; +import com.vividsolutions.jts.geom.OctagonalEnvelope; + +/** + * NDSEnvelope allows to compute an envelope for a set of NDSCoordinates. + * + * No warranties for correctness, use at own risk. + * + * @author Andreas Hessenthaler + * @since 13.02.2020 + */ +public class NDSEnvelope { + + private static OctagonalEnvelope vividEnvelope = new OctagonalEnvelope(); + + /** + * NDSEnvelope of a given set of polygon points + * + * @param polygonCoordinates + * the polygonCoordinates + */ + public NDSEnvelope(double[][] polygonCoordinates) { + for(int i = 0; i < polygonCoordinates.length; i++) { + Coordinate vividCoord = new Coordinate(polygonCoordinates[i][0], polygonCoordinates[i][1]); + vividEnvelope.expandToInclude(vividCoord); + } + } + + /** + * NDSEnvelope constructor for OctagonalEnvelope + * + * @param e + * the OctagonalEnvelope e + */ + public NDSEnvelope(OctagonalEnvelope e) { + vividEnvelope = e; + } + + /** + * Get the NDSEnvelope + * + * @return + */ + public OctagonalEnvelope getEnvelope() { + return vividEnvelope; + } + + /** + * Get the SouthWest corner of the envelope + * + * @return + */ + public NDSCoordinate getSouthWest() { + return new NDSCoordinate(vividEnvelope.getMinX(), vividEnvelope.getMinY()); + } + + /** + * Get the SouthEast corner of the envelope + * + * @return + */ + public NDSCoordinate getSouthEast() { + return new NDSCoordinate(vividEnvelope.getMaxX(), vividEnvelope.getMinY()); + } + + /** + * Get the NorthEast corner of the envelope + * + * @return + */ + public NDSCoordinate getNorthEast() { + return new NDSCoordinate(vividEnvelope.getMinX(), vividEnvelope.getMaxY()); + } + + /** + * Get the NorthWest corner of the envelope + * + * @return + */ + public NDSCoordinate getNorthWest() { + return new NDSCoordinate(vividEnvelope.getMaxX(), vividEnvelope.getMaxY()); + } + + /** + * Get the master tile for the SouthWest, SouthEast, NorthEast, NorthWest corners for a given number of levels + * + * @param maxLevels + * the maxLevels + * @return masterTileInfo + * the masterTileInfo consisting of the {masterTileLevel, masterTileNumber} + */ + public int[] getMasterTileInfo(int maxLevels) { + NDSCoordinate point0 = getSouthWest(); + NDSCoordinate point1 = getSouthEast(); + NDSCoordinate point2 = getNorthEast(); + NDSCoordinate point3 = getNorthWest(); + return getMasterTileInfo(point0, point1, point2, point3, maxLevels); + } + + /** + * Get the master tile for a set of four NDSCoordinates for a given number of levels + * + * @param maxLevels + * the maxLevels + * @return masterTileInfo + * the masterTileInfo consisting of the {masterTileLevel, masterTileNumber} + */ + public int[] getMasterTileInfo(NDSCoordinate point0, NDSCoordinate point1, NDSCoordinate point2, NDSCoordinate point3, int maxLevels) { + int masterTileLevel = -1; + int masterTileNumber = -1; + for (int li = 0; li < maxLevels; li++) { + NDSTile currTile0 = new NDSTile(li, point0); + NDSTile currTile1 = new NDSTile(li, point1); + NDSTile currTile2 = new NDSTile(li, point2); + NDSTile currTile3 = new NDSTile(li, point3); + int currTileNumber0 = currTile0.getTileNumber(); + int currTileNumber1 = currTile1.getTileNumber(); + int currTileNumber2 = currTile2.getTileNumber(); + int currTileNumber3 = currTile3.getTileNumber(); + boolean singleTileNumber = (currTileNumber0 == currTileNumber1) && (currTileNumber0 == currTileNumber2) && (currTileNumber0 == currTileNumber3); + // if at least one tile ID is different, we discard the tile IDs and keep our previously detected tile ID (i.e. on the previous level) + if(!singleTileNumber) { + break; + } + // store tile info + masterTileLevel = li; + masterTileNumber = currTileNumber0; + } + // check if valid result + if (masterTileNumber == -1) { + System.out.println(">>>ERROR: Invalid master tile ID."); + System.exit(1); + } + return new int[] {masterTileLevel, masterTileNumber}; + } +} \ No newline at end of file diff --git a/src/main/java/de/rondiplomatico/nds/NDSFloodFill.java b/src/main/java/de/rondiplomatico/nds/NDSFloodFill.java new file mode 100644 index 0000000..a3f0d4a --- /dev/null +++ b/src/main/java/de/rondiplomatico/nds/NDSFloodFill.java @@ -0,0 +1,214 @@ +package de.rondiplomatico.nds; + +import java.awt.Point; +import java.util.LinkedList; +import java.util.Queue; + +/** + * Utility class for flood fill algorithm applied to marked NDSTiles deduced from polygon points. + * + * @author Andreas Hessenthaler + * @since 13.02.2020 + */ +public class NDSFloodFill { + int[][] ff = null; + int minX = Integer.MAX_VALUE; + int maxX = Integer.MIN_VALUE; + int minY = Integer.MAX_VALUE; + int maxY = Integer.MIN_VALUE; + int dimX = -1; + int dimY = -1; + int bg = 0; + int bound = 1; + int fill = bound; + int startAtX = 0; + int startAtY = 0; + + /** + * Constructor + * + * @param f + * the image input to perform flood fill on + * assumes: + * - non-empty array + * - closed polygon without intersections with boundary color 1 and background color 0 + * @param a + * the minimum x-index + * @param b + * the maximum x-index + * @param c + * the minimum y-index + * @param d + * the maximum y-index + */ + public NDSFloodFill(int[][] f, int a, int b, int c, int d) { + ff = f; + minX = a; + maxX = b; + minY = c; + maxY = d; + dimX = maxX - minX + 1; + dimY = maxY - minY + 1; + setStartingPoint(); + } + + /* + * Iterative four-way flood fill algorithm + * + */ + public void iterativeFloodFill() { + Queue q = new LinkedList(); + q.add(new Point(startAtX, startAtY)); + while (!q.isEmpty()) { + Point p = q.remove(); + if ((p.x < 0) || (p.y < 0) || (p.x > dimX-1 || (p.y > dimY-1))) { + continue; + } + if (ff[p.x][p.y] == bg) { + ff[p.x][p.y] = fill; + q.add(new Point(p.x+1, p.y )); // right + q.add(new Point(p.x, p.y+1)); // top + q.add(new Point(p.x-1, p.y )); // left + q.add(new Point(p.x, p.y-1)); // bottom + } + } + } + + /* + * Safety net for filling left over lines after applying flood fill algorithm + * + */ + public void fillHoles() { + // first, find vertical lines + boolean vline = false; + boolean vlineEnded = false; + int currStartY = -1; + int currStopY = -1; + for (int col = 1; col < dimY-1; col++) { + vline = false; + vlineEnded = false; + for (int row = 1; row < dimX-1; row++) { + // minimum requirement: left and right neighbors are not background + if ((ff[row][col-1] != bg) && (ff[row][col+1] != bg)) { + if ((ff[row][col] == bg) + && (ff[row-1][col] != bg)) { + // beginning of vertical line + currStartY = row; + vline = true; + } else if (vline && (ff[row][col] == bg)) { + // vertical line continues + } else if (vline && (ff[row][col] != bg)) { + // end of vertical line + currStopY = row; + vlineEnded = true; + } + } else if (vline && ((ff[row][col+1] != bg) || (ff[row][col-1] != bg)) && (ff[row][col] != bg)) { + // end of horizontal line + currStopY = row; + vlineEnded = true; + } else { + vline = false; + vlineEnded = false; + } + if (vline && vlineEnded) { + // fill line and reset booleans + System.out.println(">>>>>>INFO: Filling vertical line."); + for (int idx = currStartY; idx < currStopY; idx++) { + ff[idx][col] = fill; + } + vline = false; + vlineEnded = false; + } + } + } + // second, find horizontal lines + boolean hline = false; + boolean hlineEnded = false; + for (int row = 1; row < dimX-1; row++) { + hline = false; + hlineEnded = false; + for (int col = 1; col < dimY-1; col++) { + // minimum requirement: top and bottom neighbors are not background + if ((ff[row+1][col] != bg) && (ff[row-1][col] != bg)) { + if ((ff[row][col] == bg) + && (ff[row][col-1] != bg)) { + // beginning of horizontal line + currStartY = col; + hline = true; + } else if (hline && (ff[row][col] == bg)) { + // vertical line continues + } else if (hline && (ff[row][col] != bg)) { + // end of vertical line + currStopY = col; + hlineEnded = true; + } + } else if (hline && ((ff[row+1][col] != bg) || (ff[row-1][col] != bg)) && (ff[row][col] != bg)) { + // end of vertical line + currStopY = col; + hlineEnded = true; + } else { + hline = false; + hlineEnded = false; + } + if (hline && hlineEnded) { + // fill line and reset booleans + System.out.println(">>>>>>INFO: Filling horizontal line."); + for (int idx = currStartY; idx < currStopY; idx++) { + ff[row][idx] = fill; + } + hline = false; + hlineEnded = false; + } + } + } + + return; + } + + /* + * Set initial point for flood fill algorithm + * + * @todo need to do more extensive checks if we can always find an initial point if we only check the right half of the box + */ + private void setStartingPoint() { + startAtX = (int) (0.5 * dimX); + startAtY = 0; + // we're starting on background, skip that.. + while (ff[startAtX][startAtY] == bg) { + startAtY++; + } + // we're continuing on an edge, skip that.. + while (ff[startAtX][startAtY] == bound) { + startAtY++; + } + for (int row = startAtY; row < dimY; row++, startAtY++) { + int rowSum = 0; + for (int col = startAtX; col < dimX; col++) { + if (ff[col][row-1] != ff[col][row]) { + rowSum = rowSum + ff[col][row]; + } + } + if ((rowSum % 2) == 1) { + // we've found a point that's inside + break; + } + } + return; + } + + /* + * Get initial point for flood fill algorithm + * + */ + public int[] getStartingPoint() { + return new int[] {startAtX, startAtY}; + } + + /* + * Get flood filled image + * + */ + public int[][] getFloodFill(){ + return ff; + } +} diff --git a/src/main/java/de/rondiplomatico/nds/NDSHashmap.java b/src/main/java/de/rondiplomatico/nds/NDSHashmap.java new file mode 100644 index 0000000..e124a86 --- /dev/null +++ b/src/main/java/de/rondiplomatico/nds/NDSHashmap.java @@ -0,0 +1,93 @@ +package de.rondiplomatico.nds; + +import java.util.ArrayList; +import java.util.Collections; +import java.util.HashMap; +import java.util.List; +import java.util.Map; + +/** + * NDSHashmap allows to store tileY - tileX pairs in a convenient format, + * where tileY is the key and tileX are the value. + * This is convenient because each tileY key has a (sorted) number of tileX value(s) + * and makes several tasks straighforward, such as: + * - drawing a map + * - applying flood-fill algorithm to get map tiles covered by polygon + * - etc. + * + * No warranties for correctness, use at own risk. + * + * @author Andreas Hessenthaler + * @since 28.02.2020 + */ +public class NDSHashmap { + + /** + * Compute a hash map from tile numbers. + * + * @param level + * the level + * @param tileNumbers + * the tileNumbers + * @return + */ + public Map> tileNumbersToHM(int level, int[] tileNumbers){ + // get the master tile number 0 + NDSTile masterTile = new NDSTile(0, 0); + // create hashmap of (key, value) pairs, where key is tileY and val is tileX + // note: val may contain multiple values + Map> tileHM = new HashMap>(); + for (int ti = 0; ti < tileNumbers.length; ti++) { + int[] tileXY = masterTile.getTileXYfromTileNumber(level, tileNumbers[ti]); + int key = tileXY[1]; + int newVal = tileXY[0]; + // if key already exists, add value to sorted list; else create new list + if (tileHM.containsKey(key)) { + List prevList = tileHM.get(key); + prevList.add(newVal); + Collections.sort(prevList); + tileHM.put(key, prevList); + } else { + List newList = new ArrayList(); + newList.add(newVal); + tileHM.put(key, newList); + } + } + return tileHM; + } + + /** + * Compute tile numbers from hash map. + * + * @param level + * the level + * @param hm + * the hash map + * @return + */ + public int[] hmToTileNumbers(int level, Map> hm) { + NDSTile masterTile = new NDSTile(0, 0); + int numVals = getNumberOfValuesHM(hm); + int[] filledTileIDs = new int[numVals]; + int idx = 0; + for (Map.Entry> entry : hm.entrySet()) { + int key = entry.getKey(); + List currList = entry.getValue(); + for (int idx0 = 0; idx0 < currList.size(); idx0++) { + filledTileIDs[idx] = masterTile.getTileNumberFromTileXY(level, currList.get(idx0), key); + idx++; + } + } + return filledTileIDs; + } + + + private int getNumberOfValuesHM(Map> hm) { + int numVals = 0; + for (Map.Entry> entry : hm.entrySet()) { + numVals = numVals + entry.getValue().size(); + } + return numVals; + } + +} diff --git a/src/main/java/de/rondiplomatico/nds/NDSPolyFillDemo.java b/src/main/java/de/rondiplomatico/nds/NDSPolyFillDemo.java new file mode 100644 index 0000000..0c3cf8a --- /dev/null +++ b/src/main/java/de/rondiplomatico/nds/NDSPolyFillDemo.java @@ -0,0 +1,108 @@ +package de.rondiplomatico.nds; + +import java.util.List; +import java.util.Map; + +import org.junit.Test; + +import de.rondiplomatico.nds.NDSUtils; +import de.rondiplomatico.nds.NDSHashmap; + +public class NDSPolyFillDemo { + + /** + * Main method to test selecting all tiles covered by a sample polygon. + * + * @param args + */ + @Test + public static void main(String[] args) { + + long t0 = System.currentTimeMillis(); + + // get some random polygon data for testing purposes +// double[][] polygonCoordinates = new double[][] { +// {10.5, 45.9}, +// {13.0, 50.3}, +// {15.0, 47.0}, +// {13.4, 70.0}, +// {10.5, 45.9} +// }; + + // get some random polygon data for testing purposes +// double[][] polygonCoordinates = new double[][] { +// {10.5, 45.9}, +// {13.0, 63.3}, +// {15.0, 47.0}, +// {13.4, 70.0}, +// {10.5, 45.9} +// }; + + // get some random polygon data for testing purposes +// double[][] polygonCoordinates = new double[][] { +// {10.5, 45.9}, +// {13.0, 50.3}, +// {15.0, 47.0}, +// {17.0, 50.3}, +// {20.0, 47.0}, +// {13.4, 60.0}, +// {13.4, 70.0}, +// {10.5, 45.9} +// }; + + // get some random polygon data approximating Germany + // https://www.mapsofworld.com/lat_long/germany-lat-long.html + double[][] polygonCoordinates = new double[][] { + {10.5, 45.9}, + {13.0, 45.9}, + {14.0, 49.0}, + {12.0, 50.0}, + {15.0, 51.0}, + {15.0, 54.0}, + {13.5, 54.5}, + {11.0, 54.0}, + {10.0, 55.0}, + { 8.5, 55.0}, + { 9.0, 54.0}, + { 7.0, 53.5}, + { 6.0, 52.0}, + { 6.1, 50.0}, + { 8.0, 49.0}, + { 7.5, 47.5}, + {10.5, 45.9} + }; + + // number of levels in map hierarchy + int maxLevels = 15; + // get a bounding octagonal envelope (defaults to quadrilateral in 2D case) + // for sample data corresponding to a polygon (e.g. borders of a country) + NDSEnvelope envelope = new NDSEnvelope(polygonCoordinates); + // get corresponding bounding tile ID, i.e. find level where all polygon points are on the same tile + int[] masterTileInfo = envelope.getMasterTileInfo(maxLevels); + int masterTileLevel = masterTileInfo[0]; + int masterTileID = masterTileInfo[1]; + // store the master tile + NDSTile masterTile = new NDSTile(masterTileLevel, masterTileID); + // get all tiles covered by the polygon on the tstLevel + int tstLevel = 11; + // refine polygon coordinates to avoid edges that are crossing a tile + // set refinement factor to -1 to adaptively refine + int numSamples = -1; + double[][] polygonCoordinatesRef = NDSPolygon.refinePolygon(tstLevel, polygonCoordinates, numSamples); + // let's grab all tiles with polygon points + int[] uniqueTileIDs = NDSPolygon.getUniqueTileNumbersOnLevel(tstLevel, polygonCoordinatesRef); + // dump to image file for debugging + NDSUtils.printMap(masterTile, tstLevel, uniqueTileIDs, "png", "map"); + // get x/y tile indices for tiles with polygon points + // key is tileY, val is tileX (note: val may contain multiple tileX indices) + NDSHashmap hm = new NDSHashmap(); + Map> tileHM = hm.tileNumbersToHM(tstLevel, uniqueTileIDs); + // dump to image file for debugging + Map> filledTileHM = NDSPolygon.mapFillPolygon(tileHM, uniqueTileIDs); + int[] filledTileIDs = hm.hmToTileNumbers(tstLevel, filledTileHM); + NDSUtils.printMap(masterTile, tstLevel, filledTileIDs, "png", "mapFilled"); + + long t1 = System.currentTimeMillis(); + System.out.println("\n>>>INFO: Program finished in "+(t1-t0)+" ms."); + } +} diff --git a/src/main/java/de/rondiplomatico/nds/NDSPolygon.java b/src/main/java/de/rondiplomatico/nds/NDSPolygon.java new file mode 100644 index 0000000..dfd7618 --- /dev/null +++ b/src/main/java/de/rondiplomatico/nds/NDSPolygon.java @@ -0,0 +1,146 @@ +package de.rondiplomatico.nds; + +import java.util.ArrayList; +import java.util.HashMap; +import java.util.List; +import java.util.Map; +import java.util.stream.IntStream; + +public class NDSPolygon { + + /** + * Applies a flood-fill algorithm to a polygon + * + * @param hmi + * the hash map input hmi + * @param tileNumbers + * the tileNumbers + * @return + */ + public static Map> mapFillPolygon(Map> hmi, int[] tileNumbers){ + Map> hmo = new HashMap>(); + // first get a minimum 'bounding tile box' + int minX = Integer.MAX_VALUE; + int maxX = Integer.MIN_VALUE; + int minY = Integer.MAX_VALUE; + int maxY = Integer.MIN_VALUE; + // loop over (key, values) pairs to get boundary + for (Map.Entry> entry : hmi.entrySet()) { + int key = entry.getKey(); + List currList = entry.getValue(); + // we have a sorted list, so we can just get the first and last element for the comparison + minX = Math.min(currList.get(0), minX); + maxX = Math.max(currList.get(currList.size()-1), maxX); + // key is the y-index + minY = Math.min(key, minY); + maxY = Math.max(key, maxY); + } + // get dimensions + int dimX = maxX - minX + 1; + int dimY = maxY - minY + 1; + // allocate int[][] for flood fill algorithm; default value 0 + int[][] tmp = new int[dimX][dimY]; + // loop over (key, values) pairs to label boundary + for (Map.Entry> entry : hmi.entrySet()) { + int key = entry.getKey() - minY; + List currList = entry.getValue(); + // loop over values to label boundary + for (int idx0 = 0; idx0 < currList.size(); idx0++) { + int val = currList.get(idx0) - minX; + tmp[val][key] = 1; + } + } + // perform flood fill from midpoint + NDSFloodFill ff = new NDSFloodFill(tmp, minX, maxX, minY, maxY); + ff.iterativeFloodFill(); + tmp = ff.getFloodFill(); + ff.fillHoles(); + // loop over y + for (int key = 0; key < dimY; key++) { + List currList = new ArrayList(); + // loop over x + for (int x = 0; x < dimX; x++) { + // get all elements that are on boundary or inside + if (tmp[x][key] > 0) { + currList.add(x+minX); + } + } + hmo.put(key+minY, currList); + } + return hmo; + } + + /** + * Refine polygon data by subsampling. + * + * @param level + * the level + * @param polygonCoordinates + * the polygon coordinates + * @param numSamples + * number of samples between each polygon edge (-1 requests adaptive sampling) + * @return polygonCoordinatesRef + * the refined polygon coordinates + */ + public static double[][] refinePolygon(int level, double[][] polygonCoordinates, int numSamples) { + if (numSamples == 0) { + System.out.println(">>>INFO: Polygon refinement off. Consider switching to adaptive refinement."); + return polygonCoordinates; + } + if (numSamples == -1) { + double maxDist = 0.0; + double tileSizeX = 360.0 / (double) (Math.pow(2, level+1)); + for (int idx = 0; idx < polygonCoordinates.length-1; idx++) { + double p0x = polygonCoordinates[idx ][0]; + double p0y = polygonCoordinates[idx ][1]; + double p1x = polygonCoordinates[idx+1][0]; + double p1y = polygonCoordinates[idx+1][1]; + double dist = Math.sqrt(Math.pow(p1x-p0x, 2.0) + Math.pow(p1y-p0y, 2.0)); + maxDist = Math.max(dist, maxDist); + } + double maxTargetDist = tileSizeX * 0.4; + if (maxDist <= maxTargetDist) { + numSamples = 1; + } else { + numSamples = (int) Math.ceil(maxDist / maxTargetDist); + } + System.out.println(">>>INFO: Setting adaptive refinement to "+numSamples); + } + int numCoord = polygonCoordinates.length; + int numCoordRef = (numCoord - 1) * (numSamples + 1) + 1; + double[][] polygonCoordinatesRef = new double[numCoordRef][2]; + // sample between all points except for the last and first because those coincide + for (int idx = 0; idx < numCoord-1; idx++) { + double x0 = polygonCoordinates[idx ][0]; + double x1 = polygonCoordinates[idx+1][0]; + double y0 = polygonCoordinates[idx ][1]; + double y1 = polygonCoordinates[idx+1][1]; + for (int jdx = 0; jdx < numSamples+1; jdx++) { + polygonCoordinatesRef[idx*(numSamples+1)+jdx][0] = x0 + (x1 - x0) * (double)jdx / (double)(numSamples + 1); + polygonCoordinatesRef[idx*(numSamples+1)+jdx][1] = y0 + (y1 - y0) * (double)jdx / (double)(numSamples + 1); + } + } + // final point to close polygon + polygonCoordinatesRef[numCoordRef-1][0] = polygonCoordinates[numCoord-1][0]; + polygonCoordinatesRef[numCoordRef-1][1] = polygonCoordinates[numCoord-1][1]; + // return refined polygon + return polygonCoordinatesRef; + } + + public static int[] getTileNumbersOnLevel(int level, double[][] polygonCoordinates){ + int numCoord = polygonCoordinates.length; + int[] tileNumbers = new int[numCoord]; + // loop over all coordinates to get map tile IDs + for (int idx = 0; idx < numCoord; idx++) { + NDSCoordinate currCoord = new NDSCoordinate(polygonCoordinates[idx][0], polygonCoordinates[idx][1]); + NDSTile currTile = new NDSTile(level, currCoord); + tileNumbers[idx] = currTile.getTileNumber(); + } + return tileNumbers; + } + + public static int[] getUniqueTileNumbersOnLevel(int level, double[][] polygonCoordinates){ + int[] tileNumbers = getTileNumbersOnLevel(level, polygonCoordinates); + return IntStream.of(tileNumbers).distinct().toArray(); + } +} diff --git a/src/main/java/de/rondiplomatico/nds/NDSTile.java b/src/main/java/de/rondiplomatico/nds/NDSTile.java index b3e1c2a..4cee161 100644 --- a/src/main/java/de/rondiplomatico/nds/NDSTile.java +++ b/src/main/java/de/rondiplomatico/nds/NDSTile.java @@ -1,5 +1,17 @@ package de.rondiplomatico.nds; +import java.awt.image.BufferedImage; +import java.io.File; +import java.io.IOException; +import java.util.ArrayList; +import java.util.Arrays; +import java.util.Collections; +import java.util.HashMap; +import java.util.List; +import java.util.Map; + +import javax.imageio.ImageIO; + import lombok.EqualsAndHashCode; import lombok.Getter; import lombok.ToString; @@ -11,6 +23,7 @@ * No warranties for correctness, use at own risk. * * @author Daniel Wirtz + * @author Andreas Hessenthaler * @since 20.02.2020 */ @Getter @@ -108,6 +121,121 @@ public NDSTile(int level, WGS84Coordinate coord) { this(level, new NDSCoordinate(coord.getLongitude(), coord.getLatitude())); } + /** + * Get tile indices in x / y from given tile number and level. + * + * @param level + * the level + * @param nr + * the tile number + */ + public int[] getTileXYfromTileNumber(int level, int nr) { + NDSTile tile = new NDSTile(level, nr); + // variables for storing x/y-tile indices + int tileX, tileY; + tileX = tileY = -1; + // get tile size, assuming that there are 2x1 tiles on level 0 + double tileSizeX, tileSizeY; + tileSizeX = 360.0 / (double) (Math.pow(2, level+1)); + tileSizeY = 180.0 / (double) (Math.pow(2, level )); + // get tile center and southwest corner + NDSCoordinate c = tile.getCenter(); + double lon = c.toWGS84().getLongitude() - 0.5 * tileSizeX; + double lat = c.toWGS84().getLatitude() - 0.5 * tileSizeY; + // compute tile indices + tileX = (int) Math.round((lon + 180.0) / tileSizeX); + tileY = (int) Math.round((lat + 90.0) / tileSizeY); + return new int[] {tileX, tileY}; + } + + /** + * Get tile number from indices in x / y and level. + * + * @param level + * the level + * @param tileX + * the tile x-coordinate + * @param tileY + * the tile y-coordinate + */ + public int getTileNumberFromTileXY(int level, int tileX, int tileY) { + // get tile size, assuming that there are 2x1 tiles on level 0 + double tileSizeX = 360.0 / (double) (Math.pow(2, level+1)); + double tileSizeY = 180.0 / (double) (Math.pow(2, level )); + // get center + double clon = tileX * tileSizeX + 0.5 * tileSizeX - 180.0; + double clat = tileY * tileSizeY + 0.5 * tileSizeY - 90.0; + NDSCoordinate c = new NDSCoordinate(clon, clat); + // get tile + NDSTile t = new NDSTile(level, c); + // return tile number + return t.getTileNumber(); + } + + /** + * Get tile quadkey from tile indices in x / y and level. + * + * @param level + * the level + * @param tileX + * the tile x-coordinate + * @param tileY + * the tile y-coordinate + */ + public String getTileQuadkeyFromTileXY(int level, int tileX, int tileY) { + // setup object for building quadkey + StringBuilder quadkey = new StringBuilder(); + // loop over levels + for (int i = level+1; i > 0; i--) { + char digit = '0'; + int mask = 1 << (i - 1); + if ((tileX & mask) != 0) { + digit++; + } + if ((tileY & mask) != 0) { + digit++; + digit++; + } + // add to quadkey + quadkey.append(digit); + } + // return quadkey + return quadkey.toString(); + } + + /** + * Get tile indices in x / y from tile quadkey. + * + * @param quadkey + * the quadkey + */ + public int[] getTileXYfromTileQuadkey(String quadkey) { + int[] tileXY = new int[2]; + // 0-based level index + int level = quadkey.length(); + for (int i = level; i > 0; i--) { + int mask = 1 << (i - 1); + switch (quadkey.charAt(level - i)) { + case '0': + break; + case '1': + tileXY[0] |= mask; + break; + case '2': + tileXY[1] |= mask; + break; + case '3': + tileXY[0] |= mask; + tileXY[1] |= mask; + break; + default: + System.out.println(">>>ERROR:"); + throw new IllegalArgumentException("Invalid quadkey digit sequence."); + } + } + return tileXY; + } + /** * Checks if the current Tile contains a certain coordinate. * @@ -155,7 +283,25 @@ public NDSCoordinate getCenter() { } return center; } - + + /** + * Get child tile for given number + * + * @param nr + * the child tile number (SouthWest, SouthEast, NorthEast, NorthWest) + * @return + */ + public NDSTile getChild(int nr) { + assert (0 <= nr) && (nr <= 3); + int childTileNumber = (tileNumber << 2) + nr; + if (nr == 2) { + childTileNumber = (tileNumber << 2) + 3; + } else if (nr == 3) { + childTileNumber = (tileNumber << 2) + 2; + } + return new NDSTile(level, childTileNumber); + } + /** * Creates a bounding box for the current tile. * diff --git a/src/main/java/de/rondiplomatico/nds/NDSUtils.java b/src/main/java/de/rondiplomatico/nds/NDSUtils.java new file mode 100644 index 0000000..604aa6e --- /dev/null +++ b/src/main/java/de/rondiplomatico/nds/NDSUtils.java @@ -0,0 +1,73 @@ +package de.rondiplomatico.nds; + +import java.awt.image.BufferedImage; +import java.io.File; +import java.io.IOException; +import java.util.Arrays; + +import javax.imageio.ImageIO; + +import de.rondiplomatico.nds.NDSTile; + +public class NDSUtils { + + /** + * Print map with labeled tiles. + * + * @param level + * the level + * @param tiles + * the relevant tiles + * @param printType + * the print type (e.g., png or text) + * @param fileName + * the output filename + */ + public static void printMap(NDSTile masterTile, int level, int[] tiles, String printType, String fileName) { + // get number of tiles in x/y + int numTilesX = (int) Math.pow(2, level+1); + int numTilesY = (int) Math.pow(2, level ); + // print labeled cells + if (printType.equalsIgnoreCase("text")) { + if (level > 4) { + System.out.println(">>>WARNING: Writing to console won't give readable output. Skipping.."); + return; + } + String header = new String(new char[(numTilesX+1)*2-1]).replace('\0', '-'); + for (int idxY = numTilesY-1; idxY >= 0; idxY--) { + System.out.println(header); + for (int idxX = 0; idxX < numTilesX; idxX++) { + int tileID = masterTile.getTileNumberFromTileXY(level, idxX, idxY); + // print label + // note: if check is quite inefficient.. + if (Arrays.stream(tiles).anyMatch(i -> i == tileID)) { + System.out.print("|X"); + } else { + System.out.print("| "); + } + } + System.out.println("|"); + } + System.out.println(header); + } else if (printType.equalsIgnoreCase("png")) { + if (level > 13) { + System.out.println(">>>WARNING: Required heap space is "+((long)numTilesX*(long)numTilesY*(long)4/(long)1024/(long)1024)+" MB."); + System.out.println(">>>WARNING: Writing to buffered image may fail due to limited heap space. Skipping.."); + return; + } + // create image using default value for type int: 0 + BufferedImage image = new BufferedImage(numTilesX, numTilesY, BufferedImage.TYPE_INT_RGB); + for (int idx = 0; idx < tiles.length; idx++) { + int[] tileXY = masterTile.getTileXYfromTileNumber(level, tiles[idx]); + image.setRGB(tileXY[0], numTilesY-1-tileXY[1], -16000); + } + File ImageFile = new File(System.getProperty("user.home"), fileName+".png"); + try { + ImageIO.write(image, "png", ImageFile); + } catch (IOException e) { + e.printStackTrace(); + } + } + return; + } +} diff --git a/src/test/java/de/rondiplomatico/nds/NDSCoordinateTest.java b/src/test/java/de/rondiplomatico/nds/NDSCoordinateTest.java index ed5225a..83da099 100644 --- a/src/test/java/de/rondiplomatico/nds/NDSCoordinateTest.java +++ b/src/test/java/de/rondiplomatico/nds/NDSCoordinateTest.java @@ -17,6 +17,7 @@ * Test values are taken from the NDS Format Specification, Version 2.5.4. * * @author Daniel Wirtz + * @author Andreas Hessenthaler * @since 20.02.2020 */ @Getter @@ -24,6 +25,8 @@ @EqualsAndHashCode public class NDSCoordinateTest { + private final double eps = 1E-7; + @Test public void testConstructor() { try { @@ -54,6 +57,50 @@ public void testNDSCoordinateCorners() { assertEquals(new NDSCoordinate(0.0, 0.0), new NDSCoordinate(0, 0)); } + @Test + public void testFractionBetweenWorks() { + // Tile 0 + NDSTile t = new NDSTile(0, 0); + NDSCoordinate c = t.getCenter(); + NDSBBox bb = t.getBBox(); + NDSCoordinate sw = bb.southWest(); + NDSCoordinate se = bb.southEast(); + NDSCoordinate ne = bb.northEast(); + NDSCoordinate nw = bb.northWest(); + NDSCoordinate csw = sw.fractionBetween(c, 0.5); + NDSCoordinate cse = se.fractionBetween(c, 0.5); + NDSCoordinate cne = ne.fractionBetween(c, 0.5); + NDSCoordinate cnw = nw.fractionBetween(c, 0.5); + assertEquals( 45.0, csw.toWGS84().getLongitude(), eps); + assertEquals(135.0, cse.toWGS84().getLongitude(), eps); + assertEquals( 45.0, cnw.toWGS84().getLongitude(), eps); + assertEquals(135.0, cne.toWGS84().getLongitude(), eps); + assertEquals(-45.0, csw.toWGS84().getLatitude(), eps); + assertEquals(-45.0, cse.toWGS84().getLatitude(), eps); + assertEquals( 45.0, cnw.toWGS84().getLatitude(), eps); + assertEquals( 45.0, cne.toWGS84().getLatitude(), eps); + // Tile 1 + t = new NDSTile(0, 1); + c = t.getCenter(); + bb = t.getBBox(); + sw = bb.southWest(); + se = bb.southEast(); + ne = bb.northEast(); + nw = bb.northWest(); + csw = sw.fractionBetween(c, 0.5); + cse = se.fractionBetween(c, 0.5); + cne = ne.fractionBetween(c, 0.5); + cnw = nw.fractionBetween(c, 0.5); + assertEquals(-135.0, csw.toWGS84().getLongitude(), eps); + assertEquals( -45.0, cse.toWGS84().getLongitude(), eps); + assertEquals(-135.0, cnw.toWGS84().getLongitude(), eps); + assertEquals( -45.0, cne.toWGS84().getLongitude(), eps); + assertEquals( -45.0, csw.toWGS84().getLatitude(), eps); + assertEquals( -45.0, cse.toWGS84().getLatitude(), eps); + assertEquals( 45.0, cnw.toWGS84().getLatitude(), eps); + assertEquals( 45.0, cne.toWGS84().getLatitude(), eps); + } + /** * * Verifies the values of Table 8-1 in Section "7.2.1 Coding of Coordinates", NDS Spec 2.5.4 diff --git a/src/test/java/de/rondiplomatico/nds/NDSEnvelopeTest.java b/src/test/java/de/rondiplomatico/nds/NDSEnvelopeTest.java new file mode 100644 index 0000000..baad0d0 --- /dev/null +++ b/src/test/java/de/rondiplomatico/nds/NDSEnvelopeTest.java @@ -0,0 +1,67 @@ +package de.rondiplomatico.nds; + +import static org.junit.Assert.assertEquals; + +import org.junit.Test; + +import com.vividsolutions.jts.geom.OctagonalEnvelope; + +import lombok.EqualsAndHashCode; +import lombok.Getter; +import lombok.ToString; + +/** + * Tests for the NDSEnvelope class. + * + * @author Andreas Hessenthaler + * @since 28.02.2020 + */ +@Getter +@ToString +@EqualsAndHashCode +public class NDSEnvelopeTest { + + private final double eps = 1E-7; + + @Test + public void testEnvelopeWorks() { + // rough polygon data approximating Germany + double[][] polygonCoordinates = new double[][] { + {10.5, 45.9}, + {13.0, 45.9}, + {14.0, 49.0}, + {12.0, 50.0}, + {15.0, 51.0}, + {15.0, 54.0}, + {13.5, 54.5}, + {11.0, 54.0}, + {10.0, 55.0}, + { 8.5, 55.0}, + { 9.0, 54.0}, + { 7.0, 53.5}, + { 6.0, 52.0}, + { 6.1, 50.0}, + { 8.0, 49.0}, + { 7.5, 47.5}, + {10.5, 45.9} + }; + NDSEnvelope envelope = new NDSEnvelope(polygonCoordinates); + OctagonalEnvelope vividEnvelope = envelope.getEnvelope(); + // check envelope boundaries + assertEquals( 6.0, vividEnvelope.getMinX(), eps); + assertEquals(15.0, vividEnvelope.getMaxX(), eps); + assertEquals(45.9, vividEnvelope.getMinY(), eps); + assertEquals(55.0, vividEnvelope.getMaxY(), eps); + // check bounding box corners + assertEquals(581131592357515410L, envelope.getSouthWest().getMortonCode()); + assertEquals(595825689965249734L, envelope.getSouthEast().getMortonCode()); + assertEquals(592806849050488888L, envelope.getNorthEast().getMortonCode()); + assertEquals(607500946658223212L, envelope.getNorthWest().getMortonCode()); + // check master tile details + int[] info = envelope.getMasterTileInfo(15); + int level = info[0]; + int number = info[1]; + assertEquals(3, level); + assertEquals(8, number); + } +} diff --git a/src/test/java/de/rondiplomatico/nds/NDSPolyFillDemoTest.java b/src/test/java/de/rondiplomatico/nds/NDSPolyFillDemoTest.java new file mode 100644 index 0000000..09613ab --- /dev/null +++ b/src/test/java/de/rondiplomatico/nds/NDSPolyFillDemoTest.java @@ -0,0 +1,13 @@ +package de.rondiplomatico.nds; + +import org.junit.Test; + +import de.rondiplomatico.nds.NDSPolyFillDemo; + +public class NDSPolyFillDemoTest { + @Test + public void testPolyFillDemoWorks() { + String[] test = {"Test", "Test2"}; + NDSPolyFillDemo.main(test); + } +} diff --git a/src/test/java/de/rondiplomatico/nds/NDSTileTest.java b/src/test/java/de/rondiplomatico/nds/NDSTileTest.java index 9d58ad6..4c9e9ca 100644 --- a/src/test/java/de/rondiplomatico/nds/NDSTileTest.java +++ b/src/test/java/de/rondiplomatico/nds/NDSTileTest.java @@ -6,6 +6,8 @@ import org.junit.Test; +import com.vividsolutions.jts.geom.OctagonalEnvelope; + import de.rondiplomatico.nds.NDSBBox; import de.rondiplomatico.nds.NDSCoordinate; import de.rondiplomatico.nds.NDSTile; @@ -17,6 +19,7 @@ * * * @author Daniel Wirtz + * @author Andreas Hessenthaler * @since 20.02.2020 */ public class NDSTileTest { @@ -298,6 +301,88 @@ public void testCenterWorks() { } } + @Test + public void testChildTileWorks() { + // get tile and its center + int masterTileLevel = 3; + NDSTile t = new NDSTile(masterTileLevel, 8); + NDSCoordinate c = t.getCenter(); + // store the master tile corners + NDSBBox bb = t.getBBox(); + NDSCoordinate sw = bb.southWest(); + NDSCoordinate se = bb.southEast(); + NDSCoordinate ne = bb.northEast(); + NDSCoordinate nw = bb.northWest(); + // first, let's create four points that we know are in the child tiles + int childTileLevel = masterTileLevel + 1; + NDSCoordinate csw = sw.fractionBetween(c, 0.5); + NDSCoordinate cse = se.fractionBetween(c, 0.5); + NDSCoordinate cne = ne.fractionBetween(c, 0.5); + NDSCoordinate cnw = nw.fractionBetween(c, 0.5); + // get child tiles based on master level for double-checking + NDSTile tsw = new NDSTile(childTileLevel, csw); + NDSTile tse = new NDSTile(childTileLevel, cse); + NDSTile tne = new NDSTile(childTileLevel, cne); + NDSTile tnw = new NDSTile(childTileLevel, cnw); + // now perform check + assertEquals(tsw.getTileNumber(), t.getChild(0).getTileNumber()); + assertEquals(tse.getTileNumber(), t.getChild(1).getTileNumber()); + assertEquals(tne.getTileNumber(), t.getChild(2).getTileNumber()); + assertEquals(tnw.getTileNumber(), t.getChild(3).getTileNumber()); + } + + @Test + public void testQuadkeyWorks() { + // rough polygon data approximating Germany + double[][] polygonCoordinates = new double[][] { + {10.5, 45.9}, + {13.0, 45.9}, + {14.0, 49.0}, + {12.0, 50.0}, + {15.0, 51.0}, + {15.0, 54.0}, + {13.5, 54.5}, + {11.0, 54.0}, + {10.0, 55.0}, + { 8.5, 55.0}, + { 9.0, 54.0}, + { 7.0, 53.5}, + { 6.0, 52.0}, + { 6.1, 50.0}, + { 8.0, 49.0}, + { 7.5, 47.5}, + {10.5, 45.9} + }; + // number of levels in map hierarchy + int maxLevels = 15; + // get a bounding octagonal envelope (defaults to quadrilateral in 2D case) + // for sample data corresponding to a polygon (e.g. borders of a country) + NDSEnvelope envelope = new NDSEnvelope(polygonCoordinates); + // get corresponding bounding tile ID, i.e. find level where all polygon points are on the same tile + int[] masterTileInfo = envelope.getMasterTileInfo(maxLevels); + int masterTileLevel = masterTileInfo[0]; + int masterTileID = masterTileInfo[1]; + // store the master tile + NDSTile masterTile = new NDSTile(masterTileLevel, masterTileID); + // set some random tile number and level for testing + int tstTileNumber = 7; + int tstTileLevel = 3; + int[] tstTileXY = masterTile.getTileXYfromTileNumber(tstTileLevel, tstTileNumber); + int tstTileX = tstTileXY[0]; + int tstTileY = tstTileXY[1]; + assertEquals(11, tstTileX); + assertEquals( 5, tstTileY); + String tstTileQuadkey = masterTile.getTileQuadkeyFromTileXY(tstTileLevel, tstTileX, tstTileY); + assertEquals("1213", tstTileQuadkey); + assertEquals("00313102310", masterTile.getTileQuadkeyFromTileXY(10, 486, 332)); + assertEquals("01202102332221212", masterTile.getTileQuadkeyFromTileXY(16, 35210, 21493)); + // inverse direction + tstTileXY = masterTile.getTileXYfromTileQuadkey(tstTileQuadkey); + assertEquals(tstTileX, tstTileXY[0]); + assertEquals(tstTileY, tstTileXY[1]); + assertEquals(tstTileNumber, masterTile.getTileNumberFromTileXY(tstTileLevel, tstTileX, tstTileY)); + } + @Test public void testConstructor() { try {