diff --git a/gears/src/main/java/org/hortonmachine/gears/libs/modules/GridNode.java b/gears/src/main/java/org/hortonmachine/gears/libs/modules/GridNode.java index 9b7df9da9..4c9cb6583 100644 --- a/gears/src/main/java/org/hortonmachine/gears/libs/modules/GridNode.java +++ b/gears/src/main/java/org/hortonmachine/gears/libs/modules/GridNode.java @@ -56,6 +56,11 @@ public class GridNode extends Node { private static DecimalFormat f = new DecimalFormat("0.0000000"); private List nodes; + + + public GridNode (HMRaster raster, int col, int row) { + this(raster.getIter(), raster.getCols(), raster.getRows(), raster.getXRes(), raster.getYRes(), col, row, raster.getNovalue()); + } /** * The constructor. @@ -317,13 +322,16 @@ public List getWindow( int size ) { int tmpCol = col + c; for( int r = -delta; r <= delta; r++ ) { int tmpRow = row + r; - GridNode n = new GridNode(gridIter, cols, rows, xRes, yRes, tmpCol, tmpRow, doubleNovalue); - windowNodes.add(n); + if (isInRaster(tmpCol, tmpRow)) { + GridNode n = new GridNode(gridIter, cols, rows, xRes, yRes, tmpCol, tmpRow, doubleNovalue); + windowNodes.add(n); + } } } return windowNodes; } + /** * Get the value of the elevation in one of the surrounding direction. * diff --git a/gears/src/main/java/org/hortonmachine/gears/libs/modules/HMRaster.java b/gears/src/main/java/org/hortonmachine/gears/libs/modules/HMRaster.java index c5e52188b..943d32870 100644 --- a/gears/src/main/java/org/hortonmachine/gears/libs/modules/HMRaster.java +++ b/gears/src/main/java/org/hortonmachine/gears/libs/modules/HMRaster.java @@ -152,6 +152,10 @@ public RegionMap getRegionMap() { public GridGeometry2D getGridGeometry() { return gridGeometry; } + + public RandomIter getIter() { + return iter; + } /** * @return the columns of this raster. diff --git a/gears/src/main/java/org/hortonmachine/gears/modules/r/holefiller/OmsHoleFiller.java b/gears/src/main/java/org/hortonmachine/gears/modules/r/holefiller/OmsHoleFiller.java index aa60344bd..c89f3b6a3 100644 --- a/gears/src/main/java/org/hortonmachine/gears/modules/r/holefiller/OmsHoleFiller.java +++ b/gears/src/main/java/org/hortonmachine/gears/modules/r/holefiller/OmsHoleFiller.java @@ -118,7 +118,7 @@ public class OmsHoleFiller extends HMModel { public static final String OMSHOLEFILLER_IN_RASTER_DESCRIPTION = "The raster to fill holes in."; public static final String OMSHOLEFILLER_MODE_DESCRIPTION = "The interpolation mode to use."; public static final String OMSHOLEFILLER_P_BUFFER_DESCRIPTION = "The buffer to use for interpolation."; - public static final String OMSHOLEFILLER_IN_ROI_DESCRIPTION = "The regions vector map on which to checkc for nulls."; + public static final String OMSHOLEFILLER_IN_ROI_DESCRIPTION = "The regions vector map on which to check for nulls."; public static final String OMSHOLEFILLER_OUT_RASTER_DESCRIPTION = "The new raster."; diff --git a/gears/src/main/java/org/hortonmachine/gears/modules/r/rasternull/OmsRasterMissingValuesFiller.java b/gears/src/main/java/org/hortonmachine/gears/modules/r/rasternull/OmsRasterMissingValuesFiller.java index ac368f2b2..92e5484a5 100644 --- a/gears/src/main/java/org/hortonmachine/gears/modules/r/rasternull/OmsRasterMissingValuesFiller.java +++ b/gears/src/main/java/org/hortonmachine/gears/modules/r/rasternull/OmsRasterMissingValuesFiller.java @@ -44,6 +44,10 @@ import java.util.stream.Stream; import org.geotools.coverage.grid.GridCoverage2D; +import org.geotools.data.simple.SimpleFeatureCollection; +import org.hortonmachine.gears.io.rasterreader.OmsRasterReader; +import org.hortonmachine.gears.io.rasterwriter.OmsRasterWriter; +import org.hortonmachine.gears.io.vectorreader.OmsVectorReader; import org.hortonmachine.gears.libs.modules.Direction; import org.hortonmachine.gears.libs.modules.HMConstants; import org.hortonmachine.gears.libs.modules.HMModel; @@ -54,10 +58,15 @@ import org.hortonmachine.gears.modules.r.interpolation2d.core.ISurfaceInterpolator; import org.hortonmachine.gears.modules.r.interpolation2d.core.LinearDWInterpolator; import org.hortonmachine.gears.modules.r.interpolation2d.core.TPSInterpolator; +import org.hortonmachine.gears.utils.features.FeatureUtilities; import org.locationtech.jts.geom.Coordinate; import org.locationtech.jts.geom.Envelope; +import org.locationtech.jts.geom.Geometry; import org.locationtech.jts.geom.Point; +import org.locationtech.jts.geom.prep.PreparedGeometry; +import org.locationtech.jts.geom.prep.PreparedGeometryFactory; import org.locationtech.jts.index.strtree.STRtree; +import org.locationtech.jts.operation.union.CascadedPolygonUnion; import oms3.annotations.Author; import oms3.annotations.Description; @@ -89,7 +98,11 @@ public class OmsRasterMissingValuesFiller extends HMModel { @Description(OMSRASTERNULLFILLER_IN_RASTERMASK_DESCRIPTION) @In - public GridCoverage2D inMask; + public GridCoverage2D inRasterMask; + + @Description(OMSRASTERNULLFILLER_IN_VECTORMASK_DESCRIPTION) + @In + public SimpleFeatureCollection inVectorMask; @Description(OMSRASTERNULLFILLER_pMinDistance_DESCRIPTION) @In @@ -123,6 +136,7 @@ public class OmsRasterMissingValuesFiller extends HMModel { public static final String OMSRASTERNULLFILLER_AUTHORCONTACTS = "www.hydrologis.com"; public static final String OMSRASTERNULLFILLER_IN_RASTER_DESCRIPTION = "The raster in which to fill missing values."; public static final String OMSRASTERNULLFILLER_IN_RASTERMASK_DESCRIPTION = "An optional raster mask on which to fill missing data."; + public static final String OMSRASTERNULLFILLER_IN_VECTORMASK_DESCRIPTION = "An optional vector mask on which to fill missing data."; public static final String OMSRASTERNULLFILLER_pMinDistance_DESCRIPTION = "The min distance to consider for the interpolation (in cells)."; public static final String OMSRASTERNULLFILLER_pMaxDistance_DESCRIPTION = "The max distance to consider for the interpolation (in cells)."; public static final String OMSRASTERNULLFILLER_OUT_RASTER_DESCRIPTION = "The new raster."; @@ -166,8 +180,15 @@ public void process() throws Exception { .build();) { HMRaster mask = null; - if (inMask != null) { - mask = HMRaster.fromGridCoverage(inMask); + if (inRasterMask != null) { + mask = HMRaster.fromGridCoverage(inRasterMask); + } + + PreparedGeometry vectorMask = null; + if (inVectorMask != null) { + List geoms = FeatureUtilities.featureCollectionToGeometriesList(inVectorMask, false, null); + Geometry union = CascadedPolygonUnion.union(geoms); + vectorMask = PreparedGeometryFactory.prepare(union); } List novaluePoints = new ArrayList<>(); @@ -176,17 +197,18 @@ public void process() throws Exception { novaluePoints.add(new Coordinate(col, row)); } }); + pm.message("Found holes:" + novaluePoints.size()); if (doUseOnlyBorderValues) { - processUsingOnlyBorderValues(interpolator, inData, mask, outData, novaluePoints); + processUsingOnlyBorderValues(interpolator, inData, mask, vectorMask, outData, novaluePoints); } else { - processUsingAllValues(interpolator, inData, mask, outData, novaluePoints); + processUsingAllValues(interpolator, inData, mask, vectorMask, outData, novaluePoints); } } } - private void processUsingAllValues( ISurfaceInterpolator interpolator, HMRaster inData, HMRaster mask, HMRaster outData, - List novaluePoints ) throws IOException { + private void processUsingAllValues( ISurfaceInterpolator interpolator, HMRaster inData, HMRaster mask, + PreparedGeometry vectorMask, HMRaster outData, List novaluePoints ) throws IOException { // now find touching points with values pm.beginTask("Filling holes...", novaluePoints.size()); for( Coordinate noValuePoint : novaluePoints ) { @@ -195,7 +217,15 @@ private void processUsingAllValues( ISurfaceInterpolator interpolator, HMRaster pm.worked(1); continue; } - + if (vectorMask != null) { + Coordinate w = inData.getWorld((int) noValuePoint.x, (int) noValuePoint.y); + Point nvPoint = gf.createPoint(w); + if (!vectorMask.intersects(nvPoint)) { + // do not fill when not in ROI + pm.worked(1); + continue; + } + } List result = inData.getSurroundingCells((int) noValuePoint.x, (int) noValuePoint.y, pMaxDistance, true); result = result.stream().filter(c -> !inData.isNovalue(c.z) && c.distance(noValuePoint) > pMinDistance) .collect(Collectors.toList()); @@ -215,7 +245,7 @@ private void processUsingAllValues( ISurfaceInterpolator interpolator, HMRaster } private void processUsingOnlyBorderValues( ISurfaceInterpolator interpolator, HMRaster inData, HMRaster mask, - HMRaster outData, List novaluePoints ) throws IOException { + PreparedGeometry vectorMask, HMRaster outData, List novaluePoints ) throws IOException { // now find touching points with values STRtree touchingPointsTreetree = new STRtree(); TreeSet checkSet = new TreeSet<>(); @@ -246,6 +276,14 @@ private void processUsingOnlyBorderValues( ISurfaceInterpolator interpolator, HM pm.worked(1); continue; } + if (vectorMask != null) { + Point nvPoint = gf.createPoint(noValuePoint); + if (!vectorMask.intersects(nvPoint)) { + // do not fill when not in ROI + pm.worked(1); + continue; + } + } // get points with values in range Envelope env = new Envelope(new Coordinate(noValuePoint.x, noValuePoint.y)); env.expandBy(pMaxDistance); @@ -300,4 +338,25 @@ public double getValue( List controlPoints, Coordinate interpolated } } + + public static void main( String[] args ) throws Exception { + int cells = 5; + String mode = "IDW"; + OmsRasterMissingValuesFiller transformer = new OmsRasterMissingValuesFiller(); + transformer.inRaster = OmsRasterReader + .readRaster("/storage/lavori_tmp/GEOLOGICO_TN/20240318_dati_per_Silvia/1_errori_dtm_dbm_fiumi/zscore_holes_1.tif"); +// transformer.inRasterMask = inMask; + transformer.inVectorMask = OmsVectorReader.readVector( + "/storage/lavori_tmp/GEOLOGICO_TN/20240318_dati_per_Silvia/1_errori_dtm_dbm_fiumi/fiumi_rappresentabili_buffer.shp"); + transformer.doUseOnlyBorderValues = false; + transformer.pMaxDistance = cells; + transformer.pMinDistance = 0; + transformer.pMode = mode; + transformer.process(); + GridCoverage2D outGC = transformer.outRaster; + OmsRasterWriter.writeRaster( + "/storage/lavori_tmp/GEOLOGICO_TN/20240318_dati_per_Silvia/1_errori_dtm_dbm_fiumi/zscore_holes_1_filled_" + mode + + "_" + cells + ".tif", + outGC); + } } diff --git a/gears/src/test/java/org/hortonmachine/gears/modules/TestRasterMissingFiller.java b/gears/src/test/java/org/hortonmachine/gears/modules/TestRasterMissingFiller.java index 20fc8723f..4f2c2e2d5 100644 --- a/gears/src/test/java/org/hortonmachine/gears/modules/TestRasterMissingFiller.java +++ b/gears/src/test/java/org/hortonmachine/gears/modules/TestRasterMissingFiller.java @@ -70,7 +70,7 @@ public void testRasterFiller() throws Exception { OmsRasterMissingValuesFiller transformer = new OmsRasterMissingValuesFiller(); transformer.inRaster = inCoverage; - transformer.inMask = inMask; + transformer.inRasterMask = inMask; transformer.doUseOnlyBorderValues = false; transformer.pMaxDistance = 1; transformer.pMinDistance = 0; @@ -99,7 +99,7 @@ public void testRasterFillerMasked() throws Exception { OmsRasterMissingValuesFiller transformer = new OmsRasterMissingValuesFiller(); transformer.inRaster = inCoverage; - transformer.inMask = inMask2; + transformer.inRasterMask = inMask2; transformer.doUseOnlyBorderValues = false; transformer.pMaxDistance = 1; transformer.pMinDistance = 0; @@ -128,7 +128,7 @@ public void testRasterFillerBordersOnly() throws Exception { OmsRasterMissingValuesFiller transformer = new OmsRasterMissingValuesFiller(); transformer.inRaster = inCoverage; - transformer.inMask = inMask; + transformer.inRasterMask = inMask; transformer.doUseOnlyBorderValues = false; transformer.pMaxDistance = 3; transformer.pMinDistance = 0; diff --git a/modules/src/main/java/org/hortonmachine/modules/RasterMissingValuesFiller.java b/modules/src/main/java/org/hortonmachine/modules/RasterMissingValuesFiller.java index 232425385..10855e642 100644 --- a/modules/src/main/java/org/hortonmachine/modules/RasterMissingValuesFiller.java +++ b/modules/src/main/java/org/hortonmachine/modules/RasterMissingValuesFiller.java @@ -102,7 +102,7 @@ public void process() throws Exception { OmsRasterMissingValuesFiller ormvf = new OmsRasterMissingValuesFiller(); ormvf.pm = pm; ormvf.inRaster = getRaster(inRaster); - ormvf.inMask = getRaster(inMask); + ormvf.inRasterMask = getRaster(inMask); ormvf.pMinDistance = pMinDistance; ormvf.pMaxDistance = pMaxDistance; ormvf.doUseOnlyBorderValues = doUseOnlyBorderValues;