-
Notifications
You must be signed in to change notification settings - Fork 101
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
Piecewise orbit model #1378
Piecewise orbit model #1378
Conversation
Codecov Report
@@ Coverage Diff @@
## master #1378 +/- ##
==========================================
+ Coverage 62.19% 62.38% +0.19%
==========================================
Files 89 91 +2
Lines 20215 20621 +406
Branches 3650 3760 +110
==========================================
+ Hits 12572 12865 +293
- Misses 6857 6924 +67
- Partials 786 832 +46
Help us with your feedback. Take ten seconds to tell us how you rate us. Have a feature suggestion? Share it here. |
# self.validate() | ||
|
||
def update_binary_object(self, toas=None, acc_delay=None): | ||
super().update_binary_object(toas=toas, acc_delay=acc_delay) |
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.
If all this method does is call the same method on super()
can it not be omitted completely? Several other methods here seem to be in the same situation.
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.
Updated so only super() is called when required
self.validate() | ||
|
||
def add_piecewise_param(self, param, param_unit, paramx, j): | ||
if int(j) in self.get_prefix_mapping_component("X_"): |
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.
This should be after the if j is None
, shouldn't 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.
Reordered the statements, good catch
if param is "A1": | ||
self.add_param( | ||
prefixParameter( | ||
name=param + "X_{0}".format(i), |
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.
Nicer to use f"X_{i}"
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.
Changed string formatting as suggested (here and in other places).
frozen=False, | ||
) | ||
) | ||
elif param is "T0": |
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.
For strings it's better to use ==
than is
- you could have two different strings that contain T0
that were equal but not is
.
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.
changed all instances of is to == for strings
|
||
def setup(self): | ||
# print(self.PB) | ||
super(BinaryBTPiecewise, self).setup() |
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.
Just use super().setup()
.
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.
Changed to super().setup()
# if None in T0Xs.values(): | ||
# print("Group with non-defined T0X, applying default T0 to group") | ||
# TODO set default T0 value | ||
if None not in T0Xs.values(): |
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.
Doesn't this skip defining all the parameters if any of the values is unset? That doesn't seem right?
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.
Addressed this by assigning a global A1/T0 in place of the None
it was previously being set to. The intention of this was to later re-assign None
's with the global value elsewhere, but as you spotted this would break if the values are unset
self.binary_instance.add_binary_params(xr2_name, xr2_value) | ||
self.update_binary_object() | ||
|
||
def validate(self): |
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.
Are you sure there aren't any ways the new parameters could be set that would make the model unusable? validate()
should check for invalid settings of the piecewise parameters. (Overlapping T0X intervals? T0X intervals with end before beginning? Any weirdness that won't work can be detected here and then you won't need to worry about it later.)
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.
Added additional checks. Will now raise a value error if:
- Overlapping detected
- mismatching group indexes (i.e. declaring the lower bound index as 0000 and upper bound index as 0001)
- mismatching indexes between the piecewise param and group boundaries (i.e. declaring T0X_0000 and pieces PLB/PUB_0001)
- If the end of a group boundary is chronologically before the start
self.extended_group_range = [] | ||
self.d_binarydelay_d_par_funcs = [self.d_BTdelay_d_par] | ||
if t is not None: | ||
self._t = t |
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.
It's probably a good idea to set selt._t
to None if t
is None, rather than leaving it uninitialized. That way you get a ValueError rather than an AttributeError (the latter usually means you made a typo in the name of an attribute).
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.
Changed as suggested
if input_params is not None: | ||
if self.T0X is None: | ||
self.update_input(input_params) | ||
elif self.T0X is not None: |
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.
This is just else
isn't 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.
Silly me, probably a hangover from when I thought there would be a need for extra conditions. Fixed now (along with any other iterations of this occurring)
if self.T0X is None: | ||
self.update_input(input_params) | ||
elif self.T0X is not None: | ||
self.update_input() |
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.
Isn't this doing exactly the same thing in both branches of the if
statement?
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.
Yes, removed the condition and doesn't break the tests
super().set_param_values(valDict=valDict) | ||
self.piecewise_parameter_loader(valDict=valDict) | ||
|
||
def piecewise_parameter_loader(self, valDict=None): |
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.
It's preferable if a function name is an imperative - load_piecewise_parameters
, perhaps, or setup_internal_structures
or something. This function should have a docstring so you can explain what it does and when to call 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.
Forgot to add the docstring with the update, still to do...
if ( | ||
valDict is None | ||
): # If there are no updates passed by binary_instance, sets default value (usually overwritten when reading from parfile) | ||
self.T0X_arr = self.T0 |
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.
Does this make T0X_arr
into a scalar instead of an array? Won't this break things?
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.
Apparently it didn't, confused as to why since it's been like this through testing. Changed so it lives in a list so it matches the expected form of T0X_arr after self._t is declared
else: | ||
piece_index = ( | ||
[] | ||
) # iniialise array used to count the number of pieces. Operates by seaching for "A1X_i/T0X_i" and appending i to the array. Can operate if pieces are given out of order. |
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.
Try to avoid comments on the same line as your code, it makes for nasty formatting.
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.
moved remaining comments off the same line
index | ||
) in ( | ||
piece_index | ||
): # looping through each index in order they are given (0 -> n) |
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'm not sure this is really the order you're getting. Doesn't np.unique
give you the strings in sorted (lexicographical) order, regardless of the order that they were given?
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.
So it's populating the list piece_index
with all the unique piecewise parameter keys (i.e. '0000', '0001' etc). Then looks up the corresponding values from valDict using that key. It's a matching process rather than using the iteration number of the loop.
if len(self.piecewise_parameter_information) > 0: | ||
# check = hasattr(self,"t") | ||
# print(f"Piecewise parameter loader can see t: {check}") | ||
if hasattr(self, "_t") is True: |
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.
Don't use is True
; just say if hasattr(...):
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.
Also if you set this to None in __init__
you can avoid the need tor hasattr
.
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.
removed instances of is True. There may still be some hasattr
's in other conditions that could follow the same treatment
return toas | ||
|
||
|
||
def get_toa_group_indexes(model, toas): |
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.
This isn't recognized as a test because it doesn't start with test_
. Wouldn't it be better to just use model.which_group_is_toa_in(toas)
directly?
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.
Removed redundant functions like this and similar
|
||
def get_number_of_groups(model): | ||
# returns number of groups | ||
model.setup() |
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.
Likewise it's better to use model.get_number_of_groups()
directly.
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.
Fixed as part of the previous comment
@pytest.fixture() | ||
def build_piecewise_model_with_two_pieces(build_model): | ||
# takes the basic model frame and adds 2 non-ovelerlapping pieces to it | ||
piecewise_model = deepcopy(build_model) |
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.
Fixtures always give you a fresh-created new object, so there's no need to copy 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.
Removed instances of deepcopy()
.
Just for reference, if I call for example:
m1 = build_model
then add pieces to m1
followed by: m2 = build_model
. Then attempting to add pieces results in a failure due to trying to add a piecewise parameter where one already exists. Does this look like a caching issue or is it because build_model is creating a fresh version for the test which is shared between m1 and m2?
return piecewise_model | ||
|
||
|
||
def test_add_piecewise_parameter(build_model): |
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 this test is actually checking whether the model round-trips to a par file correctly? If so, maybe call the test something that says that?
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.
There are round-trip tests elsewhere as well, it would be good to add a BTpiecewise model to those, to benefit from that testing code as well.
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.
Unfortunately I couldn't find where some of these additional tests live when updating. Still needed.
Otherwise, renamed test to something appropriate. (Likely other instances of inappropriate test names)
copy_param_dict = m3.get_params_dict(which="all") | ||
number_of_keys = 0 | ||
comparison = 0 | ||
for key, value in param_dict.items(): |
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.
Can you just ask assert param_dict == copy_param_dict
? If this fails you'll get a list of the mismatched keys/values.
comparison = 0 | ||
for key, value in param_dict.items(): | ||
number_of_keys = number_of_keys + 1 # iterates up to total number of keys | ||
for copy_key, copy_value in param_dict.items(): |
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.
Shouldn't this be copy_param_dict
?
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.
Good spot, yes it should be
|
||
|
||
def test_get_number_of_groups(build_piecewise_model_with_two_pieces): | ||
# test to make sure number of groups matches with number of added piecewise parameters |
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.
Can you put this information in the function name? We don't call test functions by hand so their names don't particularly need to be concise.
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.
Unchanged yet, fell through the cracks while updating
np.unique(index, return_counts=True)[1][1], | ||
] | ||
expected_toas_in_each_group = [10, 10] | ||
is_there_ten_toas_per_group = np.array_equal( |
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.
it's better to use numpy.testing.assert_array_equal
, it will tell you more about the mismatch
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.
Implemented numpy.testing.assert_array_equal
where appropriate. Additionally using numpy.testing.assert_allclose
where appropriate.
m_piecewise = build_piecewise_model_with_two_pieces | ||
toa = make_toas_to_go_with_two_piece_model | ||
m_piecewise.setup() | ||
rs = pint.residuals.Residuals(toa, m_piecewise) |
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.
Is this ever used?
def test_problematic_group_indexes_and_ranges(build_model): | ||
# Test to flag issues with problematic group indexes | ||
m_piecewise = build_model | ||
with pytest.raises(Exception): |
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.
Please specify what kind of exception is to be raised.
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.
Included the types of exceptions that should be raised (basically all value errors)
# Test to flag issues with problematic group indexes | ||
m_piecewise = build_model | ||
with pytest.raises(Exception): | ||
assert m_piecewise.add_group_range( |
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.
No point asserting anything if it's supposed to raise an exception. Just run the code.
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 suggested, just running the code. New function added to binary_piecewise.py
called lock_groups
which runs validate then setup. Avoids the need to run it after every piecewise parameter is added but must be run after adding the a set of/sets of piecewise parameters and date ranges for the boundaries.
assert m_piecewise.add_group_range( | ||
m_piecewise.START.value, m_piecewise.FINISH.value, j=-1 | ||
) | ||
assert m_piecewise.add_group_range( |
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.
This needs to be in a separate with
statement - the first exception exits the with
.
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.
Fixed to only include a single assertion in a with
statement. I realise most of these have been removed during the update
) | ||
|
||
|
||
def add_groups_and_make_toas(build_model, build_piecewise_model_with_two_pieces, param): |
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.
This should be three separate functions, they don't share any appreciable code.
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.
Separated out into two functions since overlapping groups are now fail validation tests.
) | ||
|
||
|
||
def convert_int_into_index(i): |
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.
Isn't there a utility function to do this in the BTpiecewise
code?
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.
It's removed now, all it was doing is making 1->0001 so I just declare them where necessary now rather than using a function
return truth_array_comparison | ||
|
||
|
||
@pytest.mark.parametrize( |
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.
You can pass a list of functions to parametrize
.
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.
Didn't know that, thanks!
"non-complete group coverage", | ||
], | ||
) | ||
def test_does_toa_lie_in_group( |
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.
The three branches here don't share any appreciable code so they should be three separate tests.
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.
Again separated into two tests, since overlapping groups are handled differently now
rs_temp = pint.residuals.Residuals( | ||
toas, m_piecewise_temp, subtract_mean=False | ||
).resids_value | ||
have_residuals_changed = np.invert(rs_temp == rs_value) |
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.
you can use rs_temp != rs_value
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.
Talk about going the long way to get the same thing! Silly oversight on my part
|
||
|
||
def get_d_delay_d_xxxx(toas, model, parameter_string): | ||
rs_temp = pint.residuals.Residuals(toas, model, subtract_mean=False) |
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.
Do you really need to construct Residuals here? Why?
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.
It appears given the updates I don't (/perhaps never needed to), now constructing residuals is reserved for tests that actually use them.
return d_delay_temp | ||
|
||
|
||
# @pytest.mark.parametrize("param, index, offset_size",[("T0X",0, 1e-5*u.d),#("T0X",-1, 1e-5*u.d),("A1X",0, 1e-5*ls),("A1X",-1, 1e-5*ls)]) |
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.
You can use pytest.mark.skip
with a reason to indicate that a test should be skipped.
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.
Thanks for the info, but this test has been removed and replaced now, informed by some of the test suggestions you make in the next comment!
This looks excellent! A little code cleanup would improve it a lot. I'd also like to suggest a few additional tests:
|
Some discussion of the
|
Hi All: What is the status of this code? Is it working? Or might it be soon? I have a new student who will be working on timing Ter5A and would love to use this, if possible. |
Adding some of my extra thoughts based on Anne's last block of suggestions.
Haven't been able to spot where these live. Seems like a good idea to implement still
Implemented as validation tests which are checked in the testing script.
I think I've made a start on testing some of these in the final test of the script. It performs some basic arithmetic and attempts to change piecewise parameter values in two independent models based on the other model. Currently passing but I doubt this is the exact form you were suggesting. Further discussion required. |
Replying to @dlakaplan's comment. Yes there does seem to be a natural way these two can interface. I would need to think about how I would implement this. Either as a |
Replying to @scottransom's comment. The code has received the much needed face-lift in a lot of places, and I'm in the process of trimming down an example notebook. It's not quite at the level to be integrated into PINT's tutorial notebooks just yet, but happy to distribute privately and run through how I'd use it. |
@poneill129 : are you able to attend the weekly PINT call to discuss this? |
Yes I will be, will this be on Thursday? If so, what time does it start? |
A completed version of this (#1545) has already been merged. So I am closing this. |
This is intended as an extension to the existing BT model currently implemented. BT Piecewise was designed to allow a user to fit for orbital phase (T0) with user-defined MJD blocks.
Preliminary work was made to allow for other piecewise orbital parameters, but since this was a tool intended for a specific problem, the development of other piecewise orbital parameters was stopped. Though the foundations for future development in this area are present.
It feels like setting the piecewise parameters from binary_piecewise to BT_piecewise is currently a jury-rig which could do with improvements. Additional concerns surround this involves when a setup() function should be called. A better understanding of this would help performance time.