#include <CSCChamberFitter.h>

Public Member Functions
	CSCChamberFitter (const edm::ParameterSet &iConfig, std::vector< CSCPairResidualsConstraint * > &residualsConstraints)

bool	fit (std::vector< CSCAlignmentCorrections * > &corrections) const

void	radiusCorrection (AlignableNavigator alignableNavigator, AlignmentParameterStore alignmentParameterStore, bool combineME11) const

virtual	~CSCChamberFitter ()

Protected Member Functions
long	alignableId (std::string alignable) const

double	chi2 (const AlgebraicVector &A, double lambda) const

double	hessian (int k, int l, double lambda) const

int	index (std::string alignable) const

bool	isFrame (int i) const

double	lhsVector (int k) const

void	walk (std::map< int, bool > &touched, int alignable) const

Protected Attributes
std::vector< std::string >	m_alignables

std::vector< CSCPairConstraint * >	m_constraints

int	m_fixed

std::vector< int >	m_frames

std::string	m_name

Detailed Description

Date:: 2010/05/06 05:57:33

Revision:: 1.2

Author: J. Pivarski - Texas A&M University pivar.nosp@m.ski@.nosp@m.physi.nosp@m.cs.t.nosp@m.amu.e.nosp@m.du

Definition at line 24 of file CSCChamberFitter.h.

Constructor & Destructor Documentation

CSCChamberFitter::CSCChamberFitter	(	const edm::ParameterSet &	iConfig,
		std::vector< CSCPairResidualsConstraint * > &	residualsConstraints
	)

Definition at line 14 of file CSCChamberFitter.cc.

