CommonAlignmentParametrization/src/BeamSpotAlignmentParameters.cc

0001 /** \file BeamSpotAlignmentParameters.cc
0002  *
0003  *  Version    : $Revision: 1.14 $
0004  *  last update: $Date: 2008/09/02 15:08:12 $
0005  *  by         : $Author: flucke $
0006  */
0007
0008 #include "FWCore/Utilities/interface/Exception.h"
0009
0010 #include "Alignment/CommonAlignment/interface/Alignable.h"
0011 #include "Alignment/CommonAlignment/interface/AlignableDetOrUnitPtr.h"
0012 #include "Alignment/CommonAlignment/interface/Utilities.h"
0013 #include "Alignment/CommonAlignmentParametrization/interface/AlignmentParametersFactory.h"
0014 #include "Alignment/CommonAlignmentParametrization/interface/BeamSpotAlignmentDerivatives.h"
0015 #include "CondFormats/Alignment/interface/Definitions.h"
0016
0017 // This class's header
0018 #include "Alignment/CommonAlignmentParametrization/interface/BeamSpotAlignmentParameters.h"
0019
0020 //__________________________________________________________________________________________________
0021 BeamSpotAlignmentParameters::BeamSpotAlignmentParameters(Alignable *ali, bool calcMis)
0022     : AlignmentParameters(ali, displacementFromAlignable(calcMis ? ali : nullptr), AlgebraicSymMatrix(N_PARAM, 0)) {}
0023
0024 //__________________________________________________________________________________________________
0025 BeamSpotAlignmentParameters::BeamSpotAlignmentParameters(Alignable *alignable,
0026                                                          const AlgebraicVector &parameters,
0027                                                          const AlgebraicSymMatrix &covMatrix)
0028     : AlignmentParameters(alignable, parameters, covMatrix) {
0029   if (parameters.num_row() != N_PARAM) {
0030     throw cms::Exception("BadParameters")
0031         << "in BeamSpotAlignmentParameters(): " << parameters.num_row() << " instead of " << N_PARAM << " parameters.";
0032   }
0033 }
0034
0035 //__________________________________________________________________________________________________
0036 BeamSpotAlignmentParameters::BeamSpotAlignmentParameters(Alignable *alignable,
0037                                                          const AlgebraicVector &parameters,
0038                                                          const AlgebraicSymMatrix &covMatrix,
0039                                                          const std::vector<bool> &selection)
0040     : AlignmentParameters(alignable, parameters, covMatrix, selection) {
0041   if (parameters.num_row() != N_PARAM) {
0042     throw cms::Exception("BadParameters")
0043         << "in BeamSpotAlignmentParameters(): " << parameters.num_row() << " instead of " << N_PARAM << " parameters.";
0044   }
0045 }
0046
0047 //__________________________________________________________________________________________________
0048 BeamSpotAlignmentParameters::~BeamSpotAlignmentParameters() {}
0049
0050 //__________________________________________________________________________________________________
0051 BeamSpotAlignmentParameters *BeamSpotAlignmentParameters::clone(const AlgebraicVector &parameters,
0052                                                                 const AlgebraicSymMatrix &covMatrix) const {
0053   BeamSpotAlignmentParameters *rbap = new BeamSpotAlignmentParameters(alignable(), parameters, covMatrix, selector());
0054
0055   if (userVariables())
0056     rbap->setUserVariables(userVariables()->clone());
0057   rbap->setValid(isValid());
0058
0059   return rbap;
0060 }
0061
0062 //__________________________________________________________________________________________________
0063 BeamSpotAlignmentParameters *BeamSpotAlignmentParameters::cloneFromSelected(const AlgebraicVector &parameters,
0064                                                                             const AlgebraicSymMatrix &covMatrix) const {
0065   BeamSpotAlignmentParameters *rbap = new BeamSpotAlignmentParameters(
0066       alignable(), expandVector(parameters, selector()), expandSymMatrix(covMatrix, selector()), selector());
0067
0068   if (userVariables())
0069     rbap->setUserVariables(userVariables()->clone());
0070   rbap->setValid(isValid());
0071
0072   return rbap;
0073 }
0074
0075 //__________________________________________________________________________________________________
0076 AlgebraicMatrix BeamSpotAlignmentParameters::derivatives(const TrajectoryStateOnSurface &tsos,
0077                                                          const AlignableDetOrUnitPtr &alidet) const {
0078   const Alignable *ali = this->alignable();  // Alignable of these parameters
0079
0080   if (ali == alidet) {  // same alignable => same frame
0081     return BeamSpotAlignmentDerivatives()(tsos);
0082   } else {
0083     throw cms::Exception("MisMatch") << "BeamSpotAlignmentParameters::derivatives: The hit alignable must "
0084                                         "match the "
0085                                      << "aligned one, i.e. these parameters make only sense for "
0086                                         "AlignableBeamSpot.\n";
0087     return AlgebraicMatrix(N_PARAM, 2);  // please compiler
0088   }
0089 }
0090
0091 //__________________________________________________________________________________________________
0092 AlgebraicMatrix BeamSpotAlignmentParameters::selectedDerivatives(const TrajectoryStateOnSurface &tsos,
0093                                                                  const AlignableDetOrUnitPtr &alignableDet) const {
0094   const AlgebraicMatrix dev = this->derivatives(tsos, alignableDet);
0095
0096   int ncols = dev.num_col();
0097   int nrows = dev.num_row();
0098   int nsel = numSelected();
0099
0100   AlgebraicMatrix seldev(nsel, ncols);
0101
0102   int ir2 = 0;
0103   for (int irow = 0; irow < nrows; ++irow) {
0104     if (selector()[irow]) {
0105       for (int icol = 0; icol < ncols; ++icol)
0106         seldev[ir2][icol] = dev[irow][icol];
0107       ++ir2;
0108     }
0109   }
0110
0111   return seldev;
0112 }
0113
0114 //__________________________________________________________________________________________________
0115 AlgebraicVector BeamSpotAlignmentParameters::translation(void) const {
0116   AlgebraicVector shift(3);
0117   for (int i = 0; i < 2; ++i)
0118     shift[i] = theData->parameters()[i];
0119   shift[2] = 0.0;
0120
0121   return shift;
0122 }
0123
0124 //__________________________________________________________________________________________________
0125 AlgebraicVector BeamSpotAlignmentParameters::rotation(void) const {
0126   AlgebraicVector rot(3);
0127
0128   double dxdz = theData->parameters()[2];
0129   double dydz = theData->parameters()[3];
0130   double angleY = std::atan(dxdz);
0131   double angleX = -std::atan(dydz);
0132
0133   align::RotationType rotY(std::cos(angleY), 0., -std::sin(angleY), 0., 1., 0., std::sin(angleY), 0., std::cos(angleY));
0134
0135   align::RotationType rotX(1., 0., 0., 0., std::cos(angleX), std::sin(angleX), 0., -std::sin(angleX), std::cos(angleX));
0136
0137   align::EulerAngles angles = align::toAngles(rotY * rotX);
0138
0139   rot[0] = angles(1);
0140   rot[1] = angles(2);
0141   rot[2] = angles(3);
0142
0143   return rot;
0144 }
0145
0146 //__________________________________________________________________________________________________
0147 void BeamSpotAlignmentParameters::apply() {
0148   Alignable *alignable = this->alignable();
0149   if (!alignable) {
0150     throw cms::Exception("BadParameters") << "BeamSpotAlignmentParameters::apply: parameters without alignable";
0151   }
0152
0153   // Translation in local frame
0154   AlgebraicVector shift = this->translation();  // fixme: should be LocalVector
0155
0156   // Translation local->global
0157   align::GlobalVector gv(shift[0], shift[1], shift[2]);
0158   alignable->move(gv);
0159
0160   // Rotation in local frame
0161   align::EulerAngles angles = this->rotation();
0162   // original code:
0163   //  alignable->rotateInLocalFrame( align::toMatrix(angles) );
0164   // correct for rounding errors:
0165   align::RotationType rot = alignable->surface().toGlobal(align::toMatrix(angles));
0166   align::rectify(rot);
0167   alignable->rotateInGlobalFrame(rot);
0168 }
0169
0170 //__________________________________________________________________________________________________
0171 int BeamSpotAlignmentParameters::type() const { return AlignmentParametersFactory::kBeamSpot; }
0172
0173 //__________________________________________________________________________________________________
0174 AlgebraicVector BeamSpotAlignmentParameters::globalParameters(void) const {
0175   AlgebraicVector m_GlobalParameters(N_PARAM, 0);
0176
0177   const AlgebraicVector shift = translation();  // fixme: should return LocalVector
0178
0179   const align::GlobalVector dg(shift[0], shift[1], 0);
0180
0181   m_GlobalParameters[0] = dg.x();
0182   m_GlobalParameters[1] = dg.y();
0183
0184   align::LocalVector lv(0.0, 0.0, 1.0);
0185   align::GlobalVector gv = theAlignable->surface().toGlobal(lv);
0186
0187   double dxdz = gv.x() / gv.z();
0188   double dydz = gv.x() / gv.z();
0189
0190   m_GlobalParameters[2] = dxdz;
0191   m_GlobalParameters[3] = dydz;
0192
0193   return m_GlobalParameters;
0194 }
0195
0196 //__________________________________________________________________________________________________
0197 void BeamSpotAlignmentParameters::print(void) const {
0198   std::cout << "Contents of BeamSpotAlignmentParameters:"
0199             << "\nParameters: " << theData->parameters() << "\nCovariance: " << theData->covariance() << std::endl;
0200 }
0201
0202 //__________________________________________________________________________________________________
0203 AlgebraicVector BeamSpotAlignmentParameters::displacementFromAlignable(const Alignable *ali) {
0204   AlgebraicVector displacement(N_PARAM);
0205
0206   if (ali) {
0207     const align::RotationType &dR = ali->rotation();
0208
0209     const align::LocalVector shifts(ali->globalRotation() * (dR.transposed() * ali->displacement().basicVector()));
0210
0211     align::GlobalVector gv(0.0, 0.0, 1.0);
0212     align::LocalVector lv(dR.transposed() * gv.basicVector());
0213
0214     displacement[0] = shifts.x();
0215     displacement[1] = shifts.y();
0216     displacement[2] = lv.x() / lv.z();
0217     displacement[3] = lv.y() / lv.z();
0218   }
0219
0220   return displacement;
0221 }