CMS 3D CMS Logo

List of all members | Public Member Functions | Protected Member Functions | Protected Attributes
pftools::LinearCalibrator Class Reference

#include <LinearCalibrator.h>

Inheritance diagram for pftools::LinearCalibrator:
pftools::Calibrator

Public Member Functions

LinearCalibratorclone () const
 
LinearCalibratorcreate () const
 
 LinearCalibrator ()
 
virtual ~LinearCalibrator ()
 
- Public Member Functions inherited from pftools::Calibrator
void addDetectorElement (DetectorElementPtr const de)
 
void addParticleDeposit (ParticleDepositPtr pd)
 
 Calibrator ()
 
std::map< DetectorElementPtr, double > getCalibrationCoefficients () noexcept(false)
 
std::vector< ParticleDepositPtrgetParticles ()
 
int hasParticles () const
 
virtual ~Calibrator ()
 

Protected Member Functions

virtual std::map< DetectorElementPtr, double > getCalibrationCoefficientsCore () noexcept(false)
 
virtual TMatrixD & getHessian (const TMatrixD &eij, TMatrixD &hess, const TVectorD &truthE) const
 
virtual TVectorD & getProjections (const TMatrixD &eij, TVectorD &proj, const TVectorD &truthE) const
 
virtual void initEijMatrix (TMatrixD &eij, TVectorD &truthE)
 
 LinearCalibrator (const LinearCalibrator &lc)
 
virtual void populateDetElIndex ()
 

Protected Attributes

std::map< DetectorElementPtr, unsigned > myDetElIndex
 
- Protected Attributes inherited from pftools::Calibrator
std::vector< DetectorElementPtrmyDetectorElements
 
std::vector< ParticleDepositPtrmyParticleDeposits
 

Detailed Description

Definition at line 24 of file LinearCalibrator.h.

Constructor & Destructor Documentation

LinearCalibrator::LinearCalibrator ( )

Definition at line 42 of file LinearCalibrator.cc.

Referenced by clone(), and create().

42  {
43 
44 }
LinearCalibrator::~LinearCalibrator ( )
virtual

Definition at line 46 of file LinearCalibrator.cc.

46  {
47 }
LinearCalibrator::LinearCalibrator ( const LinearCalibrator lc)
protected

Definition at line 52 of file LinearCalibrator.cc.

References pftools::Calibrator::myDetectorElements, and pftools::Calibrator::myParticleDeposits.

52  {
55 }
std::vector< ParticleDepositPtr > myParticleDeposits
Definition: Calibrator.h:61
std::vector< DetectorElementPtr > myDetectorElements
Definition: Calibrator.h:60

Member Function Documentation

LinearCalibrator * LinearCalibrator::clone ( void  ) const
virtual

Implements pftools::Calibrator.

Definition at line 48 of file LinearCalibrator.cc.

References LinearCalibrator().

48  {
49  return new LinearCalibrator(*this);
50 }
LinearCalibrator * LinearCalibrator::create ( ) const
virtual

Implements pftools::Calibrator.

Definition at line 57 of file LinearCalibrator.cc.

References LinearCalibrator().

57  {
58  return new LinearCalibrator();
59 }
std::map< DetectorElementPtr, double > LinearCalibrator::getCalibrationCoefficientsCore ( )
protectedvirtualnoexcept

Reimplemented from pftools::Calibrator.

Definition at line 61 of file LinearCalibrator.cc.

References gather_cfg::cout, MillePedeFileConverter_cfg::e, getHessian(), getProjections(), pftools::Calibrator::hasParticles(), initEijMatrix(), myDetElIndex, and convertSQLiteXML::ok.

