ValidationFormats/interface/MaterialAccountingStep.h

0001 #ifndef MaterialAccountingStep_h
0002 #define MaterialAccountingStep_h
0003
0004 #include <utility>
0005 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0006
0007 // struct to keep material accounting informations on a per-step basis
0008 // TODO split segment info (in, out) into separate child class
0009 class MaterialAccountingStep {
0010 public:
0011   MaterialAccountingStep(void) : m_length(0.), m_radiationLengths(0.), m_energyLoss(0.), m_in(), m_out() {}
0012
0013   MaterialAccountingStep(double position, double radlen, double loss, const GlobalPoint& in, const GlobalPoint& out)
0014       : m_length(position), m_radiationLengths(radlen), m_energyLoss(loss), m_in(in), m_out(out) {}
0015
0016   void clear(void) {
0017     m_length = 0.;
0018     m_radiationLengths = 0.;
0019     m_energyLoss = 0.;
0020     m_in = GlobalPoint();
0021     m_out = GlobalPoint();
0022   }
0023
0024 private:
0025   double m_length;
0026   double m_radiationLengths;
0027   double m_energyLoss;
0028   GlobalPoint m_in;
0029   GlobalPoint m_out;
0030
0031 public:
0032   double length(void) const { return m_length; }
0033
0034   double radiationLengths(void) const { return m_radiationLengths; }
0035
0036   double energyLoss(void) const { return m_energyLoss; }
0037
0038   const GlobalPoint& in(void) const { return m_in; }
0039
0040   const GlobalPoint& out(void) const { return m_out; }
0041
0042   /// split the step (0..1) in (0..f) + (f..1) using linear interpolation
0043   std::pair<MaterialAccountingStep, MaterialAccountingStep> split(double fraction) const {
0044     // no check is done to ensure that 0 <= f <= 1 !
0045     GlobalPoint limit(m_in.x() * fraction + m_out.x() * (1. - fraction),
0046                       m_in.y() * fraction + m_out.y() * (1. - fraction),
0047                       m_in.z() * fraction + m_out.z() * (1. - fraction));
0048
0049     MaterialAccountingStep part1(
0050         fraction * m_length, fraction * m_radiationLengths, fraction * m_energyLoss, m_in, limit);
0051
0052     MaterialAccountingStep part2(
0053         (1 - fraction) * m_length, (1 - fraction) * m_radiationLengths, (1 - fraction) * m_energyLoss, limit, m_out);
0054     return std::make_pair(part1, part2);
0055   }
0056
0057   /// assignement operator
0058   MaterialAccountingStep& operator=(const MaterialAccountingStep& step) {
0059     m_length = step.m_length;
0060     m_radiationLengths = step.m_radiationLengths;
0061     m_energyLoss = step.m_energyLoss;
0062     m_in = step.m_in;
0063     m_out = step.m_out;
0064     return *this;
0065   }
0066
0067   /// add a step
0068   MaterialAccountingStep& operator+=(const MaterialAccountingStep& step) {
0069     m_length += step.m_length;
0070     m_radiationLengths += step.m_radiationLengths;
0071     m_energyLoss += step.m_energyLoss;
0072
0073     // assume that perp2 is 0 only for uninitialized steps
0074     if ((m_in.perp2() == 0.0) or (step.m_in.perp2() < m_in.perp2()))
0075       m_in = step.m_in;
0076
0077     if ((m_out.perp2() == 0.0) or (step.m_out.perp2() > m_out.perp2()))
0078       m_out = step.m_out;
0079
0080     return *this;
0081   }
0082
0083   /// subtract a step
0084   MaterialAccountingStep& operator-=(const MaterialAccountingStep& step) {
0085     m_length -= step.m_length;
0086     m_radiationLengths -= step.m_radiationLengths;
0087     m_energyLoss -= step.m_energyLoss;
0088
0089     // can anything more sensible be done for m_in and/or m_out ?
0090     if ((step.m_in.perp2() <= m_in.perp2()) and (step.m_out.perp2() >= m_in.perp2()))
0091       m_in = step.m_out;
0092
0093     if ((step.m_out.perp2() >= m_out.perp2()) and (step.m_in.perp2() <= m_out.perp2()))
0094       m_out = step.m_in;
0095
0096     return *this;
0097   }
0098
0099   /// multiply two steps, usefull to implement (co)variance
0100   MaterialAccountingStep& operator*=(const MaterialAccountingStep& step) {
0101     m_length *= step.m_length;
0102     m_radiationLengths *= step.m_radiationLengths;
0103     m_energyLoss *= step.m_energyLoss;
0104     return *this;
0105   }
0106
0107   /// multiply by a scalar
0108   MaterialAccountingStep& operator*=(double x) {
0109     m_length *= x;
0110     m_radiationLengths *= x;
0111     m_energyLoss *= x;
0112     return *this;
0113   }
0114
0115   /// divide by a scalar
0116   MaterialAccountingStep& operator/=(double x) {
0117     m_length /= x;
0118     m_radiationLengths /= x;
0119     m_energyLoss /= x;
0120     return *this;
0121   }
0122 };
0123
0124 inline MaterialAccountingStep operator+(const MaterialAccountingStep& x, const MaterialAccountingStep& y) {
0125   MaterialAccountingStep step(x);
0126   step += y;
0127   return step;
0128 }
0129
0130 inline MaterialAccountingStep operator-(const MaterialAccountingStep& x, const MaterialAccountingStep& y) {
0131   MaterialAccountingStep step(x);
0132   step -= y;
0133   return step;
0134 }
0135
0136 inline MaterialAccountingStep operator*(const MaterialAccountingStep& x, const MaterialAccountingStep& y) {
0137   MaterialAccountingStep step(x);
0138   step *= y;
0139   return step;
0140 }
0141
0142 inline MaterialAccountingStep operator*(const MaterialAccountingStep& x, double y) {
0143   MaterialAccountingStep step(x);
0144   step *= y;
0145   return step;
0146 }
0147
0148 inline MaterialAccountingStep operator*(double y, const MaterialAccountingStep& x) {
0149   MaterialAccountingStep step(x);
0150   step *= y;
0151   return step;
0152 }
0153
0154 inline MaterialAccountingStep operator/(const MaterialAccountingStep& x, double y) {
0155   MaterialAccountingStep step(x);
0156   step /= y;
0157   return step;
0158 }
0159
0160 #endif  // MaterialAccountingStep_h