manaakiwhenua / rhealpixdggs-py

rHEALPixDGGS-py implements code to both define an rHEALPix DGGS Reference System and to perform topological queries on its Identfiers. Our roadmap is for v1.0 to be fully compliant with OGC Topic 21 v2.0 / ISO 19170-1:2020.
GNU Lesser General Public License v3.0
18 stars 4 forks source link

R88676887 #15

Open ChocopieKewpie opened 11 months ago

ChocopieKewpie commented 11 months ago

Trying to draw a geometry for resolution 8 rhealpix for the entirety of New Zealand, getting issues with this particular cell ID.

R88676887

It tries to return coordinates that are nowhere near New Zealand. Error: input coordinates (-1.570796, 3.778075) are out of bounds


File ~\AppData\Local\miniconda3\envs\h3panda_exp\Lib\site-packages\rhealpixdggs\dggs.py:2245 in vertices
    result = [self.rdggs.rhealpix(*p, inverse=True) for p in result]

  File ~\AppData\Local\miniconda3\envs\h3panda_exp\Lib\site-packages\rhealpixdggs\dggs.py:2245 in <listcomp>
    result = [self.rdggs.rhealpix(*p, inverse=True) for p in result]

  File ~\AppData\Local\miniconda3\envs\h3panda_exp\Lib\site-packages\rhealpixdggs\dggs.py:453 in rhealpix
    return f(u, v, inverse=inverse)

  File ~\AppData\Local\miniconda3\envs\h3panda_exp\Lib\site-packages\rhealpixdggs\projection_wrapper.py:121 in __call__
    lam, phi = f(u, v, radians=radians, inverse=True)

  File ~\AppData\Local\miniconda3\envs\h3panda_exp\Lib\site-packages\rhealpixdggs\pj_rhealpix.py:507 in f
    lam, phi = array(

TypeError: iteration over a 0-d array
ChocopieKewpie commented 11 months ago

@ndemaio

rggibb commented 11 months ago

Unfortunately this is a recurrence of Issue #6 which wwe thought was resolved by this commit. The commentary for the commit describes the issue and the approach used to resolve it. But also note that the accceptance testing was premised by: 'I can confirm those two cells compute properly now. But only after the source code is updated as described in https://github.com/manaakiwhenua/rhealpixdggs-py/issues/10'. So it might be wiorth double checking that the version of code yopu are using has both of those changes.

ChocopieKewpie commented 11 months ago

Hi Robert,

I've been testing out @ndemaio 's attempt to bring the library to python 3.11. (https://github.com/manaakiwhenua/rhealpixdggs-py/tree/demaion/raster2dggs)

The problem seems to be back, perhaps the solution to #6 is somehow missed?

Attempting to compute cell boundaries for "R88888", which is supposed to be solved by the commit (https://github.com/manaakiwhenua/rhealpixdggs-py/commit/fd4e75a2bcc6830b7d82d7f1ef69603c3d7d6ed8), throws an error. Makes me think that we are missing something in our attempt to bring the code to 3.11.

Will see what went wrong.

rggibb commented 11 months ago

Ok .. the problem arises (ie prior to the #6 solution) because at the point in the code where the projection code is called to determine the lat/long it tests the unprojected coord to test whether to use the the Colignon or the Equal Area Lambert cylindrical projection. For any coord on the boundary between the two projections, rounding errors can result in the wrong projection being chosen, so the solution, was to add the CellID to the calling parameters, so that the projection test can be based on the CellID instead of coords. The specific code changes are tracked in tne commit I mention in my previous comment, and it should be fairly quickly evident if the code you were using has got that change in it or not. Hasppy to have a chat if necessary

alpha-beta-soup commented 2 days ago

@ChocopieKewpie can you please provide a full example to demonstrate the issue? I attempted a basic replication but not currently seeing any issues. (I'm using the bleeding edge.)

Attempting to compute cell boundaries for "R88888", which is supposed to be solved by the commit (https://github.com/manaakiwhenua/rhealpixdggs-py/commit/fd4e75a2bcc6830b7d82d7f1ef69603c3d7d6ed8), throws an error.

>>> from rhealpixdggs.dggs import RHEALPixDGGS
>>> rdggs = RHEALPixDGGS()
>>> c = rdggs.cell(('R', 8,8,8,8,8))
>>> c.boundary()
[(np.float64(19973926.00295788), np.float64(-4962593.986301907)), (np.float64(20015109.3555413), np.float64(-4962593.986301907)), (np.float64(20015109.3555413), np.float64(-5003777.338885326)), (np.float64(19973926.00295788), np.float64(-5003777.338885326))]
>>> c.boundary(plane=False)
[(np.float64(179.62962962962962), np.float64(-41.51722626880009)), (np.float64(-180.0), np.float64(-41.51722626880009)), (np.float64(-180.0), np.float64(-41.93785391016015)), (np.float64(179.62962962962962), np.float64(-41.93785391016015))]

Trying to draw a geometry for resolution 8 rhealpix for the entirety of New Zealand, getting issues with this particular cell ID. R88676887

>>> c = rdggs.cell(('R', 8,8,6,7,6,8,8,7))
>>> c.boundary(plane=False)
[(np.float64(171.45404663923182), np.float64(-41.92222652426085)), (np.float64(171.46776406035664), np.float64(-41.92222652426085)), (np.float64(171.46776406035664), np.float64(-41.93785391016015)), (np.float64(171.45404663923182), np.float64(-41.93785391016015))]