References alignableId(), CSCDetId::chamber(), HLT_FULL_cff::constraint, createBeamHaloJobs::constraints, CSCDetId::endcap(), relativeConstraints::error, Exception, edm::ParameterSet::getParameter(), mps_fire::i, index(), dqmiolumiharvest::j, m_alignables, m_constraints, m_fixed, m_frames, m_name, CSCDetId::ring(), CSCDetId::station(), AlCaHLTBitMon_QueryRunRegistry::string, relativeConstraints::value, and walk().

                                                                                                   {
   m_name = iConfig.getParameter<std::string>("name");
   m_alignables = iConfig.getParameter<std::vector<std::string> >("alignables");
   if (m_alignables.empty()) {
     throw cms::Exception("BadConfig") << "Fitter " << m_name << " has no alignables!" << std::endl;
   }
 
   int i = 0;
   for (std::vector<std::string>::const_iterator alignable = m_alignables.begin(); alignable != m_alignables.end();
        ++alignable) {
     if (alignableId(*alignable) == -1)
       m_frames.push_back(i);
     i++;
   }
 
   m_fixed = -1;
   std::string fixed = iConfig.getParameter<std::string>("fixed");
   if (!fixed.empty()) {
     int i = 0;
     for (std::vector<std::string>::const_iterator alignable = m_alignables.begin(); alignable != m_alignables.end();
          ++alignable) {
       if (fixed == *alignable) {
         m_fixed = i;
       }
       i++;
     }
     if (m_fixed == -1)
       throw cms::Exception("BadConfig") << "Cannot fix unrecognized alignable " << fixed << std::endl;
   }
 
   int numConstraints = 0;
   std::vector<edm::ParameterSet> constraints = iConfig.getParameter<std::vector<edm::ParameterSet> >("constraints");
   for (std::vector<edm::ParameterSet>::const_iterator constraint = constraints.begin(); constraint != constraints.end();
        ++constraint) {
     int i = index(constraint->getParameter<std::string>("i"));
     int j = index(constraint->getParameter<std::string>("j"));
     double value = constraint->getParameter<double>("value");
     double error = constraint->getParameter<double>("error");
 
     if (i < 0)
       throw cms::Exception("BadConfig") << "Unrecognized alignable " << constraint->getParameter<std::string>("i")
                                         << " in constraint " << numConstraints << " of fitter " << m_name << std::endl;
     if (j < 0)
       throw cms::Exception("BadConfig") << "Unrecognized alignable " << constraint->getParameter<std::string>("j")
                                         << " in constraint " << numConstraints << " of fitter " << m_name << std::endl;
     if (error <= 0.)
       throw cms::Exception("BadConfig") << "Non-positive uncertainty in constraint " << numConstraints << " of fitter "
                                         << m_name << std::endl;
     if (i == j)
       throw cms::Exception("BadConfig") << "Self-connection from " << constraint->getParameter<std::string>("i")
                                         << " to " << constraint->getParameter<std::string>("j")
                                         << " is not allowed in constraint " << numConstraints << " of fitter " << m_name
                                         << std::endl;
 
     m_constraints.push_back(new CSCPairConstraint(i, j, value, error));
     numConstraints++;
   }
 
   // insert CSCPairResidualsConstraints
   for (unsigned int i = 0; i < m_alignables.size(); i++) {
     std::string alignable_i = m_alignables[i];
     long id_i = alignableId(alignable_i);
     if (id_i != -1) {
       CSCDetId cscid_i(id_i);
 
       for (unsigned int j = 0; j < m_alignables.size(); j++) {
         std::string alignable_j = m_alignables[j];
         long id_j = alignableId(alignable_j);
         if (i != j && id_j != -1) {
           CSCDetId cscid_j(id_j);
 
           if (!(cscid_i.station() == 1 && cscid_i.ring() == 3 && cscid_j.station() == 1 && cscid_j.ring() == 3)) {
             int next_chamber = cscid_i.chamber() + 1;
             if (cscid_i.station() > 1 && cscid_i.ring() == 1 && next_chamber == 19)
               next_chamber = 1;
             else if (!(cscid_i.station() > 1 && cscid_i.ring() == 1) && next_chamber == 37)
               next_chamber = 1;
             if (cscid_i.endcap() == cscid_j.endcap() && cscid_i.station() == cscid_j.station() &&
                 cscid_i.ring() == cscid_j.ring() && next_chamber == cscid_j.chamber()) {
               CSCPairResidualsConstraint *residualsConstraint =
                   new CSCPairResidualsConstraint(residualsConstraints.size(), i, j, cscid_i, cscid_j);
               m_constraints.push_back(residualsConstraint);
               residualsConstraints.push_back(residualsConstraint);
               numConstraints++;
             }
           }
         }
       }
     }
   }
 
   std::map<int, bool> touched;
   for (unsigned int i = 0; i < m_alignables.size(); i++)
     touched[i] = false;
   walk(touched, 0);
   for (unsigned int i = 0; i < m_alignables.size(); i++) {
     if (!touched[i])
       throw cms::Exception("BadConfig") << "Fitter " << m_name << " is not a connected graph (no way to get to "
                                         << m_alignables[i] << " from " << m_alignables[0] << ", for instance)"
                                         << std::endl;
   }
 }

virtual CSCChamberFitter::~CSCChamberFitter ( )

inlinevirtual

Definition at line 27 of file CSCChamberFitter.h.

27 {};

Member Function Documentation

long CSCChamberFitter::alignableId ( std::string alignable ) const

protected

Definition at line 143 of file CSCChamberFitter.cc.

References CSCDetId, Reference_intrackfit_cff::endcap, relativeConstraints::ring, and relativeConstraints::station.

