C2SM / iconarray

R&D Python package for processing ICON data with xarray
https://c2sm.github.io/iconarray/
BSD 3-Clause "New" or "Revised" License
10 stars 2 forks source link

Implement ind_from_latlon using a BallTree #59

Closed leuty closed 1 year ago

leuty commented 1 year ago

I tried using ind_from_latlon on with the new ICON-1 domain of MCH and it was kinda slow (20 Min). I implemented a new version using a BallTree.

I also added functionality to return the n closest points. That is super useful if one want to compare against instruments and avoid double penalties by slight location offsets.

leuty commented 1 year ago

Note, an version for spherical earth can be found here:

https://github.com/Unidata/python-workshop/blob/fall-2016/notebooks/netcdf-by-coordinates.ipynb

AnnikaLau commented 1 year ago

Thanks for the update/improvement! I'm still going through it but one thing I'm not sure of yet, would ind[0] be the closest of the closest points?

leuty commented 1 year ago

@AnnikaLau I think the distinction is made by either passing k=n or k=[n] to cKDTree.query, see here.

Would it help if I re-implemented the verbose option for ind_from_latlon ? Then you can test.

AnnikaLau commented 1 year ago

@AnnikaLau I think the distinction is made by either passing k=n or k=[n] to cKDTree.query, see here.

Would it help if I re-implemented the verbose option for ind_from_latlon ? Then you can test.

Yes I think the verbose option is useful for the case k=1, which we use for the icon-vis examples. So that the users are aware of that the location of the nearest point might be quite different. So it would be great if you could re-implement it for that case.

leuty commented 1 year ago

@AnnikaLau So that the users are aware of that the location of the nearest point might be quite different.

I dont understand. Does it yield a different location for your case?

AnnikaLau commented 1 year ago

@AnnikaLau So that the users are aware of that the location of the nearest point might be quite different.

I dont understand. Does it yield a different location for your case?

It takes the nearest point to the coordinates given. So in the case the grid cells are big, this can be quite a bit off.

AnnikaLau commented 1 year ago

I created a corresponding merge request for icon-vis: https://github.com/C2SM/icon-vis/pull/39

leuty commented 1 year ago

@AnnikaLau I realized that my approach did not work with global domains. I now use another Tree Data Structure and Haversines to compute the distances. The downside is that it adds another dependency: skikit-learn.

AnnikaLau commented 1 year ago

Radian and degree are now mixed up for input and output. I would only give degree as input and output as it is much more intuitive.

AnnikaLau commented 1 year ago

Apart from one last thing, which should be fixed, I think it looks good and should be merged together with https://github.com/C2SM/icon-vis/pull/39. However, the tests are currently failing on icon-vis, so we should wait until this is fixed before we can merge it...