Invalid pos/index within RasterLayer (geographic spatial location not exposed to GeoSpace / RasterLayer)?

[GondekNP]

Thanks for all your work on this library - I am looking forward to using it for our vegetation modelling efforts!

Describe the bug
The pos attribute of a Cell (within a RasterLayer) does not contain valid x, y coordinate in physical space in my specified CRS. This is despite the fact that the RasterLayer is displayed in the proper location in physical space via the accompanying Solara Visualization.

Expected behavior
Based on the mesa-geo docs of Cell, I expect the 'pos' attribute to represent the x, y location of the cell in physical space, within the CRS specified when the RasterLayer is instantiated (in my case, from the GeoTIFF created by previous instantiation of the RasterLayer and exporting via the to_filemethod.

More generally, I expected some way to query the RasterLayer object to grab relevant Cells based on their location in space. I only noticed this (potential) bug when trying out the various methods for grabbing neighboring cells, and noticing than none of them seemed to use a spatial argument, only indexes of the underlying RasterLayer.

To Reproduce
There is a devcontainer setup within this repo, which was adopted from mesa-geo's rainfall model example. There are a few artifacts of the original rainfall model's behavior here, as I haven't successfully gotten through a complete step of the model yet.

I am using the following, on an Ubuntu-based image running Python 3.12.7 :

  - mesa >= 3.0.3
  - mesa-geo >= 0.9.0a0
  - solara >= 1.41.0

To run, with Docker, VSCode and Extensions - Dev Containers installed, one can run:

Cmd + P + Reopen in Dev Container

From there, one can run the solara app using the provided launch configuration in the VSCode debugger, or run Solara as normal from the terminal (using the already installed mesa conda environment.

Additional context

As alluded to above - I am having some trouble understanding the relationship between the agents (which exist in vector / continuous space) and the underlying cells (which exist in raster / discrete space).

It seems like a very natural need for these modelling approaches to provide an interface between the two (for example, in our vegetation case, the vegetation agents must know something about the aridity of the underlying cells, but the aridity is a feature of the cell, not the agent itself) - which leads me to believe I am misunderstanding something about the design of mesa / mesa-geo in general. Since this seems like such a foundational concept, I am pretty tempted to say this is a misunderstanding of the tools you've built on my part - but, since I am on an alpha release and so much seems to have been updated for mesa 3.0+, I figured I would see if this is truly a bug.

Thank you again for all your efforts and I am excited to give mesa / mesa-geo a shot for our academic applications!

[wang-boyu]

Hi, thanks for the detailed writeup and sorry about the confusion. The rasterlayer was implemented some time ago and it deserves some improvements.

Regarding how it currently works, there are two ways of indexing a cell:

• using cell.pos, which is (x, y), with an origin at the lower left corner of the raster grid.

  This is mainly because the current implementation of RasterLayer was largely copied over from Mesa's Grid space. As such, cells can be indexed with rasterlayer.cells[x][y] where x, y = cell.pos, similar to how to index agents in grid space in Mesa. It can also be used in neighborhood queries such as raster_layer.get_neighboring_cells(cell.pos)

• using cell.indices, which is in (row, col) format with an origin at the upper left corner of the raster grid.

  This is usually how matrix indexing works with for example numpy array. So raster_layer.get_raster("attr_name")[row, col] will return the cell at [row, col] of the numpy matrix.

Based on the mesa-geo docs of Cell, I expect the 'pos' attribute to represent the x, y location of the cell in physical space

This absolutely makes sense. Unfortunately till now there's no built-in method or attribute to retrieve the physical location of a cell in longitude & latitude. This will not work with the rasterlayer.cells[x][y] syntax if its pos is the real physical location.

In Mesa-Geo the physical locations are handled by the geometry attribute of GeoAgent, whereas Cell does not have geometry at the moment (should it be a Point or a square Polygon?). But, the physical locations of cells should be able to be computed by raster_layer.transform * cell.indices (using Rasterio's Transforms). The resulting (x, y) is the physical location (i.e., lon & lat), not to be confused with the cell.pos attribute mentioned above.

This is despite the fact that the RasterLayer is displayed in the proper location in physical space via the accompanying Solara Visualization.

I'm guessing this works because the geometry attribute is set correctly: https://github.com/SchmidtDSE/mesa_abm_poc/blob/600b7495e466b7879b9d6df0974193808211c0a4/vegetation/patch/model.py#L61-L63

I hope this clarifies! We wish to refactor raster layer in the near future (as mentioned in #201). This is the last one or two things that's blocking a v1.0 release of Mesa-Geo.

[GondekNP]

@wang-boyu - Thanks so much for the detailed markup! I appreciate it greatly and this clarifies a lot for me.

Regarding the specifics task of getting the underlying cells, I addressed it using the technique you suggest:

        self.float_indices = ~self.model.space.raster_layer._transform * \
            (np.float64(geometry.x), np.float64(geometry.y))

        self.indices = (int(self.float_indices[0]), int(self.float_indices[1]))

Where self in this case references the GeoAgent. This is performed during initiation and is fine for our use case, as our agents don't move spatial location throughout the simulation (they are plants 😉)

I obviously don't have the specifics for mesa-geo's case, but from my perspective, some primitive linking the GeoAgent to its underlying Cell would be greatly appreciated - if the Cell had a geometry, I could use some operation to find the intersecting cell by iterating over potential cells, but in my case at least my agent doesn't really care about its spatial location or the spatial location of the intersecting cell - only the attributes that are inherent to that cell. Just my two cents and would be happy to be a guinea pig for future abstractions that apply to our use case!

Thanks again - since the bug in question is not a bug, I will close.