Referenced by CSCChamberFitter(), and fit().

                                                             {
   if (alignable.size() != 9)
     return -1;
 
   if (alignable[0] == 'M' && alignable[1] == 'E') {
     int endcap = -1;
     if (alignable[2] == '+')
       endcap = 1;
     else if (alignable[2] == '-')
       endcap = 2;
 
     if (endcap != -1) {
       int station = -1;
       if (alignable[3] == '1')
         station = 1;
       else if (alignable[3] == '2')
         station = 2;
       else if (alignable[3] == '3')
         station = 3;
       else if (alignable[3] == '4')
         station = 4;
 
       if (alignable[4] == '/' && station != -1) {
         int ring = -1;
         if (alignable[5] == '1')
           ring = 1;
         else if (alignable[5] == '2')
           ring = 2;
         else if (alignable[5] == '3')
           ring = 3;
         else if (alignable[5] == '4')
           ring = 4;
         if (station > 1 && ring > 2)
           return -1;
 
         if (alignable[6] == '/' && ring != -1) {
           int chamber = -1;
           if (alignable[7] == '0' && alignable[8] == '1')
             chamber = 1;
           else if (alignable[7] == '0' && alignable[8] == '2')
             chamber = 2;
           else if (alignable[7] == '0' && alignable[8] == '3')
             chamber = 3;
           else if (alignable[7] == '0' && alignable[8] == '4')
             chamber = 4;
           else if (alignable[7] == '0' && alignable[8] == '5')
             chamber = 5;
           else if (alignable[7] == '0' && alignable[8] == '6')
             chamber = 6;
           else if (alignable[7] == '0' && alignable[8] == '7')
             chamber = 7;
           else if (alignable[7] == '0' && alignable[8] == '8')
             chamber = 8;
           else if (alignable[7] == '0' && alignable[8] == '9')
             chamber = 9;
           else if (alignable[7] == '1' && alignable[8] == '0')
             chamber = 10;
           else if (alignable[7] == '1' && alignable[8] == '1')
             chamber = 11;
           else if (alignable[7] == '1' && alignable[8] == '2')
             chamber = 12;
           else if (alignable[7] == '1' && alignable[8] == '3')
             chamber = 13;
           else if (alignable[7] == '1' && alignable[8] == '4')
             chamber = 14;
           else if (alignable[7] == '1' && alignable[8] == '5')
             chamber = 15;
           else if (alignable[7] == '1' && alignable[8] == '6')
             chamber = 16;
           else if (alignable[7] == '1' && alignable[8] == '7')
             chamber = 17;
           else if (alignable[7] == '1' && alignable[8] == '8')
             chamber = 18;
           else if (alignable[7] == '1' && alignable[8] == '9')
             chamber = 19;
           else if (alignable[7] == '2' && alignable[8] == '0')
             chamber = 20;
           else if (alignable[7] == '2' && alignable[8] == '1')
             chamber = 21;
           else if (alignable[7] == '2' && alignable[8] == '2')
             chamber = 22;
           else if (alignable[7] == '2' && alignable[8] == '3')
             chamber = 23;
           else if (alignable[7] == '2' && alignable[8] == '4')
             chamber = 24;
           else if (alignable[7] == '2' && alignable[8] == '5')
             chamber = 25;
           else if (alignable[7] == '2' && alignable[8] == '6')
             chamber = 26;
           else if (alignable[7] == '2' && alignable[8] == '7')
             chamber = 27;
           else if (alignable[7] == '2' && alignable[8] == '8')
             chamber = 28;
           else if (alignable[7] == '2' && alignable[8] == '9')
             chamber = 29;
           else if (alignable[7] == '3' && alignable[8] == '0')
             chamber = 30;
           else if (alignable[7] == '3' && alignable[8] == '1')
             chamber = 31;
           else if (alignable[7] == '3' && alignable[8] == '2')
             chamber = 32;
           else if (alignable[7] == '3' && alignable[8] == '3')
             chamber = 33;
           else if (alignable[7] == '3' && alignable[8] == '4')
             chamber = 34;
           else if (alignable[7] == '3' && alignable[8] == '5')
             chamber = 35;
           else if (alignable[7] == '3' && alignable[8] == '6')
             chamber = 36;
 
           if (station > 1 && ring == 1 && chamber > 18)
             return -1;
 
           if (chamber != -1) {
             return CSCDetId(endcap, station, ring, chamber, 0).rawId();
           }
         }
       }
     }
   }
 
   return -1;
 }

double CSCChamberFitter::chi2	(	const AlgebraicVector &	A,
		double	lambda
	)		const

protected

Definition at line 275 of file CSCChamberFitter.cc.

References HLT_FULL_cff::constraint, mps_fire::i, isFrame(), m_alignables, m_constraints, m_fixed, funct::pow(), and alignCSCRings::s.

