rafaqz / Rasters.jl

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

Geometric lookups and "vector data cubes" #748

Open asinghvi17 opened 1 month ago

asinghvi17 commented 1 month ago

It's super easy to create dimension axes that have vector lookups:

julia> using CairoMakie, DimensionalData
julia> points = Point2f.(1, 1:10)
10-element Vector{Point{2, Float32}}:
 [1.0, 1.0]
 [1.0, 2.0]
 [1.0, 3.0]
 [1.0, 4.0]
 [1.0, 5.0]
 [1.0, 6.0]
 [1.0, 7.0]
 [1.0, 8.0]
 [1.0, 9.0]
 [1.0, 10.0]

julia> rand(X(DimensionalData.Categorical(points)))
╭────────────────────────────────╮
│ 10-element DimArray{Float64,1} │
├────────────────────────────────┴─────────────────────────────────────── dims ┐
  ↓ X Categorical{Point{2, Float32}} [Float32[1.0, 1.0], Float32[1.0, 2.0], …, Float32[1.0, 9.0], Float32[1.0, 10.0]] ForwardOrdered
└──────────────────────────────────────────────────────────────────────────────┘
 Float32[1.0, 1.0]   0.122667
 Float32[1.0, 2.0]   0.670468
 Float32[1.0, 3.0]   0.398075
 Float32[1.0, 4.0]   0.275494
 Float32[1.0, 5.0]   0.437735
 Float32[1.0, 6.0]   0.333305
 Float32[1.0, 7.0]   0.43346
 Float32[1.0, 8.0]   0.740607
 Float32[1.0, 9.0]   0.48268
 Float32[1.0, 10.0]  0.221809

but some questions arise here.

  1. How do we plot them? For scatter I could see just using specapi to forward color as the values of the dimvector. What about 2d arrays?
  2. For e.g. weather station data, you would have a time axis and a position or geometry axis. How would people use data from here? What are common access patterns?
  3. Can we make selectors like Contains, Touches, .., At, etc work, potentially using GeometryOps? Do we want an extension for this? I think this will involve geometric predicates at least.
  4. Do we want to construct spatial indices for such geometric lookups?
  5. Can we set GI.geometrycolumns of a raster to point to geometric lookups? Do we want to stash this in the metadata?
rafaqz commented 1 month ago

Dimensional data.jl already handles point lookups a bit, see mergedims. But polygons etc need more work. Probably we make a GeometryLookups.jl package so everyone can use them?

rafaqz commented 1 month ago
  1. scatter default is good for points with color as the values. if there is another dimension like time we plot points in 3d? and the same for lines/polygons etc
  2. This already works for points with MergedLookup, we can use that syntax for everything else too
  3. See merged.jl for the point implementation, for polygons we would need special lookups.
  4. That could be better than what MergedLookup does currently. Its slow.
  5. Yep we can do that. One hitch is that geometry can be either values in the array or in a lookup, but we can detect that. (a confusing part of the "vector data cube" talk is people use and discuss both but theyre totally different things)
asinghvi17 commented 1 month ago

I just realized that geometric lookups are a great way to represent discretized global grids, since each grid cell is just an integer - a DGGSLookup could be a subtype of some AbstractGeometricLookup, and with a nice API we could have super easy plotting / subsetting / rasterization to a Cartesian grid!

Screenshot 2024-09-20 at 4 09 04 PM

(image from the DGGS talk at SDSL 2024)

rafaqz commented 1 month ago

For sure. Except for the nesting, and there should be optimisations for finding polygons on these grids?