-
Notifications
You must be signed in to change notification settings - Fork 4
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
Initial Orbit Determination #31
Conversation
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.
Many thanks for this contribution @LuEdRaMo! This is an awesome contribution to the preliminary OD capabilities in NEOs.jl! Not sure if this is ready for review, but left a couple of comments around; hope this helps!
function NEOParameters(params::NEOParameters{T}; kwargs...) where {T <: AbstractFloat} | ||
fields = fieldnames(NEOParameters{T}) | ||
vals = Vector{Any}(undef, length(fields)) | ||
for i in eachindex(vals) | ||
if fields[i] in keys(kwargs) | ||
vals[i] = kwargs[fields[i]] | ||
else | ||
vals[i] = getfield(params, i) | ||
end | ||
end | ||
|
||
return NEOParameters{T}(vals...) | ||
end |
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.
Parameters.jl and Base.@kwdef
implement similar functionality, have you considered using either of those to implement that functionality here?
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 implemented this manually to have more flexibility with NEOParameters
, specially wrt the debiasing tables. The user specifies the table as a string (e.g. debias_table = "2018"
) but the struct saves the corresponding variables (a string, a Healpix
varible and a matrix), and I do not know if that is possible with the options you mention.
# Scalig factors | ||
@test all(sol1.scalings .< 1e-6) | ||
# Compatibility with JPL | ||
@test all(abs.(sol1() - JPL) ./ sigmas(sol1) .< 1.6) |
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.
IIUC, in this test you update sol
with more observations, but adding the new observations worsens the comparison to the JPL reference. I would expect the comparison to improve, any idea why this is happening?
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.
JPL gives no information about which observations they discard as outliers or which weights they use; therefore, a strict comparison with JPL is not possible. Furthermore, the first orbit is well constrained because it is computed over the best quality astrometry, the observations added for the second orbit do not span a significant amount of extra time and may only introduce noise and/or outliers. Yet, the comparison stays within a 3 sigma threshold wrt JPL. This particular NEO requires a detailed filtering of the observations.
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.
Okay, if comparison wrt a JPL solution is not straightforward, then probably it makes sense to test for improvement of the post-fit statistics for each solution, e.g.:
@test nrms(sol1) < nrms(sol)
@test all.(sigmas(sol1) .< sigmas(sol))
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 NRMS increases from 0.29 to 0.37 from sol
to sol1
but sigmas(sol1)
are smaller than sigmas(sol)
by a factor of 4. I've added the test
@test all(sigmas(sol1) .< sigmas(sol))
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.
Alright, then I'd only suggest to add the NRMS test commented for now and leave as a TODO to try to understand a bit better the differences wrt to the JPL solution.
@LuEdRaMo thank you for adding the suggestions we discussed offline; I left a couple of comments but I aside from those I think this PR is pretty much ready to be merged; any further improvements we can always address in another PR. This has been really an awesome contribution, truly many thanks! |
Just one minor detail: since the DataFrames extension is not a weakdep anymore, I'd suggest to rename |
Done. |
Just pushed a couple of minor fixes, merging! |
Excellent! Thanks a lot @LuEdRaMo |
This PR introduces a new function
gaussinitcond
which transforms a vector of optical observations into a vector of initial conditions viagauss_method
.