Referenced by fit().

                                                                            {
   double sumFixed = 0.;
 
   if (m_fixed == -1) {
     for (unsigned int i = 0; i < m_alignables.size(); i++) {
       if (!isFrame(i)) {
         sumFixed += A[i];
       }
     }
   } else {
     sumFixed = A[m_fixed];
   }
 
   double s = lambda * sumFixed * sumFixed;
   for (std::vector<CSCPairConstraint *>::const_iterator constraint = m_constraints.begin();
        constraint != m_constraints.end();
        ++constraint) {
     if ((*constraint)->valid()) {
       s += pow((*constraint)->value() - A[(*constraint)->i()] + A[(*constraint)->j()], 2) / (*constraint)->error() /
            (*constraint)->error();
     }
   }
   return s;
 }

bool CSCChamberFitter::fit ( std::vector< CSCAlignmentCorrections * > & corrections ) const

Definition at line 343 of file CSCChamberFitter.cc.

References A, alignableId(), chi2(), HLT_FULL_cff::constraint, pixelCPEforGPU::correction(), CSCDetId, hessian(), mps_fire::i, infinity, CSCAlignmentCorrections::insertCorrection(), CSCAlignmentCorrections::insertMode(), CSCAlignmentCorrections::insertResidual(), isFrame(), dqmiolumiharvest::j, isotrackApplyRegressor::k, cmsLHEtoEOSManager::l, lhsVector(), m_alignables, m_constraints, m_name, mathSSE::sqrt(), and cms::cuda::V.

Referenced by trackingPlots.Iteration::modules().

                                                                                   {
   double lambda = 1. / infinity / infinity;
 
   AlgebraicVector A(m_alignables.size());
   AlgebraicVector V(m_alignables.size());
   AlgebraicMatrix M(m_alignables.size(), m_alignables.size());
 
   for (unsigned int k = 0; k < m_alignables.size(); k++) {
     A[k] = 0.;
     V[k] = lhsVector(k);
     for (unsigned int l = 0; l < m_alignables.size(); l++) {
       M[k][l] = hessian(k, l, lambda);
     }
   }
 
   double oldchi2 = chi2(A, lambda);
 
   int ierr;
   M.invert(ierr);
   if (ierr != 0) {
     edm::LogError("CSCOverlapsAlignmentAlgorithm")
         << "Matrix inversion failed for fitter " << m_name << " matrix is " << M << std::endl;
     return false;
   }
 
   A = M * V;  // that's the alignment step
 
   CSCAlignmentCorrections *correction = new CSCAlignmentCorrections(m_name, oldchi2, chi2(A, lambda));
 
   for (unsigned int i = 0; i < m_alignables.size(); i++) {
     if (!isFrame(i)) {
       correction->insertCorrection(m_alignables[i], CSCDetId(alignableId(m_alignables[i])), A[i]);
     }
   }
 
   // we have to switch to a completely different linear algebrea
   // package because CLHEP doesn't compute
   // eigenvectors/diagonalization (?!?)
   TMatrixD tmatrix(m_alignables.size(), m_alignables.size());
   for (unsigned int i = 0; i < m_alignables.size(); i++) {
     for (unsigned int j = 0; j < m_alignables.size(); j++) {
       tmatrix[i][j] = M[i][j];
     }
   }
   TMatrixDEigen tmatrixdeigen(tmatrix);
   const TMatrixD &basis = tmatrixdeigen.GetEigenVectors();
   TMatrixD invbasis = tmatrixdeigen.GetEigenVectors();
   invbasis.Invert();
   TMatrixD diagonalized = invbasis * (tmatrix * basis);
 
   for (unsigned int i = 0; i < m_alignables.size(); i++) {
     std::vector<double> coefficient;
     std::vector<std::string> modename;
     std::vector<long> modeid;
     for (unsigned int j = 0; j < m_alignables.size(); j++) {
       coefficient.push_back(invbasis[i][j]);
       modename.push_back(m_alignables[j]);
       modeid.push_back(alignableId(m_alignables[j]));
     }
 
     correction->insertMode(
         coefficient, modename, modeid, sqrt(2. * fabs(diagonalized[i][i])) * (diagonalized[i][i] >= 0. ? 1. : -1.));
   }
 
   for (std::vector<CSCPairConstraint *>::const_iterator constraint = m_constraints.begin();
        constraint != m_constraints.end();
        ++constraint) {
     if ((*constraint)->valid()) {
       double residual = (*constraint)->value() - A[(*constraint)->i()] + A[(*constraint)->j()];
       correction->insertResidual(m_alignables[(*constraint)->i()],
                                  m_alignables[(*constraint)->j()],
                                  (*constraint)->value(),
                                  (*constraint)->error(),
                                  residual,
                                  residual / (*constraint)->error());
     }
   }
 
   corrections.push_back(correction);
   return true;
 }

