precompute finite value mask for kdtree query #2201

chrishavlin · 2023-06-27T15:43:22Z

Rationale

scipy>=1.11 now raises an error when querying a kdtree object with coordinate arrays that contain nan or inf values. This PR adds a mask for finite values before that query occurs.

Implications

Closes #2199

chrishavlin · 2023-06-27T15:44:30Z

lib/cartopy/img_transform.py

@@ -268,16 +268,24 @@ def regrid(array, source_x_coords, source_y_coords, source_proj, target_proj,
                                              target_x_points.flatten(),
                                              target_y_points.flatten())

+    # Find mask of valid points before querying kdtree: scipy >= 1.11 errors
+    # when querying nan points, might as well use for pykdtree too.
+    indices = np.zeros(target_xyz.shape[0], dtype=int)


was I right in having the indices array default to 0? I based this on the indices[mask]=0 line after this, wasn't entirely certain.

I think this is fine to use 0 here as you noted because I think it is just settting it to some valid index. But, I think we still need to add these non-finite locations you found here into the mask below as well? Something like mask = ~finite_xyz | (indices >= len(xyz))

lib/cartopy/img_transform.py

lib/cartopy/tests/test_img_transform.py

…prefilter

chrishavlin · 2023-07-03T14:59:44Z

@greglucas I think the actual code changes are ready for (final?) review.

The main outstanding question I have is how to handle tests. To actually test these changes with scipy, I'd need to refactor how the scipy/pykdtree choice happens. I'm happy to do that -- my initial inclination is to add a function to set which package to use (defaulting to pykdtree) so that there is a way to toggle between the two in the test environment. Is it OK to add that to this PR or would it be better as a followup PR?

Edit: it's also hard for me to tell if these test failures are from this PR, so any tips on checking that would be useful...

rcomer · 2023-07-05T20:31:48Z

I may have misunderstood but, for testing, could you monkeypatch img_transform._is_pykdtree to False?

chrishavlin · 2023-07-05T22:11:22Z

I may have misunderstood but, for testing, could you monkeypatch img_transform._is_pykdtree to False?

not exactly -- because the scipy kdtree import only happens on initial module import. So if you monkeypatch _is_pykdtree to False, you get an error cause the scipy kdtree hasn't been imported yet.

chrishavlin · 2023-07-05T22:20:44Z

One option would be to store the kdtree function handle in a variable to allow monkey patching during tests.

In pseudo-code, I was thinking along the lines of:

def _get_kdtree_function(use_pykdtree = True):
    if use_pykdtree and pykdtree is available:
         return True, pykdtree.kdtree.KDTree
    elif use_pykdtree is False and scipy is available:
        return False, scipy.spatial.cKDTree
    else:
          error 

_is_pydktree, _kdtree_func = _get_kdtree_function()

and then using that _kdtree_func handle down in regrid. That would allow you to monkeypatch both the boolean and the function handle in the test environment.

chrishavlin · 2023-07-07T20:03:14Z

OK, found a simple-ish solution. Realized that if the scipy and pykdtree kdtree object handles were imported into a module attribute similar to _is_pykdtree, the kdtree object could be monkeypatched in a simple way during testing (in addition to the the _is_pykdtree flag). I updated the new test to switch between using scipy and pykdtree.

So I don't see anything else to add here at the moment, let me know if anything more is needed (or if there's a better approach to testing both scipy and pykdtree here).

And from what I can tell, the failing tests look to be the same tests that are failing in other PRs, so I don't think this PR is causing any new failures.

chrishavlin · 2023-07-07T20:05:23Z

oops, forgot to setup pre-commit locally... will fix that in a moment.

chrishavlin · 2023-07-07T20:06:37Z

pre-commit.ci autofix

for more information, see https://pre-commit.ci

greglucas

I think this is good to go from a code perspective, I've just made a few minor nits that you can take or leave.

greglucas · 2023-07-12T02:25:56Z

lib/cartopy/tests/test_img_transform.py

+        _ = img_trans.regrid(data, lons, lats, data_trans, target_prj,
+                             target_x, target_y)
+    else:
+        _ = img_trans.regrid(data, lons, lats, data_trans, target_prj,
+                             target_x, target_y)


Can you dedent this and only keep one of them and only do the monkeypatching in the if-block, but the regridding outside?

Yes! Will do

greglucas · 2023-07-12T02:28:38Z

lib/cartopy/img_transform.py

@@ -13,11 +13,11 @@


 try:
-    import pykdtree.kdtree
+    from pykdtree.kdtree import KDTree as _kdtree_handle


Suggested change

from pykdtree.kdtree import KDTree as _kdtree_handle

