/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_5_3_14/src/Alignment/SurveyAnalysis/src/SurveyAlignment.cc

Go to the documentation of this file.

#include "Alignment/CommonAlignment/interface/Alignable.h"
#include "Alignment/CommonAlignment/interface/AlignmentParameters.h"
#include "Alignment/SurveyAnalysis/interface/SurveyOutput.h"

#include "Alignment/SurveyAnalysis/interface/SurveyAlignment.h"

using namespace align;

SurveyAlignment::SurveyAlignment(const Alignables& sensors,
                                 const std::vector<StructureType>& levels):
  theSensors(sensors),
  theLevels(levels)
{
}

void SurveyAlignment::shiftSensors()
{
  unsigned int nSensor = theSensors.size();

  for (unsigned int i = 0; i < nSensor; ++i)
  {
    Alignable* ali = theSensors[i];

    const AlignableSurface& surf = ali->surface();
    const AlgebraicVector&  pars = ali->alignmentParameters()->parameters();

    EulerAngles angles(3);

    angles(1) = pars[3]; angles(2) = pars[4]; angles(3) = pars[5];

    RotationType rot = surf.toGlobal( toMatrix(angles) );

    rectify(rot); // correct for rounding errors

    ali->move( surf.toGlobal( align::LocalVector(pars[0], pars[1], pars[2]) ) );
    ali->rotateInGlobalFrame(rot);
  }
}

void SurveyAlignment::iterate(unsigned int nIteration,
                              const std::string& fileName,
                              bool bias)
{
  static const double tolerance = 1e-4; // convergence criteria

  SurveyOutput out(theSensors, fileName);

  out.write(0);

  for (unsigned int i = 1; i <= nIteration; ++i)
  {
    std::cout << "***** Iteration " << i << " *****\n";
    findAlignPars(bias);
    shiftSensors();
    out.write(i);

  // Check convergence

    double parChi2 = 0.;

    unsigned int nSensor = theSensors.size();

    for (unsigned int j = 0; j < nSensor; ++j)
    {
      AlignmentParameters* alignPar = theSensors[j]->alignmentParameters();

      const AlgebraicVector&    par = alignPar->parameters();
      const AlgebraicSymMatrix& cov = alignPar->covariance();

      int dummy;

      parChi2 += cov.inverse(dummy).similarity(par);
    }

    parChi2 /= static_cast<double>(nSensor);
    std::cout << "chi2 = " << parChi2 << std::endl;
    if (parChi2 < tolerance) break; // converges, so exit loop
  }
}