double CSCChamberFitter::hessian	(	int	k,
		int	l,
		double	lambda
	)		const

protected

Definition at line 316 of file CSCChamberFitter.cc.

References HLT_FULL_cff::constraint, ztail::d, infinity, isFrame(), isotrackApplyRegressor::k, cmsLHEtoEOSManager::l, m_constraints, m_fixed, and alignCSCRings::s.

Referenced by fit().

                                                                   {
   double s = 0.;
 
   if (m_fixed == -1) {
     if (!isFrame(k) && !isFrame(l))
       s += 2. * lambda;
   } else {
     if (k == l && l == m_fixed)
       s += 2. * lambda;
   }
 
   for (std::vector<CSCPairConstraint *>::const_iterator constraint = m_constraints.begin();
        constraint != m_constraints.end();
        ++constraint) {
     double d = 2. / infinity / infinity;
     if ((*constraint)->valid()) {
       d = 2. / (*constraint)->error() / (*constraint)->error();
     }
 
     if (k == l && ((*constraint)->i() == k || (*constraint)->j() == k))
       s += d;
     if (((*constraint)->i() == k && (*constraint)->j() == l) || ((*constraint)->j() == k && (*constraint)->i() == l))
       s -= d;
   }
   return s;
 }

int CSCChamberFitter::index ( std::string alignable ) const

protected

Definition at line 118 of file CSCChamberFitter.cc.

References a, mps_fire::i, and m_alignables.

Referenced by BeautifulSoup.PageElement::_invert(), and CSCChamberFitter().

                                                      {
   int i = 0;
   for (std::vector<std::string>::const_iterator a = m_alignables.begin(); a != m_alignables.end(); ++a) {
     if (*a == alignable)
       return i;
     i++;
   }
   return -1;
 }

bool CSCChamberFitter::isFrame ( int i ) const

protected

Definition at line 267 of file CSCChamberFitter.cc.

References m_frames.

Referenced by chi2(), fit(), and hessian().

                                           {
   for (std::vector<int>::const_iterator frame = m_frames.begin(); frame != m_frames.end(); ++frame) {
     if (i == *frame)
       return true;
   }
   return false;
 }

double CSCChamberFitter::lhsVector ( int k ) const

protected

Definition at line 300 of file CSCChamberFitter.cc.

References HLT_FULL_cff::constraint, ztail::d, isotrackApplyRegressor::k, m_constraints, and alignCSCRings::s.

Referenced by fit().

                                               {
   double s = 0.;
   for (std::vector<CSCPairConstraint *>::const_iterator constraint = m_constraints.begin();
        constraint != m_constraints.end();
        ++constraint) {
     if ((*constraint)->valid()) {
       double d = 2. * (*constraint)->value() / (*constraint)->error() / (*constraint)->error();
       if ((*constraint)->i() == k)
         s += d;
       if ((*constraint)->j() == k)
         s -= d;
     }
   }
   return s;
 }

void CSCChamberFitter::radiusCorrection	(	AlignableNavigator *	alignableNavigator,
		AlignmentParameterStore *	alignmentParameterStore,
		bool	combineME11
	)		const

Definition at line 426 of file CSCChamberFitter.cc.

