astromfit.h

#ifndef ASTROMFIT__H
#define ASTROMFIT__H

#include <string>
#include <iostream>
#include <map>

#include "ccdimage.h"

#include "eigenstuff.h"
#include "tripletlist.h"
#include "distortionmodel.h"

class Associations;
class DistortionModel;

class AstromFit {
  private : 

  Associations &_assoc; 
  std::string _WhatToFit;
  bool _fittingDistortions, _fittingPos, _fittingRefrac, _fittingPM;
  DistortionModel * _distortionModel;
  int _LastNTrip; // last triplet count, used to speed up allocation
  double _referenceColor; // average color
  unsigned _nRefrac;
  std::vector<double> _refracCoefficient; // fit parameter
  unsigned int _refracPosInMatrix; // where it stands
  double _JDRef; // average Julian date

  // counts in parameter subsets.
  unsigned int _nParDistortions;
  // unsigned int _nParStars; unused presently
  unsigned int _nParTot;
  unsigned _nMeasuredStars;
  
 public :

  AstromFit (Associations &A, DistortionModel *D);
  

  bool Minimize(const std::string &WhatToFit);

  void LSDerivatives(const CcdImage &Ccd, 
             TripletList &TList, Eigen::VectorXd &Rhs) const;

  void LSDerivatives(const CcdImageList  &L, 
             TripletList &TList, Eigen::VectorXd &Rhs);


  void LSDerivatives(TripletList &TList, Eigen::VectorXd &Rhs);

  void AssignIndices(const std::string &WhatToFit);



  void OffsetParams(const Eigen::VectorXd &Delta);

  struct Chi2
  {
    double chi2;
    unsigned ndof;
    Chi2() : chi2(0), ndof(0) {};
    friend std::ostream& operator << (std::ostream& s, const Chi2 &C)
    {
      s << "Chi2/ndof : " << C.chi2 << '/' << C.ndof << '=' <<  C.chi2/C.ndof; return s;
    }

    // Addentry has a third argument in order to make it compatible with an 
    //other stat accumulator.
    void AddEntry(double Inc, unsigned Dof, const MeasuredStar *M) 
    {chi2+= Inc; ndof += Dof;}

    void operator += (const Chi2 &R) 
    {chi2 += R.chi2; ndof += R.ndof;}

  }; // end of struct Chi2

  Chi2 ComputeChi2() const;

  unsigned RemoveOutliers(const double &NSigCut);

  void MakeResTuple(const std::string &TupleName) const;

  std::vector<double> RefractionCoefficients() const 
    { return _refracCoefficient;}

 private : 

  Point TransformFittedStar(const FittedStar &F,
                const Gtransfo * Sky2TP,
                const Point &RefractionVector,
                const double RefractionCoeff,
                const double Jd) const;

  template <class ListType, class Accum> 
    void AccumulateStatImageList(ListType &L, Accum &A) const;

  template <class ImType, class Accum> 
    void AccumulateStatImage(ImType &I, Accum &A) const;
  
  void GetMeasuredStarIndices(const MeasuredStar &Ms, 
                  std::vector<unsigned> &Indices) const;

};


#endif /* ASTROMFIT__H */