feat: Add regrid functionality based on match geometry#1597
feat: Add regrid functionality based on match geometry#1597eliascapriles-NOAA wants to merge 52 commits into
Conversation
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## main #1597 +/- ##
==========================================
+ Coverage 85.58% 85.61% +0.03%
==========================================
Files 79 78 -1
Lines 6998 7169 +171
==========================================
+ Hits 5989 6138 +149
- Misses 1009 1031 +22
Flags with carried forward coverage won't be shown. Click here to find out more. ☔ View full report in Codecov by Sentry. 🚀 New features to boost your workflow:
|
|
Hey @eliascapriles-NOAA : Thanks for the PR! Below are some comments based on quick look:
|
|
Sounds good. I will get started on implementing the changes ! |
|
Hi Wu-Jung sorry for the delay ! The tests uncovered a couple of bugs that I wanted to fix before submitting for a new PR. I have added unit tests for my helper function, and integration test using SV data for the regridding function. Let me know if there are any changes ! |
Merged changes in main to regrid for PR
for more information, see https://pre-commit.ci
|
Thanks @eliascapriles-NOAA ! I'm going to ask @LOCEANlloydizard to help review this since you had most of the discussions with him. @LOCEANlloydizard - Could you please take a look? feel free to ping me for discussions. Thanks! |
|
Sounds good ! Thanks @leewujung |
LOCEANlloydizard
left a comment
There was a problem hiding this comment.
Hey @eliascapriles-NOAA, thx for the PR!
Following our discussion, a small recap (you probably noted more!):
- the errors from CI are now gone with the merge of pandas < 3, just need to address the assertionError
- remove the +20 padding at the end of the function
- double-check that pings are actually aligned before regridding (there’s already a helper in echopype that does this, see getting_started notebook)
- update the notebook to call the function from your PR, so I can run it on my side too
- we can go from there for other points !
One thing I wanted to flag more generally: right now the regridding is done by looping over channel in Python, and inside that running an apply_ufunc that vectorises over ping_time × range_sample.
This is more a design question than a bug: are we happy with looping over channels in Python and vectorising over pings with apply_ufunc, or should we try to make the whole regridding run as one channel-aware vectorised operation? I know you looked into it.., and i could have a look as well! and maybe @leewujung would have recommendations for this?
Cheers!
|
Hey @eliascapriles-NOAA @LOCEANlloydizard : Not sure if I am really making a good suggestion, but since |
|
Hi @leewujung LLoyd brought this up yesterday. The reason I currently have the channels in a loop is because my apply_ufunc is parallelizing the function across ping_time as specificed by the Echoview algorithm. However, I will try to rework my function to parallelize across the channel dimension as well |
| ) | ||
|
|
||
| if (target_channel is not None) == (target_grid is not None): | ||
| raise ValueError("Provide only one of target_channel or target_grid, not both.") |
There was a problem hiding this comment.
Shouldn't the == be and ?
Otherwise, if target_channel is None and target_grid is also None, these 2 checks would be equal and then the value error would not make sense.
This is not caught by the current tests also, so would be good to add a test for this.
There was a problem hiding this comment.
I've actually removed this line as it was a duplicate of the one below
if (target_channel is None) == (target_grid is None):
So I think it's OK like this? we're raising if both are None (True == True) or if both are provided (False == False)? And we have a small test (test_resample_requires_exactly_one_target) but I've added the match for the error message, thanks!
| if (target_channel is None) == (target_grid is None): | ||
| raise ValueError("Provide exactly one of target_channel or target_grid.") |
There was a problem hiding this comment.
| if exist_reversed_time(ds_Sv, "ping_time"): | ||
| coerce_increasing_time(ds_Sv) |
There was a problem hiding this comment.
Suggest adding a user warning here to let them know that there are reversed ping times and they have been corrected.
There was a problem hiding this comment.
i've added
warnings.warn(
"Reversed ping_time values detected. The dataset has been reordered "
"to increasing ping_time before resampling.",
UserWarning,
)
but we don't have a test for this.. should we?
| def regrid(): | ||
| return 1 | ||
| def resample_to_geometry( | ||
| ds_Sv, target_variable: str, target_channel: str | None = None, target_grid: xr.DataArray = None |
There was a problem hiding this comment.
target_variable is missing in the Parameters list.
Also, should we set its default to "Sv"?
There was a problem hiding this comment.
thx! i've added target_variable and set it to default!
| if target_variable in LOG_VARIABLES: | ||
| source_linear = _log2lin(ds_source) | ||
| else: | ||
| source_linear = ds_source |
There was a problem hiding this comment.
This assumes that unless the user wants to regrid a variable in LOG_VARIABLES, the regridding will just be done on the provided values of the variable. Since we allow users to use this function on arbitrary target_variable, I feel we should expose a boolean input argument to give user control of whether the variable they want to regrid should be done in the linear or log domain.
There was a problem hiding this comment.
yes agreed! I've replaced LOG_VARIABLES with an is_log boolean argument, giving users control over the resampling domain!
| ) | ||
|
|
||
| # Convert back to log domain | ||
| result_linear = result_linear.where(result_linear > 0) |
There was a problem hiding this comment.
I am curious -- what are the conditions when result_linear would be =0 or <0? Given that this is a weight interpolation, the result_linear would always be >=0, right? So is this just to remove =0 results? and if so, why?
There was a problem hiding this comment.
thx! i now moved result_linear = result_linear.where(result_linear > 0) only inside the is_log condition to prevents zero linear values from becoming -inf when converted back to the log domain.. what do you think?
| new_ds = xr.Dataset( | ||
| data_vars={ | ||
| target_variable: ds_combined, | ||
| "echo_range": echo_range_aligned, | ||
| "frequency_nominal": ds_Sv["frequency_nominal"], | ||
| } | ||
| ) |
There was a problem hiding this comment.
We should also have channel as a dimension? Other parts of echopype anchors on using channel as the main coordinate (since there could be 2 different channels of the same frequency, say in a echosounder configuration containing 2 transducers of the same frequencies pointing in 2 different directions), so it'd be good to do that here as well.
There was a problem hiding this comment.
thx! i think channel is preserved as a dimension through xr.concat? (ds_combined = xr.concat(aligned_arrays, dim="channel")) but I’ve now also made it explicit as a coordinate in the output dataset!
| range_var_attrs = np.round(target_range_da.diff("range_sample").mean().values, decimals=4) | ||
| if target_channel: | ||
| new_ds[target_variable].attrs["resampling_mode"] = "target_channel" | ||
| new_ds[target_variable].attrs["target_channel"] = target_channel | ||
| new_ds[target_variable].attrs["grid_spacing"] = str(range_var_attrs) + "m" | ||
| else: | ||
| new_ds[target_variable].attrs["resampling_mode"] = "target_grid" | ||
| new_ds[target_variable].attrs["target_grid"] = target_grid | ||
| new_ds[target_variable].attrs["grid_spacing"] = str(range_var_attrs) + "m" |
There was a problem hiding this comment.
This works right now because echopype currently does not produce range vectors that are of unequal spacing (because right now we only allow constant temperature and salinity when computing sound speed to generate Sv). This will change soon (hopefully) as we incorporate the capability to compute Sv based on a sound speed profile (ie varying along range or depth).
I am not sure what is the goal for including the grid_spacing attribute, but renaming it to mean_grid_spacing would be more accurate.
There was a problem hiding this comment.
yep thanks! I’m investigating how the regridding handles irregular grids across pings, and I agree that in this case the mean spacing is not informative.. so I’ve removed the attribute!
| if "echo_range" in ds_Sv: | ||
| new_ds["echo_range"].attrs = ds_Sv["echo_range"].attrs |
There was a problem hiding this comment.
I don't understand the logic here -- the regridding computation is entirely based on echo_range, so is there any cases in which echo_range does not exist in ds_Sv?
There was a problem hiding this comment.
true! removed the condition
| return ping_time_bin_resvalue, ping_time_bin_resunit_label | ||
|
|
||
|
|
||
| def get_valid_max_depth_ping_index( |
There was a problem hiding this comment.
this was a leftover from an earlier implementation before we simplified the logic! I've removed it, thx!
leewujung
left a comment
There was a problem hiding this comment.
@eliascapriles-NOAA @LOCEANlloydizard : Thanks for this PR, awesome to have this in echopype! I think it's in great shape with just a few things needing attention. Please take a look on the inline comments. Some are very minor textual changes and some would require some (small) code changes.
Co-authored-by: Wu-Jung Lee <leewujung@gmail.com>
for more information, see https://pre-commit.ci
|
Hey @leewujung @eliascapriles-NOAA, i've applied the changes and fixed the incoherence/errors! I've also updated the notebook (echostack-org/echopype-examples#101). We should be good, let me know if something seems off! cheers! |
Implements regrid_all_channel function that would allow for users to match the sampling rates of a channels in a ds_Sv dataset to a specific channel