References AlignableDetOrUnitPtr::alignable(), AlignableNavigator::alignableFromDetId(), Alignable::alignmentParameters(), AlignmentParameterStore::applyParameters(), CSCDetId::chamber(), AlignmentParameters::cloneFromSelected(), HLT_FULL_cff::constraint, alignCSCRings::e, CSCDetId::endcap(), gpuClustering::id, CSCPairResidualsConstraint::id_i(), m_constraints, M_PI, submitPVValidationJobs::params, applyRadialCorrections::radial_correction, CSCPairResidualsConstraint::radius(), CSCDetId::ring(), Alignable::setAlignmentParameters(), AlignmentParameters::setValid(), CSCDetId::station(), CSCPairResidualsConstraint::valid(), and CSCPairResidualsConstraint::value().

                                                                 {
   double sum_phipos_residuals = 0.;
   double num_valid = 0.;
   double sum_radius = 0.;
   double num_total = 0.;
   for (std::vector<CSCPairConstraint *>::const_iterator constraint = m_constraints.begin();
        constraint != m_constraints.end();
        ++constraint) {
     CSCPairResidualsConstraint *residualsConstraint = dynamic_cast<CSCPairResidualsConstraint *>(*constraint);
     if (residualsConstraint != nullptr) {
       if (residualsConstraint->valid()) {
         sum_phipos_residuals += residualsConstraint->value();
         num_valid += 1.;
       }
 
       sum_radius += residualsConstraint->radius(true);
       num_total += 1.;
     }
   }
   if (num_valid == 0. || num_total == 0.)
     return;
   double average_phi_residual = sum_phipos_residuals / num_valid;
   double average_radius = sum_radius / num_total;
 
   double radial_correction = average_phi_residual * average_radius * num_total / (2. * M_PI);
 
   for (std::vector<CSCPairConstraint *>::const_iterator constraint = m_constraints.begin();
        constraint != m_constraints.end();
        ++constraint) {
     CSCPairResidualsConstraint *residualsConstraint = dynamic_cast<CSCPairResidualsConstraint *>(*constraint);
     if (residualsConstraint != nullptr) {
       const DetId id(residualsConstraint->id_i());
       Alignable *alignable = alignableNavigator->alignableFromDetId(id).alignable();
       Alignable *also = nullptr;
       if (combineME11 && residualsConstraint->id_i().station() == 1 && residualsConstraint->id_i().ring() == 1) {
         CSCDetId alsoid(residualsConstraint->id_i().endcap(), 1, 4, residualsConstraint->id_i().chamber(), 0);
         const DetId alsoid2(alsoid);
         also = alignableNavigator->alignableFromDetId(alsoid2).alignable();
       }
 
       AlgebraicVector params(6);
       AlgebraicSymMatrix cov(6);
 
       params[1] = radial_correction;
       cov[1][1] = 1e-6;
 
       AlignmentParameters *parnew = alignable->alignmentParameters()->cloneFromSelected(params, cov);
       alignable->setAlignmentParameters(parnew);
       alignmentParameterStore->applyParameters(alignable);
       alignable->alignmentParameters()->setValid(true);
       if (also != nullptr) {
         AlignmentParameters *parnew2 = also->alignmentParameters()->cloneFromSelected(params, cov);
         also->setAlignmentParameters(parnew2);
         alignmentParameterStore->applyParameters(also);
         also->alignmentParameters()->setValid(true);
       }
     }
   }
 }

void CSCChamberFitter::walk	(	std::map< int, bool > &	touched,
		int	alignable
	)		const

protected

Definition at line 128 of file CSCChamberFitter.cc.

References HLT_FULL_cff::constraint, and m_constraints.

Referenced by CSCChamberFitter().

                                                                            {
   touched[alignable] = true;
 
   for (std::vector<CSCPairConstraint *>::const_iterator constraint = m_constraints.begin();
        constraint != m_constraints.end();
        ++constraint) {
     if (alignable == (*constraint)->i() || alignable == (*constraint)->j()) {
       if (!touched[(*constraint)->i()])
         walk(touched, (*constraint)->i());
       if (!touched[(*constraint)->j()])
         walk(touched, (*constraint)->j());
     }
   }
 }