Add tests for unitful cellarea
#799
Conversation
|
Hmm with units, probably we need to check they are the same units as the crs and otherwise convert. You can (We will need a Unitful.jl extension) I think if the crs units are "US survey foot" and the lookup has meter units we can just error. I guess a dumb first round check is is it in "meters" and then do the |
|
That sounds like a way bigger effort :D, but doable at some point via Proj. It might be nice for reproject, actually, since we could discover the CRS's units (maybe stash that in lookup metadata) and this would all optimize away at compile-time. Then indexing with arbitrary units could just work TM via Unitful's automatic conversion. Most of the unit finding work would happen in Proj. |
| # The tests below fail because Unitful apparently can't convert to Float64... | ||
| # (see https://github.com/PainterQubits/Unitful.jl/issues/742) | ||
| # But we also have to re-convert back to the original unit type, which cellarea currently | ||
| # doesn't do. | ||
|
|
There was a problem hiding this comment.
| # The tests below fail because Unitful apparently can't convert to Float64... | |
| # (see https://github.com/PainterQubits/Unitful.jl/issues/742) | |
| # But we also have to re-convert back to the original unit type, which cellarea currently | |
| # doesn't do. | |
| # This needs fixes for unit handling, probably using `ustrip(u"m", x)` in a Unitful.jl extension. That will also require checking the units of the crs. | |
Unitful is like a free test suite for physical models, not converting is by design!
If we are going to make units work here we need to do ustrip(u"m", x) to ensure we are actually getting metres coming in. If people pass in kilometers and we convert back and forth with Proj.jl treating them as meters in the middle, we are really breaking the contract.
There was a problem hiding this comment.
That makes sense. Maybe we have a maybe_ustrip(unit, x) defined in Rasters, its definition for Unitful unit in RastersUnitfulExt, and its definition for Proj units in RastersProjUnitfulExt?
(This is getting to be a mouthful..)
There was a problem hiding this comment.
Yeah those names! But makes sense. Somehow we need some mapping from Proj units to Unitful.
There was a problem hiding this comment.
Why does conversion from Proj to Unitful have to live in Rasters.jl? Shouldn't this be in a ProjUnitfulExt?
There was a problem hiding this comment.
Hmm, we could introduce a nice interface to batch transform an array, but do it out of place in Proj.jl. Then the Unitful stuff could also live there.
There was a problem hiding this comment.
I was initially going to put that here because this calls the C-API directly, but if it doesn't have to then there's no reason to stop this from calling Proj.
| ux_3857 = rebuild(x_3857; val = rebuild(val(x_3857); data = val(x_3857) .* Unitful.u"m", span = Regular(val(x_3857).span.step * Unitful.u"m"))) | ||
| uy_3857 = rebuild(y_3857; val = rebuild(val(y_3857); data = val(y_3857) .* Unitful.u"m", span = Regular(val(y_3857).span.step * Unitful.u"m"))) | ||
| unitful_dimz_3857 = (ux_3857, uy_3857) |
There was a problem hiding this comment.
This is also super complicated, should we have DD overloads to make constructing lookups with unitful easier?
There was a problem hiding this comment.
I think set will just work this out for us, no need to rebuild manually?
Does what it says on the tin. Unfortunately it looks like Unitful can't
convertto Float64 (see PainterQubits/Unitful.jl#742) so going to Proj is where this fails.It also occurs to me that we could try to promote and rebuild e.g. bignumbers or Unitful quantities when going into
reproject. That would allow unitful axes to pass through reproject unchanged. But that's a different PR.