dted reader implementation in Java

DtedReader.java

package whee;

import java.io.Closeable;
import java.io.IOException;
import java.nio.ByteBuffer;
import java.nio.ByteOrder;
import java.nio.ShortBuffer;
import java.nio.channels.FileChannel;
import java.nio.file.Files;
import java.nio.file.Path;
import java.nio.file.Paths;
import java.nio.file.StandardOpenOption;
import java.util.ArrayList;
import java.util.List;
import java.util.Objects;
import java.util.stream.IntStream;
import java.util.stream.Stream;

import org.locationtech.jts.geom.Coordinate;
import org.locationtech.jts.geom.Envelope;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;

import whee.util.Geometries;

/**
 * Barebones DTED reader class.
 * <p>
 * References:
 * <ul>
 * <li>MIL-PRF-89020B (https://dds.cr.usgs.gov/srtm/version2_1/Documentation/MIL-PDF-89020B.pdf)</li>
 * <li>NASA WorldWind (https://github.com/nasa/World-Wind-Java/blob/93ac47f196692a674eab3c8cc33f824062f7aa59/WorldWind/src/gov/nasa/worldwind/formats/dted/DTED.java)</li>
 * </ul>
 */
public final class DtedReader implements Closeable {
    private static final int SIZE_UHL = 80;
    private static final int SIZE_DSI = 648;
    private static final int SIZE_ACC = 2700;
    private static final int SIZE_RECORD_HEADER = 8;
    private static final int SIZE_RECORD_CHECKSUM = 4;

    private static final int OFFSET_UHL = 0;
    private static final int OFFSET_DSI = OFFSET_UHL + SIZE_UHL;
    private static final int OFFSET_ACC = OFFSET_DSI + SIZE_DSI;
    private static final int OFFSET_DATA = OFFSET_ACC + SIZE_ACC;

    private static final short NODATA = -32767;

    private static final Logger LOG = LoggerFactory.getLogger(DtedReader.class);

    private final FileChannel channel;
    private final Path path;
    private final int level;
    private final double originX;
    private final double originY;
    private final double pixelWidth;
    private final double pixelHeight;
    private final int xCount;
    private final int yCount;
    private final int recordSize;

    public DtedReader(String path) throws IOException {
        this(Paths.get(path));
    }

    public DtedReader(Path path) throws IOException {
        this.path = path;

        LOG.debug("Open DTED file (path={}, size={})",
                this.path,
                String.format("%,d", Files.size(this.path)));

        this.channel = FileChannel.open(this.path, StandardOpenOption.READ);
        try {
            // Read UHL header
            final byte[] uhlHeader = readHeader(channel, "UHL", OFFSET_UHL, SIZE_UHL);
            this.xCount = parseInt("xCount", uhlHeader, 47, 4);
            this.yCount = parseInt("yCount", uhlHeader, 51, 4);
            this.pixelWidth = 1D / (xCount - 1);
            this.pixelHeight = 1D / (yCount - 1);
            this.originX = parseDms("longitude", uhlHeader, 4) - pixelWidth / 2;   // Account for pixel overlap
            this.originY = parseDms("latitude", uhlHeader, 12) + pixelHeight / 2;  // Account for pixel overlap

            // Read DSI header
            final byte[] dsiHeader = readHeader(channel, "DSI", OFFSET_DSI, SIZE_DSI);
            this.level = parseLevel(dsiHeader);

            // Misc
            this.recordSize = SIZE_RECORD_HEADER + xCount * Short.SIZE / Byte.SIZE + SIZE_RECORD_CHECKSUM;

            // Dump diagnostic info
            LOG.debug("DTED Properties:\n" +
                      "---\n" +
                      "        Path : {}\n" +
                      "       Level : {}\n" +
                      "      Origin : ({}, {})\n" +
                      "  Pixel Size : {} x {}\n" +
                      " Pixel Count : {} (estimated)\n" +
                      " Record Size : {} bytes\n" +
                      "---",
                    path,
                    level,
                    originX,
                    originY,
                    pixelWidth,
                    pixelHeight,
                    String.format("%,d", xCount * yCount),
                    recordSize);
        }
        catch (Exception e) {
            // Prevent leaking open file descriptors if init fails
            this.channel.close();

            // Rethrow exception
            throw e;
        }
    }

    @Override
    public void close() throws IOException {
        LOG.debug("Closing DTED file (file={})", path);
        channel.close();
    }

    public Envelope getExtent() {
        return Geometries.envelope(
                originX,
                originY,
                originX + (pixelWidth * xCount),
                originY + (pixelHeight * yCount));
    }

    public int getLevel() {
        return level;
    }

    public double getOriginX() {
        return originX;
    }

    public double getOriginY() {
        return originY;
    }

    public double getPixelHeight() {
        return pixelHeight;
    }

    public double getPixelWidth() {
        return pixelWidth;
    }

    public int getXCount() {
        return xCount;
    }

    public int getYCount() {
        return yCount;
    }

