Skip to content

[GH-3214] Deprecate RS_MapAlgebra in favour of Python UDFs - #3212

Draft
jiayuasu wants to merge 9 commits into
apache:masterfrom
jiayuasu:docs/raster-udf-guide
Draft

[GH-3214] Deprecate RS_MapAlgebra in favour of Python UDFs#3212
jiayuasu wants to merge 9 commits into
apache:masterfrom
jiayuasu:docs/raster-udf-guide

Conversation

@jiayuasu

@jiayuasu jiayuasu commented Jul 29, 2026

Copy link
Copy Markdown
Member

Did you read the Contributor Guide?

Is this PR related to a ticket?

What changes were proposed in this PR?

Raster-returning Python UDFs landed in #2956 and ship in 1.9.1, but they are documented only as a ~55-line subsection of docs/tutorial/raster.md — no page under docs/api/sql/, no nav entry, and no link to them from the map algebra docs. Map algebra by contrast has a concept page in nav, a per-function reference page, an index row, and a worked example in the tutorial walkthrough, so a reader comparing the two options sees them very unevenly.

This adds docs/api/sql/Raster-UDF.md (plus the zh translation and both nav entries) as a concept page beside Raster-map-algebra.md, presenting UDFs as the recommended way to process rasters and covering both ways of reading pixels:

  • raster-to-scalar and raster-to-raster, in both the DataFrame form and the register + expr form (the DataFrame form is used throughout test_flexible_bands.py but appeared in no doc)
  • NDVI written as a Jiffle script and as NumPy, side by side
  • two-raster UDFs — the equivalent of RS_MapAlgebra's five-argument form
  • using rasterio inside a UDF: as_rasterio() → compute → with_bands() round-trips for grid-preserving operations such as fillnodata. This works today but was not documented, and the tutorial's only rasterio mention showed ds.read(1) with no return path, which made it look impossible.
  • setting NODATA on the output with the new nodata= argument, and why the inherited default is usually wrong for a derived raster
  • Limits, stated plainly: the output stays on the input's grid, and uint32, int8, int64/uint64 don't survive the round trip cleanly
  • Scala UDFs over GridCoverage2D, which work in both directions and — unlike the Python path — have no restriction on the output grid

Also cross-links Raster-map-algebra.md (+ zh), RS_MapAlgebra.md, Raster-Functions.md, and docs/tutorial/raster.md (+ zh), and fixes the RS_SetBandNoDataValue link path in the new tutorial text — that page lives under Raster-Operators/, not Raster-Band-Accessors/.

