Add brigtness calc to comp all orbits#384
Conversation
…b.com/sblunt/orbitize into add-brigtness-calc-to-comp-all-orbits
…b.com/sblunt/orbitize into add-brigtness-calc-to-comp-all-orbits
…b.com/sblunt/orbitize into add-brigtness-calc-to-comp-all-orbits
…b.com/sblunt/orbitize into add-brigtness-calc-to-comp-all-orbits
…b.com/sblunt/orbitize into add-brigtness-calc-to-comp-all-orbits
This reverts commit d1f285c.
semaphoreP
left a comment
There was a problem hiding this comment.
Nice, I like it. I think we have some documentation updates (especially since this is a breaking change to compute_all_orbits(). I also had some questions below about how this will work in the future.
|
|
||
| return mean_anom | ||
|
|
||
| def times2trueanom_and_eccanom( |
|
|
||
| alpha = (1/np.pi)*(np.sin(B)+(np.pi-B)*np.cos(B)) | ||
|
|
||
| albedo = 0.5 # NOTE: we're only fitting relative changes in brightness, so the actual value of albedo doesn't matter |
There was a problem hiding this comment.
if we're fitting relative changes, how will that work? Are we going to analytically marginalize over some linear scale factor term in the future in lnlike?
There was a problem hiding this comment.
Great question. Noting for posterity here our offline conversation that it indeed makes sense to fit a linear scale factor. I just raised issue #416 in response to this comment.
| tanom, eanom = kepler.times2trueanom_and_eccanom(sma, epochs, mtot, ecc, tau, tau_ref_epoch=self.tau_ref_epoch) | ||
|
|
||
|
|
||
| R = (sma*(1-ecc**2))/(1+ecc*np.cos(tanom)) | ||
|
|
||
| z = (R)*(-np.cos(argp)*np.sin(inc)*np.sin(tanom)-np.cos(tanom)*np.sin(inc)*np.sin(argp)) | ||
|
|
||
| B = np.arctan2(-R, z)+ np.pi | ||
|
|
||
| alpha = (1/np.pi)*(np.sin(B)+(np.pi-B)*np.cos(B)) | ||
|
|
||
| albedo = 0.5 # NOTE: we're only fitting relative changes in brightness, so the actual value of albedo doesn't matter | ||
| brightness = albedo*alpha/R**2 |
There was a problem hiding this comment.
maybe not for this PR, but we should make this part of code customizable: we can replace it with a user-provided function or class, so that we can handle custom and non-Lambertian models. We could refactor this into a separate function here or make an issue to do it.
There was a problem hiding this comment.
Agreed! I raised issue #417 to keep track. Let's do it in a future PR.
| vz (np.array of float): N_epochs x N_bodies x N_orbits array of | ||
| radial velocities at each epoch. |
There was a problem hiding this comment.
update docstring to add fourth return element. I think we also (maybe separately) need to update the orbitize release notes about this breaking function signature change for compute_all_orbits(). It may break people's code that use it.
There was a problem hiding this comment.
Good call. I'll update the release notes after I pull into main, and will be sure to note this is a breaking change.
|
|
||
| def test_compute_posteriors(): | ||
| """ | ||
| Test that a short mcmc runs to completion |
There was a problem hiding this comment.
is this test just to show that adding brightness to the csv doesn't break anything? are we actually fitting for it right now?
There was a problem hiding this comment.
Update the docstring to explain the purpose of the test
|
Thanks for the review @semaphoreP! I'm going to pull the changes to main and release a new version now. The CI tests on this branch are broken because they're using an old version of the coveralls action. They should be fixed when pulled into main. |
No description provided.