LSSTApplications
19.0.0-14-gb0260a2+72efe9b372,20.0.0+7927753e06,20.0.0+8829bf0056,20.0.0+995114c5d2,20.0.0+b6f4b2abd1,20.0.0+bddc4f4cbe,20.0.0-1-g253301a+8829bf0056,20.0.0-1-g2b7511a+0d71a2d77f,20.0.0-1-g5b95a8c+7461dd0434,20.0.0-12-g321c96ea+23efe4bbff,20.0.0-16-gfab17e72e+fdf35455f6,20.0.0-2-g0070d88+ba3ffc8f0b,20.0.0-2-g4dae9ad+ee58a624b3,20.0.0-2-g61b8584+5d3db074ba,20.0.0-2-gb780d76+d529cf1a41,20.0.0-2-ged6426c+226a441f5f,20.0.0-2-gf072044+8829bf0056,20.0.0-2-gf1f7952+ee58a624b3,20.0.0-20-geae50cf+e37fec0aee,20.0.0-25-g3dcad98+544a109665,20.0.0-25-g5eafb0f+ee58a624b3,20.0.0-27-g64178ef+f1f297b00a,20.0.0-3-g4cc78c6+e0676b0dc8,20.0.0-3-g8f21e14+4fd2c12c9a,20.0.0-3-gbd60e8c+187b78b4b8,20.0.0-3-gbecbe05+48431fa087,20.0.0-38-ge4adf513+a12e1f8e37,20.0.0-4-g97dc21a+544a109665,20.0.0-4-gb4befbc+087873070b,20.0.0-4-gf910f65+5d3db074ba,20.0.0-5-gdfe0fee+199202a608,20.0.0-5-gfbfe500+d529cf1a41,20.0.0-6-g64f541c+d529cf1a41,20.0.0-6-g9a5b7a1+a1cd37312e,20.0.0-68-ga3f3dda+5fca18c6a4,20.0.0-9-g4aef684+e18322736b,w.2020.45
LSSTDataManagementBasePackage
|
Class that handles the astrometric least squares problem.
More...
#include <AstrometryFit.h>
Class that handles the astrometric least squares problem.
This is the class that actually computes the quantities required to carry out a LS astrometric fit wrt distortion mappings and coordinates of common objects. Namely it computes the Jacobian and gradient of the chi2 (w.r.t. parameters), and the Chi2 itself. It interfaces with the actual modelling of distortions via a mimimum virtual interface AstrometryModel, and the actual mappings via an other virtual interface : Mapping.
In short AstrometryFit aims at computing derivatives of least quares. The terms of the chi2 are of two kinds:
kind 1 -> (T(X_M) - p(F))^T W (T(X_M) - p(F))
with X_M is a measured (2d) position in CCD coordinates, F refers to the position of the object in some space, defined in practise by p. There is one such term per measurement. The default setup would be that p is the projection from sky to some tangent plane and hence T maps the CCD coordinates onto this TP. p is obtained via the DistorsionModel and can be different for all CcdImage's. Depending on what is beeing fitted, one could imagine cases where the projector p is the same for all CcdImages.
Kind 2 -> (p'(F)-p'(R))^T W_R (p'(F)-p'(R)) R refers to some externally-provided reference object position, and p' to some projector from sky to some plane. The reference objects define the overall coordinate frame, which is required when all T and all F are fitted simultaneously. There is one such term per external reference object. There can be more F (fitted) objects than R (reference) objects.
In the same framework, one can fit relative transforms between images by setting p = Identity for all input CcdImages and not fitting T for one of the CcdImage's. One does not need reference object and would then naturally not have any Kind 2 terms.
Definition at line 78 of file AstrometryFit.h.
◆ AstrometryFit() [1/3]
this is the only constructor
Definition at line 46 of file AstrometryFit.cc.
49 _astrometryModel(astrometryModel),
50 _refractionCoefficient(0),
64 _referenceColor += i->color;
69 _referenceColor /= double(count);
70 if (_sigCol > 0) _sigCol =
sqrt(_sigCol / count -
std::pow(_referenceColor, 2));
72 LOGLS_INFO(
_log,
"Reference Color: " << _referenceColor <<
" sig " << _sigCol);
◆ AstrometryFit() [2/3]
lsst::jointcal::AstrometryFit::AstrometryFit |
( |
AstrometryFit const & |
| ) |
|
|
delete |
No copy or move: there is only ever one fitter of a given type.
◆ AstrometryFit() [3/3]
lsst::jointcal::AstrometryFit::AstrometryFit |
( |
AstrometryFit && |
| ) |
|
|
delete |
◆ assignIndices()
void lsst::jointcal::AstrometryFit::assignIndices |
( |
std::string const & |
whatToFit | ) |
|
|
overridevirtual |
Set parameters to fit and assign indices in the big matrix.
- Parameters
-
[in] | whatToFit | Valid strings: zero or more of "Distortions", "Positions", "Refrac", "PM" which define which parameter set is going to be variable when computing derivatives (leastSquareDerivatives) and minimizing (minimize()). whatToFit="Positions Distortions" will minimize w.r.t mappings and objects positions, and not w.r.t proper motions and refraction modeling. However if proper motions and/or refraction parameters have already been set, then they are accounted for when computing residuals. The string is forwarded to the AstrometryModel, and it can then be used to turn subsets of distortion parameter on or off, if the AstrometryModel implements such a thing. |
Implements lsst::jointcal::FitterBase.
Definition at line 451 of file AstrometryFit.cc.
454 _fittingDistortions = (
_whatToFit.
find(
"Distortions") != std::string::npos);
455 _fittingPos = (
_whatToFit.
find(
"Positions") != std::string::npos);
456 _fittingRefrac = (
_whatToFit.
find(
"Refrac") != std::string::npos);
457 if (_sigCol == 0 && _fittingRefrac) {
459 "Cannot fit refraction coefficients without a color lever arm. Ignoring refraction.");
460 _fittingRefrac =
false;
465 _nParDistortions = 0;
466 if (_fittingDistortions) _nParDistortions = _astrometryModel->assignIndices(
_whatToFit, 0);
470 FittedStarList &fittedStarList =
_associations->fittedStarList;
471 for (
auto &fittedStar : fittedStarList) {
476 fittedStar->setIndexInMatrix(ipar);
478 if ((_fittingPM)&fittedStar->mightMove) ipar +=
NPAR_PM;
481 _nParPositions = ipar - _nParDistortions;
482 if (_fittingRefrac) {
483 _refracPosInMatrix = ipar;
◆ checkStuff()
void lsst::jointcal::AstrometryFit::checkStuff |
( |
| ) |
|
DEBUGGING routine.
Definition at line 517 of file AstrometryFit.cc.
519 const char *what2fit[] = {
"Positions",
522 "Positions Distortions",
524 "Distortions Refrac",
525 "Positions Distortions Refrac"};
527 const char *what2fit[] = {
"Positions",
"Distortions",
"Positions Distortions"};
529 for (
unsigned k = 0; k <
sizeof(what2fit) /
sizeof(what2fit[0]); ++k) {
531 TripletList tripletList(10000);
536 jacobian.setFromTriplets(tripletList.begin(), tripletList.end());
◆ computeChi2()
Returns the chi2 for the current state.
Definition at line 42 of file FitterBase.cc.
◆ findOutliers()
Find Measurements and references contributing more than a cut, computed as.
\[ <chi2> + nSigmaCut + rms(chi2). \]
The outliers are NOT removed, and no refit is done.
After returning from here, there are still measurements that contribute above the cut, but their contribution should be evaluated after a refit before discarding them.
- Parameters
-
[in] | nSigmaCut | Number of sigma to select on. |
[out] | msOutliers | list of MeasuredStar outliers to populate |
[out] | fsOutliers | list of FittedStar outliers to populate |
- Returns
- Total number of outliers that were removed.
Definition at line 52 of file FitterBase.cc.
63 size_t nval = chi2List.size();
64 if (nval == 0)
return 0;
65 sort(chi2List.begin(), chi2List.end());
66 double median = (nval & 1) ? chi2List[nval / 2].chi2
67 : 0.5 * (chi2List[nval / 2 - 1].chi2 + chi2List[nval / 2].chi2);
68 auto averageAndSigma = chi2List.computeAverageAndSigma();
69 LOGLS_DEBUG(
_log,
"findOutliers chi2 stat: mean/median/sigma " << averageAndSigma.first <<
'/' << median
70 <<
'/' << averageAndSigma.second);
71 double cut = averageAndSigma.first + nSigmaCut * averageAndSigma.second;
77 Eigen::VectorXi affectedParams(
_nParTot);
78 affectedParams.setZero();
82 for (
auto chi2 = chi2List.rbegin(); chi2 != chi2List.rend(); ++chi2) {
83 if (chi2->chi2 < cut)
break;
88 auto measuredStar = std::dynamic_pointer_cast<MeasuredStar>(chi2->star);
90 if (measuredStar ==
nullptr) {
92 fittedStar = std::dynamic_pointer_cast<FittedStar>(chi2->star);
93 if (fittedStar->getMeasurementCount() == 0) {
94 LOGLS_WARN(
_log,
"FittedStar with no measuredStars found as an outlier: " << *fittedStar);
99 indices.
push_back(fittedStar->getIndexInMatrix());
100 LOGLS_TRACE(
_log,
"Removing refStar " << *(fittedStar->getRefStar()) <<
" chi2: " << chi2->chi2);
105 auto tempFittedStar = measuredStar->getFittedStar();
106 if (tempFittedStar->getMeasurementCount() == 1 && tempFittedStar->getRefStar() ==
nullptr) {
107 LOGLS_WARN(
_log,
"FittedStar with 1 measuredStar and no refStar found as an outlier: "
112 LOGLS_TRACE(
_log,
"Removing measStar " << *measuredStar <<
" chi2: " << chi2->chi2);
120 for (
auto const &i : indices) {
121 if (affectedParams(i) != 0) {
128 if (measuredStar ==
nullptr) {
130 fsOutliers.push_back(fittedStar);
133 msOutliers.push_back(measuredStar);
136 for (
auto const &i : indices) {
142 LOGLS_INFO(
_log,
"findOutliers: found " << msOutliers.size() <<
" meas outliers and " << fsOutliers.size()
143 <<
" ref outliers ");
◆ freezeErrorTransform()
void lsst::jointcal::AstrometryFit::freezeErrorTransform |
( |
| ) |
|
|
inline |
The transformations used to propagate errors are freezed to the current state.
The routine can be called when the mappings are roughly in place. After the call, the transformations used to propage errors are no longer affected when updating the mappings. This allows to have an exactly linear fit, which can be useful.
Definition at line 117 of file AstrometryFit.h.
117 { _astrometryModel->freezeErrorTransform(); }
◆ getModel()
Return the model being fit.
Definition at line 122 of file AstrometryFit.h.
122 {
return _astrometryModel; }
◆ leastSquareDerivatives()
void lsst::jointcal::FitterBase::leastSquareDerivatives |
( |
TripletList & |
tripletList, |
|
|
Eigen::VectorXd & |
grad |
|
) |
| const |
|
inherited |
Evaluates the chI^2 derivatives (Jacobian and gradient) for the current whatToFit setting.
The Jacobian is given as triplets in a sparse matrix, the gradient as a dense vector. The parameters which vary, and their indices, are to be set using assignIndices.
- Parameters
-
tripletList | tripletList of (row,col,value) representing the Jacobian of the chi2. |
grad | The gradient of the chi2. |
Definition at line 317 of file FitterBase.cc.
319 for (
auto const &ccdImage : ccdImageList) {
◆ minimize()
MinimizeResult lsst::jointcal::FitterBase::minimize |
( |
std::string const & |
whatToFit, |
|
|
double const |
nSigmaCut = 0 , |
|
|
bool const |
doRankUpdate = true , |
|
|
bool const |
doLineSearch = false , |
|
|
std::string const & |
dumpMatrixFile = "" |
|
) |
| |
|
inherited |
Does a 1 step minimization, assuming a linear model.
This is a complete Newton Raphson step. Compute first and second derivatives, solve for the step and apply it, with an optional line search.
It calls assignIndices, leastSquareDerivatives, solves the linear system and calls offsetParams, then removes outliers in a loop if requested. Relies on sparse linear algebra via Eigen's CholmodSupport package.
- Parameters
-
[in] | whatToFit | See child method assignIndices for valid string values. |
[in] | nSigmaCut | How many sigma to reject outliers at. Outlier rejection ignored for nSigmaCut=0. |
[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. |
[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 in the domain [-1, 2], but if the scale factor is far from 1.0, then the problem is likely in a significantly non-linear regime. |
[in] | dumpMatrixFile | Write the pre-fit Hessian matrix and gradient to the files with "-mat.txt" and "-grad.txt". Be aware, this requires a large increase in memory usage to create a dense matrix before writing it; the output file may be large. Writing the matrix can be helpful for debugging bad fits. Read it and compute the real eigenvalues (recall that the Hessian is symmetric by construction) with numpy: hessian = np.matrix(np.loadtxt("dumpMatrixFile-mat.txt"))
values, vectors = np.linalg.eigh(hessian)
|
- Returns
- Return code describing success/failure of fit.
- Note
- When fitting one parameter set by itself (e.g. "Model"), the system is purely linear (assuming there are no cross-terms in the derivatives, e.g. the SimpleAstrometryModel), which should result in the optimal chi2 after a single step. This can be used to debug the fitter by fitting that parameter set twice in a row: the second run with the same "whatToFit" will produce no change in the fitted parameters, if the calculations and indices are defined correctly.
Definition at line 171 of file FitterBase.cc.
179 TripletList tripletList(nTrip);
188 LOGLS_DEBUG(
_log,
"End of triplet filling, ntrip = " << tripletList.size());
194 << hessian.rows() <<
" non-zeros=" << hessian.nonZeros()
195 <<
" filling-frac = " << hessian.nonZeros() /
std::pow(hessian.rows(), 2));
197 if (dumpMatrixFile !=
"") {
198 if (hessian.rows() * hessian.cols() > 2e8) {
199 LOGLS_WARN(
_log,
"Hessian matrix is too big to dump to file, with rows, columns: "
200 << hessian.rows() <<
", " << hessian.cols());
202 dumpMatrixAndGradient(hessian, grad, dumpMatrixFile,
_log);
207 if (chol.info() != Eigen::Success) {
217 Eigen::VectorXd delta = chol.solve(grad);
219 scale = _lineSearch(delta);
225 LOGL_ERROR(
_log,
"chi2 is not finite. Aborting outlier rejection.");
229 if (currentChi2.chi2 > oldChi2 && totalMeasOutliers + totalRefOutliers != 0) {
230 LOGL_WARN(
_log,
"chi2 went up, skipping outlier rejection loop");
234 oldChi2 = currentChi2.chi2;
236 if (nSigmaCut == 0)
break;
237 MeasuredStarList msOutliers;
238 FittedStarList fsOutliers;
241 totalMeasOutliers += msOutliers.size();
242 totalRefOutliers += fsOutliers.size();
243 if (nOutliers == 0)
break;
244 TripletList outlierTriplets(nOutliers);
254 H.setFromTriplets(outlierTriplets.begin(), outlierTriplets.end());
255 chol.update(H,
false );
266 LOGLS_DEBUG(
_log,
"Triplets recomputed, ntrip = " << nextTripletList.size());
268 hessian = createHessian(
_nParTot, nextTripletList);
269 nextTripletList.clear();
272 "Restarting factorization, hessian: dim="
273 << hessian.rows() <<
" non-zeros=" << hessian.nonZeros()
274 <<
" filling-frac = " << hessian.nonZeros() /
std::pow(hessian.rows(), 2));
275 chol.compute(hessian);
276 if (chol.info() != Eigen::Success) {
284 if (nSigmaCut != 0) {
285 LOGLS_INFO(
_log,
"Number of outliers (Measured + Reference = Total): "
286 << totalMeasOutliers <<
" + " << totalRefOutliers <<
" = "
287 << totalMeasOutliers + totalRefOutliers);
◆ offsetParams()
void lsst::jointcal::AstrometryFit::offsetParams |
( |
Eigen::VectorXd const & |
delta | ) |
|
|
overridevirtual |
Offset the parameters by the requested quantities.
The used parameter layout is the one from the last call to assignIndices or minimize(). There is no easy way to check that the current setting of whatToFit and the provided Delta vector are compatible: we can only test the size.
- Parameters
-
[in] | delta | vector of offsets to apply |
Implements lsst::jointcal::FitterBase.
Definition at line 489 of file AstrometryFit.cc.
491 throw LSST_EXCEPT(pex::exceptions::InvalidParameterError,
492 "AstrometryFit::offsetParams : the provided vector length is not compatible with "
493 "the current whatToFit setting");
494 if (_fittingDistortions) _astrometryModel->offsetParams(delta);
497 FittedStarList &fittedStarList =
_associations->fittedStarList;
498 for (
auto const &i : fittedStarList) {
503 Eigen::Index index = fs.getIndexInMatrix();
504 fs.x += delta(index);
505 fs.y += delta(index + 1);
506 if ((_fittingPM)&fs.mightMove) {
507 fs.pmx += delta(index + 2);
508 fs.pmy += delta(index + 3);
512 if (_fittingRefrac) {
513 _refractionCoefficient += delta(_refracPosInMatrix);
◆ operator=() [1/2]
◆ operator=() [2/2]
◆ outliersContributions()
Contributions to derivatives from (presumably) outlier terms.
No discarding done.
Definition at line 292 of file FitterBase.cc.
294 for (
auto &outlier : msOutliers) {
295 MeasuredStarList tmp;
296 tmp.push_back(outlier);
297 const CcdImage &ccdImage = outlier->getCcdImage();
◆ removeMeasOutliers()
void lsst::jointcal::FitterBase::removeMeasOutliers |
( |
MeasuredStarList & |
outliers | ) |
|
|
protectedinherited |
Remove measuredStar outliers from the fit. No Refit done.
Definition at line 303 of file FitterBase.cc.
304 for (
auto &measuredStar : outliers) {
305 auto fittedStar = measuredStar->getFittedStar();
306 measuredStar->setValid(
false);
307 fittedStar->getMeasurementCount()--;
◆ removeRefOutliers()
void lsst::jointcal::FitterBase::removeRefOutliers |
( |
FittedStarList & |
outliers | ) |
|
|
protectedinherited |
Remove refStar outliers from the fit. No Refit done.
Definition at line 311 of file FitterBase.cc.
312 for (
auto &fittedStar : outliers) {
313 fittedStar->setRefStar(
nullptr);
◆ saveChi2Contributions()
void lsst::jointcal::FitterBase::saveChi2Contributions |
( |
std::string const & |
baseName | ) |
const |
|
virtualinherited |
Save the full chi2 term per star that was used in the minimization, for debugging.
Saves results to text files "baseName-meas.csv" and "baseName-ref.csv" for the MeasuredStar and RefStar contributions, respectively. This method is mostly useful for debugging: we will probably want to create a better persistence system for jointcal's internal representations in the future (see DM-12446).
Definition at line 325 of file FitterBase.cc.
327 auto pos = baseName.
find(replaceStr);
329 measFilename.replace(pos, replaceStr.
size(),
"-meas.csv");
331 refFilename.replace(pos, replaceStr.
size(),
"-ref.csv");
◆ saveChi2MeasContributions()
void lsst::jointcal::AstrometryFit::saveChi2MeasContributions |
( |
std::string const & |
filename | ) |
const |
|
overrideprotectedvirtual |
Save a CSV file containing residuals of measurement terms.
Implements lsst::jointcal::FitterBase.
Definition at line 542 of file AstrometryFit.cc.
546 ofile <<
"#id" << separator <<
"xccd" << separator <<
"yccd " << separator;
547 ofile <<
"rx" << separator <<
"ry" << separator;
548 ofile <<
"xtp" << separator <<
"ytp" << separator;
549 ofile <<
"mag" << separator <<
"mjd" << separator;
550 ofile <<
"xErr" << separator <<
"yErr" << separator <<
"xyCov" << separator;
551 ofile <<
"xtpi" << separator <<
"ytpi" << separator;
552 ofile <<
"rxi" << separator <<
"ryi" << separator;
553 ofile <<
"color" << separator <<
"fsindex" << separator;
554 ofile <<
"ra" << separator <<
"dec" << separator;
555 ofile <<
"chi2" << separator <<
"nm" << separator;
556 ofile <<
"chip" << separator <<
"visit" <<
std::endl;
558 ofile <<
"#id in source catalog" << separator <<
"coordinates in CCD (pixels)" << separator << separator;
559 ofile <<
"residual on TP (degrees)" << separator << separator;
560 ofile <<
"transformed coordinate in TP (degrees)" << separator << separator;
561 ofile <<
"rough magnitude" << separator <<
"Modified Julian Date of the measurement" << separator;
562 ofile <<
"transformed measurement uncertainty (degrees)" << separator << separator << separator;
563 ofile <<
"as-read position on TP (degrees)" << separator << separator;
564 ofile <<
"as-read residual on TP (degrees)" << separator << separator;
565 ofile <<
"currently unused" << separator <<
"unique index of the fittedStar" << separator;
566 ofile <<
"on sky position of fittedStar" << separator << separator;
567 ofile <<
"contribution to Chi2 (2D dofs)" << separator <<
"number of measurements of this fittedStar"
569 ofile <<
"chip id" << separator <<
"visit id" <<
std::endl;
572 for (
auto const &ccdImage : ccdImageList) {
573 const MeasuredStarList &cat = ccdImage->getCatalogForFit();
574 const AstrometryMapping *mapping = _astrometryModel->getMapping(*ccdImage);
575 const auto readTransform = ccdImage->getReadWcs();
576 const Point &refractionVector = ccdImage->getRefractionVector();
577 double mjd = ccdImage->getMjd() - _JDRef;
578 for (
auto const &ms : cat) {
579 if (!ms->isValid())
continue;
581 FatPoint inPos = *ms;
582 tweakAstromMeasurementErrors(inPos, *ms, _posError);
583 mapping->transformPosAndErrors(inPos, tpPos);
584 auto sky2TP = _astrometryModel->getSkyToTangentPlane(*ccdImage);
586 compose(*sky2TP, *readTransform);
587 FatPoint inputTpPos = readPixToTangentPlane->apply(inPos);
590 Point fittedStarInTP =
591 transformFittedStar(*fs, *sky2TP, refractionVector, _refractionCoefficient, mjd);
592 Point res = tpPos - fittedStarInTP;
593 Point inputRes = inputTpPos - fittedStarInTP;
594 double det = tpPos.vx * tpPos.vy -
std::pow(tpPos.vxy, 2);
595 double wxx = tpPos.vy /
det;
596 double wyy = tpPos.vx /
det;
597 double wxy = -tpPos.vxy /
det;
598 double chi2 = wxx * res.x * res.x + wyy * res.y * res.y + 2 * wxy * res.x * res.y;
601 ofile << ms->getId() << separator << ms->x << separator << ms->y << separator;
602 ofile << res.x << separator << res.y << separator;
603 ofile << tpPos.x << separator << tpPos.y << separator;
604 ofile << fs->getMag() << separator << mjd << separator;
605 ofile << tpPos.vx << separator << tpPos.vy << separator << tpPos.vxy << separator;
606 ofile << inputTpPos.x << separator << inputTpPos.y << separator;
607 ofile << inputRes.x << separator << inputRes.y << separator;
608 ofile << fs->color << separator << fs->getIndexInMatrix() << separator;
609 ofile << fs->x << separator << fs->y << separator;
610 ofile << chi2 << separator << fs->getMeasurementCount() << separator;
611 ofile << ccdImage->getCcdId() << separator << ccdImage->getVisit() <<
std::endl;
◆ saveChi2RefContributions()
void lsst::jointcal::AstrometryFit::saveChi2RefContributions |
( |
std::string const & |
filename | ) |
const |
|
overrideprotectedvirtual |
Save a CSV file containing residuals of reference terms.
Implements lsst::jointcal::FitterBase.
Definition at line 616 of file AstrometryFit.cc.
620 ofile <<
"#ra" << separator <<
"dec " << separator;
621 ofile <<
"rx" << separator <<
"ry" << separator;
622 ofile <<
"mag" << separator;
623 ofile <<
"xErr" << separator <<
"yErr" << separator <<
"xyCov" << separator;
624 ofile <<
"color" << separator <<
"fsindex" << separator;
625 ofile <<
"chi2" << separator <<
"nm" <<
std::endl;
627 ofile <<
"#coordinates of fittedStar (degrees)" << separator << separator;
628 ofile <<
"residual on TP (degrees)" << separator << separator;
629 ofile <<
"magnitude" << separator;
630 ofile <<
"refStar transformed measurement uncertainty (degrees)" << separator << separator << separator;
631 ofile <<
"currently unused" << separator <<
"unique index of the fittedStar" << separator;
632 ofile <<
"refStar contribution to Chi2 (2D dofs)" << separator
633 <<
"number of measurements of this FittedStar" <<
std::endl;
636 const FittedStarList &fittedStarList =
_associations->fittedStarList;
637 TanRaDecToPixel proj(AstrometryTransformLinear(), Point(0., 0.));
638 for (
auto const &i : fittedStarList) {
639 const FittedStar &fs = *i;
640 const RefStar *rs = fs.getRefStar();
641 if (rs ==
nullptr)
continue;
642 proj.setTangentPoint(fs);
645 proj.transformPosAndErrors(*rs, rsProj);
646 double rx = rsProj.x;
647 double ry = rsProj.y;
648 double det = rsProj.vx * rsProj.vy -
std::pow(rsProj.vxy, 2);
649 double wxx = rsProj.vy /
det;
650 double wyy = rsProj.vx /
det;
651 double wxy = -rsProj.vxy /
det;
652 double chi2 = wxx *
std::pow(rx, 2) + 2 * wxy * rx * ry + wyy *
std::pow(ry, 2);
655 ofile << fs.x << separator << fs.y << separator;
656 ofile << rx << separator << ry << separator;
657 ofile << fs.getMag() << separator;
658 ofile << rsProj.vx << separator << rsProj.vy << separator << rsProj.vxy << separator;
659 ofile << fs.color << separator << fs.getIndexInMatrix() << separator;
660 ofile << chi2 << separator << fs.getMeasurementCount() <<
std::endl;
◆ _associations
◆ _lastNTrip
Eigen::Index lsst::jointcal::FitterBase::_lastNTrip |
|
protectedinherited |
◆ _log
◆ _nMeasuredStars
Eigen::Index lsst::jointcal::FitterBase::_nMeasuredStars |
|
protectedinherited |
◆ _nParTot
Eigen::Index lsst::jointcal::FitterBase::_nParTot |
|
protectedinherited |
◆ _whatToFit
The documentation for this class was generated from the following files:
- /j/snowflake/release/lsstsw/stack/cb4e2dc/Linux64/jointcal/19.0.0-14-gb0260a2+72efe9b372/include/lsst/jointcal/AstrometryFit.h
- /j/snowflake/release/lsstsw/stack/cb4e2dc/Linux64/jointcal/19.0.0-14-gb0260a2+72efe9b372/src/AstrometryFit.cc
T setprecision(T... args)
virtual void leastSquareDerivativesReference(FittedStarList const &fittedStarList, TripletList &tripletList, Eigen::VectorXd &grad) const =0
Compute the derivatives of the reference terms.
virtual void leastSquareDerivativesMeasurement(CcdImage const &ccdImage, TripletList &tripletList, Eigen::VectorXd &grad, MeasuredStarList const *measuredStarList=nullptr) const =0
Compute the derivatives of the measured stars and model for one CcdImage.
std::list< std::shared_ptr< CcdImage > > CcdImageList
#define LOGL_WARN(logger, message...)
virtual void accumulateStatImageList(CcdImageList const &ccdImageList, Chi2Accumulator &accum) const =0
Compute the chi2 (per star or total, depending on which Chi2Accumulator is used) for measurements.
void assignIndices(std::string const &whatToFit) override
Set parameters to fit and assign indices in the big matrix.
#define LOGLS_INFO(logger, message)
#define LOGLS_ERROR(logger, message)
#define LOGLS_TRACE(logger, message)
virtual void accumulateStatRefStars(Chi2Accumulator &accum) const =0
Compute the chi2 (per star or total, depending on which Chi2Accumulator is used) for RefStars.
#define LOGLS_WARN(logger, message)
#define LOGL_ERROR(logger, message...)
std::size_t findOutliers(double nSigmaCut, MeasuredStarList &msOutliers, FittedStarList &fsOutliers) const
Find Measurements and references contributing more than a cut, computed as.
std::shared_ptr< Associations > _associations
void leastSquareDerivatives(TripletList &tripletList, Eigen::VectorXd &grad) const
Evaluates the chI^2 derivatives (Jacobian and gradient) for the current whatToFit setting.
FitterBase(std::shared_ptr< Associations > associations)
virtual void getIndicesOfMeasuredStar(MeasuredStar const &measuredStar, IndexVector &indices) const =0
Set the indices of a measured star from the full matrix, for outlier removal.
MinimizeResult
Return value of minimize()
#define LSST_EXCEPT(type,...)
Create an exception with a given type.
virtual void saveChi2MeasContributions(std::string const &filename) const =0
Save a CSV file containing residuals of measurement terms.
#define LOGLS_DEBUG(logger, message)
void outliersContributions(MeasuredStarList &msOutliers, FittedStarList &fsOutliers, TripletList &tripletList, Eigen::VectorXd &grad)
Contributions to derivatives from (presumably) outlier terms.
Chi2Statistic computeChi2() const
Returns the chi2 for the current state.
virtual void assignIndices(std::string const &whatToFit)=0
Set parameters to fit and assign indices in the big matrix.
void removeMeasOutliers(MeasuredStarList &outliers)
Remove measuredStar outliers from the fit. No Refit done.
Eigen::SparseMatrix< double, 0, Eigen::Index > SparseMatrixD
virtual void offsetParams(Eigen::VectorXd const &delta)=0
Offset the parameters by the requested quantities.
def scale(algorithm, min, max=None, frame=None)
std::unique_ptr< AstrometryTransform > compose(AstrometryTransform const &left, AstrometryTransform const &right)
Returns a pointer to a composition of transforms, representing left(right()).
void removeRefOutliers(FittedStarList &outliers)
Remove refStar outliers from the fit. No Refit done.
Eigen::Index _nMeasuredStars
virtual void saveChi2RefContributions(std::string const &filename) const =0
Save a CSV file containing residuals of reference terms.