    public List<Coordinate> read() throws IOException {
        channel.position(OFFSET_DATA);

        final int recordSize = SIZE_RECORD_HEADER + xCount * Short.SIZE / Byte.SIZE + SIZE_RECORD_CHECKSUM;

        final List<Coordinate> coordinates = new ArrayList<>();

        final ByteBuffer byteBuffer = ByteBuffer.allocate(recordSize).order(ByteOrder.BIG_ENDIAN);

        /*
         * Pixel Iteration
         *
         * According to MIL-PRF-89020B, cell origin is the pixel in the
         * Southwestern corner.  The iteration order here differs from the way
         * GDAL iterates.
         *
         * Assuming a 3x3 pixel grid, the iteration order is as follows:
         *
         *     3 6 9
         *     2 5 8
         *     1 4 7
         */
        for (int xOffset = 0; xOffset < xCount; xOffset++) {
            channel.read(byteBuffer);
            byteBuffer.flip();

            final ShortBuffer shortBuffer = byteBuffer.asShortBuffer();
            for (int yOffset = yCount - 1; yOffset >= 0; yOffset--) {
                double x = originX + pixelWidth * xOffset;
                double y = originY + pixelHeight * yOffset;
                short z = shortBuffer.get(yOffset + 4);  // skip 4 shorts of header

                // Discard NULL data
                if (z == NODATA) {
                    continue;
                }

                coordinates.add(new Coordinate(x, y, z));
            }
        }

        return coordinates;
    }

    public Stream<Coordinate> stream() {
        final ByteBuffer byteBuffer = ByteBuffer.allocate(recordSize).order(ByteOrder.BIG_ENDIAN);

        /*
         * Pixel Iteration
         *
         * According to MIL-PRF-89020B, cell origin is the pixel in the
         * Southwestern corner.  The iteration order here differs from
         * the way GDAL iterates.
         *
         * Assuming a 3x3 pixel grid, the iteration order is as follows:
         *
         *     3 6 9
         *     2 5 8
         *     1 4 7
         */
        return IntStream.range(0, xCount * yCount)
                .mapToObj(i -> {
                    final int xOffset = i / xCount;
                    final int yOffset = yCount - (i % yCount) - 1;
                    final int recordOffset = xOffset * recordSize;

                    // Read new data into the buffer whenever xOffset changes
                    if (i % xCount == 0) {
                        try {
                            channel.read(byteBuffer, OFFSET_DATA + recordOffset);
                            byteBuffer.flip();
                        }
                        catch (IOException e) {
                            throw new RuntimeException(e);
                        }
                    }

                    double x = originX + pixelWidth * xOffset;
                    double y = originY + pixelHeight * yOffset;
                    short z = byteBuffer.asShortBuffer().get(yOffset + 4);  // skip 4 shorts of header

                    // Discard NULL data
                    if (z == NODATA) {
                        return null;
                    }

                    return new Coordinate(x, y, z);
                })
                .filter(Objects::nonNull);
    }

    //
    // Helpers
    //

    private static double parseDouble(String context, byte[] bytes, int offset, int length) throws IOException {
        final String value = new String(bytes, offset, length);
        try {
            return Double.parseDouble(value);
        }
        catch (NumberFormatException e) {
            throw new IOException("unreadable double value for '" + context + "' (value='" + value + "')", e);
        }
    }

    private static double parseDms(String context, byte[] bytes, int offset) throws IOException {
        final String raw = new String(bytes, offset, 8).trim().toUpperCase();

        if (LOG.isTraceEnabled()) {
            LOG.trace("parseDms: \033[34m{}\033[0m", raw);
        }

        // Validate format
        if (!raw.matches("^\\d{3}\\d{2}\\d{2}[NSEW]$")) {
            throw new IOException("UHL contains malformed " + context + " (value='" + raw + "')");
        }

        // Extract
        double value = parseDouble(context + " degree", raw.getBytes(), 0, 3);    // Degrees
        value += parseDouble(context + " minutes", raw.getBytes(), 3, 2) / 60;    // Minutes
        value += parseDouble(context + " seconds", raw.getBytes(), 5, 2) / 3600;  // Seconds
        value *= (raw.endsWith("W") || raw.endsWith("S")) ? -1 : 1;               // Hemisphere

        // Validate value
        if ((raw.endsWith("E") || raw.endsWith("W")) && (value < -180 || value > 180)) {
            throw new IOException("out of bounds longitude '" + raw + "'");
        }
        else if ((raw.endsWith("N") || raw.endsWith("S")) && (value < -90 || value > 90)) {
            throw new IOException("out of bounds latitude '" + raw + "'");
        }

        return value;
    }

    private static int parseInt(String context, byte[] bytes, int offset, int length) throws IOException {
        final String value = new String(bytes, offset, length);
        try {
            return Integer.parseInt(value);
        }
        catch (NumberFormatException e) {
            throw new IOException("unreadable integer value for '" + context + "' (value='" + value + "')", e);
        }
    }

    private static int parseLevel(byte[] dsiHeader) throws IOException {
        final String raw = new String(dsiHeader, 59, 5);

        // Validate format
        if (!raw.matches("^DTED\\d$")) {
            throw new IOException("DSI contains unknown DTED level '" + raw + "'");
        }

        return parseInt("dted level", dsiHeader, 63, 1);
    }

    private static byte[] readHeader(FileChannel channel, String id, int offset, int size) throws IOException {
        final byte[] value = new byte[size];

        channel.position(offset);
        channel.read(ByteBuffer.wrap(value));

        // Validate basic format
        if (!id.equals(new String(value, 0, id.length()))) {
            throw new IOException("malformed " + id + " header");
        }

        if (LOG.isTraceEnabled()) {
            LOG.trace("readHeader: \033[34m{}\033[0m", new String(value));
        }

        return value;
    }
}