From b7f3492c1a3e242d8121a9c7bdc2a36abb05a8f5 Mon Sep 17 00:00:00 2001 From: Ben Steffensmeier Date: Thu, 18 Oct 2012 12:39:01 -0500 Subject: [PATCH] Issue #1279 change GridCoverage to choose more reasonable longitude ranges, use a 0 central meridian for huge grids to avoid normalization in geotools, fix intersections in GridGeometryProjections to work for unnormalized areas, fix data wrapper world wrap detection to work for any crs which doesn't normalize. Former-commit-id: 9bef14fdb1639aad28cbb0298920053eabc56180 [formerly aac3144f518ac21b69b98d3b21830efb94bcd0c8] [formerly 2f7a04f8efb8eade6cf9abb627ba34fac4def73a] [formerly e3cc30e8a8c3d4930a0dd6581976d42f12828130 [formerly 2f7a04f8efb8eade6cf9abb627ba34fac4def73a [formerly 908b9ec666e5dfb50ef4fe3899bf2288bf4cecc5]]] Former-commit-id: e3cc30e8a8c3d4930a0dd6581976d42f12828130 Former-commit-id: f7f14001910f26e749c9069418cd703088ba92a9 [formerly a67d5e2b6c1750e962391094fe05e41639a40f98] Former-commit-id: 6748ebb4ecc063b14e726d669c5882bd1c17d152 --- .../spatial/projections/GridCoverage.java | 21 ++++++++++++++--- .../projections/LatLonGridCoverage.java | 11 ++++++++- .../uf/common/geospatial/MapUtil.java | 23 ++++++++++++++++--- .../data/AbstractDataWrapper.java | 6 ++--- 4 files changed, 50 insertions(+), 11 deletions(-) diff --git a/edexOsgi/com.raytheon.uf.common.dataplugin.grib/src/com/raytheon/uf/common/dataplugin/grib/spatial/projections/GridCoverage.java b/edexOsgi/com.raytheon.uf.common.dataplugin.grib/src/com/raytheon/uf/common/dataplugin/grib/spatial/projections/GridCoverage.java index dc17fb6cd2..e6fc0f1667 100644 --- a/edexOsgi/com.raytheon.uf.common.dataplugin.grib/src/com/raytheon/uf/common/dataplugin/grib/spatial/projections/GridCoverage.java +++ b/edexOsgi/com.raytheon.uf.common.dataplugin.grib/src/com/raytheon/uf/common/dataplugin/grib/spatial/projections/GridCoverage.java @@ -377,8 +377,8 @@ public abstract class GridCoverage extends PersistableDataObject implements if (gridGeometry == null) { gridGeometry = MapUtil.getGridGeometry(this); } - - return gridGeometry; + + return gridGeometry; } public void setGridGeometry(GridGeometry2D gridGeometry) { @@ -595,6 +595,20 @@ public abstract class GridCoverage extends PersistableDataObject implements maxLon -= 1.0E-12; minLon += 1.0E-12; } + + // Normalize the range by shifting 360 degrees to bring the range + // within +/-360 degree as much as possible. For example the + // Canadian-NH model gets calculated as 179.7 to 540.3 but it + // works better to use -180.3 to 180.3. + while (minLon > 0 && maxLon > 360) { + minLon -= 360; + maxLon -= 360; + } + // Normalize the low end. + while (minLon < -360 && maxLon < 0) { + minLon += 360; + maxLon += 360; + } try { geometry = MapUtil.createGeometry(minLat, minLon, maxLat, maxLon); @@ -743,7 +757,8 @@ public abstract class GridCoverage extends PersistableDataObject implements if (isSubGridded()) { String subGridName = getName(); int index = subGridName.lastIndexOf(SUBGRID_TOKEN); - if (index >= 0 && index + SUBGRID_TOKEN.length() < subGridName.length()) { + if (index >= 0 + && index + SUBGRID_TOKEN.length() < subGridName.length()) { model = subGridName.substring(index + SUBGRID_TOKEN.length()); } } diff --git a/edexOsgi/com.raytheon.uf.common.dataplugin.grib/src/com/raytheon/uf/common/dataplugin/grib/spatial/projections/LatLonGridCoverage.java b/edexOsgi/com.raytheon.uf.common.dataplugin.grib/src/com/raytheon/uf/common/dataplugin/grib/spatial/projections/LatLonGridCoverage.java index f384ad2476..9cc4dffd6b 100644 --- a/edexOsgi/com.raytheon.uf.common.dataplugin.grib/src/com/raytheon/uf/common/dataplugin/grib/spatial/projections/LatLonGridCoverage.java +++ b/edexOsgi/com.raytheon.uf.common.dataplugin.grib/src/com/raytheon/uf/common/dataplugin/grib/spatial/projections/LatLonGridCoverage.java @@ -97,7 +97,16 @@ public class LatLonGridCoverage extends GridCoverage { double maxLon = minLon + dx * nx; double centralMeridian = (minLon + maxLon) / 2.0; - centralMeridian = MapUtil.correctLon(centralMeridian); + if (dx * nx <= 360) { + centralMeridian = MapUtil.correctLon(centralMeridian); + } else { + // For almost all map projections geotools will clip all math + // transforms to be within +-180 of the central meridian. For grids + // that wrap around the world more than once this is a problem. When + // the central Meridian is 0.0 then geotools does not do this + // clipping, which works much better. + centralMeridian = 0.0; + } crs = MapUtil.constructEquidistantCylindrical( MapUtil.AWIPS_EARTH_RADIUS, MapUtil.AWIPS_EARTH_RADIUS, centralMeridian, 0); diff --git a/edexOsgi/com.raytheon.uf.common.geospatial/src/com/raytheon/uf/common/geospatial/MapUtil.java b/edexOsgi/com.raytheon.uf.common.geospatial/src/com/raytheon/uf/common/geospatial/MapUtil.java index 8bfba07b62..f71f344d93 100644 --- a/edexOsgi/com.raytheon.uf.common.geospatial/src/com/raytheon/uf/common/geospatial/MapUtil.java +++ b/edexOsgi/com.raytheon.uf.common.geospatial/src/com/raytheon/uf/common/geospatial/MapUtil.java @@ -385,11 +385,28 @@ public class MapUtil { LATLON_PROJECTION, false); ReferencedEnvelope newTargetREnv = targetREnv.transform( LATLON_PROJECTION, false); - + com.vividsolutions.jts.geom.Envelope intersection = newTargetREnv + .intersection(newSourceEnv); + // Its possible to get two envelopes that don't intersect in a common + // space, for example one could have longitude from -200 to -160 and + // another could have longitude from 160 to 200. Even though these are + // the same range, they don't intersect. These two loops will shift the + // data 360 degrees in the x direction until all intersections are + // found. + while (newSourceEnv.getMaxX() > newTargetREnv.getMinX()) { + newSourceEnv.translate(-360, 0); + intersection.expandToInclude(newTargetREnv + .intersection(newSourceEnv)); + } + while (newSourceEnv.getMinX() < newTargetREnv.getMaxX()) { + newSourceEnv.translate(360, 0); + intersection.expandToInclude(newTargetREnv + .intersection(newSourceEnv)); + } // Get the newEnvelope ReferencedEnvelope newEnv = new ReferencedEnvelope(JTS.getEnvelope2D( - newTargetREnv.intersection(newSourceEnv), LATLON_PROJECTION), - LATLON_PROJECTION); + intersection, LATLON_PROJECTION), LATLON_PROJECTION); + newEnv = newEnv.transform(targetCRS, false, 500); // Calculate nx and ny, start with the number of original grid // points in the intersection and then adjust to the new aspect diff --git a/edexOsgi/com.raytheon.uf.common.geospatial/src/com/raytheon/uf/common/geospatial/interpolation/data/AbstractDataWrapper.java b/edexOsgi/com.raytheon.uf.common.geospatial/src/com/raytheon/uf/common/geospatial/interpolation/data/AbstractDataWrapper.java index a61e87f9fd..54fdb2e2d8 100644 --- a/edexOsgi/com.raytheon.uf.common.geospatial/src/com/raytheon/uf/common/geospatial/interpolation/data/AbstractDataWrapper.java +++ b/edexOsgi/com.raytheon.uf.common.geospatial/src/com/raytheon/uf/common/geospatial/interpolation/data/AbstractDataWrapper.java @@ -27,8 +27,6 @@ import org.opengis.referencing.crs.CoordinateReferenceSystem; import org.opengis.referencing.datum.PixelInCell; import org.opengis.referencing.operation.MathTransform; -import com.raytheon.uf.common.geospatial.MapUtil; - /** * * Abstract class for any data implementation that can act as both a source and @@ -138,8 +136,8 @@ public abstract class AbstractDataWrapper implements DataSource, grid2crs.transform(corner2, corner2); crs2LatLon.transform(corner1, corner1); crs2LatLon.transform(corner2, corner2); - corner1.x = MapUtil.correctLon(corner1.x); - corner2.x = MapUtil.correctLon(corner2.x); + corner1.x = corner1.x - 360; + corner2.x = corner2.x - 360; crs2LatLon.inverse().transform(corner1, corner1); crs2LatLon.inverse().transform(corner2, corner2); grid2crs.inverse().transform(corner1, corner1);