/* This program is free software: you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public License
as published by the Free Software Foundation, either version 3 of
the License, or (at your option) any later version.
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with this program. If not, see <http://www.gnu.org/licenses/>. */
package org.opentripplanner.graph_builder.impl.ned;
import java.io.BufferedInputStream;
import java.io.DataInputStream;
import java.io.IOException;
import java.io.InputStream;
/**
* A Vertical datum specified as a grid of offsets from NAD83
* http://vdatum.noaa.gov/dev/gtx_info.html
*
* @author novalis
*
*/
public class VerticalDatum {
double lowerLeftLatitude;
double lowerLeftLongitude;
double deltaLatitude;
double deltaLongitude;
float[][] datum;
public VerticalDatum(double lowerLeftLongitude, double lowerLeftLatitude, double width,
double height, float[][] datum) {
this.lowerLeftLongitude = lowerLeftLongitude;
this.lowerLeftLatitude = lowerLeftLatitude;
this.deltaLongitude = width;
this.deltaLatitude = height;
this.datum = datum;
}
double interpolatedHeight(double longitude, double latitude) {
// because VDatums can cross -180 longitude, there's some complication in interpolating
// longitude
double lowerLeftAdjusted = lowerLeftLongitude;
if (lowerLeftLongitude + deltaLongitude > 180) {
// then lowerLeftLongitude must be left of 180
if (longitude < lowerLeftLongitude) {
// lowerLeftLongitude is right of 180 (probablY)
lowerLeftAdjusted -= 360; // this transforms lowerLeft so that interpolation will
// work
if (longitude < lowerLeftLongitude) {
throw new RuntimeException("longitude out of range");
}
}
}
if (longitude > lowerLeftAdjusted + deltaLongitude || longitude < lowerLeftAdjusted) {
throw new RuntimeException("longitude out of range");
}
if (latitude < lowerLeftLatitude || latitude > lowerLeftLatitude + deltaLatitude) {
throw new RuntimeException("latitude out of range");
}
double longitudeNormalized = (longitude - lowerLeftAdjusted) / deltaLongitude;
double latitudeNormalized = (latitude - lowerLeftLatitude) / deltaLatitude;
int rows = datum.length;
int columns = datum[0].length;
int x1 = (int) Math.floor(longitudeNormalized * columns);
int y1 = (int) Math.floor(latitudeNormalized * rows);
double gridXFraction = longitudeNormalized * columns - x1;
double gridYFraction = latitudeNormalized * rows - y1;
return datum[y1][x1] * gridXFraction * gridYFraction + datum[y1][x1 + 1]
* (1 - gridXFraction) * gridYFraction + datum[y1 + 1][x1] * gridXFraction
* (1 - gridYFraction) + datum[y1 + 1][x1 + 1] * (1 - gridXFraction)
* (1 - gridYFraction);
}
boolean covers(double longitude, double latitude) {
double lowerLeftAdjusted = lowerLeftLongitude;
if (lowerLeftLongitude + deltaLongitude > 180) {
// then lowerLeftLongitude must be left of 180
if (longitude < lowerLeftLongitude) {
// lowerLeftLongitude is right of 180 (probablY)
lowerLeftAdjusted -= 360; // this transforms lowerLeft so that interpolation will
// work
if (longitude < lowerLeftLongitude) {
return false;
}
}
}
if (longitude > lowerLeftAdjusted + deltaLongitude || longitude < lowerLeftAdjusted) {
return false;
}
if (latitude < lowerLeftLatitude || latitude > lowerLeftLatitude + deltaLatitude) {
return false;
}
return true;
}
public static VerticalDatum fromGTX (InputStream inputStream) throws IOException {
DataInputStream stream = new DataInputStream(new BufferedInputStream(inputStream));
double lowerLeftLatitude = stream.readDouble();
double lowerLeftLongitude = stream.readDouble();
if (lowerLeftLongitude > 180) {
lowerLeftLongitude -= 360; //convert to standard coordinates
}
double deltaLatitude = stream.readDouble();
double deltaLongitude = stream.readDouble();
int nRows = stream.readInt();
int nColumns = stream.readInt();
float[][] data = new float[nRows][nColumns];
for (int y = 0; y < nRows; ++y) {
for (int x = 0; x < nColumns; ++x) {
data[y][x] = stream.readFloat();
}
}
return new VerticalDatum(lowerLeftLongitude, lowerLeftLatitude, deltaLatitude * nColumns, deltaLongitude * nRows, data);
}
}