Add fitting functionality using Ceres-Solver #189
Conversation
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
added 11 commits
October 28, 2024 12:52
mmccrackan
changed the title
mhasself
requested changes
Jan 31, 2025
| int get_dtype(const bp::object &); | ||
|
|
||
| template <typename T> | ||
| T _calculate_median(const T*, const int); |
Comment on lines
+100
to
+103
| // Model independent Negative Log Likelihood for generalized | ||
| // unconstrained minimization | ||
| template <typename Model> | ||
| struct NegLogLikelihood |
There was a problem hiding this comment.
I think it's important to point out, in the comment and/or name of the class, that this is what you should use for fitting power spectra. One should use this only on data whose residuals are expected to be Chi2(1)-distributed. Is it accurate to call it PowerSpectrumCostFunction?
| } | ||
|
|
||
| template <typename T> | ||
| bool _invert_matrix(const T* matrix, T* inverse, const int n) { |
Comment on lines
+451
to
+452
| " p: parameter array (float32/64) with dimensions (ndets, nparams). This is modified in place.\n" | ||
| " c: uncertaintiy array (float32/64) with dimensions (ndets, nparams). This is modified in place.\n" |
There was a problem hiding this comment.
For p and c -- do the values going in matter? Or can it just be zeros/empty?
| "method that minimizes a log likelihood. OMP is used to parallelize across dets (rows)." | ||
| "\n" | ||
| "Args:\n" | ||
| " f: frequency array (float32/64) with dimensions (nsamps,).\n" |
There was a problem hiding this comment.
Are there requirements on f? Must it be positive definite; must it be strictly increasing?
Comment on lines
+234
to
+250
| // index for f > lower fwhite | ||
| for (int i = 0; i < nsamps; ++i) { | ||
| if (f[i] > fwhite_lower) { | ||
| fwhite_i.push_back(i); | ||
| break; | ||
| } | ||
| } | ||
|
|
||
| // index for f < upper fwhite | ||
| for (int i = nsamps - 1; i >= 0; --i) { | ||
| if (f[i] < fwhite_upper) { | ||
| fwhite_i.push_back(i); | ||
| break; | ||
| } | ||
| } | ||
|
|
||
| int fwhite_size = fwhite_i[1] - fwhite_i[0] + 1; |
There was a problem hiding this comment.
This looks segfaulty. Make it robust against fwhite_lower / fwhite_upper pathologically dodging your tests. A good way is like
int fwhite_lo = 0;
while (fwhite_lo < nsamps && f[fwhite_lo] < fwhite_lower)
fwhite_lo++;
fwhite_hi = fwhite_lo;
whilte (...
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
This branch is an attempt at adding fitting functionality to
so3gusing theCeres-Solverlibrary to support both non-linear least squares and general optimization problems. It is used to duplicate the existing noise-fitting function insotodlib. The idea of this function is to:The cost function classes and fitting functions have been written in a model-independent manner, so they should be easy to adapt for other data and model types.
Ceresincludes options for bounds/constraints and supports auto differentiation so derivatives and gradients don't need be calculated manually. The general minimization routine forCeresdoesn't seem to have uncertainty calculations (its least-squares routines do) so I have manually calculated the inverted Hessian matrix to get the parameter uncertainties.I put most of the code into new
fitting_ops.cxxandfitting_ops.hfiles instead of putting it all inarray_ops.cxx. I added anarray_ops.hfile to hold shared functions declarations.I've added a cmake directory with cmake files for
Ceres,Eigen,Gflag, andGlogwith the last 3 being dependencies ofCeres-Solver. We could move these intospt3gat some point, but this also seems to work in my tests.Eigenis a nice optimized vector, matrix, and linear algebra library that might be useful on its own, but it is not required to useEigenwhen working withCeres.Building the Docker image for tests will take a few minutes longer now as I found that Ceres needed to be built from source as the necessary version doesn't appear to be available through
apt-get libceres-dev.