DM-14510 implement line search in jointcal #91
Conversation
parejkoj
force-pushed
the
tickets/DM-14510
branch
4 times, most recently
from
db1beec
to
6acc78f
Compare
include/lsst/jointcal/FitterBase.h
Outdated
| * | ||
| * @param[in] whatToFit See child method assignIndices for valid string values. | ||
| * @param[in] nSigmaCut How many sigma to reject outliers at. Outlier | ||
| * rejection ignored for nSigmaCut=0. | ||
| * @param[in] doRankUpdate Use CholmodSimplicialLDLT2.update() to do a fast rank update after outlier | ||
| * removal; otherwise do a slower full recomputation of the matrix. | ||
| * Only matters if nSigmaCut != 0. | ||
| * @param[in] doLineSearch Use boost's brent_find_minima to perform a line search after the gradient | ||
| * solution is found, and apply the scale factor to the computed offsets. | ||
| * The line search is done for scale factors in [-1, 2], but if scale factor |
python/lsst/jointcal/jointcal.py
Outdated
| @@ -594,8 +600,11 @@ def _fit_photometry(self, associations, dataName=None): | |||
| model = lsst.jointcal.ConstrainedPhotometryModel(associations.getCcdImageList(), | |||
| self.focalPlaneBBox, | |||
| visitOrder=self.config.photometryVisitOrder) | |||
| # potentially nonlinear problem, so we may need a line search to converge. | |||
| doLineSearch = True if self.config.allowLineSearch else False | |||
There was a problem hiding this comment.
Why not just doLineSearch = self.config.allowLineSearch?
There was a problem hiding this comment.
Because I'm an idiot. Thanks.
src/FitterBase.cc
Outdated
| @@ -299,5 +305,18 @@ void FitterBase::saveChi2Contributions(std::string const &baseName) const { | |||
| saveChi2RefContributions(refTuple); | |||
| } | |||
|
|
|||
| float FitterBase::_lineSearch(Eigen::VectorXd const &delta) { | |||
There was a problem hiding this comment.
Are you sure this should be float? Especially since you are only assigning it to scale (which is double) above.
There was a problem hiding this comment.
I have no idea why I made that a float.
src/FitterBase.cc
Outdated
| offsetParams(-scale * delta); | ||
| return chi2.chi2; | ||
| }; | ||
| auto result = boost::math::tools::brent_find_minima(func, -1.0, 2.0, std::numeric_limits<double>::digits); |
There was a problem hiding this comment.
Note that in principle, the minima can not be located to greater accuracy than the square root of machine epsilon (for 64-bit double, sqrt(1e-16)≅1e-8), therefore the value of bits will be ignored if it's greater than half the number of bits in the mantissa of T.
Since ::digits gets you the latter this will be ignored. Did you intend that? If so please add a note.
There was a problem hiding this comment.
Their example does exactly this:
For simplicity, we set the precision parameter bits to std::numeric_limits::digits, which is effectively the maximum possible i.e. std::numeric_limits::digits/2.
Given that, and the lack of clarity on what "ignored" means in this context (I think the implication is "we will use the maximum possible precision"), I changed it to digits/2 and also added a note.
src/FitterBase.cc
Outdated
| offsetParams(scale * delta); | ||
| auto chi2 = computeChi2(); | ||
| // reset the system to where it was before offsetting. | ||
| offsetParams(-scale * delta); |
There was a problem hiding this comment.
I'm not a great fan of this setting/resetting, or the state changing in general for something that could be a const function, but I assume it is unavoidable?
There was a problem hiding this comment.
It is necessary to change the state in order to compute the chi2 at the new position.
Another option would be to make a copy of the model and offset the copy, but the model is very large, so allocate/copy/add is probably of comparable speed to add/subtract (in place) as implemented, but uses more memory and requires implementing copy for the models (which have many shared pointers to the underlying transforms, each of which would have to be deep copied). Also, any numerical imprecision induced by the state change should be well below the overall precision desired by the model.
| @@ -172,6 +172,16 @@ def test_jointcalTask_2_visits_constrainedPhotometry_no_rank_update(self): | |||
|
|
|||
| self._testJointcalTask(2, None, None, pa1, metrics=metrics) | |||
|
|
|||
| def test_jointcalTask_2_visits_constrainedPhotometry_lineSearch(self): | |||
| """Activating the line search should only slightly change the chi2 in this case.""" | |||
There was a problem hiding this comment.
Hmm, can this expected change be formalized? Checking exact values doesn't seem particularly robust.
There was a problem hiding this comment.
I don't know how I would go about formalizing the difference here: the values in these tests are almost all empirical (see the note in the readme). I did add a hopefully clarifying comment.
typedefing std::pair meant that anything that could be cast to `<int,int>` was becoming a CcdImageKey, so e.g. `std::out << somepair` was using our overload!
Uses boost::brent_find_minima() to perform the line search. Add allowLineSearch config option to turn on lineSearch when it's relevant.
parejkoj
force-pushed
the
tickets/DM-14510
branch
from
0a1ec8d
to
a16da5b
Compare
No description provided.