Allow single point queries with Opendatacube

[whisperingpixel]

Description

Currently, the intersection of the feature geometries with the pixels are based on the centroids. This is okay and best for Polygons but points are not loaded at all and for lines the centroid points need to be hit like here:

• Why is the feature needed?

Only polygons can be currently used.

• What should the feature do?
  The make_geocube function (

  semantique/semantique/extent.py

  Lines 184 to 189 in 2244625

  	raster_obj = geocube.api.core.make_geocube(
  	vector_data = vector_obj,
  	measurements = ["index"],
  	output_crs = crs,
  	resolution = tuple(resolution),
  	)["index"]

  ) allows a parameter centroid (default) and all_touch via the rasterize function (https://corteva.github.io/geocube/stable/geocube.html#geocube.rasterize.rasterize_image).

There are three situations:

1. Only Polygons

Use centroid intersection (all_touch = false)

2. Only Lines and/or only Points

Use touch (all_touch = true)

3. Polygons and (Lines or Points)
   This is a bit difficult. The perfect solution would be to use centroid intersection for polygons, otherwise touch (e.g. splitting the collection and merging the raster results). However, we need this feature rather soon (for points), therefore the acceptable solution would be to use centroid intersection in this case as well.

Additional context

[whisperingpixel]

Related: #14

[whisperingpixel]

Unfortunately, #21 changed the loading behaviour of the geometries but did not solve the problem with point-based queries.

[luukvdmeer]

Hmm, could this be an issue with Opendatacube? When I do point queries using the semantique test datacube it works fine (also saving and reading the GeoTIFF).

Do you have the layout of the current EODC datacube such that I can test it there?

[luukvdmeer]

Source of the problem found by @whisperingpixel, thanks! It turns out that deep in the OpenDataCube code there is the assumption of a regular grid and the cell spacing of the first two cells is used to calculate the resolution for the entire grid. Obviously, this is not possible for an AOI consisting of a single point (since the grid only has a single cell in that case).

There is the option to hand over a fallback_resolution. We need to set a reasonable default here, e.g. 1. Since this is only used in the case of a point it will be a single pixel anyways corresponding on the pixel size. However, setting this hard-coded could lead to unwanted effects later. We need further testing and understanding what the effect of different sizes is.

@whisperingpixel just to confirm, this is not an issue when the AOI consists of multiple points, right? In that case the grid should be regular (except when the points are so close to each other they all fall in the same cell).

[whisperingpixel]

@luukvdmeer - Thanks! Yes, I confirm: If there are two or more points it works.

[luukvdmeer]

@whisperingpixel This should work now, as the spatial resolution is now always stored in the arrays in the way ODC expects it. This means the "fallback resolution" is equal to the resolution the user defined when executing the query. Unfortunately I cannot test good if it works since I don't have a working ODC example. So we should confirm this after updating on the workstation.

[luukvdmeer]

@whisperingpixel Just as a reminder, we should confirm this on the workstation/server whenever the latest semantique version is installed (no urgency from my side of course ;)

[whisperingpixel]

Of course! @augustinh22 is currently setting this up for testing.