Rasters.jl
Rasters.jl copied to clipboard
writing a Regular Intervals to netcdf produces Explicit Intervals
Writing a Raster
to a netcdf and reading it back in again yields a Raster
with Explicit Intervals
, even if the raster originally had Regular Intervals
.
E.g.
using Rasters, NCDatasets
x = X(-50:1:50; sampling=Rasters.Intervals());
y = Y(-50:1:50; sampling=Rasters.Intervals());
ras = Raster(rand(x, y))
write("myraster.nc", ras; force = true)
ras_nc = Raster("myraster.nc")
This isn't always trivial. For example, rasterizing a polygon to ras
works, but to ras_nc
does not work.
using ArchGDAL
pointvec = [
(-20.0, 30.0),
(-20.0, 10.0),
(0.0, 10.0),
(0.0, 30.0),
(-20.0, 30.0),
]
poly = ArchGDAL.createpolygon(pointvec)
rasterize(poly; to = ras, fill = true)
rasterize(poly; to = ras_nc, fill = true) # errors
Yes we should try to minimise how much we get Explicit
lookups. But also we should try to convert Explicit
to Regular
before we rasterize as sometimes that will be possible.
Actually rasterizing into Explicit
and Irregular intervals is harder and slow, because e.g. line burning is assuming a straight line accross the grid so it does the angle calculation once.
Maybe in practice most lines only cross one or a few pixels, so it wont make much difference?
I think for normal center
rasterization it wont matter as much, just a level of indirection to check the current position from a vector as we step along the rows. We have very good performance so it might even still be relatively fast.
So:
- [ ] dont write
Explicit
intervals if we dont have to - [x] ignore bounds matrices on read when we don't need them
- [ ] try to convert
Explicit
to regular in rasterize/mask/zonal etc where this matters - [ ] write rasterization and line-burning algs that can handle irregular grids
Not entirely sure if this is totally fixed, but the MWE I gave no longer errors
I think getting Regular
instead of Explicit
on load will fix 99% of the problem