62  {
63 // std::cout << __PRETTY_FUNCTION__
64 // << ": determining linear calibration coefficients...\n";
65  if (!hasParticles()) {
66  //I have no particles to calibrate to - throw exception.
67  PFToolsException me("Calibrator has no particles for calibration!");
68  throw me;
69  }
70  //std::cout << "\tGetting eij matrix...\n";
71  TMatrixD eij;
72  TVectorD truthE;
73  initEijMatrix(eij, truthE);
74 
75 // std::cout << "\tEij matrix:\n";
76 // printMat(std::cout, eij);
77 
78  //std::cout << "\tGetting projections...\n";
79  TVectorD proj;
80  TMatrixD hess;
81 
82  getProjections(eij, proj, truthE);
83  //std::cout << "\tProjections:\n";
84  //printVec(std::cout, proj);
85  getHessian(eij, hess, truthE);
86 
87  TDecompLU lu;
88  lu.SetMatrix(hess);
89  //std::cout << "\tHessian:\n";
90  //printMat(std::cout, hess);
91 
92  lu.SetTol(1e-25);
93  TMatrixD hessInv = lu.Invert();
94  //std::cout <<"\tInverse Hessian:\n";
95  //printMat(std::cout, hessInv);
96  TVectorD calibsSolved(eij.GetNcols());
97 
98  bool ok(true);
99  calibsSolved = lu.Solve(proj, ok);
100  if (ok) {
101  //std::cout << "\tLU reports ok.\n";
102  }
103  else {
104  std::cout << "\tWARNING: LU reports NOT ok.\n";
105  //This usually happens when you've asked the calibrator to solve for the 'a' term, without including
106  //a dummy 'a' term in each particle deposit.
107  //Make sure you do
108  //Deposition dOffset(offset, eta, phi, 1.0);
109  //particle->addRecDeposition(dOffset);
110  //particle->addTruthDeposition(dOffset);
111  std::cout
112  << "\tDid you forget to add a dummy offset deposition to each particle?\n";
113  std::cout << "\tThrowing an exception!"<< std::endl;
115  me("TDecompLU did not converge successfully when finding calibrations. Did you forget to add a dummy offset deposition to each particle?");
116  throw me;
117  }
118 
119  //std::cout << "\tCalibrations: \n";
120  //printVec(std::cout, calibsSolved);
121 
122  std::map<DetectorElementPtr, double> answers;
123  for (std::map<DetectorElementPtr, unsigned>::iterator it =
124  myDetElIndex.begin(); it != myDetElIndex.end(); ++it) {
125  DetectorElementPtr de = (*it).first;
126  answers[de] = calibsSolved[(*it).second];
127  }
128  return answers;
129 }
virtual TVectorD & getProjections(const TMatrixD &eij, TVectorD &proj, const TVectorD &truthE) const
int hasParticles() const
Definition: Calibrator.h:46
boost::shared_ptr< DetectorElement > DetectorElementPtr
General purpose exception class for use by classes in the pftools namespace.
virtual TMatrixD & getHessian(const TMatrixD &eij, TMatrixD &hess, const TVectorD &truthE) const
virtual void initEijMatrix(TMatrixD &eij, TVectorD &truthE)
std::map< DetectorElementPtr, unsigned > myDetElIndex
TMatrixD & LinearCalibrator::getHessian ( const TMatrixD &  eij,
TMatrixD &  hess,
const TVectorD &  truthE 
) const
protectedvirtual

Definition at line 179 of file LinearCalibrator.cc.

References i, j, gen::n, and funct::pow().

Referenced by getCalibrationCoefficientsCore().

180  {
181  //std::cout << __PRETTY_FUNCTION__ << "\n";
182  unsigned nCalibs(eij.GetNcols());
183  hess.ResizeTo(nCalibs, nCalibs);
184  hess.Zero();
185 
186  for (unsigned i(0); i < nCalibs; ++i) {
187  for (unsigned j(0); j < nCalibs; ++j) {
188  for (int n(0); n < eij.GetNrows(); ++n) {
189  hess[i][j] += eij[n][i] * eij[n][j]/ pow(truthE[n], 2.0);
190  }
191  }
192  }
193 
194  return hess;
195 }
int i
Definition: DBlmapReader.cc:9
int j
Definition: DBlmapReader.cc:9
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
TVectorD & LinearCalibrator::getProjections ( const TMatrixD &  eij,
TVectorD &  proj,
const TVectorD &  truthE 
) const
protectedvirtual

