Newer
Older
TB_Chris / TbAlignment / src / .svn / text-base / TbAlignmentMinuitBase.cpp.svn-base
  1. #include "TbAlignmentMinuitBase.h"
  2.  
  3. namespace Tb {
  4. TbAlignmentMinuitBase* gTbAlignmentMinuitBase(nullptr);
  5.  
  6. //=============================================================================
  7. // Function to be minimised
  8. //=============================================================================
  9. void fcn(int&, double* g, double& f, double* par, int) {
  10. gTbAlignmentMinuitBase->chi2(f, par, g);
  11. }
  12. }
  13.  
  14. //=============================================================================
  15. // Standard constructor, initializes variables
  16. //=============================================================================
  17. TbAlignmentMinuitBase::TbAlignmentMinuitBase(const std::string& type,
  18. const std::string& name,
  19. const IInterface* parent)
  20. : TbAlignmentBase(type, name, parent) {
  21.  
  22. declareProperty("NIterations", m_nIterations = 3);
  23. declareProperty("FitStrategy", m_fitStrategy = 2);
  24. }
  25.  
  26. //=============================================================================
  27. // Destructor
  28. //=============================================================================
  29. TbAlignmentMinuitBase::~TbAlignmentMinuitBase() {}
  30.  
  31. //=============================================================================
  32. // Update the module positions and orientations
  33. //=============================================================================
  34. void TbAlignmentMinuitBase::updateGeometry() {
  35.  
  36. info() << "Updating alignment constants:" << endmsg;
  37. for (auto im = m_modules.begin(), end = m_modules.end(); im != end; ++im) {
  38. const unsigned int i = im - m_modules.begin();
  39. const double x = (*im)->x() + m_fitter->GetParameter(6 * i + 0);
  40. const double y = (*im)->y() + m_fitter->GetParameter(6 * i + 1);
  41. const double z = (*im)->z() + m_fitter->GetParameter(6 * i + 2);
  42. const double rotx = (*im)->rotX() + m_fitter->GetParameter(6 * i + 3);
  43. const double roty = (*im)->rotY() + m_fitter->GetParameter(6 * i + 4);
  44. const double rotz = (*im)->rotZ() + m_fitter->GetParameter(6 * i + 5);
  45. info() << format("%2i", i) << format(" %-15s", (*im)->id().c_str())
  46. << format(" %10.3f", (*im)->x()) << format(" %10.3f", (*im)->y())
  47. << format(" %10.3f", (*im)->z()) << format(" %10.3f", (*im)->rotX())
  48. << format(" %10.3f", (*im)->rotY())
  49. << format(" %10.3f", (*im)->rotZ()) << endmsg;
  50. (*im)->setAlignment(x, y, z, rotx, roty, rotz, 0., 0., 0., 0., 0., 0.);
  51. }
  52. }
  53.  
  54. //=============================================================================
  55. // Main alignment function
  56. //=============================================================================
  57. void TbAlignmentMinuitBase::align(
  58. std::vector<TbAlignmentTrack*>& alignmentTracks) {
  59.  
  60. m_tracks = alignmentTracks;
  61. // Set up the fitter.
  62. m_fitter = new TFitter(60);
  63. Tb::gTbAlignmentMinuitBase = this;
  64. m_fitter->SetFCN(Tb::fcn);
  65. // Set the output level (-1: none, 0: minimum: 1: normal).
  66. double arglist = -1.;
  67. if (msgLevel(MSG::DEBUG)) arglist = 1.;
  68. m_fitter->ExecuteCommand("SET PRINTOUT", &arglist, 1);
  69. // Set the value defining parameter errors.
  70. arglist = 1.;
  71. m_fitter->ExecuteCommand("SET ERR", &arglist, 1);
  72. // Set the strategy for calculating derivatives.
  73. arglist = m_fitStrategy;
  74. m_fitter->ExecuteCommand("SET STRATEGY", &arglist, 1);
  75. // Set the floating point accuracy.
  76. arglist = 1.e-12;
  77. m_fitter->ExecuteCommand("SET EPS", &arglist, 1);
  78. setParameters();
  79. }
  80.  
  81. //=============================================================================
  82. // Set the initial values of the fit parameters.
  83. //=============================================================================
  84. void TbAlignmentMinuitBase::setParameters() {
  85.  
  86. info() << "Setting initial values of alignment parameters" << endmsg;
  87. for (auto im = m_modules.begin(), end = m_modules.end(); im != end; ++im) {
  88. const unsigned int i = im - m_modules.begin();
  89. const double dx = (*im)->dX();
  90. const double dy = (*im)->dY();
  91. const double dz = (*im)->dZ();
  92. const double drx = (*im)->dRotX();
  93. const double dry = (*im)->dRotY();
  94. const double drz = (*im)->dRotZ();
  95. const std::string id = (*im)->id();
  96. m_fitter->SetParameter(6 * i + 0, (id + " dx ").c_str(), dx, 0.001, 0., 0.);
  97. m_fitter->SetParameter(6 * i + 1, (id + " dy ").c_str(), dy, 0.001, 0., 0.);
  98. m_fitter->SetParameter(6 * i + 2, (id + " dz ").c_str(), dz, 1., 0., 0.);
  99. m_fitter->SetParameter(6 * i + 3, (id + " dRx").c_str(), drx, 0.001, 0.,
  100. 0.);
  101. m_fitter->SetParameter(6 * i + 4, (id + " dRy").c_str(), dry, 0.001, 0.,
  102. 0.);
  103. m_fitter->SetParameter(6 * i + 5, (id + " dRz").c_str(), drz, 0.001, 0.,
  104. 0.);
  105. }
  106. }