rafaqz / Rasters.jl

Raster manipulation for the Julia language
MIT License
213 stars 36 forks source link

writing a Regular Intervals to netcdf produces Explicit Intervals #632

Open tiemvanderdeure opened 6 months ago

tiemvanderdeure commented 6 months ago

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
rafaqz commented 6 months ago

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:

  1. [ ] dont write Explicit intervals if we dont have to
  2. [x] ignore bounds matrices on read when we don't need them
  3. [ ] try to convert Explicit to regular in rasterize/mask/zonal etc where this matters
  4. [ ] write rasterization and line-burning algs that can handle irregular grids
tiemvanderdeure commented 6 months ago

Not entirely sure if this is totally fixed, but the MWE I gave no longer errors

rafaqz commented 6 months ago

I think getting Regular instead of Explicit on load will fix 99% of the problem