Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,9 @@
*/
package org.apache.sedona.common.raster;

import static org.apache.sedona.common.utils.RasterUtils.getRasterScaleY;
import static org.apache.sedona.common.utils.RasterUtils.getRasterUpperLeftY;

import java.awt.image.WritableRaster;
import java.math.BigDecimal;
import java.math.RoundingMode;
Expand Down Expand Up @@ -69,6 +72,11 @@ protected static List<Object> rasterize(
rasterized = coverageFactory.create("rasterized", params.writableRaster, rasterExtent);
}

// If original raster is bottom-up, flip the rasterized result vertically
if (params.bottomUp) {
rasterized = RasterUtils.flipVerticallyGridSpace(rasterized);
}

// Return results compatible with the original function
List<Object> objects = new ArrayList<>();
objects.add(params.writableRaster);
Expand Down Expand Up @@ -136,16 +144,22 @@ private static void rasterizePoint(
int x = startX;
int y = startY;

for (double worldY = geomExtent.getMinY();
worldY < geomExtent.getMaxY();
worldY += params.scaleY, y++) {
for (double worldY = geomExtent.getMaxY();
worldY > geomExtent.getMinY();
worldY += params.scaleY, y--) {
x = startX;
for (double worldX = geomExtent.getMinX();
worldX < geomExtent.getMaxX();
worldX += params.scaleX, x++) {

// Flip y-axis (since raster Y starts from top-left)
int yIndex = -y - 1;
// Make zero indexed
int yIndex = y - 1;

// Adjust for bottom-up rasters
// Reverse the y index
if (params.bottomUp) {
yIndex = params.writableRaster.getHeight() - 1 - yIndex;
}

// Create envelope for this pixel
double cellMaxX = worldX + params.scaleX;
Expand Down Expand Up @@ -182,9 +196,9 @@ private static void rasterizeLineString(
}

double x0 = (start.x - params.upperLeftX) / params.scaleX;
double y0 = (params.upperLeftY - start.y) / params.scaleY;
double y0 = (start.y - params.upperLeftY) / params.scaleY;
double x1 = (end.x - params.upperLeftX) / params.scaleX;
double y1 = (params.upperLeftY - end.y) / params.scaleY;
double y1 = (end.y - params.upperLeftY) / params.scaleY;

// Apply Bresenham for this segment
drawLineBresenham(params, x0, y0, x1, y1, value, 0.2);
Expand Down Expand Up @@ -221,6 +235,12 @@ private static void drawLineBresenham(
int rasterX = (int) (Math.floor(x));
int rasterY = (int) (Math.floor(y));

// Adjust for bottom-up rasters
// Reverse the y index
if (params.bottomUp) {
rasterY = params.writableRaster.getHeight() - 1 - rasterY;
}

// Only write if within raster bounds
if (rasterX >= 0
&& rasterX < params.writableRaster.getWidth()
Expand Down Expand Up @@ -338,9 +358,9 @@ static ReferencedEnvelope rasterizeGeomExtent(

// Using BigDecimal to avoid floating point errors
double upperLeftX = metadata[0];
double upperLeftY = metadata[1];
double upperLeftY = getRasterUpperLeftY(metadata);
double scaleX = metadata[4];
double scaleY = metadata[5];
double scaleY = getRasterScaleY(metadata);

// Compute the aligned min/max values
double alignedMinX =
Expand Down Expand Up @@ -464,13 +484,14 @@ private static RasterizationParams calculateRasterizationParams(
upperLeftY = geomExtent.getMaxY();
} else {
upperLeftX = metadata[0];
upperLeftY = metadata[1];
upperLeftY = getRasterUpperLeftY(metadata);
}

WritableRaster writableRaster;
if (useGeometryExtent) {
int geomExtentWidth = (int) (Math.round(geomExtent.getWidth() / metadata[4]));
int geomExtentHeight = (int) (Math.round(geomExtent.getHeight() / -metadata[5]));
int geomExtentHeight =
(int) (Math.round(geomExtent.getHeight() / -getRasterScaleY(metadata)));

writableRaster =
RasterFactory.createBandedRaster(
Expand All @@ -483,19 +504,30 @@ private static RasterizationParams calculateRasterizationParams(
RasterUtils.getDataTypeCode(pixelType), rasterWidth, rasterHeight, 1, null);
}

boolean bottomUp = metadata[5] > 0;

return new RasterizationParams(
writableRaster, pixelType, metadata[4], -metadata[5], upperLeftX, upperLeftY);
writableRaster,
pixelType,
metadata[4],
getRasterScaleY(metadata),
upperLeftX,
upperLeftY,
bottomUp);
}

private static void validateRasterMetadata(double[] rasterMetadata) {
if (rasterMetadata[4] < 0) {
throw new IllegalArgumentException("ScaleX cannot be negative");
}
if (rasterMetadata[5] > 0) {
throw new IllegalArgumentException("ScaleY must be negative");
throw new IllegalArgumentException(
String.format(
"Invalid raster metadata: scaleX is negative (%.2f). Right-to-left rasters are not supported.",
rasterMetadata[4]));
}
if (rasterMetadata[6] != 0 || rasterMetadata[7] != 0) {
throw new IllegalArgumentException("SkewX and SkewY must be zero");
throw new IllegalArgumentException(
String.format(
"Invalid raster metadata: skewX is %.2f and skewY is %.2f. Both values must be zero.",
rasterMetadata[6], rasterMetadata[7]));
}
}

Expand All @@ -506,20 +538,23 @@ private static class RasterizationParams {
double scaleY;
double upperLeftX;
double upperLeftY;
boolean bottomUp;

RasterizationParams(
WritableRaster writableRaster,
String pixelType,
double scaleX,
double scaleY,
double upperLeftX,
double upperLeftY) {
double upperLeftY,
boolean bottomUp) {
this.writableRaster = writableRaster;
this.pixelType = pixelType;
this.scaleX = scaleX;
this.scaleY = scaleY;
this.upperLeftX = upperLeftX;
this.upperLeftY = upperLeftY;
this.bottomUp = bottomUp;
}
}

Expand Down Expand Up @@ -602,15 +637,15 @@ private static Map<Double, TreeSet<Double>> computeScanlineIntersections(
// Using BigDecimal to avoid floating point errors
double yStart =
Math.round(
(BigDecimal.valueOf(params.upperLeftY)
.subtract(BigDecimal.valueOf(worldP1.y))
(BigDecimal.valueOf(worldP1.y)
.subtract(BigDecimal.valueOf(params.upperLeftY))
.divide(BigDecimal.valueOf(params.scaleY), RoundingMode.CEILING))
.doubleValue());

double yEnd =
Math.round(
(BigDecimal.valueOf(params.upperLeftY)
.subtract(BigDecimal.valueOf(worldP2.y))
(BigDecimal.valueOf(worldP2.y)
.subtract(BigDecimal.valueOf(params.upperLeftY))
.divide(BigDecimal.valueOf(params.scaleY), RoundingMode.FLOOR))
.doubleValue());

Expand All @@ -619,7 +654,7 @@ private static Map<Double, TreeSet<Double>> computeScanlineIntersections(
yStart = Math.min((params.writableRaster.getHeight() - 0.5), Math.abs(yStart) - 0.5);

double p1X = (worldP1.x - params.upperLeftX) / params.scaleX;
double p1Y = (params.upperLeftY - worldP1.y) / params.scaleY;
double p1Y = (worldP1.y - params.upperLeftY) / params.scaleY;

if (worldP1.x == worldP2.x) {
// Vertical line case: directly set xIntercept. Avoid divide by zero error when
Expand Down Expand Up @@ -695,6 +730,13 @@ private static void fillPolygon(

for (Map.Entry<Integer, List<int[]>> entry : scanlineFillRanges.entrySet()) {
int y = entry.getKey();

// Adjust for bottom-up rasters
// Reverse the y index
if (params.bottomUp) {
y = params.writableRaster.getHeight() - 1 - y;
}

for (int[] range : entry.getValue()) {
if (range.length == 1) {
params.writableRaster.setSample(range[0], y, 0, value);
Expand Down
166 changes: 166 additions & 0 deletions common/src/main/java/org/apache/sedona/common/utils/RasterUtils.java
Original file line number Diff line number Diff line change
Expand Up @@ -984,4 +984,170 @@ public static GridCoverage2D applyRasterMask(GridCoverage2D raster, GridCoverage
false);
return modifiedRaster;
}

/**
* This function calculates the upper-left Y world coordinate of the raster. It accounts for both
* bottom-up and top-down raster conventions.
*
* <p>For bottom-up rasters (metadata[5] > 0), the upper-left Y world coordinate is directly taken
* from metadata[1]. - For top-down rasters (metadata[5] < 0), the upper-left Y coordinate is
* adjusted by adding the product of the raster height and the scale factor (metadata[5]).
*
* @param metadata The metadata array containing raster properties.
* @return The upper-left Y world coordinate of the raster.
*/
public static double getRasterUpperLeftY(double[] metadata) {
if (metadata[5] < 0) {
return metadata[1];
} else {
int rasterHeight = (int) metadata[3];
return metadata[1] + (rasterHeight * metadata[5]);
}
}

/**
* This function calculates the Y scale factor of the raster. It ensures that the scale factor is
* always negative, regardless of the raster convention.
*
* <p>For bottom-up rasters (metadata[5] > 0), the scale factor is returned as-is. - For top-down
* rasters (metadata[5] < 0), the scale factor is negated to maintain consistency.
*
* @param metadata The metadata array containing raster properties.
* @return The Y scale factor of the raster.
*/
public static double getRasterScaleY(double[] metadata) {
if (metadata[5] < 0) {
return metadata[5];
} else {
return -metadata[5];
}
}

/**
* GRID SPACE: Grid space refers to the raster's index space (col, row). Pixel values remain at
* the same array indices, but the grid-to-CRS affine transform is rewritten.
*
* <p>EFFECT: - The raster is reinterpreted spatially (top-down ↔ bottom-up). - Pixel values are
* NOT moved in memory. - Pixel world coordinates DO change. - The world envelope (bounding box)
* is preserved.
*
* <p>Use this when you want to change how the raster is oriented in CRS space without touching
* the pixel data.
*/
public static GridCoverage2D flipVerticallyGridSpace(GridCoverage2D raster) {

GridGeometry2D gridGeometry = raster.getGridGeometry();

MathTransform mt = gridGeometry.getGridToCRS2D();
if (!(mt instanceof AffineTransform2D)) {
throw new UnsupportedOperationException(
"flipVertically only supports AffineTransform2D grid-to-CRS transforms");
}

AffineTransform2D affine = (AffineTransform2D) mt;

RenderedImage image = raster.getRenderedImage();
int height = image.getHeight();

/*
* Affine transform definition (GeoTools / GDAL compatible):
*
* worldX = scaleX * col + shearX * row + translateX
* worldY = shearY * col + scaleY * row + translateY
*/
double scaleX = affine.getScaleX();
double shearX = affine.getShearX();
double shearY = affine.getShearY();
double scaleY = affine.getScaleY();
double translateX = affine.getTranslateX();
double translateY = affine.getTranslateY();

/*
* Apply grid-space vertical flip:
*
* row = (H - 1) - row'
*/
double newShearX = -shearX;
double newScaleY = -scaleY;

double newTranslateX = translateX + shearX * (height - 1);
double newTranslateY = translateY + scaleY * (height - 1);

AffineTransform2D flippedAffine =
new AffineTransform2D(scaleX, shearY, newShearX, newScaleY, newTranslateX, newTranslateY);

GridGeometry2D newGridGeometry =
new GridGeometry2D(
gridGeometry.getGridRange(),
flippedAffine,
gridGeometry.getCoordinateReferenceSystem());

return RasterUtils.clone(
raster.getRenderedImage(),
newGridGeometry,
raster.getSampleDimensions(),
raster,
null,
true);
}

/**
* PIXEL SPACE: Pixel space refers to pixel values anchored at their world coordinates. Pixel
* values are physically rearranged so that each value remains at the same CRS location.
*
* <p>EFFECT: - Raster rows are reversed in memory. - The affine transform is updated accordingly.
* - Pixel world coordinates are preserved. - The world envelope (bounding box) is preserved.
*
* <p>Use this when pixel values must remain fixed at their geographic locations.
*/
public static GridCoverage2D flipVerticallyPixelSpace(GridCoverage2D raster) {
GridGeometry2D gridGeometry = raster.getGridGeometry();

MathTransform mt = gridGeometry.getGridToCRS2D();
if (!(mt instanceof AffineTransform2D)) {
throw new UnsupportedOperationException(
"flipVertically only supports AffineTransform2D grid-to-CRS transforms");
}

AffineTransform2D affine = (AffineTransform2D) mt;

RenderedImage image = raster.getRenderedImage();
int height = image.getHeight();

// Flip the affine transformation
double scaleX = affine.getScaleX();
double shearX = affine.getShearX();
double shearY = affine.getShearY();
double scaleY = affine.getScaleY();
double translateX = affine.getTranslateX();
double translateY = affine.getTranslateY();

double newShearX = -shearX;
double newScaleY = -scaleY;
double newTranslateX = translateX + shearX * (height - 1);
double newTranslateY = translateY + scaleY * (height - 1);

AffineTransform2D flippedAffine =
new AffineTransform2D(scaleX, shearY, newShearX, newScaleY, newTranslateX, newTranslateY);

// Create a new WritableRaster to preserve pixel values
Raster originalRaster = RasterUtils.getRaster(image);
WritableRaster flippedRaster = originalRaster.createCompatibleWritableRaster();

// Copy pixel values to the flipped raster
for (int y = 0; y < height; y++) {
for (int x = 0; x < image.getWidth(); x++) {
flippedRaster.setPixel(x, height - y - 1, originalRaster.getPixel(x, y, (double[]) null));
}
}

GridGeometry2D newGridGeometry =
new GridGeometry2D(
gridGeometry.getGridRange(),
flippedAffine,
gridGeometry.getCoordinateReferenceSystem());

return RasterUtils.clone(
flippedRaster, newGridGeometry, raster.getSampleDimensions(), raster, null, true);
}
}
Loading
Loading