V1_1_1.java [src/java/m/rse/wepot] Revision: Date:
package m.rse.wepot;
import csip.*;
import csip.annotations.*;
import static csip.annotations.ResourceType.*;
import oms3.annotations.Description;
import oms3.annotations.Name;
import java.io.*;
import java.sql.*;
import java.util.*;
import java.util.logging.Level;
import javax.ws.rs.Path;
import org.codehaus.jettison.json.JSONArray;
import org.codehaus.jettison.json.JSONException;
import org.codehaus.jettison.json.JSONObject;
import m.rse.cfactor.utils.Const;
/**
*
* @author Shaun Case
* @author wl
* @author od
*/
//@Name("wepot")
@Description("This service computes water and wind erodibility potentials for "
+ "an area of analysis (AoA). The service clips SSURGO soil mapunits "
+ "with AoA geometry, determines the dominant soil component in the "
+ "AoA, gets parameters from the SSURGO component table, including a "
+ "climate factor from the C Factor layer, and computes the following "
+ "equations: Wind Erosion Potential = C*I/T ; "
+ "Water Erosion Potential = K*(LS)/T")
//@Path("m/wepot/1.1.1")
//@Polling(first = 5000, next = 2000)
// gdal
@Resource(type = REFERENCE, file = "gdallocationinfo", id = "gdal-lin")
@Resource(type = EXECUTABLE, file = "/bin/win/gdallocationinfo.exe", id = "gdal-win")
@Resource(type = ARCHIVE, file = "/bin/win/gdal_dlls.zip")
@Resource(type = FILE, file = "/data/us_cvalues_topo2ras_masked.tfw")
@Resource(type = FILE, file = "/data/us_cvalues_topo2ras_masked.tif", id = "tif")
@Resource(type = FILE, file = "/data/us_cvalues_topo2ras_masked.tif.aux.xml")
// db
@Resource(type = JDBC, file = "${cfactor.db}", id = "cfactor", env = {
"removeAbandonedTimeout=10",
"jdbcInterceptors=org.apache.tomcat.jdbc.pool.interceptor.ConnectionState;"
+ "org.apache.tomcat.jdbc.pool.interceptor.StatementFinalizer;"
+ "org.apache.tomcat.jdbc.pool.interceptor.ResetAbandonedTimer"
})
public class V1_1_1 extends ModelDataService {
static final int JSON_LATITUDE = 1;
static final int JSON_LONGITUDE = 0;
static final String RSE_AOAID = "AoAId";
static final String AOA_GEOMETRY = "aoa_geometry";
static final double[][] LightleWeesiesSlopeLength = {
{0, .75, 100},
{.75, 1.5, 200},
{1.5, 2.5, 300},
{2.5, 3.5, 200},
{3.5, 4.5, 180},
{4.5, 5.5, 160},
{5.5, 6.5, 150},
{6.5, 7.5, 140},
{7.5, 8.5, 130},
{8.5, 9.5, 125},
{9.5, 10.5, 120},
{10.5, 11.5, 110},
{11.5, 12.5, 100},
{12.5, 13.5, 90},
{13.5, 14.5, 80},
{14.5, 15.5, 70},
{15.5, 17.5, 60},
{17.5, -1, 50}
};
static final double[][] PalouseSlopeLength = {
{0, 5.5, 350},
{5.5, 10.5, 275},
{10.5, 15.5, 225},
{15.5, 20.5, 175},
{20.5, 25.5, 150},
{25.5, 35.5, 125},
{35.5, -1, 100}
};
private String AoAId = "Error";
private JSONObject aoaGeometry;
private PolygonLatLon polygonData;
private String error_msg = "";
//Result variables
private double aoa_cfactor = Const.UNKNOWN_CFACTOR;
private String aoa_dom_comp = "";
private String aoa_dom_compname = "";
private String aoa_dom_wind_comp = "";
private String aoa_dom_wind_compname = "";
private double aoa_ifactor = Const.UNKNOWN_CFACTOR;
private double aoa_tfactor = Const.UNKNOWN_CFACTOR;
private double aoa_wind_tfactor = Const.UNKNOWN_CFACTOR;
private double aoa_lsfactor = Const.UNKNOWN_CFACTOR;
private double aoa_kfactor = Const.UNKNOWN_CFACTOR;
private double aoa_wind_ep = Const.UNKNOWN_CFACTOR;
private double aoa_water_ep = Const.UNKNOWN_CFACTOR;
private boolean isPalouse = false;
@Override
protected void preProcess() throws ServiceException {
//These functions will throw a usuable user message if they cannot find these values.
// Just catch it and return it to the user.
AoAId = getStringParam(RSE_AOAID);
aoaGeometry = getJSONParam(AOA_GEOMETRY);
}
@Override
protected void doProcess() throws Exception {
try {
JSONArray features = aoaGeometry.optJSONArray("features");
if (buildPolygon(features)) {
polygonData.makeESRIRotation(PolygonLatLon.ESRIContourType.outer);
//Got a valid Polygon now.
//HINT: You can get the WKT from the PolygonLatLon object to build an SQL query via the "toWKT" funciton,
// and then encasing that WKT in the typical "geometry::STPolyFromText(<KWT goes here>,4326)", etc.
LOG.info("THE POLYGON AS WE KNOW IT IS=" + polygonData.toWKT());
String cfactor_polygon = "geometry::STPolyFromText(" + polygonData.toWKT() + ",4326)";
LOG.info("getting centroid for polygon now");
Centroid centroid = getCentroid(cfactor_polygon);
LOG.info("About to get C Factor Raster value for lat=" + centroid.lat + " long=" + centroid.lon + " ");
aoa_cfactor = getCFactorRaster(centroid.lat, centroid.lon) / 100.0;
ArrayList<SsurgoAttributes> intersectionAttributes = ssurgoIntersect(polygonData.toWKT());
isPalouse = isPalouseArea(polygonData.toWKT());
if (null != intersectionAttributes) {
if (windCriticalComponentLookup(intersectionAttributes)) {
aoa_wind_ep = aoa_cfactor * aoa_ifactor / aoa_wind_tfactor;
if (waterCriticalComponentLookup(intersectionAttributes)) {
aoa_water_ep = aoa_kfactor * aoa_lsfactor / aoa_tfactor;
} else {
throw new ServiceException(" Could not calculate the water erodibility potenial. Missing component information. ");
}
} else {
error_msg += " Could not calculate the wind erodibility potenial. Missing component information. ";
}
}
}
} catch (JSONException ex) {
error_msg += "Cannot proceed with processing the request: " + ex.getMessage();
throw new ServiceException("Cannot proceed with processing the request", ex);
} catch (SQLException se) {
error_msg += "SQL exception getting cfactor: " + se.getMessage();
} catch (ServiceException sx) {
error_msg += "Service exception getting cfactor: " + sx.getMessage();
} catch (IOException ioe) {
error_msg += "IO Exception getting cFactor value from raster geotiff layer: " + ioe.getMessage();
}
}
@Override
protected void postProcess() throws Exception {
putResult("aoa_id", AoAId);
putResult("aoa_dom_water_comp", aoa_dom_comp);
putResult("aoa_dom_water_compname", aoa_dom_compname);
putResult("aoa_dom_wind_comp", aoa_dom_wind_comp);
putResult("aoa_dom_wind_compname", aoa_dom_wind_compname);
putResult("aoa_cfactor", aoa_cfactor);
putResult("aoa_ifactor", aoa_ifactor);
putResult("aoa_tfactor", aoa_tfactor);
putResult("aoa_wind_tfactor", aoa_wind_tfactor);
putResult("aoa_lsfactor", aoa_lsfactor);
putResult("aoa_kfactor", aoa_kfactor);
putResult("aoa_wind_ep", aoa_wind_ep);
putResult("aoa_water_ep", aoa_water_ep);
}
//Private functions
private boolean buildPolygon(JSONArray features) throws JSONException {
boolean ret_val = false;
if (null != features) {
if (features.length() > 0) {
JSONObject feature = features.getJSONObject(0).optJSONObject("geometry");
if (isGeometryPolygonType(feature)) {
JSONArray polygon = feature.optJSONArray("coordinates");
if (null != polygon) {
if (polygon.length() > 0) {
polygonData = readPolygonCoordinates(polygon.getJSONArray(0));
if (null != polygonData) {
if (!polygonData.isValid()) {
error_msg += "Invalid latitude and/or longitude data contained in this polygon. Cannot proceed with processing of this request. ";
} else {
ret_val = true;
}
}//No else needed here, the error message will be built in readPolygonCoordinates()
} else {
error_msg = "No coordinates found associated with the polygon specified in feature collection number: 1 . ";
}
} else {
error_msg = "No coordinates found associated with the polygon specified in feature collection number: 1 . ";
}
}//No else needed here, the error message will be built in isGeometryPolygonType()
} else {
error_msg = "No geometry found associated with this feature. ";
}
} else {
error_msg = "Cannot process request JSON, missing features. ";
}
return ret_val;
}
private boolean isGeometryPolygonType(JSONObject geometry) {
boolean ret_val = false;
if (null != geometry) {
String geometryType = geometry.optString("type");
if (!geometryType.isEmpty()) {
if (geometryType.toLowerCase().equals("polygon")) {
ret_val = true;
} else {
error_msg = "No valid geometry type found in the feature collection number: 1 . Looking for 'Polygon'. ";
}
} else {
error_msg = "No geometry type specified in the feature collection number: 1 . ";
}
} else {
error_msg = "No geometry found in the feature collection number: 1 . ";
}
return ret_val;
}
private PolygonLatLon readPolygonCoordinates(JSONArray shape) throws JSONException {
PolygonLatLon tPolygon = null;
if (null != shape) {
tPolygon = new PolygonLatLon();
for (int k = 0; k < shape.length(); k++) {
JSONArray jPoint = shape.getJSONArray(k);
tPolygon.add(jPoint.getDouble(JSON_LATITUDE), jPoint.getDouble(JSON_LONGITUDE));
}
} else {
error_msg = "Cannot find the latitude and longitude values for this polygon. ";
}
return tPolygon;
}
// If the WKTPolygon of interest intersects any polygons in the Palouse, MLRA-9, region, this funciton will return true.
private boolean isPalouseArea(String wktPolygon) {
boolean ret_val = false;
if (null != wktPolygon) {
String query = "SELECT id"
+ " FROM ssurgo.mlra_9 "
+ " WHERE geometry::STGeomFromText(component_boundary.STAsText(),4326).STIntersects(geometry::STGeomFromText(" + wktPolygon + ",4326))=1"
+ " AND geometry::STGeomFromText(" + wktPolygon + ",4326).STIsValid()=1;";
LOG.log(Level.INFO, "ssurgo Palouse query={0}", query);
try (Connection conn = getResourceJDBC("cfactor");
Statement statement = conn.createStatement();) {
try (ResultSet results = statement.executeQuery(query)) {
LOG.info(" past the query execute");
if (null != results) {
if (results.next()) {
ret_val = true;
}
}
}
} catch (SQLException | ServiceException ex) {
error_msg += "Cannot continue processing this request in the find Palouse region procedure: " + ex.getMessage();
LOG.log(Level.SEVERE, "error", ex);
}
}
return ret_val;
}
// Returns the intersected mukeys
private ArrayList<SsurgoAttributes> ssurgoIntersect(String wktPolygon) {
ArrayList<SsurgoAttributes> ret_val = null;
if (null != wktPolygon) {
ArrayList<SsurgoAttributes> attributeList = new ArrayList<>();
HashMap<String, SsurgoAttributes> mukeyAttributes = new HashMap<>();
String query = "SELECT areasymbol, musym, mukey, geography::STGeomFromText(intersectPoly.STAsText(), 4326).MakeValid().STArea() / 4046.86 as sizeIntersectionAcres, "
+ " geography::STGeomFromText(" + wktPolygon + ", 4326).MakeValid().STArea() / 4046.86 as aoaArea"
+ " FROM"
+ " (SELECT m.areasymbol, m.musym, m.mukey, "
+ " m.the_geom.STIntersection(geometry:: STGeomFromText (" + wktPolygon + ",0)).MakeValid() as intersectPoly"
+ " FROM ssurgo.soilmu_a as m "
+ " WITH (index(geom_sidx))"
+ " WHERE m.the_geom.STIntersects(geometry:: STGeomFromText (" + wktPolygon + ",0)) = 1"
+ " AND m.the_geom.STIsValid()=1"
+ " AND (geometry:: STGeomFromText (" + wktPolygon + ",0).STIsValid())=1)"
+ " as a;";
LOG.info("ssurgo query=" + query);
try (Connection conn = getResourceJDBC("cfactor");
Statement statement = conn.createStatement()) {
ResultSet results = statement.executeQuery(query);
LOG.info(" past the query execute");
if (null != results) {
while (results.next()) {
String tKey = results.getString("areasymbol") + ":" + results.getString("mukey");
SsurgoAttributes tAttributes = mukeyAttributes.get(tKey);
//Merge, i.e. sum, areas of duplicate AoA:mukey combinations
if (null != tAttributes) {
tAttributes.gid_area += results.getDouble("sizeIntersectionAcres");
} else {
tAttributes = new SsurgoAttributes();
tAttributes.gid = results.getString("musym");
tAttributes.AoA_Id = results.getString("areasymbol");
tAttributes.mukey = results.getString("mukey");
tAttributes.gid_area = results.getDouble("sizeIntersectionAcres");
tAttributes.aoa_area = results.getDouble("aoaArea");
mukeyAttributes.put(tKey, tAttributes);
attributeList.add(tAttributes);
}
}
ret_val = attributeList;
} else {
error_msg += "No results from the intersect query for this geometry. ";
}
} catch (SQLException | csip.ServiceException ex) {
error_msg += "Cannot continue processing this request in the intersect procedure: " + ex.getMessage();
LOG.log(Level.SEVERE, "error", ex);
ret_val = null;
}
}
return ret_val;
}
private boolean windCriticalComponentLookup(ArrayList<SsurgoAttributes> intersectionAttributeList) {
boolean ret_val = true;
String dom_comp = "";
String dom_compname = "";
double ifact = 0.0;
double tfact = 0.0;
double sandtotal = 0.0;
double area = 0.0;
for (SsurgoAttributes attributes : intersectionAttributeList) {
if (attributes.gid_area >= (attributes.aoa_area * 0.10)) {
componentLookup(attributes, false);
if (attributes.horizons.size() > 0) {
for (ComponentAttributes tComponent : attributes.horizons) {
if ((tComponent.sand_pct >= sandtotal) && (tComponent.sand_pct > 0)) {
if (sandtotal == tComponent.sand_pct) {
if (tComponent.area_pct > area) {
sandtotal = tComponent.sand_pct;
ifact = tComponent.ifact;
tfact = tComponent.tfact;
dom_comp = tComponent.cokey;
dom_compname = tComponent.compname;
area = tComponent.area_pct;
}
} else {
sandtotal = tComponent.sand_pct;
ifact = tComponent.ifact;
tfact = tComponent.tfact;
dom_comp = tComponent.cokey;
dom_compname = tComponent.compname;
area = tComponent.area_pct;
}
}
}
}
}
}
aoa_dom_wind_comp = dom_comp;
aoa_dom_wind_compname = dom_compname;
aoa_ifactor = ifact;
aoa_wind_tfactor = tfact;
if (dom_comp.isEmpty()) {
error_msg += " Could not find any wind critical components in the database. No sand found or no components were greater than 10% of the AOI. ";
ret_val = false;
}
return ret_val;
}
private boolean waterCriticalComponentLookup(ArrayList<SsurgoAttributes> intersectionAttributeList) {
boolean ret_val = true;
String dom_comp = "";
String dom_compname = "";
double ifact = 0.0;
double lsfactor = 0.0;
double tfactor = 0.0;
double kffact = 0.0;
double area = 0.0;
for (SsurgoAttributes attributes : intersectionAttributeList) {
if (attributes.gid_area >= (attributes.aoa_area * 0.10)) {
componentLookup(attributes, true);
if (attributes.horizons.size() > 0) {
for (ComponentAttributes tComponent : attributes.horizons) {
if ((tComponent.kffact >= kffact) && (tComponent.kffact > 0)) {
if (kffact == tComponent.kffact) {
if (tComponent.area_pct > area) {
kffact = tComponent.kffact;
lsfactor = tComponent.lsfact;
tfactor = tComponent.tfact;
ifact = tComponent.ifact;
dom_comp = tComponent.cokey;
dom_compname = tComponent.compname;
area = tComponent.area_pct;
}
} else {
kffact = tComponent.kffact;
lsfactor = tComponent.lsfact;
tfactor = tComponent.tfact;
ifact = tComponent.ifact;
dom_comp = tComponent.cokey;
dom_compname = tComponent.compname;
area = tComponent.area_pct;
}
}
}
}
}
}
aoa_dom_comp = dom_comp;
aoa_dom_compname = dom_compname;
aoa_lsfactor = lsfactor;
aoa_tfactor = tfactor;
aoa_ifactor = ifact;
aoa_kfactor = kffact;
if (dom_comp.isEmpty()) {
error_msg += " Could not find any water critical components in the database. No kffact found or no components were greater than 10% of the AOI. ";
ret_val = false;
}
return ret_val;
}
private boolean componentLookup(SsurgoAttributes attribute, boolean useKFact) {
boolean ret_val = true;
double lambda = 0.0;
if ((null != attribute.mukey) && !attribute.mukey.isEmpty()) {
try (Connection conn = getResourceJDBC("cfactor");
Statement statement = conn.createStatement()) {
String query = "SELECT component.cokey, component.compname, component.wei as ifact, component.slope_r, component.slopelenusle_r as lsfact, component.tfact, component.comppct_r, "
+ " chorizon.cokey, chorizon.kffact, chorizon.kwfact, chorizon.sandtotal_r "
+ " FROM ssurgo.component "
+ " INNER JOIN ssurgo.chorizon on component.cokey=chorizon.cokey "
+ " WHERE component.mukey='" + attribute.mukey + "' "
+ " AND component.slope_r is not null"
+ " AND UPPER(component.majcompflag)='YES' "
+ " AND chorizon.hzdept_r=0 "
+ ((useKFact) ? " AND (chorizon.kffact is not null OR chorizon.kwfact is not null )" : " AND chorizon.sandtotal_r is not null ") + ";";
LOG.info(" ssurgo query=" + query);
ResultSet results = statement.executeQuery(query);
LOG.info(" past the query execute");
attribute.horizons = new ArrayList<>();
while (results.next()) {
ComponentAttributes tHorizon = new ComponentAttributes();
double slope_r = results.getDouble("slope_r");
tHorizon.cokey = results.getString("cokey");
tHorizon.compname = results.getString("compname");
tHorizon.ifact = results.getDouble("ifact");
tHorizon.tfact = results.getDouble("tfact");
tHorizon.slope_r = slope_r;
tHorizon.area_pct = (results.getDouble("comppct_r") / 100.0) * (attribute.gid_area / attribute.aoa_area);
if (useKFact) { // Water Erodibility needed.
// Look up tables for slope length
double[][] slopeLength;
if (isPalouse) {
slopeLength = PalouseSlopeLength;
} else {
slopeLength = LightleWeesiesSlopeLength;
}
for (int i = 0; i < slopeLength.length; i++) {
if (slopeLength[i][1] != -1) {
if ((slope_r >= slopeLength[i][0]) && (slope_r < slopeLength[i][1])) {
lambda = slopeLength[i][2];
break;
} else if (slope_r >= slopeLength[i][0]) {
lambda = slopeLength[i][2];
break;
}
}
}
double length_m = ((slope_r >= 5.0) ? 0.5 : (slope_r >= 3.5) ? 0.4 : (slope_r >= 1.0) ? 0.3 : 0.2);
double theta = Math.atan(slope_r / 100.0);
tHorizon.lsfact = Math.pow((lambda / 72.6), length_m) * (65.41 * Math.pow(Math.sin(theta), 2.0) + 4.56 * Math.sin(theta) + 0.065); //results.getDouble("lsfact");
tHorizon.kffact = results.getDouble("kffact");
if (results.wasNull()) {
tHorizon.kffact = results.getDouble("kwfact");
if (results.wasNull()) {
error_msg = "Invalid kffact selected, empty columns in both kffact and kwfact";
ret_val = false;
} else {
attribute.horizons.add(tHorizon);
ret_val = true;
}
} else {
attribute.horizons.add(tHorizon);
ret_val = true;
}
} else {
tHorizon.lsfact = 0.0;
tHorizon.sand_pct = results.getDouble("sandtotal_r");
attribute.horizons.add(tHorizon);
ret_val = true;
}
}
if (attribute.horizons.size() <= 0) {
ret_val = false;
error_msg += "No soil components and/or horizon data found for this mukey: " + attribute.mukey + " .";
}
} catch (SQLException | csip.ServiceException ex) {
error_msg += "Could not continue processing the component lookups for this mapunit intersect: " + ex.getMessage();
}
}
return ret_val;
}
/**
* @author Shaun Case
* @param polygon
* @return
* @throws SQLException
*/
Centroid getCentroid(String polygon) throws SQLException, ServiceException {
Centroid cent = new Centroid();
try (Connection c = getResourceJDBC("cfactor")) {
try (Statement s = c.createStatement()) {
String find_centroid = "select (" + polygon + ".STCentroid().STY) as lat, (" + polygon + ".STCentroid().STX) as long;";
try (ResultSet r = s.executeQuery(find_centroid)) {
LOG.info("query sql=" + find_centroid);
if (r.getMetaData().getColumnCount() != 2) {
LOG.severe("invalid columns in getCentroid query");
throw new ServiceException("invalid columns in getCentroid query");
}
LOG.info("getting centroid for polygon");
if (r.next()) {
cent.lat = r.getString(1);
cent.lon = r.getString(2);
}
}
}
}
return cent;
}
/**
* @author Wes Lloyd
* @param slat
* @param slon
* @return
* @throws SQLException
* @throws ServiceException
* @throws IOException Synchronized in an attempt to reduce random shell
* failures...
*/
double getCFactorRaster(String slat, String slon) throws SQLException, ServiceException, IOException {
double c_fact = Const.UNKNOWN_CFACTOR;
//File tiffFile = new File("/tmp/us_cvalues_topo2ras_masked.tif");
String gdal = Config.getString("gdal.deployment", Const.RSE_DEPLOYMENT_LINUX);
Executable gdal_exe = getResourceExe(gdal.equals(Const.RSE_DEPLOYMENT_WINDOWS) ? "gdal-win" : "gdal-lin");
// to do
//pc.exe = (gdal.equals(Const.RSE_DEPLOYMENT_WINDOWS) ? gdalwinpath : "") + "gdallocationinfo";
LOG.info("gdal is equal to=" + gdal);
File tifptr = getResourceFile("tif");
gdal_exe.setArguments(tifptr.toString(), slon, slat, "-wgs84");
StringWriter stdout = new StringWriter();
StringWriter stderr = new StringWriter();
gdal_exe.redirectOutput(stdout);
gdal_exe.redirectError(stderr);
int res = gdal_exe.exec();
LOG.info("RET=" + res);
LOG.info("STDOUT=" + stdout.toString());
LOG.info("STDERR=" + stderr.toString());
String output[] = stdout.toString().split("\\r?\\n");
if (res == 0 && output.length > 0) {
String val[];
val = output[output.length - 1].split(" ");
if ((null != val) && (val.length > 0) && !val[val.length - 1].isEmpty()) {
c_fact = Double.parseDouble(val[val.length - 1]);
} else {
LOG.severe("gdallocationinfo returned no parseable output. Cannot proceed. ");
throw new csip.ServiceException("Cannot proceed with cfactor calculation; gdallocationinfo returned no parseable output: " + stdout.toString() + " err " + stderr.toString());
}
} else {
LOG.severe("gdallocationinfo returned no output on stdout. Cannot proceed. ");
throw new csip.ServiceException("Cannot proceed with cfactor calculation: " + stderr.toString());
}
return c_fact;
}
///// inner classes
static class SsurgoAttributes {
public String gid = "";
public String AoA_Id = "";
public String mukey = "";
public double gid_area;
public double aoa_area = 1.0;
public ArrayList<ComponentAttributes> horizons;
}
static class ComponentAttributes {
public String cokey = "";
public String compname = "";
public double ifact;
public double slope_r;
public double tfact;
public double kffact;
public double kwfact;
public double area_pct;
public double lsfact;
public double sand_pct;
}
static class Centroid {
public String lon = "0.0";
public String lat = "0.0";
}
public static class PolygonLatLon {
private ArrayList<PointLatLon> points = new ArrayList<>();
private boolean bValid = true;
public enum ESRIContourType {
outer, inner, unknown
};
protected ESRIContourType myRotation;
public PolygonLatLon() {
}
public PolygonLatLon(ArrayList<PointLatLon> tPoints) {
if (null != tPoints) {
for (PointLatLon tPoint : tPoints) {
if (!tPoint.isValid()) {
bValid = false;
break;
}
}
}
}
public void add(PointLatLon tPoint) {
points.add(tPoint);
if (!tPoint.isValid()) {
bValid = false;
}
}
public void add(double Lat, double Lon) {
PointLatLon tPoint = new PointLatLon(Lat, Lon);
points.add(tPoint);
if (!tPoint.isValid()) {
bValid = false;
}
}
public String toWKT() {
String ret_val = "'POLYGON((";
int currentPoint = 0;
if (bValid) {
for (PointLatLon tPoint : points) {
if (currentPoint > 0) {
ret_val += ", ";
}
// Remember Lon is X and Lat is Y...
ret_val += tPoint.getLon() + " " + tPoint.getLat();
currentPoint++;
}
ret_val += "))'";
} else {
ret_val = "";
}
return ret_val;
}
public boolean isValid() {
getRotation();
return bValid;
}
protected void swapRotation() {
ArrayList<PointLatLon> newPoints = new ArrayList<>();
for (int j = points.size() - 1; j >= 0; j--) {
PointLatLon tPoint = new PointLatLon(points.get(j).getLat(), points.get(j).getLon());
newPoints.add(tPoint);
}
points.clear();
points = newPoints;
}
public boolean makeESRIRotation(ESRIContourType contourType) {
boolean ret_val = false;
if (isValid()) {
if (((contourType == ESRIContourType.outer) && (myRotation == ESRIContourType.inner))
|| ((contourType == ESRIContourType.inner) && (myRotation == ESRIContourType.outer))) {
swapRotation();
}
ret_val = true;
}
return ret_val;
}
protected void getRotation() {
int j;
myRotation = ESRIContourType.unknown;
if (bValid) {
// Calculate cartesian area...test for positive or negative value. This tells the rotation of points CW versus CCW.
// NOTE: Don't use this method to calc "real" areas of Lat/Lon! This is just a fast way of finding out the order of points.
// Lat/Lon must be projected first to get actual area using this type of cross-product method.
if (points.size() > 2) {
double contour_area = 0.0;
for (j = 0; j < points.size() - 2; j++) {
PointLatLon a = points.get(j);
PointLatLon b = points.get(j + 1);
contour_area += a.getLon() * b.getLat() - b.getLon() * a.getLat();
}
PointLatLon a = points.get(j);
PointLatLon b = points.get(0);
contour_area += a.getLon() * b.getLat() - b.getLon() * a.getLat();
if (contour_area < 0) {
myRotation = ESRIContourType.inner; //CCW
} else {
myRotation = ESRIContourType.outer; //CW
}
}
}
}
public static class PointLatLon {
private double Lat;
private double Lon;
private boolean bValid;
public PointLatLon(double Lat, double Lon) {
this.Lat = Lat;
this.Lon = Lon;
bValid = true;
if (Math.abs(Lon) > 180) {
bValid = false;
}
if (Math.abs(Lat) > 90) {
bValid = false;
}
}
public boolean isValid() {
return bValid;
}
public double getLat() {
return Lat;
}
public double getLon() {
return Lon;
}
public void setLon(double Lon) {
this.Lon = Lon;
if (Math.abs(Lon) > 180) {
bValid = false;
}
}
public void setLat(double Lat) {
this.Lat = Lat;
if (Math.abs(Lat) > 90) {
bValid = false;
}
}
}
}
}