Added snap_segments - Githubissues

jploveless commented 2 years ago

Called as first line in process_segment, so that newly created segments have additional properties added to them Other modifications here:

japan_segment_new.csv includes a correction to which segments are replaced by the Nankai mesh (one is missing in japan_segment.csv)
japan.msh includes a single edited coordinate, so that its southeast corner is perfectly coincident with the east corner of sagami.msh
notebooks/celeri_dense_new.ipynb runs through using snap_segments. This NB needed to be modified to recalculate elastic partials, because the number of segments changes after running snap_segments. Results here look nearly identical to original, which is good.

review-notebook-app[bot] commented 2 years ago

Check out this pull request on

See visual diffs & provide feedback on Jupyter Notebooks.

Powered by ReviewNB

jploveless commented 2 years ago

Hmm…tested on the Japan model, but there are big problems with the global model, so I'm working on it now

brendanjmeade commented 2 years ago

This is exciting and it's cool that the global closure test caught something!

jploveless commented 2 years ago

I'm a little confused why it's still failing. I thought the key issue was a hanging segment that resulting from using the full global_segment.csv file, which references all of the global meshes, but only reading in the single Cascadia mesh as specified in the global_mesh_parameters.json file. Previously, snap_segments would connect Nankai segments (the first mesh referenced in global_segment.csv) with Cascadia segments (the only, and therefore the first, mesh in global_mesh_parameters.json). This actually made sense from an indexing perspective but is clearly wrong. Now, I can create an adjusted segment structure that conforms to the Cascadia mesh shape, but the tests are still failing.

If I run the tests locally, they pass.

brendanjmeade commented 2 years ago

@jploveless Thanks for this. I've looked at the error message a bit. The error appears to be that the call to np.allclose is comparing arrays of two different lengths:

    assert np.allclose(all_edge_idxs, all_edge_idxs_stored)
...
E   ValueError: operands could not be broadcast together with shapes (17284,) (17252,)

The second number is: 2 x (# number of segments in the segment file)

The first number is 32 larger which suggests that somehow there are an additional 16 segments in all_edge_idxs vs. all_edge_idxs_stored?

Having said that I haven't looked at the code yet. I will and happy to discuss.

jploveless commented 2 years ago

keep_segment has the original segments deleted (those with patch_flag = 1 and patch_file_name not equal to 1) (lines 557-559 of celeri.py) and then has the mesh edge segments appended (line 590). The mesh edges therefore replace the original segments.

But it does seem that the segment count is wrong as your previous message indicates, and I could be doing this incorrectly, although in the Japan case the data frame sizes seem correct. I’m confused as to why closure can complete in the Japan case (a full dense run is successful) and why just running closure within a notebook works fine for this global case but this test fails. I’ll look into it more tomorrow.

On Thu, May 19, 2022 at 10:15 PM Brendan Meade @.***> wrote:

@.**** commented on this pull request.

In celeri/celeri.py https://github.com/brendanjmeade/celeri/pull/101#discussion_r877675484:

closest_edge_idx = np.argmin(hang_to_mesh_dist, axis=1)

Replace segment coordinates with closest mesh coordinate

Using a loop because we need to evaluate whether to replace endpoint 1 or 2

for i in range(len(closest_edge_idx)):

if hanging_idx[i] < len(keep_segment.lon1):

keep_segment.loc[hanging_idx[i], "lon1"] = ecoords[0, closest_edge_idx[i]]

keep_segment.loc[hanging_idx[i], "lat1"] = ecoords[1, closest_edge_idx[i]]

else:

keep_segment.loc[hanging_idx[i] - len(keep_segment.lon1), "lon2"] = ecoords[

0, closest_edge_idx[i]

]

keep_segment.loc[hanging_idx[i] - len(keep_segment.lon1), "lat2"] = ecoords[

1, closest_edge_idx[i]

]

Merge with mesh edge segments

new_segment = keep_segment.append(all_edge_segment)

@jploveless https://github.com/jploveless Reading through this now. Should appending happen for the Cascadia case? If so do we need to delete the list elements with the previous coordinates?

— Reply to this email directly, view it on GitHub https://github.com/brendanjmeade/celeri/pull/101#pullrequestreview-979441872, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABZQE4RITH32I2KZRGZDPWDVK3YSZANCNFSM5WMO3UMA . You are receiving this because you were mentioned.Message ID: @.***>

brendanjmeade commented 2 years ago

keep_segment has the original segments deleted (those with patch_flag = 1 and patch_file_name not equal to 1) (lines 557-559 of celeri.py) and then has the mesh edge segments appended (line 590). The mesh edges therefore replace the original segments.

Understood. Is it 16 segments that should be replaced?

But it does seem that the segment count is wrong as your previous message indicates, and I could be doing this incorrectly, although in the Japan case the data frame sizes seem correct. I’m confused as to why closure can complete in the Japan case (a full dense run is successful) and why just running closure within a notebook works fine for this global case but this test fails. I’ll look into it more tomorrow.

That is super strange! Does the closure work when it is all run through celeri_solve? Does it past pytest locally? Definitely super weird...and sounds like most of my days

jploveless commented 2 years ago

Understood. Is it 16 segments that should be replaced?

Yes, it looks like it's net 16 segments that are getting replaced (39 removed, 55 added). I'm going to try to figure out where all_edge_idxs and all_edge_idxs_stored are generated. The new segments are added very early in the code, so I'm not sure how the closure algorithm is "remembering" the original length, but I'm still a bit fuzzy on mutable variables, etc. and so maybe I'm doing something wrong in trying to update the segment structure with these new segments.

That is super strange! Does the closure work when it is all run through celeri_solve? Does it past pytest locally? Definitely super weird...and sounds like most of my days

Shoot, I was still running celeri_dense.ipynb. I'll try celeri_solve. And I know I had said that it ran through pytest locally, but that's actually not true: it's throwing the same error as here.

jploveless commented 2 years ago

I'm going to try to figure out where all_edge_idxs and all_edge_idxs_stored are generated.

Oh — looks like all_edge_idxs_stored is loaded from a file, /tests/global_closure_test_data.npy and so it makes sense that these are not the same length. This explains why the new segment dataframe, reflecting changes from snap_segments, can pass the actual closure algorithm but fails test_basic.py. I'm not sure of the best path forward here: save a new .npy file with the updated segment dataframe?

jploveless commented 2 years ago

@brendanjmeade The current tests/global_closure_test_data.npy reflects the 16 additional segments; tests/global_closure_test_data_old.npy is the previous file. I'm not sure if this is the right path forward (we'd need to change it again once all meshes are integrated in the global model files) but it does confirm where the issue with failing checks was occurring.

brendanjmeade commented 2 years ago

@jploveless Cool! In terms of testing, maybe we can develop a set of input files that are specifically designed for the Cascadia and other cases. For example, a command file like global_command_cascadia_for_pytest.json and an associated .npy file? That's a down the road thing.

brendanjmeade commented 2 years ago

@jploveless This passes tests locally for me as well!

brendanjmeade / celeri

Added snap_segments #101

@.**** commented on this pull request.

Replace segment coordinates with closest mesh coordinate

Using a loop because we need to evaluate whether to replace endpoint 1 or 2

Merge with mesh edge segments