Deprecating RS_MapAlgebra (#3214)

With the replacement documented, RS_MapAlgebra is marked deprecated as of v1.9.1, following the precedent set by ST_Envelope_Aggr in 1.8.1. Notices go on the map algebra concept page (+ zh), the RS_MapAlgebra reference page, the Raster-Functions.md index row and the tutorial walkthrough, and @Deprecated goes on the three MapAlgebra.mapAlgebra overloads in sedona-common.

The notices state that the function will be removed in a future version and that existing scripts should be migrated. It stays registered and working for now — rasteralgebraTest still passes at 166 tests — and this PR does not remove anything. The array-based functions on the map algebra page (RS_BandAsArray, RS_Add, RS_NormalizedDifference, …) are not deprecated. The removal itself, which would also drop the Jiffle dependency and the janino shading workaround in common/pom.xml that exists for #1945, is a separate change for a later version.

The docs no longer present the two as alternatives to weigh. The tutorial's map-algebra reference section is now "Per-pixel calculations" and shows the two-raster change-detection case as a UDF, and Raster-Functions.md points at Raster UDFs as the recommended path. The NDVI section keeps both forms, since that one is a migration aid rather than a choice.

The tutorial walkthrough still computes NDVI with RS_MapAlgebra and now carries a deprecation warning pointing at the UDF version, rather than being rewritten — the walkthrough's ndvi-math.svg illustrates the Jiffle script itself, so switching it over is a larger change worth doing on its own.

Fixing NODATA on UDF output (#3213)

A UDF could not set the NODATA of the raster it returned, so recommending UDFs as the replacement for a function that takes noDataValue as an argument would have left a real gap.

GridSampleDimensionSerializer.write() emits a nodata double "for interoperability with Python RasterType", but read() discarded it and rebuilt the sample dimension from the Kryo category blob alone. with_bands() replays those blobs from its source raster, so the value Python declared never reached the JVM. Two symptoms: a mask carved out of a scene silently inherited that scene's NODATA, and for bands added beyond the source's band count Python and the JVM reported different values for the same band.

read() now reconciles the two through the existing RasterUtils.createSampleDimensionWithNoDataValue(). The JVM writer emits RasterUtils.getNoDataValue(), derived from the categories it also writes, so the two always agree on JVM-to-JVM round trips and the helper returns its argument unchanged there — the change is a no-op on that path. A declared NaN means "nothing declared" and leaves the categories alone, which is also the path taken by rasters serialized before this change. This mirrors the reconcileColorModel() that #2956 added for the colour model.

with_bands() gains nodata=, taking a scalar or one value per band, and widened bands now report the value they will actually receive instead of NaN. The NODATA section of the new Raster UDF page is therefore a documented parameter rather than a documented limitation.

How was this patch tested?

Every code example was executed against a local 1.9.1-SNAPSHOT build (-Dspark=3.4 -Dscala=2.12 -Dgeotools) rather than written from reading the source. The stated limits and error messages are observed behaviour:

  • the NDVI UDF, the two-raster UDF, and the fillnodata rasterio round trip all produce single-band rasters with SRID and scale carried over from the input
  • with_bands() on a mismatched shape raises exactly the ValueError quoted on the page
  • int8 -2 reads back as 254; int64/uint64 raise ValueError
  • all three udf.register forms shown work
  • the Scala examples were run as a ScalaTest suite, including one that returns a 7×5 EPSG:4326 raster from a 4×3 EPSG:3857 input, which is what backs the claim that Scala UDFs aren't pinned to the input grid

The two NODATA symptoms in #3213 were reproduced against the unfixed build before the fix was written: a 4→1 band UDF over an input with NODATA 0 on band 1 produced an output with NODATA 0, giving RS_Count(rast, 1, true) = 6 against RS_Count(rast, 1, false) = 12 on a 12-pixel mask; and widening a 4-band raster whose band 4 has NODATA -1 to 8 bands gave Python [nan, nan, nan, -1.0, nan, nan, nan, nan] while SQL reported -1.0 for bands 4–8. Both are covered by the new tests below and no longer reproduce.

Rather than leaving those runs as a one-off, every example on the page is now a committed test, so the documentation cannot drift from what runs:

  • python/tests/raster/test_raster_udf_docs.py — raster to scalar in both the DataFrame and register-for-SQL forms, raster to raster, nodata= as a scalar and as a per-band sequence, the two-raster UDF, the rasterio fillnodata round trip, the grid-mismatch ValueError, and the three dtype claims including int8 -2 reading back as 254. The NDVI test asserts the UDF and RS_MapAlgebra forms agree to 1e-6, so the migration path the page recommends is checked rather than asserted — they differ only by the 1e-10 the UDF adds to the denominator to avoid 0/0.
  • spark/common/src/test/scala/org/apache/sedona/sql/rasterUdfTest.scala — the Scala examples, including the one returning a 7×5 EPSG:4326 raster from a 4×3 EPSG:3857 input.

New tests for the NODATA fix:

  • test_flexible_bands.py gains five cases — inherited default, scalar override, per-band sequence override, Python and the JVM agreeing about widened bands, and sequence-length validation
  • SerdeTest gains two — a declared nodata that disagrees with the categories (patched into the serialized bytes, which is what Python produces) and a raster with no NODATA at all, confirming a declared NaN is not turned into a nodata category

Suites run green: sedona-common at 1196 tests, the Python raster suite at 59, rasteralgebraTest at 166, rasterUdfTest at 3.

uv run mkdocs build produces 15 warnings, identical to the count on unmodified master — the change introduces none. Both locales render, the new nav entry appears in each, and the #limits / #scala-and-java anchors resolve in both. The zh page needed its explicit anchors wrapped in {% raw %} (the convention already used in raster.zh.md) because {#...} otherwise opens an unterminated Jinja comment and the macros plugin replaces the page with a syntax error. Full pre-commit suite passes.

Did this PR include necessary documentation updates?

  • Yes, I am adding a new API. I am using the current SNAPSHOT version number in vX.Y.Z format.
  • Yes, I have updated the documentation.

Raster-returning Python UDFs landed in apache#2956 but are documented only as a
~55-line subsection of the raster tutorial, with no API page, no nav entry, and
no path to them from the map algebra docs. Map algebra meanwhile has a concept
page in nav, a per-function reference page and an index row, so the two options
are presented very unevenly.

Adds docs/api/sql/Raster-UDF.md (plus the zh translation and nav entries)
covering both read paths:

- when a UDF is a better fit than RS_MapAlgebra, and when it isn't
- raster-to-scalar and raster-to-raster, in both the DataFrame and the
  register-plus-expr forms
- NDVI written as Jiffle and as NumPy, side by side
- two-raster UDFs, the equivalent of RS_MapAlgebra's five-argument form
- using rasterio inside a UDF, which round-trips through with_bands() for
  grid-preserving operations and was not documented anywhere
- the limits: output stays on the input's grid, NODATA is inherited rather
  than settable, and three NumPy dtypes that don't survive the round trip
- Scala UDFs over GridCoverage2D, which have no restriction on the output grid

Every example was executed against a local 1.9.1-SNAPSHOT build, including the
Scala ones; the stated limits and error messages are the observed behaviour.
Also cross-links Raster-map-algebra.md, RS_MapAlgebra.md, Raster-Functions.md
and the tutorial, and corrects the RS_SetBandNoDataValue path used in the new
tutorial text (it lives under Raster-Operators, not Raster-Band-Accessors).
Two changes that finish making Python UDFs a complete replacement for
RS_MapAlgebra.

Fix NODATA on rasters returned from Python UDFs (closes apache#3213).
GridSampleDimensionSerializer.write() emits a nodata double "for
interoperability with Python RasterType", but read() discarded it and rebuilt
the sample dimension from the Kryo category blob alone. with_bands() replays
those blobs from its source raster, so the value Python declared never reached
the JVM: a mask carved out of a scene silently inherited that scene's NODATA,
and for bands added beyond the source's band count Python and the JVM reported
different values for the same band.

read() now reconciles the two through the existing
RasterUtils.createSampleDimensionWithNoDataValue(). The JVM writer emits
RasterUtils.getNoDataValue(), derived from the categories it also writes, so the
two always agree on JVM-to-JVM round trips and the helper returns its argument
unchanged there. A declared NaN means "nothing declared" and leaves the
categories alone, which is also the path taken by rasters serialized before this
change. This mirrors the reconcileColorModel() that apache#2956 added for the colour
model.

with_bands() gains nodata=, taking a scalar or one value per band, so a UDF can
state what its output means instead of inheriting. Widened bands now also report
the value they will actually receive rather than NaN, removing the Python/JVM
disagreement.

Deprecate RS_MapAlgebra in favour of Raster UDFs (closes apache#3214), following the
ST_Envelope_Aggr precedent: notices on the concept page, the function reference,
the function index and the tutorial, plus @deprecated on the three
MapAlgebra.mapAlgebra overloads. The function stays registered and working, and
nothing is scheduled for removal. The array-based functions on the map algebra
page are not deprecated. Removing Jiffle, and with it the janino shading
workaround in common/pom.xml, would be a separate decision for a later major
version.

The NODATA section of the new Raster UDF page changes from a documented
limitation to a documented parameter.

Tests: five new cases in test_flexible_bands.py covering inherited default,
scalar and per-band override, Python/JVM agreement when widening, and sequence
length validation; two in SerdeTest covering a declared value that disagrees
with the categories and a raster with no NODATA at all. sedona-common is green
at 1196 tests, the Python raster suite at 47, and rasteralgebraTest at 166.
@jiayuasu jiayuasu changed the title [DOCS] Add a Raster UDF page and link it from map algebra [GH-3214] Deprecate RS_MapAlgebra in favour of Python UDFs Jul 29, 2026
@jiayuasu
jiayuasu marked this pull request as draft July 30, 2026 02:09
@jiayuasu jiayuasu added this to the sedona-1.9.1 milestone Jul 30, 2026
jiayuasu added 7 commits July 29, 2026 19:20
…with tests

Three follow-ups on the review.

The deprecation notices now say RS_MapAlgebra will be removed in a future
version rather than that nothing is scheduled for removal, and tell readers to
migrate. Updated on the concept page (+ zh), the function reference, the
function index row, the tutorial walkthrough (+ zh), and the @deprecated javadoc
on all three MapAlgebra.mapAlgebra overloads.

Dropped the pandas_udf section from the Raster UDF page and its zh translation.

Added tests that mirror every example on the page, so the documentation cannot
drift from what runs:

- python/tests/raster/test_raster_udf_docs.py — raster to scalar in both the
  DataFrame and register-for-SQL forms, raster to raster, nodata= as a scalar
  and as a per-band sequence, the two-raster UDF, the rasterio fillnodata round
  trip, the grid-mismatch ValueError, and the three dtype claims including int8
  -2 reading back as 254
- spark/common/.../rasterUdfTest.scala — the Scala examples, including the one
  that returns a 7x5 EPSG:4326 raster from a 4x3 EPSG:3857 input, which is what
  backs the claim that Scala UDFs are not pinned to the input grid

The NDVI test asserts the UDF and RS_MapAlgebra forms agree to 1e-6, so the
migration path the page recommends is checked rather than asserted. They differ
only by the 1e-10 the UDF adds to the denominator to avoid 0/0.

Python raster suite is now 59 tests, sedona-common 1196, rasterUdfTest 3.
The page opened with a table weighing RS_MapAlgebra against Python UDFs and told
readers to "reach for RS_MapAlgebra when the operation is per-pixel arithmetic".
That framing contradicts the deprecation: UDFs are the way forward, not one of
two options. Removed the section from both locales and replaced it with a short
statement in the intro that UDFs are recommended and RS_MapAlgebra is going
away, pointing at the NDVI comparison for migration.

Same correction elsewhere:

- Raster-map-algebra.md (+ zh) loses the "a Raster UDF is the better fit if you
  need a library" paragraph, which the deprecation warning above it supersedes.
- Raster-Functions.md points at Raster UDFs as the recommended path rather than
  offering the two side by side, while noting the array-based functions are not
  deprecated.
- The tutorial's "Map algebra" reference section becomes "Per-pixel
  calculations" and shows the two-raster change detection as a UDF instead of an
  RS_MapAlgebra query.

The NDVI section keeps both forms — that one is a migration aid, not a choice.
Picks up the proj4sedona 0.1.6 bump (apache#3215), which fixes the CRS metadata
failures in geoparquetIOTests seen on this PR's CI run.
Review follow-ups on the NODATA work.

Encoding the declared nodata needs the sample dimension type of the pixels it
will describe, and read() ran before those pixels were deserialized, so it used
the source raster's type. A byte raster widened to float64 by a UDF with
nodata=-99999 therefore failed deserialization with IllegalArgumentException,
even though the output dtype holds that value fine. Reconciliation moves to
Serde.SerializableState.read(), after the image is available, and takes the type
from its SampleModel. GridSampleDimensionSerializer now returns the declared
value alongside the sample dimension instead of applying it, and RasterUtils
gains a createSampleDimensionWithNoDataValue overload that takes the type
explicitly.

A declared NaN previously left the replayed categories alone, which made the
documented nodata=float("nan") a no-op whenever the source had a NODATA: the
source's value survived into the output and silently skewed RS_Count and zonal
statistics. NaN now drops the category via the existing
RasterUtils.removeNoDataValue, but only when there is a real value to drop, so
categories carrying a NaN nodata — which Sedona already reports as null — are
untouched.

with_bands(nodata=...) accepts NumPy scalars. np.float32(-9999) is not a Python
float, so it fell through to the sequence branch and raised TypeError:
'numpy.float32' object is not iterable.

The fillnodata recipe on the Raster UDF page could not work as written.
as_rasterio() does not propagate NODATA into the MEM dataset, so src.nodata is
None and read_masks() reports every pixel valid, leaving holes unfilled. The
example now builds the mask from as_numpy_masked(), and a warning documents the
gap. The old test only exercised a raster with no NODATA, so it passed on a
no-op; it now marks a value as NODATA and asserts the hole is interpolated while
every other pixel survives.

The tutorial's two-raster block gains the imports it was missing and a comment
explaining the two DataFrames it operates on.

New tests: nodata wider than the source dtype, NaN clearing an inherited value,
NumPy scalar acceptance, and the as_rasterio NODATA gap pinned as a known
limitation. Python raster suite 63, sedona-common 1198, geoparquetIOTests 223.
Review round two.

The two-raster example subtracted raw NODATA sentinels. as_numpy() returns them
verbatim, so a hole in one input became a large bogus difference and a hole in
both cancelled into a plausible zero, and nodata= only labelled the output
without marking which pixels were invalid. The example now reads through
as_numpy_masked(), which substitutes NaN, and writes the sentinel back with
np.where before declaring it. A warning spells out why, since this is the trap
anyone writing multi-raster arithmetic will hit. Fixed on the Raster UDF page,
its zh translation, and both tutorials.

The fillnodata example returned without declaring a NODATA, so the output
inherited the input's marker even though every hole was gone — any filled pixel
landing on that value would have read as invalid. It now returns
nodata=float("nan").

Its test was also too weak to catch that: asserting approx(5.0, abs=2.0) around
a hole whose original value was 5.0 passed whether or not the fill happened. It
now asserts the hole moved off 5.0, that every other pixel is untouched, that
the output declares no NODATA, and that RS_Count skips nothing.

createSampleDimensionWithNoDataValue split range-valued no data categories
around the new value while preserving their name, leaving no data-named
fragments behind. getNoDataValue() reads the first of those, so replacing a
NODATA of [0..10] with 5 still reported 0. Existing no data categories are now
dropped whole before the split.

Rejecting a nodata the output dtype cannot hold moves to with_bands(), where the
dtype is known and the error can name the mistake. A fractional or out-of-range
value on an integral band previously produced a sample dimension that failed on
read with GeoTools' "Illegal argument: value=5.5".

The scalar check narrows from np.number, which accepted np.complex64 and
silently discarded the imaginary part, to real integer and floating types, and
now unwraps 0-d arrays instead of failing on "iteration over a 0-d array".

New tests: holes surviving two-raster arithmetic, range-valued no data
replacement, non-real and 0-d nodata, and integral-dtype rejection.
sedona-common 1200, Python raster 67, Scala raster 169.
Review round three.

Inherited NODATA now goes through the same validation as an explicit nodata=.
Narrowing a float64 scene whose NODATA is -9999 to uint8 used to keep declaring
-9999 — a value no uint8 pixel can hold — so every hole silently read as valid.
with_bands() now raises with guidance instead.

Validation itself now targets the JVM storage type rather than the nominal
NumPy dtype, and reinterprets the two documented storage mismatches the same
way the pixel data is reinterpreted: int8 nodata=-2 is stored and reported as
254, uint32 values above 2^31-1 as their signed equivalent, so the declared
value keeps matching the pixels it marks. float32 values are rounded to the
nearest float32 — declaring 0.1 verbatim could never equal the float32 0.1 the
pixels hold. Non-finite values other than NaN are rejected for every dtype.

Two problems in createSampleDimensionWithNoDataValue:

- Splitting a real-valued data category around an interior no data value built
  exclusive bounds touching the single-valued no data category, which
  GridSampleDimension rejects ("Ranges [v..MAX] and [v..v] overlap"), and a
  fractional value could not be cast to an integral range class at all. Real
  splits now use closed ranges bounded by the adjacent doubles. The split point
  and the category value are both the TypeMap-wrapped value — on a float band
  99.2 encodes as 99.19999694824219, and boundaries computed from the raw
  double land on the wrong side of the category they should bracket.
- The early return kept a range-valued no data category whenever its minimum
  equalled the requested value, so declaring 5 over NODATA [0..10] still
  reported 0. Existing no data categories are dropped whole before the split,
  and the serde reconciliation strips them first, keeping only an exact
  single-valued match untouched.

The fillnodata recipe no longer claims the filling is complete. fillnodata only
interpolates within max_search_distance (100 pixels by default) of valid data,
so a hole larger than that keeps its sentinel cells — clearing the NODATA
declaration, as the recipe used to, turned those cells into valid-looking data.
The recipe now keeps the declaration, with the trade-off explained inline, and
a new test pins unreachable holes staying invalid.

The raw-sentinel warning is no longer confined to the two-raster section: the
"Reading pixel data" section warns that as_numpy() feeds holes into arithmetic
and comparisons as ordinary numbers, the mask example in the NODATA section
reads through as_numpy_masked(), and the NDVI example notes it reads raw values
exactly like the Jiffle script it mirrors. The zh page gains explicit
reading-pixel-data and two-rasters anchors, fixing its broken fragment link.

Tests: one-sided and shared holes through two-raster arithmetic, unreachable
fill holes staying flagged, inherited-narrowing rejection, int8/uint32
reinterpretation observed through RS_BandNoDataValue and RS_Count, float32
coercion, non-finite rejection, range-min reconcile collapse, and real and
widened fractional splits. sedona-common 1204, Python raster 73, Scala raster
169.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

Deprecate RS_MapAlgebra in favour of Python UDFs Python UDF rasters cannot set NODATA, and Python and the JVM disagree about it

1 participant