jiayuasu commented on code in PR #1150:
URL: https://github.com/apache/sedona/pull/1150#discussion_r1432974466
##########
common/src/main/java/org/apache/sedona/common/raster/MapAlgebra.java:
##########
@@ -180,6 +174,86 @@ public static GridCoverage2D mapAlgebra(GridCoverage2D
gridCoverage2D, String pi
}
}
+ public static GridCoverage2D mapAlgebra(GridCoverage2D gridCoverage2D,
String pixelType, String script) {
+ return mapAlgebra(gridCoverage2D, pixelType, script, null);
+ }
+
+ private static ColorModel fetchColorModel(ColorModel originalColorModel,
WritableRaster resultRaster) {
+ if (originalColorModel.isCompatibleRaster(resultRaster)) {
+ return originalColorModel;
+ }else {
+ return PlanarImage.createColorModel(resultRaster.getSampleModel());
+ }
+ }
+
+ public static GridCoverage2D mapAlgebra(GridCoverage2D rast0,
GridCoverage2D rast1, String pixelType, String script, Double noDataValue) {
+ if (rast0 == null || rast1 == null || script == null) {
+ return null;
+ }
+ RasterUtils.isRasterSameShape(rast0, rast1);
+
+ RenderedImage renderedImageRast0 = rast0.getRenderedImage();
+ int rasterDataType = pixelType != null ?
RasterUtils.getDataTypeCode(pixelType) :
renderedImageRast0.getSampleModel().getDataType();
+ int width = renderedImageRast0.getWidth();
+ int height = renderedImageRast0.getHeight();
+ // ImageUtils.createConstantImage is slow, manually constructing a
buffered image proved to be faster.
+ // It also eliminates the data-copying overhead when converting raster
data types after running jiffle script.
+ WritableRaster resultRaster =
RasterFactory.createBandedRaster(DataBuffer.TYPE_DOUBLE, width, height, 1,
null);
+
+ ColorModel cmRast0 =
fetchColorModel(renderedImageRast0.getColorModel(), resultRaster);
+ RenderedImage renderedImageRast1 = rast1.getRenderedImage();
+ ColorModel cmRast1 =
fetchColorModel(renderedImageRast1.getColorModel(), resultRaster);
+
+ // todo temp, figure out if needed
+ if (!cmRast0.equals(cmRast1)) {
+ throw new IllegalArgumentException("Color Model did not match.
Provide rasters that has the same properties.");
+ }
+
+ WritableRenderedImage resultImage = new BufferedImage(cmRast0,
resultRaster, false, null);
+ try {
+ String prevScript = previousScript.get();
+ JiffleDirectRuntime prevRuntime = previousRuntime.get();
+ JiffleDirectRuntime runtime;
+ if (prevRuntime != null && script.equals(prevScript)) {
+ // Reuse the runtime to avoid recompiling the script
+ runtime = prevRuntime;
+ runtime.setSourceImage("rast0", renderedImageRast0);
+ runtime.setSourceImage("rast1", renderedImageRast1);
+ runtime.setDestinationImage("out", resultImage);
+ runtime.setDefaultBounds();
+ } else {
+ JiffleBuilder builder = new JiffleBuilder();
+ runtime = builder.script(script)
+ .source("rast0", renderedImageRast0)
+ .source("rast1", renderedImageRast1)
+ .dest("out", resultImage)
+ .getRuntime();
+ previousScript.set(script);
+ previousRuntime.set(runtime);
+ }
+
+ // TODO is it worth to deduplicate code from 236-246 in this
function
Review Comment:
Please fix TODOs
##########
spark/common/src/test/scala/org/apache/sedona/sql/rasteralgebraTest.scala:
##########
@@ -1296,6 +1296,58 @@ class rasteralgebraTest extends TestBaseScala with
BeforeAndAfter with GivenWhen
assertEquals("UNSIGNED_8BITS", result)
}
+ it("Passed RS_MapAlgebra with two raster columns") {
+ val df1 = sparkSession.read.format("binaryFile")
+ .option("recursiveFileLookup", "true")
+ .option("pathGlobFilter", "*.tif*")
+ .load(resourceFolder + "raster")
+ .selectExpr("path", "RS_FromGeoTiff(content) as rast0")
+ val df2 = sparkSession.read.format("binaryFile")
+ .option("recursiveFileLookup", "true")
+ .option("pathGlobFilter", "*.tif*")
+ .load(resourceFolder + "raster")
+ .selectExpr("path", "RS_FromGeoTiff(content) as rast1")
+ .withColumnRenamed("path", "path1")
+ val df = df1.join(df2, df1("path") === df2("path1"),
"inner").select("rast0", "rast1", "path")
+ Seq(null, "b", "s", "i", "f", "d").foreach { pixelType =>
+ val pixelTypeExpr = if (pixelType == null) null else s"'$pixelType'"
+ val dfResult = df.withColumn("rast_2", expr(s"RS_MapAlgebra(rast0,
rast1, $pixelTypeExpr, 'out[0] = rast0[0] * 0.5 + rast1[0] * 0.5;', null)"))
Review Comment:
Can the script have `out[0]`? What is the difference between `out` and
`out[0]`?
##########
common/src/test/java/org/apache/sedona/common/raster/MapAlgebraTest.java:
##########
@@ -436,6 +436,126 @@ public void testCountValue() {
assertEquals(expected, actual);
}
+ @Test
+ public void testMapAlgebra2Rasters() throws FactoryException {
+ Random random = new Random();
+ String[] pixelTypes = {null, "b", "i", "s", "us", "f", "d"};
+ for (String pixelType : pixelTypes) {
+ int width = random.nextInt(100) + 10;
+ int height = random.nextInt(100) + 10;
+ testMapAlgebra2Rasters(width, height, pixelType, null);
+ testMapAlgebra2Rasters(width, height, pixelType, 100.0);
+ testMapAlgebra2RastersMultiBand(width, height, pixelType, null);
+ testMapAlgebra2RastersMultiBand(width, height, pixelType, 100.0);
+ }
+ }
+
+ private void testMapAlgebra2RastersMultiBand(int width, int height, String
pixelType, Double noDataValue) throws FactoryException {
+ GridCoverage2D rast0 = RasterConstructors.makeEmptyRaster(2, "b",
width, height, 10, 20, 1);
+ GridCoverage2D rast1 = RasterConstructors.makeEmptyRaster(2, "b",
width, height, 10, 20, 1);
+ double[] band1 = new double[width * height];
+ double[] band2 = new double[width * height];
+ double[] band3 = new double[width * height];
+ double[] band4 = new double[width * height];
+ for (int i = 0; i < band1.length; i++) {
+ band1[i] = Math.random() * 10;
+ band2[i] = Math.random() * 10;
+ band3[i] = Math.random() * 10;
+ band4[i] = Math.random() * 10;
+ }
+ rast0 = MapAlgebra.addBandFromArray(rast0, band1, 1);
+ rast0 = MapAlgebra.addBandFromArray(rast0, band2, 2);
+ rast1 = MapAlgebra.addBandFromArray(rast1, band3, 1);
+ rast1 = MapAlgebra.addBandFromArray(rast1, band4, 2);
+ GridCoverage2D result = MapAlgebra.mapAlgebra(rast0, rast1, pixelType,
"out = (rast0[0] + rast0[1] + rast1[0] + rast1[1]) * 0.4;", noDataValue);
+ double actualNoDataValue =
RasterUtils.getNoDataValue(result.getSampleDimension(0));
+ if (noDataValue != null) {
+ Assert.assertEquals(noDataValue, actualNoDataValue, 1e-9);
+ } else {
+ Assert.assertTrue(Double.isNaN(actualNoDataValue));
+ }
+
+ int resultDataType =
result.getRenderedImage().getSampleModel().getDataType();
+ int expectedDataType;
+ if (pixelType != null) {
+ expectedDataType = RasterUtils.getDataTypeCode(pixelType);
+ } else {
+ expectedDataType =
rast0.getRenderedImage().getSampleModel().getDataType();
+ }
+ Assert.assertEquals(expectedDataType, resultDataType);
+
+ Assert.assertEquals(rast0.getGridGeometry().getGridToCRS2D(),
result.getGridGeometry().getGridToCRS2D());
+ band1 = MapAlgebra.bandAsArray(rast0, 1);
+ band2 = MapAlgebra.bandAsArray(rast0, 2);
+ band3 = MapAlgebra.bandAsArray(rast1, 1);
+ band4 = MapAlgebra.bandAsArray(rast1, 2);
+ double[] bandResult = MapAlgebra.bandAsArray(result, 1);
+ Assert.assertEquals(band1.length, bandResult.length);
+ for (int i = 0; i < band1.length; i++) {
+ double expected = (band1[i] + band2[i] + band3[i] + band4[i]) *
0.4;
+ double actual = bandResult[i];
+ switch (resultDataType) {
+ case DataBuffer.TYPE_BYTE:
+ case DataBuffer.TYPE_SHORT:
+ case DataBuffer.TYPE_USHORT:
+ case DataBuffer.TYPE_INT:
+ Assert.assertEquals((int) expected, (int) actual);
+ break;
+ default:
+ Assert.assertEquals(expected, actual, 1e-3);
+ }
+ }
+ }
+
+ private void testMapAlgebra2Rasters(int width, int height, String
pixelType, Double noDataValue) throws FactoryException {
+ GridCoverage2D rast0 = RasterConstructors.makeEmptyRaster(1, "b",
width, height, 10, 20, 1);
+ GridCoverage2D rast1 = RasterConstructors.makeEmptyRaster(1, "b",
width, height, 10, 20, 1);
+ double[] band1 = new double[width * height];
+ double[] band2 = new double[width * height];
+ for (int i = 0; i < band1.length; i++) {
+ band1[i] = Math.random() * 10;
+ band2[i] = Math.random() * 10;
+ }
+ rast0 = MapAlgebra.addBandFromArray(rast0, band1, 1);
+ rast1 = MapAlgebra.addBandFromArray(rast1, band2, 1);
+ GridCoverage2D result = MapAlgebra.mapAlgebra(rast0, rast1, pixelType,
"out = (rast0[0] + rast1[0]) * 0.4;", noDataValue);
+ double actualNoDataValue =
RasterUtils.getNoDataValue(result.getSampleDimension(0));
+ if (noDataValue != null) {
+ Assert.assertEquals(noDataValue, actualNoDataValue, 1e-9);
+ } else {
+ Assert.assertTrue(Double.isNaN(actualNoDataValue));
+ }
+
+ int resultDataType =
result.getRenderedImage().getSampleModel().getDataType();
+ int expectedDataType;
+ if (pixelType != null) {
+ expectedDataType = RasterUtils.getDataTypeCode(pixelType);
+ } else {
+ expectedDataType =
rast0.getRenderedImage().getSampleModel().getDataType();
+ }
+ Assert.assertEquals(expectedDataType, resultDataType);
+
+ Assert.assertEquals(rast0.getGridGeometry().getGridToCRS2D(),
result.getGridGeometry().getGridToCRS2D());
+ band1 = MapAlgebra.bandAsArray(rast0, 1);
+ band2 = MapAlgebra.bandAsArray(rast1, 1);
+ double[] bandResult = MapAlgebra.bandAsArray(result, 1);
+ Assert.assertEquals(band1.length, bandResult.length);
+ for (int i = 0; i < band1.length; i++) {
+ double expected = (band1[i] + band2[i]) * 0.4;
+ double actual = bandResult[i];
+ switch (resultDataType) {
+ case DataBuffer.TYPE_BYTE:
+ case DataBuffer.TYPE_SHORT:
+ case DataBuffer.TYPE_USHORT:
+ case DataBuffer.TYPE_INT:
+ Assert.assertEquals((int) expected, (int) actual);
+ break;
+ default:
+ Assert.assertEquals(expected, actual, 1e-3);
Review Comment:
I think we already define a variable in the TestBase as the allowed precision
##########
common/src/main/java/org/apache/sedona/common/raster/MapAlgebra.java:
##########
@@ -180,6 +174,86 @@ public static GridCoverage2D mapAlgebra(GridCoverage2D
gridCoverage2D, String pi
}
}
+ public static GridCoverage2D mapAlgebra(GridCoverage2D gridCoverage2D,
String pixelType, String script) {
+ return mapAlgebra(gridCoverage2D, pixelType, script, null);
+ }
+
+ private static ColorModel fetchColorModel(ColorModel originalColorModel,
WritableRaster resultRaster) {
+ if (originalColorModel.isCompatibleRaster(resultRaster)) {
+ return originalColorModel;
+ }else {
+ return PlanarImage.createColorModel(resultRaster.getSampleModel());
+ }
+ }
+
+ public static GridCoverage2D mapAlgebra(GridCoverage2D rast0,
GridCoverage2D rast1, String pixelType, String script, Double noDataValue) {
+ if (rast0 == null || rast1 == null || script == null) {
+ return null;
+ }
+ RasterUtils.isRasterSameShape(rast0, rast1);
+
+ RenderedImage renderedImageRast0 = rast0.getRenderedImage();
+ int rasterDataType = pixelType != null ?
RasterUtils.getDataTypeCode(pixelType) :
renderedImageRast0.getSampleModel().getDataType();
+ int width = renderedImageRast0.getWidth();
+ int height = renderedImageRast0.getHeight();
+ // ImageUtils.createConstantImage is slow, manually constructing a
buffered image proved to be faster.
+ // It also eliminates the data-copying overhead when converting raster
data types after running jiffle script.
+ WritableRaster resultRaster =
RasterFactory.createBandedRaster(DataBuffer.TYPE_DOUBLE, width, height, 1,
null);
+
+ ColorModel cmRast0 =
fetchColorModel(renderedImageRast0.getColorModel(), resultRaster);
+ RenderedImage renderedImageRast1 = rast1.getRenderedImage();
+ ColorModel cmRast1 =
fetchColorModel(renderedImageRast1.getColorModel(), resultRaster);
+
+ // todo temp, figure out if needed
+ if (!cmRast0.equals(cmRast1)) {
+ throw new IllegalArgumentException("Color Model did not match.
Provide rasters that has the same properties.");
+ }
+
+ WritableRenderedImage resultImage = new BufferedImage(cmRast0,
resultRaster, false, null);
+ try {
+ String prevScript = previousScript.get();
+ JiffleDirectRuntime prevRuntime = previousRuntime.get();
+ JiffleDirectRuntime runtime;
+ if (prevRuntime != null && script.equals(prevScript)) {
+ // Reuse the runtime to avoid recompiling the script
+ runtime = prevRuntime;
+ runtime.setSourceImage("rast0", renderedImageRast0);
+ runtime.setSourceImage("rast1", renderedImageRast1);
+ runtime.setDestinationImage("out", resultImage);
+ runtime.setDefaultBounds();
+ } else {
+ JiffleBuilder builder = new JiffleBuilder();
+ runtime = builder.script(script)
+ .source("rast0", renderedImageRast0)
+ .source("rast1", renderedImageRast1)
+ .dest("out", resultImage)
+ .getRuntime();
+ previousScript.set(script);
+ previousRuntime.set(runtime);
+ }
+
+ // TODO is it worth to deduplicate code from 236-246 in this
function
+ // TODO and code from 160-171?
+
+ runtime.evaluateAll(null);
+
+ // If pixelType does not match with the data type of the result
image (should be double since Jiffle only supports
+ // double destination image), we need to convert the resultImage
to the specified pixel type.
+ if (rasterDataType != resultImage.getSampleModel().getDataType()) {
Review Comment:
What is the band data type after the map algebra? If the map algebra result
is a double, and the original data type is Int, how will the final value be?
Will it round the values?
##########
spark/common/src/test/scala/org/apache/sedona/sql/rasteralgebraTest.scala:
##########
@@ -1296,6 +1296,58 @@ class rasteralgebraTest extends TestBaseScala with
BeforeAndAfter with GivenWhen
assertEquals("UNSIGNED_8BITS", result)
}
+ it("Passed RS_MapAlgebra with two raster columns") {
+ val df1 = sparkSession.read.format("binaryFile")
+ .option("recursiveFileLookup", "true")
+ .option("pathGlobFilter", "*.tif*")
+ .load(resourceFolder + "raster")
+ .selectExpr("path", "RS_FromGeoTiff(content) as rast0")
+ val df2 = sparkSession.read.format("binaryFile")
Review Comment:
No need to read the same dataset twice. Please use self-join:
https://www.sparkcodehub.com/spark-self-join-in-spark-sql-and-dataframe
--
This is an automated message from the Apache Git Service.
To respond to the message, please log on to GitHub and use the
URL above to go to the specific comment.
To unsubscribe, e-mail: issues-unsubscr...@sedona.apache.org
For queries about this service, please contact Infrastructure at:
us...@infra.apache.org