from pykdtree.kdtree import KDTree as _kdtreeClass

Thoughts on making this explicit that it is a Class you're using below? I wasn't entirely clear what handle meant at first glance.

lib/cartopy/img_transform.py

chrishavlin · 2023-07-12T20:01:40Z

Thanks for the review! Will address the changes in the next few days. (edit: done!)

chrishavlin · 2023-07-13T22:19:32Z

Totally OK with keeping this open til #2209 is merged then merging here to make extra sure that I'm not actually breaking any tests here.

greglucas · 2023-07-14T14:12:37Z

Thank you @chrishavlin! This was a very nice fix.

precompute finite value mask for kdtree query

0236e6b

chrishavlin commented Jun 27, 2023

View reviewed changes

lib/cartopy/img_transform.py Outdated Show resolved Hide resolved

add a test

04b63f6

chrishavlin commented Jun 27, 2023

View reviewed changes

lib/cartopy/tests/test_img_transform.py Outdated Show resolved Hide resolved

chrishavlin added 3 commits June 27, 2023 11:19

linting fixes

e8919e2

add the non-finite points to the mask

7d0eacc

Merge remote-tracking branch 'upstream/main' into scipy_1pt11_kdtree_…

fc0f1f3

…prefilter

valeriupredoi mentioned this pull request Jul 6, 2023

Recipe carvalhais is failing due to incompatibility of cartopy with scipy>=1.11.0 ESMValGroup/ESMValTool#3281

Closed

trexfeathers mentioned this pull request Jul 7, 2023

[CI Bot] environment lockfiles auto-update SciTools/iris#5366

Merged

add _kdtree_handle

175a27d

pre-commit-ci bot and others added 2 commits July 7, 2023 20:06

[pre-commit.ci] auto fixes from pre-commit.com hooks

fd08137

for more information, see https://pre-commit.ci

shorten test name

68b438c

greglucas approved these changes Jul 12, 2023

View reviewed changes

review suggestions

58c5fe5

greglucas merged commit 41880a8 into SciTools:main Jul 14, 2023

greglucas added this to the 0.22 milestone Jul 14, 2023

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

precompute finite value mask for kdtree query #2201

precompute finite value mask for kdtree query #2201

chrishavlin commented Jun 27, 2023

chrishavlin Jun 27, 2023

greglucas Jun 28, 2023

chrishavlin commented Jul 3, 2023 •

edited

Loading

rcomer commented Jul 5, 2023

chrishavlin commented Jul 5, 2023

chrishavlin commented Jul 5, 2023

chrishavlin commented Jul 7, 2023 •

edited

Loading

chrishavlin commented Jul 7, 2023

chrishavlin commented Jul 7, 2023

greglucas left a comment

greglucas Jul 12, 2023

chrishavlin Jul 12, 2023

chrishavlin Jul 13, 2023

greglucas Jul 12, 2023

chrishavlin Jul 13, 2023

chrishavlin commented Jul 12, 2023 •

edited

Loading

chrishavlin commented Jul 13, 2023

greglucas commented Jul 14, 2023

	from pykdtree.kdtree import KDTree as _kdtree_handle
	from pykdtree.kdtree import KDTree as _kdtreeClass

precompute finite value mask for kdtree query #2201

precompute finite value mask for kdtree query #2201

Conversation

chrishavlin commented Jun 27, 2023

Rationale

Implications

chrishavlin Jun 27, 2023

Choose a reason for hiding this comment

greglucas Jun 28, 2023

Choose a reason for hiding this comment

chrishavlin commented Jul 3, 2023 • edited Loading

rcomer commented Jul 5, 2023

chrishavlin commented Jul 5, 2023

chrishavlin commented Jul 5, 2023

chrishavlin commented Jul 7, 2023 • edited Loading

chrishavlin commented Jul 7, 2023

chrishavlin commented Jul 7, 2023

greglucas left a comment

Choose a reason for hiding this comment

greglucas Jul 12, 2023

Choose a reason for hiding this comment

chrishavlin Jul 12, 2023

Choose a reason for hiding this comment

chrishavlin Jul 13, 2023

Choose a reason for hiding this comment

greglucas Jul 12, 2023

Choose a reason for hiding this comment

chrishavlin Jul 13, 2023

Choose a reason for hiding this comment

chrishavlin commented Jul 12, 2023 • edited Loading

chrishavlin commented Jul 13, 2023

greglucas commented Jul 14, 2023

chrishavlin commented Jul 3, 2023 •

edited

Loading

chrishavlin commented Jul 7, 2023 •

edited

Loading

chrishavlin commented Jul 12, 2023 •

edited

Loading