TrackingTools/interface/MuonErrorMatrix.h

0001 #ifndef MUONERRORMATRIX_H
0002 #define MUONERRORMATRIX_H
0003
0004 /** \class MuonErrorMatrix
0005  *
0006  * This class holds an error matrix paramertization. Either the error matrix value or scale factors.
0007  * The implementation uses TProfile3D with pt, eta, phi axis.
0008  * Error/Scale factor matrix is obtained using get(GlobalVector momentum)
0009  *
0010  *
0011  * \author Jean-Roch Vlimant  UCSB
0012  * \author Finn Rebassoo      UCSB
0013  */
0014
0015 #include <TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h>
0016 #include <TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h>
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018
0019 #include <TFile.h>
0020 #include <TProfile3D.h>
0021 #include <TString.h>
0022 #include <TAxis.h>
0023
0024 class MuonErrorMatrix {
0025 public:
0026   /// enum type to define if the class is used as a tool or to be created
0027   enum action { use, constructor };
0028   /// constructor from a parameter set
0029   MuonErrorMatrix(const edm::ParameterSet& pset);
0030
0031   /// destructor
0032   ~MuonErrorMatrix();
0033
0034   /// close the root file attached to the class
0035   void close();
0036
0037   /// main method to be used. Retrieve a 5x5 symetrical matrix according to parametrization of error or scale factor
0038   CurvilinearTrajectoryError get(GlobalVector momentum, bool convolute = true);
0039   CurvilinearTrajectoryError getFast(GlobalVector momentum);
0040
0041   /// multiply term by term the two matrix
0042   static void multiply(CurvilinearTrajectoryError& initial_error, const CurvilinearTrajectoryError& scale_error);
0043
0044   /// divide term by term the two matrix
0045   static bool divide(CurvilinearTrajectoryError& num_error, const CurvilinearTrajectoryError& denom_error);
0046
0047   /// actually get access to the TProfile3D used for the parametrization
0048   inline TProfile3D* get(int i, int j) { return Index(i, j); }
0049   inline unsigned int index(int i, int j) { return Pindex(i, j); }
0050
0051   /// names of the variables of the 5x5 error matrix
0052   static const TString vars[5];
0053
0054   /// provide the numerical value used. sigma or correlation factor
0055   static double Term(const AlgebraicSymMatrix55& curv, int i, int j);
0056
0057   ///method to get the bin index, taking care of under/overlow: first(1)/last(GetNbins())returned
0058   int findBin(TAxis* axis, double value);
0059
0060   /// convert sigma2/COV -> sigma/rho
0061   void simpleTerm(const AlgebraicSymMatrix55& input, AlgebraicSymMatrix55& output);
0062
0063   ///convert sigma/rho -> sigma2/COV
0064   void complicatedTerm(const AlgebraicSymMatrix55& input, AlgebraicSymMatrix55& output);
0065
0066   /// adjust the error matrix on the state
0067   void adjust(FreeTrajectoryState& state);
0068
0069   /// adjust the error matrix on the state
0070   void adjust(TrajectoryStateOnSurface& state);
0071
0072 private:
0073   /// log category: "MuonErrorMatrix"
0074   std::string theCategory;
0075
0076   /// the attached root file, where the parametrization is saved
0077   TDirectory* theD;
0078   /// 15 TProfile, each holding he parametrization of each term of the 5x5
0079   TProfile3D* theData[15];
0080   TProfile3D* theData_fast[5][5];
0081
0082   /// decide whether to scale of to assigne terms
0083   enum TermAction { error, scale, assign };
0084   TermAction theTermAction[15];
0085
0086   /// internal methods to get the index of a matrix term.
0087   inline int Pindex(int i, int j) {
0088     static const int offset[5] = {0, 5, 5 + 4, 5 + 4 + 3, 5 + 4 + 3 + 2};
0089     return offset[i] + abs(j - i);
0090   }
0091   /// internal method to get access to the profiles
0092   inline TProfile3D* Index(int i, int j) { return theData[Pindex(i, j)]; }
0093
0094   /// internal method that retreives the value of the parametrization for term i,j
0095   double Value(GlobalVector& momentum, int i, int j, bool convolute = true);
0096   /// internal method that retreives the error on the value of the parametrization for term i,j
0097   double Rms(GlobalVector& momentum, int i, int j);
0098 };
0099
0100 #endif