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