diff --git a/cave/com.raytheon.viz.hydrobase/src/com/raytheon/viz/hydrobase/data/HydroGeoProcessor.java b/cave/com.raytheon.viz.hydrobase/src/com/raytheon/viz/hydrobase/data/HydroGeoProcessor.java index 0b7a57f6db..b4fe992fb0 100644 --- a/cave/com.raytheon.viz.hydrobase/src/com/raytheon/viz/hydrobase/data/HydroGeoProcessor.java +++ b/cave/com.raytheon.viz.hydrobase/src/com/raytheon/viz/hydrobase/data/HydroGeoProcessor.java @@ -22,11 +22,11 @@ package com.raytheon.viz.hydrobase.data; import java.util.ArrayList; import java.util.List; +import com.raytheon.uf.common.geospatial.MapUtil; +import com.raytheon.uf.common.hydro.spatial.HRAP; import com.raytheon.viz.hydrocommon.util.HrapUtil; import com.vividsolutions.jts.geom.Coordinate; -import com.vividsolutions.jts.geom.GeometryFactory; -import com.vividsolutions.jts.geom.LinearRing; -import com.vividsolutions.jts.geom.Point; +import com.vividsolutions.jts.geom.Geometry; import com.vividsolutions.jts.geom.Polygon; /** @@ -47,18 +47,36 @@ import com.vividsolutions.jts.geom.Polygon; */ public class HydroGeoProcessor { - public static GeometryFactory factory = new GeometryFactory(); + private final Geometry[][] hrapGeometries; - public HydroGeoProcessor() { + /** + * Constructor. Set up the hrapGeometry cache. + * + * @throws Exception + */ + public HydroGeoProcessor() throws Exception { + HRAP hrap = HRAP.getInstance(); + Coordinate ll = hrap.getLatLonLL(); + Coordinate ur = hrap.getLatLonUR(); + + Coordinate hrapLL = HrapUtil.latLonToHrap(ll); + Coordinate hrapUR = HrapUtil.latLonToHrap(ur); + int cols = (int) Math.floor(hrapUR.x - hrapLL.x); + int rows = (int) Math.floor(hrapUR.y - hrapLL.y); + hrapGeometries = new Geometry[rows][cols]; } - public HrapBinList getHrapBinList(GeoAreaData geoData) { + /** + * Get the HrapBinList. + * + * @param GeoAreaData + * data for the feature to check + */ + public HrapBinList getHrapBinList(GeoAreaData geoData) throws Exception { List coords = getPointsFromArea(geoData); Coordinate[] minMaxXY = getMinMaxXY(coords); - LinearRing lr = factory.createLinearRing(coords - .toArray(new Coordinate[0])); - Polygon poly = factory.createPolygon(lr, null); + Polygon poly = MapUtil.getPolygon(coords.toArray(new Coordinate[0])); Coordinate minC = minMaxXY[0]; Coordinate maxC = minMaxXY[1]; @@ -66,15 +84,10 @@ public class HydroGeoProcessor { Coordinate hrapMin = HrapUtil.latLonToHrap(minC); Coordinate hrapMax = HrapUtil.latLonToHrap(maxC); - double wfoMinX = hrapMin.x; - double wfoMinY = hrapMin.y; - double wfoMaxX = hrapMax.x; - double wfoMaxY = hrapMax.y; - - double maxRow = Math.floor(wfoMaxY); - double maxCol = Math.floor(wfoMaxX); - double minRow = Math.floor(wfoMinY); - double minCol = Math.floor(wfoMinX); + int maxRow = (int) Math.floor(hrapMax.y); + int maxCol = (int) Math.floor(hrapMax.x); + int minRow = (int) Math.floor(hrapMin.y); + int minCol = (int) Math.floor(hrapMin.x); /* expand the box to make sure polygon has been covered */ minRow -= 2; @@ -82,29 +95,36 @@ public class HydroGeoProcessor { maxRow += 2; maxCol += 2; + int rows = maxRow - minRow; + int cols = maxCol - minCol; + int rowCtr = 0; - double rowNum = 0; - double startCol = 0; - double endCol = 0; + int rowNum = 0; + int colNum = 0; + int startCol = 0; + int endCol = 0; int binCtr = 0; double area = 0; HrapBinList binList = new HrapBinList(); - for (double r = minRow + 0.5; r <= maxRow; r++) { // row - rowNum = r; + for (int r = 0; r < rows; r++) { + rowNum = r + minRow; startCol = -1; - - for (double c = minCol + 0.5; c <= maxCol; c++) { - Coordinate coord = new Coordinate(c, r); - Coordinate gridCell = HrapUtil.hrapToLatLon(coord); - Point p = factory.createPoint(gridCell); - if (poly.intersects(p)) { // inside - endCol = c; + colNum = 0; + for (int c = 0; c < cols; c++) { + colNum = c + minCol; + Coordinate coord = new Coordinate(colNum, rowNum); + if (hrapGeometries[rowNum][colNum] == null) { + hrapGeometries[rowNum][colNum] = HrapUtil + .getGridCellPolygon(coord); + } + if (poly.intersects(hrapGeometries[rowNum][colNum])) { + endCol = c + cols; binCtr++; if (startCol == -1) { // First cell in the row - startCol = c; + startCol = c + cols; rowCtr++; } area += HrapUtil.getHrapBinArea(coord); @@ -139,13 +159,15 @@ public class HydroGeoProcessor { for (Coordinate c : coords) { if (c.x > maxX) { maxX = c.x; - } else if (c.x < minX) { + } + if (c.x < minX) { minX = c.x; } if (c.y > maxY) { maxY = c.y; - } else if (c.y < minY) { + } + if (c.y < minY) { minY = c.y; } } diff --git a/cave/com.raytheon.viz.hydrobase/src/com/raytheon/viz/hydrobase/dialogs/ArealDefinitionsDlg.java b/cave/com.raytheon.viz.hydrobase/src/com/raytheon/viz/hydrobase/dialogs/ArealDefinitionsDlg.java index d51cfb6116..ec6c29e1f6 100644 --- a/cave/com.raytheon.viz.hydrobase/src/com/raytheon/viz/hydrobase/dialogs/ArealDefinitionsDlg.java +++ b/cave/com.raytheon.viz.hydrobase/src/com/raytheon/viz/hydrobase/dialogs/ArealDefinitionsDlg.java @@ -938,12 +938,29 @@ public class ArealDefinitionsDlg extends CaveSWTDialog { } // Load the linesegs table - HydroGeoProcessor proc = new HydroGeoProcessor(); + HydroGeoProcessor proc; + try { + proc = new HydroGeoProcessor(); + } catch (Exception e) { + statusHandler.error( + "Error initializing the Hydro GeoProcessor", e); + log("Error initializing the Hydro GeoProcessor"); + return; + } if (selectedType != ArealTypeSelection.RESERVOIRS) { for (GeoAreaData data : geoDataList) { + /* do the main processing */ - HrapBinList binList = proc.getHrapBinList(data); + HrapBinList binList; + try { + binList = proc.getHrapBinList(data); + } catch (Exception e) { + statusHandler.error("Error processing input data", e); + log("Error processing input data"); + return; + } + log("Processing area " + data.getAreaId() + ":" + " Writing " + binList.getNumRows() + " rows"); dman.putLineSegs(data.getAreaId(), binList); diff --git a/cave/com.raytheon.viz.hydrocommon/src/com/raytheon/viz/hydrocommon/util/HrapUtil.java b/cave/com.raytheon.viz.hydrocommon/src/com/raytheon/viz/hydrocommon/util/HrapUtil.java index 3eaed33aa7..d4ec62f6c2 100644 --- a/cave/com.raytheon.viz.hydrocommon/src/com/raytheon/viz/hydrocommon/util/HrapUtil.java +++ b/cave/com.raytheon.viz.hydrocommon/src/com/raytheon/viz/hydrocommon/util/HrapUtil.java @@ -19,74 +19,118 @@ **/ package com.raytheon.viz.hydrocommon.util; +import java.awt.Rectangle; + import org.opengis.metadata.spatial.PixelOrientation; +import com.raytheon.uf.common.geospatial.MapUtil; +import com.raytheon.uf.common.geospatial.ReferencedCoordinate; +import com.raytheon.uf.common.geospatial.ReferencedObject.Type; import com.raytheon.uf.common.hydro.spatial.HRAP; +import com.raytheon.uf.common.hydro.spatial.HRAPCoordinates; +import com.raytheon.uf.common.hydro.spatial.HRAPSubGrid; +import com.raytheon.uf.common.status.IUFStatusHandler; +import com.raytheon.uf.common.status.UFStatus; +import com.raytheon.uf.common.status.UFStatus.Priority; import com.vividsolutions.jts.geom.Coordinate; +import com.vividsolutions.jts.geom.Polygon; /** * Lat/Lon to HRAP and HRAP to Lat/Lon conversion for HydroBase. * * 
- *
+ * 
  * SOFTWARE HISTORY
- *
+ * 
  * Date         Ticket#    Engineer    Description
  * ------------ ---------- ----------- --------------------------
- * Sep 15, 2009 2772       mpduff     Initial creation
- *
+ * Sep 15, 2009 2772       mpduff      Initial creation
+ * Jan 07, 2016 5217       mpduff      Added getGridCellPolygon and HRAPSubGrid processing.
+ * 
  *
- * + * * @author mpduff - * @version 1.0 + * @version 1.0 */ public class HrapUtil { + private static final IUFStatusHandler statusHandler = UFStatus + .getHandler(HrapUtil.class); + private static final double NMILE_PER_DEG_LAT = 60; + private static final double PI = 3.14159; - private static final double RAD_PER_DEG = PI/180; + + private static final double RAD_PER_DEG = PI / 180; + private static final double KM_PER_NMILE = 1.852; + private static final double KM_PER_MILE = 1.609344; - private static final double MILES_PER_NMILE = KM_PER_NMILE/KM_PER_MILE; - + + private static final double MILES_PER_NMILE = KM_PER_NMILE / KM_PER_MILE; + + private static HRAPSubGrid subGrid; + /** * Convert a lat/lon set to HRAP coordinates. * * @param latlon - * The latitude/longitude Coordinate - * @return - * HRAP Coordinate + * The latitude/longitude Coordinate + * @return HRAP Coordinate */ public static Coordinate latLonToHrap(Coordinate latlon) { HRAP hrap = HRAP.getInstance(); Coordinate gridCoord = null; try { - gridCoord = hrap.latLonToGridCoordinate(latlon, PixelOrientation.LOWER_LEFT); + gridCoord = hrap.latLonToGridCoordinate(latlon, + PixelOrientation.LOWER_LEFT); } catch (Exception e) { e.printStackTrace(); } return gridCoord; } - + /** * Convert an HRAP coordinate to lat/lon. * * @param gridCoord - * the HRAP Coordinate - * @return - * Coordinate - the lat/lon Coordinate + * the HRAP Coordinate + * @return Coordinate - the lat/lon Coordinate */ public static Coordinate hrapToLatLon(Coordinate gridCoord) { HRAP hrap = HRAP.getInstance(); - Coordinate latlon = new Coordinate(0,0); + Coordinate latlon = new Coordinate(0, 0); try { - latlon = hrap.gridCoordinateToLatLon(gridCoord, PixelOrientation.LOWER_LEFT); + latlon = hrap.gridCoordinateToLatLon(gridCoord, + PixelOrientation.LOWER_LEFT); } catch (Exception e) { e.printStackTrace(); } - + return latlon; } - + + /** + * Get the polygon that covers the HRAP grid cell. + * + * @param gridCoord + * The coordinate of the grid cell + * @return + */ + public static Polygon getGridCellPolygon(Coordinate gridCoord) { + Coordinate[] coors = new Coordinate[5]; + + coors[0] = getHRAPLatLon(new Coordinate(gridCoord.x, gridCoord.y)); + coors[1] = getHRAPLatLon(new Coordinate(gridCoord.x + 1, gridCoord.y)); + coors[2] = getHRAPLatLon(new Coordinate(gridCoord.x + 1, + gridCoord.y + 1)); + coors[3] = getHRAPLatLon(new Coordinate(gridCoord.x, gridCoord.y + 1)); + + // complete the square + coors[4] = coors[0]; + + return MapUtil.getPolygon(coors); + } + /** * compute area of an Hrap bin. * @@ -100,47 +144,94 @@ public class HrapUtil { * conversion factor to square miles. * * @param row - * The row + * The row * @param col - * The column - * @return - * The area + * The column + * @return The area */ public static double getHrapBinArea(Coordinate coord) { - double area; - - double latLength; - double lonLength; - - double nmileXLength; - double nmileYLength; - - double lat1; - double lat2; - double lon1; - double lon2; - - /* get lat and lon positions */ - Coordinate ll = hrapToLatLon(coord); - Coordinate ll2 = hrapToLatLon(new Coordinate(coord.x + 1, coord.y + 1)); - - lat1 = ll.y; - lat2 = ll2.y; - lon1 = ll.x; - lon2 = ll2.x; - - /* get differences in lat and lon */ - latLength = Math.abs(lat1 - lat2); - lonLength = Math.abs(lon1 - lon2); - - /* get the nautical miles from lat and lon lengths */ - nmileYLength = NMILE_PER_DEG_LAT * latLength; - nmileXLength = lonLength * NMILE_PER_DEG_LAT*(Math.cos(RAD_PER_DEG*((lat1+lat2)/2.0))); - - /* convert to miles and get the square */ - area = (MILES_PER_NMILE * nmileYLength) * - (MILES_PER_NMILE * nmileXLength); - - return area; + double area; + + double latLength; + double lonLength; + + double nmileXLength; + double nmileYLength; + + double lat1; + double lat2; + double lon1; + double lon2; + + /* get lat and lon positions */ + Coordinate ll = hrapToLatLon(coord); + Coordinate ll2 = hrapToLatLon(new Coordinate(coord.x + 1, coord.y + 1)); + + lat1 = ll.y; + lat2 = ll2.y; + lon1 = ll.x; + lon2 = ll2.x; + + /* get differences in lat and lon */ + latLength = Math.abs(lat1 - lat2); + lonLength = Math.abs(lon1 - lon2); + + /* get the nautical miles from lat and lon lengths */ + nmileYLength = NMILE_PER_DEG_LAT * latLength; + nmileXLength = lonLength * NMILE_PER_DEG_LAT + * (Math.cos(RAD_PER_DEG * ((lat1 + lat2) / 2.0))); + + /* convert to miles and get the square */ + area = (MILES_PER_NMILE * nmileYLength) + * (MILES_PER_NMILE * nmileXLength); + + return area; + } + + /** + * Get the HRAP subgrid + * + * @param grib + * factor + * @throws Exception + */ + public static HRAPSubGrid getHrapSubGrid() throws Exception { + return getHrapSubGrid(1); + } + + /** + * Get the HRAP subgrid + * + * @param grib + * factor + * @throws Exception + */ + public static HRAPSubGrid getHrapSubGrid(int hrapGribFactor) + throws Exception { + if (subGrid == null) { + Rectangle rect = HRAPCoordinates.getHRAPCoordinates(); + subGrid = new HRAPSubGrid(rect, hrapGribFactor); + } + + return subGrid; + } + + /** + * Gets the gribPoint + * + * @param gridPoint + * @return + */ + private static Coordinate getHRAPLatLon(Coordinate gridPoint) { + try { + ReferencedCoordinate rc = new ReferencedCoordinate(gridPoint, + getHrapSubGrid().getHRAP().getGridGeometry(), + Type.GRID_CORNER); + gridPoint = rc.asLatLon(); + } catch (Exception e) { + statusHandler.handle(Priority.ERROR, + "Unable translate grid coordinate: " + gridPoint); + } + return gridPoint; } }