001/* 002 * Copyright (c) 2003 Objectix Pty Ltd All rights reserved. 003 * 004 * This library is free software; you can redistribute it and/or 005 * modify it under the terms of the GNU Lesser General Public 006 * License as published by the Free Software Foundation. 007 * 008 * THIS SOFTWARE IS PROVIDED "AS IS" AND ANY EXPRESSED OR IMPLIED 009 * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES 010 * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE 011 * DISCLAIMED. IN NO EVENT SHALL OBJECTIX PTY LTD BE LIABLE FOR ANY 012 * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 013 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE 014 * GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS 015 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, 016 * WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE 017 * OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, 018 * EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 019 */ 020package org.openstreetmap.josm.data.projection.datum; 021 022import java.io.IOException; 023import java.io.InputStream; 024import java.io.Serializable; 025import java.nio.charset.StandardCharsets; 026import java.util.Arrays; 027import java.util.Objects; 028 029import org.openstreetmap.josm.tools.Logging; 030import org.openstreetmap.josm.tools.Utils; 031 032/** 033 * Models the NTv2 Sub Grid within a Grid Shift File. 034 * 035 * @author Peter Yuill 036 * Modified for JOSM : 037 * - removed the RandomAccessFile mode (Pieren) 038 * - read grid file by single bytes. Workaround for a bug in some VM not supporting 039 * file reading by group of 4 bytes from a jar file. 040 * - removed the Cloneable interface 041 * @since 2507 042 */ 043public class NTV2SubGrid implements Serializable { 044 045 private static final long serialVersionUID = 1L; 046 047 private final String subGridName; 048 private final String parentSubGridName; 049 private final String created; 050 private final String updated; 051 private final double minLat; 052 private final double maxLat; 053 private final double minLon; 054 private final double maxLon; 055 private final double latInterval; 056 private final double lonInterval; 057 private final int nodeCount; 058 059 private final int lonColumnCount; 060 private final int latRowCount; 061 private final float[] latShift; 062 private final float[] lonShift; 063 private float[] latAccuracy; 064 private float[] lonAccuracy; 065 066 private NTV2SubGrid[] subGrid; 067 068 /** 069 * Construct a Sub Grid from an InputStream, loading the node data into 070 * arrays in this object. 071 * 072 * @param in GridShiftFile InputStream 073 * @param bigEndian is the file bigEndian? 074 * @param loadAccuracy is the node Accuracy data to be loaded? 075 * @throws IOException if any I/O error occurs 076 */ 077 public NTV2SubGrid(InputStream in, boolean bigEndian, boolean loadAccuracy) throws IOException { 078 byte[] b8 = new byte[8]; 079 byte[] b4 = new byte[4]; 080 byte[] b1 = new byte[1]; 081 readBytes(in, b8); 082 readBytes(in, b8); 083 subGridName = new String(b8, StandardCharsets.UTF_8).trim(); 084 readBytes(in, b8); 085 readBytes(in, b8); 086 parentSubGridName = new String(b8, StandardCharsets.UTF_8).trim(); 087 readBytes(in, b8); 088 readBytes(in, b8); 089 created = new String(b8, StandardCharsets.UTF_8); 090 readBytes(in, b8); 091 readBytes(in, b8); 092 updated = new String(b8, StandardCharsets.UTF_8); 093 readBytes(in, b8); 094 readBytes(in, b8); 095 minLat = NTV2Util.getDouble(b8, bigEndian); 096 readBytes(in, b8); 097 readBytes(in, b8); 098 maxLat = NTV2Util.getDouble(b8, bigEndian); 099 readBytes(in, b8); 100 readBytes(in, b8); 101 minLon = NTV2Util.getDouble(b8, bigEndian); 102 readBytes(in, b8); 103 readBytes(in, b8); 104 maxLon = NTV2Util.getDouble(b8, bigEndian); 105 readBytes(in, b8); 106 readBytes(in, b8); 107 latInterval = NTV2Util.getDouble(b8, bigEndian); 108 readBytes(in, b8); 109 readBytes(in, b8); 110 lonInterval = NTV2Util.getDouble(b8, bigEndian); 111 lonColumnCount = 1 + (int) ((maxLon - minLon) / lonInterval); 112 latRowCount = 1 + (int) ((maxLat - minLat) / latInterval); 113 readBytes(in, b8); 114 readBytes(in, b8); 115 nodeCount = NTV2Util.getInt(b8, bigEndian); 116 if (nodeCount != lonColumnCount * latRowCount) 117 throw new IllegalStateException("SubGrid " + subGridName + " has inconsistent grid dimesions"); 118 latShift = new float[nodeCount]; 119 lonShift = new float[nodeCount]; 120 if (loadAccuracy) { 121 latAccuracy = new float[nodeCount]; 122 lonAccuracy = new float[nodeCount]; 123 } 124 125 for (int i = 0; i < nodeCount; i++) { 126 // Read the grid file byte after byte. This is a workaround about a bug in 127 // certain VM which are not able to read byte blocks when the resource file is in a .jar file (Pieren) 128 readBytes(in, b1); b4[0] = b1[0]; 129 readBytes(in, b1); b4[1] = b1[0]; 130 readBytes(in, b1); b4[2] = b1[0]; 131 readBytes(in, b1); b4[3] = b1[0]; 132 latShift[i] = NTV2Util.getFloat(b4, bigEndian); 133 readBytes(in, b1); b4[0] = b1[0]; 134 readBytes(in, b1); b4[1] = b1[0]; 135 readBytes(in, b1); b4[2] = b1[0]; 136 readBytes(in, b1); b4[3] = b1[0]; 137 lonShift[i] = NTV2Util.getFloat(b4, bigEndian); 138 readBytes(in, b1); b4[0] = b1[0]; 139 readBytes(in, b1); b4[1] = b1[0]; 140 readBytes(in, b1); b4[2] = b1[0]; 141 readBytes(in, b1); b4[3] = b1[0]; 142 if (loadAccuracy) { 143 latAccuracy[i] = NTV2Util.getFloat(b4, bigEndian); 144 } 145 readBytes(in, b1); b4[0] = b1[0]; 146 readBytes(in, b1); b4[1] = b1[0]; 147 readBytes(in, b1); b4[2] = b1[0]; 148 readBytes(in, b1); b4[3] = b1[0]; 149 if (loadAccuracy) { 150 lonAccuracy[i] = NTV2Util.getFloat(b4, bigEndian); 151 } 152 } 153 } 154 155 private static void readBytes(InputStream in, byte[] b) throws IOException { 156 if (in.read(b) < b.length) { 157 Logging.error("Failed to read expected amount of bytes ("+ b.length +") from stream"); 158 } 159 } 160 161 /** 162 * Tests if a specified coordinate is within this Sub Grid 163 * or one of its Sub Grids. If the coordinate is outside 164 * this Sub Grid, null is returned. If the coordinate is 165 * within this Sub Grid, but not within any of its Sub Grids, 166 * this Sub Grid is returned. If the coordinate is within 167 * one of this Sub Grid's Sub Grids, the method is called 168 * recursively on the child Sub Grid. 169 * 170 * @param lon Longitude in Positive West Seconds 171 * @param lat Latitude in Seconds 172 * @return the Sub Grid containing the Coordinate or null 173 */ 174 public NTV2SubGrid getSubGridForCoord(double lon, double lat) { 175 return !isCoordWithin(lon, lat) 176 ? null 177 : subGrid == null 178 ? this 179 : Arrays.stream(subGrid) 180 .filter(aSubGrid -> aSubGrid.isCoordWithin(lon, lat)) 181 .map(aSubGrid -> aSubGrid.getSubGridForCoord(lon, lat)) 182 .filter(Objects::nonNull) 183 .findFirst().orElse(this); 184 } 185 186 /** 187 * Tests if a specified coordinate is within this Sub Grid. 188 * A coordinate on either outer edge (maximum Latitude or 189 * maximum Longitude) is deemed to be outside the grid. 190 * 191 * @param lon Longitude in Positive West Seconds 192 * @param lat Latitude in Seconds 193 * @return true or false 194 */ 195 private boolean isCoordWithin(double lon, double lat) { 196 return (lon >= minLon) && (lon < maxLon) && (lat >= minLat) && (lat < maxLat); 197 } 198 199 /** 200 * Bi-Linear interpolation of four nearest node values as described in 201 * 'GDAit Software Architecture Manual' produced by the <a 202 * href='http://www.dtpli.vic.gov.au/property-and-land-titles/geodesy/geocentric-datum-of-australia-1994-gda94/gda94-useful-tools'> 203 * Geomatics Department of the University of Melbourne</a> 204 * @param a value at the A node 205 * @param b value at the B node 206 * @param c value at the C node 207 * @param d value at the D node 208 * @param x Longitude factor 209 * @param y Latitude factor 210 * @return interpolated value 211 */ 212 private static double interpolate(float a, float b, float c, float d, double x, double y) { 213 return a + (((double) b - (double) a) * x) + (((double) c - (double) a) * y) + 214 (((double) a + (double) d - b - c) * x * y); 215 } 216 217 /** 218 * Interpolate shift and accuracy values for a coordinate in the 'from' datum 219 * of the GridShiftFile. The algorithm is described in 220 * 'GDAit Software Architecture Manual' produced by the <a 221 * href='http://www.dtpli.vic.gov.au/property-and-land-titles/geodesy/geocentric-datum-of-australia-1994-gda94/gda94-useful-tools'> 222 * Geomatics Department of the University of Melbourne</a> 223 * <p>This method is thread safe for both memory based and file based node data. 224 * @param gs GridShift object containing the coordinate to shift and the shift values 225 */ 226 public void interpolateGridShift(NTV2GridShift gs) { 227 int lonIndex = (int) ((gs.getLonPositiveWestSeconds() - minLon) / lonInterval); 228 int latIndex = (int) ((gs.getLatSeconds() - minLat) / latInterval); 229 230 double x = (gs.getLonPositiveWestSeconds() - (minLon + (lonInterval * lonIndex))) / lonInterval; 231 double y = (gs.getLatSeconds() - (minLat + (latInterval * latIndex))) / latInterval; 232 233 // Find the nodes at the four corners of the cell 234 235 int indexA = lonIndex + (latIndex * lonColumnCount); 236 int indexB = indexA + 1; 237 int indexC = indexA + lonColumnCount; 238 int indexD = indexC + 1; 239 240 gs.setLonShiftPositiveWestSeconds(interpolate( 241 lonShift[indexA], lonShift[indexB], lonShift[indexC], lonShift[indexD], x, y)); 242 243 gs.setLatShiftSeconds(interpolate( 244 latShift[indexA], latShift[indexB], latShift[indexC], latShift[indexD], x, y)); 245 246 if (lonAccuracy == null) { 247 gs.setLonAccuracyAvailable(false); 248 } else { 249 gs.setLonAccuracyAvailable(true); 250 gs.setLonAccuracySeconds(interpolate( 251 lonAccuracy[indexA], lonAccuracy[indexB], lonAccuracy[indexC], lonAccuracy[indexD], x, y)); 252 } 253 254 if (latAccuracy == null) { 255 gs.setLatAccuracyAvailable(false); 256 } else { 257 gs.setLatAccuracyAvailable(true); 258 gs.setLatAccuracySeconds(interpolate( 259 latAccuracy[indexA], latAccuracy[indexB], latAccuracy[indexC], latAccuracy[indexD], x, y)); 260 } 261 } 262 263 /** 264 * Returns the parent sub grid name. 265 * @return the parent sub grid name 266 */ 267 public String getParentSubGridName() { 268 return parentSubGridName; 269 } 270 271 /** 272 * Returns the sub grid name. 273 * @return the sub grid name 274 */ 275 public String getSubGridName() { 276 return subGridName; 277 } 278 279 /** 280 * Returns the node count. 281 * @return the node count 282 */ 283 public int getNodeCount() { 284 return nodeCount; 285 } 286 287 /** 288 * Returns the sub grid count. 289 * @return the sub grid count 290 */ 291 public int getSubGridCount() { 292 return subGrid == null ? 0 : subGrid.length; 293 } 294 295 /** 296 * Set an array of Sub Grids of this sub grid 297 * @param subGrid subgrids 298 */ 299 public void setSubGridArray(NTV2SubGrid... subGrid) { 300 this.subGrid = Utils.copyArray(subGrid); 301 } 302 303 @Override 304 public String toString() { 305 return subGridName; 306 } 307 308 /** 309 * Returns textual details about the sub grid. 310 * @return textual details about the sub grid 311 */ 312 public String getDetails() { 313 return new StringBuilder(256) 314 .append("Sub Grid : ") 315 .append(subGridName) 316 .append("\nParent : ") 317 .append(parentSubGridName) 318 .append("\nCreated : ") 319 .append(created) 320 .append("\nUpdated : ") 321 .append(updated) 322 .append("\nMin Lat : ") 323 .append(minLat) 324 .append("\nMax Lat : ") 325 .append(maxLat) 326 .append("\nMin Lon : ") 327 .append(minLon) 328 .append("\nMax Lon : ") 329 .append(maxLon) 330 .append("\nLat Intvl: ") 331 .append(latInterval) 332 .append("\nLon Intvl: ") 333 .append(lonInterval) 334 .append("\nNode Cnt : ") 335 .append(nodeCount) 336 .toString(); 337 } 338 339 /** 340 * Get maximum latitude value 341 * @return maximum latitude 342 */ 343 public double getMaxLat() { 344 return maxLat; 345 } 346 347 /** 348 * Get maximum longitude value 349 * @return maximum longitude 350 */ 351 public double getMaxLon() { 352 return maxLon; 353 } 354 355 /** 356 * Get minimum latitude value 357 * @return minimum latitude 358 */ 359 public double getMinLat() { 360 return minLat; 361 } 362 363 /** 364 * Get minimum longitude value 365 * @return minimum longitude 366 */ 367 public double getMinLon() { 368 return minLon; 369 } 370}