Definition at line 164 of file LinearCalibrator.cc.

References i, and j.

Referenced by getCalibrationCoefficientsCore().

165  {
166  //std::cout << __PRETTY_FUNCTION__ << "\n";
167  proj.ResizeTo(eij.GetNcols());
168  proj.Zero();
169 
170  for (int j(0); j < eij.GetNcols(); ++j) {
171  for (int i(0); i < eij.GetNrows(); ++i) {
172  proj[j] += eij[i][j] / truthE[i];
173  }
174  }
175 
176  return proj;
177 }
int i
Definition: DBlmapReader.cc:9
int j
Definition: DBlmapReader.cc:9
void LinearCalibrator::initEijMatrix ( TMatrixD &  eij,
TVectorD &  truthE 
)
protectedvirtual

Definition at line 131 of file LinearCalibrator.cc.

References diffTreeTool::index, pftools::Calibrator::myDetectorElements, myDetElIndex, pftools::Calibrator::myParticleDeposits, AlCaHLTBitMon_ParallelJobs::p, and populateDetElIndex().

Referenced by getCalibrationCoefficientsCore().

131  {
132  //std::cout << __PRETTY_FUNCTION__ << "\n";
133  //std::cout << "\tGetting detector element indices...\n";
135  eij.Clear();
136  eij.Zero();
137 
138  truthE.ResizeTo(myParticleDeposits.size());
139  truthE.Zero();
140 
141  //First determine number of calibration constants.
142 
143  eij.ResizeTo(myParticleDeposits.size(), myDetElIndex.size());
144 
145  //loop over all particle deposits
146  unsigned index(0);
147  for (std::vector<ParticleDepositPtr>::const_iterator cit =
148  myParticleDeposits.begin(); cit != myParticleDeposits.end(); ++cit) {
149  ParticleDepositPtr p = *cit;
150  //get each of the relevant detector elements
151 
152  for (std::vector<DetectorElementPtr>::const_iterator detElIt =
153  myDetectorElements.begin(); detElIt != myDetectorElements.end(); ++detElIt) {
154  DetectorElementPtr de = *detElIt;
155  eij[index][myDetElIndex[de]] = p->getRecEnergy(de);
156  //truthE[p->getId()] += p->getTruthEnergy(de);
157  }
158  truthE[index] = p->getTruthEnergy();
159  ++index;
160  }
161 
162 }
std::vector< ParticleDepositPtr > myParticleDeposits
Definition: Calibrator.h:61
boost::shared_ptr< ParticleDeposit > ParticleDepositPtr
boost::shared_ptr< DetectorElement > DetectorElementPtr
std::vector< DetectorElementPtr > myDetectorElements
Definition: Calibrator.h:60
std::map< DetectorElementPtr, unsigned > myDetElIndex
void LinearCalibrator::populateDetElIndex ( )
protectedvirtual

Definition at line 197 of file LinearCalibrator.cc.

References diffTreeTool::index, pftools::Calibrator::myDetectorElements, and myDetElIndex.

Referenced by initEijMatrix().

197  {
198  //reserve index = 0 for the constant term, if we're told to compute it
199  unsigned index(0);
200 
201  myDetElIndex.clear();
202  //loop over known detector elements, and assign a unique row/column index to each
203  for (std::vector<DetectorElementPtr>::const_iterator cit =
204  myDetectorElements.begin(); cit != myDetectorElements.end(); ++cit) {
205  DetectorElementPtr de = *cit;
206  //std::cout << "\tGot element: "<< *de;
207  //check we don't have duplicate detector elements
208  if (myDetElIndex.count(de) == 0) {
209  myDetElIndex[de] = index;
210  ++index;
211  }
212  }
213 
214 }
boost::shared_ptr< DetectorElement > DetectorElementPtr
std::vector< DetectorElementPtr > myDetectorElements
Definition: Calibrator.h:60
std::map< DetectorElementPtr, unsigned > myDetElIndex

Member Data Documentation

std::map<DetectorElementPtr, unsigned> pftools::LinearCalibrator::myDetElIndex
protected