-
Notifications
You must be signed in to change notification settings - Fork 2
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Prism forward models return nan on singular points #30
Conversation
Return nan when computing the diagonal tensor components on singular points (vertices or edges perpendicular to the direction of the tensor component). When computing on the greatest face (eastern,northern or top face) make the respective diagonal component to return the limit approaching from the outside. There's no need to change the implementation, just add a term of 4 pi to the result.
Start creating special cases for singular points on the magnetic forward modelling for prisms. Returns nans when observation point falls in prism vertices and edges. Solves the issue with the limit on east, north and top faces. Lacks of tests and decide what to do with internal points (where the magnetic field integral is not valid anymore).
Add test functions that test if magnetic forward modelling functions for prisms return nans when the observation points fall in a singular point. Add test function that checks the symmetry of the magnetic field on prism faces.
When forward modelling the magnetic field of prisms, return nan if the computation point is inside the prism. Test it with a new function.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Looks good to me! Left a very minor style suggestion.
choclo/prism/_gravity.py
Outdated
if is_point_on_northing_edge( | ||
easting, northing, upward, prism | ||
) or is_point_on_upward_edge(easting, northing, upward, prism): |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
As a style suggestion, it would probably look nicer if both conditions were wrapped in ()
.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Even if I add ()
around them, Black reformats it as it is right now.
I can disable black on that if statement:
if (
is_point_on_northing_edge(easting, northing, upward, prism)
or is_point_on_upward_edge(easting, northing, upward, prism)
): # fmt: skip
return np.nan
Does it look better?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think it does... so I made the changes and pushed them. If you disagree with disabling Black on this, we can revert it.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Looks good to me 👍🏽
Strange that black doesn't respect the parenthesis.
Improve the style of the if statements that check if observation points fall in a singular point for the diagonal tensor components. Disable black on those lines to avoid reformatting them to the original style.
Sorry, I missed that 😢 That looks great 👍 |
Thanks @LL-Geo! Don't worry about the timing, we can always make changes afterwards. |
Gravity and magnetic forward modelling functions for prisms return
np.nan
when the observation points falls in a singular point. These singular points vary for the field that is being computed. For the gravity tensor components, these are all the vertices and edges perpendicular to the direction(s) of the tensor component. For the magnetic fields, these are all the vertices, all the edges and all the interior points. This PR also fixes the values of the fields on the faces: diagonal gravity tensor components and magnetic field components on prism faces return the limit of the field approaching from outside. To achieve that, a 4π constant must be added when computing on the east, north or top faces. Add tests for these new features and include explanations in the docstrings.Relevant issues
Fixes #36