MaterialAccountingStep.h

Go to the documentation of this file.

#ifndef MaterialAccountingStep_h
#define MaterialAccountingStep_h

#include <utility>
#include "DataFormats/GeometryVector/interface/GlobalPoint.h"

// struct to keep material accounting informations on a per-step basis
// TODO split segment info (in, out) into separate child class
class MaterialAccountingStep {
public:
  MaterialAccountingStep(void) : m_length(0.), m_radiationLengths(0.), m_energyLoss(0.), m_in(), m_out() {}

  MaterialAccountingStep(double position, double radlen, double loss, const GlobalPoint& in, const GlobalPoint& out)
      : m_length(position), m_radiationLengths(radlen), m_energyLoss(loss), m_in(in), m_out(out) {}

  void clear(void) {
    m_length = 0.;
    m_radiationLengths = 0.;
    m_energyLoss = 0.;
    m_in = GlobalPoint();
    m_out = GlobalPoint();
  }

private:
  double m_length;
  double m_radiationLengths;
  double m_energyLoss;
  GlobalPoint m_in;
  GlobalPoint m_out;

public:
  double length(void) const { return m_length; }

  double radiationLengths(void) const { return m_radiationLengths; }

  double energyLoss(void) const { return m_energyLoss; }

  const GlobalPoint& in(void) const { return m_in; }

  const GlobalPoint& out(void) const { return m_out; }

  std::pair<MaterialAccountingStep, MaterialAccountingStep> split(double fraction) const {
    // no check is done to ensure that 0 <= f <= 1 !
    GlobalPoint limit(m_in.x() * fraction + m_out.x() * (1. - fraction),
                      m_in.y() * fraction + m_out.y() * (1. - fraction),
                      m_in.z() * fraction + m_out.z() * (1. - fraction));

    MaterialAccountingStep part1(
        fraction * m_length, fraction * m_radiationLengths, fraction * m_energyLoss, m_in, limit);

    MaterialAccountingStep part2(
        (1 - fraction) * m_length, (1 - fraction) * m_radiationLengths, (1 - fraction) * m_energyLoss, limit, m_out);
    return std::make_pair(part1, part2);
  }

  MaterialAccountingStep& operator+=(const MaterialAccountingStep& step) {
    m_length += step.m_length;
    m_radiationLengths += step.m_radiationLengths;
    m_energyLoss += step.m_energyLoss;

    // assume that perp2 is 0 only for uninitialized steps
    if ((m_in.perp2() == 0.0) or (step.m_in.perp2() < m_in.perp2()))
      m_in = step.m_in;

    if ((m_out.perp2() == 0.0) or (step.m_out.perp2() > m_out.perp2()))
      m_out = step.m_out;

    return *this;
  }

  MaterialAccountingStep& operator-=(const MaterialAccountingStep& step) {
    m_length -= step.m_length;
    m_radiationLengths -= step.m_radiationLengths;
    m_energyLoss -= step.m_energyLoss;

    // can anything more sensible be done for m_in and/or m_out ?
    if ((step.m_in.perp2() <= m_in.perp2()) and (step.m_out.perp2() >= m_in.perp2()))
      m_in = step.m_out;

    if ((step.m_out.perp2() >= m_out.perp2()) and (step.m_in.perp2() <= m_out.perp2()))
      m_out = step.m_in;

    return *this;
  }

  MaterialAccountingStep& operator*=(const MaterialAccountingStep& step) {
    m_length *= step.m_length;
    m_radiationLengths *= step.m_radiationLengths;
    m_energyLoss *= step.m_energyLoss;
    return *this;
  }

  MaterialAccountingStep& operator*=(double x) {
    m_length *= x;
    m_radiationLengths *= x;
    m_energyLoss *= x;
    return *this;
  }

  MaterialAccountingStep& operator/=(double x) {
    m_length /= x;
    m_radiationLengths /= x;
    m_energyLoss /= x;
    return *this;
  }
};

inline MaterialAccountingStep operator+(const MaterialAccountingStep& x, const MaterialAccountingStep& y) {
  MaterialAccountingStep step(x);
  step += y;
  return step;
}

inline MaterialAccountingStep operator-(const MaterialAccountingStep& x, const MaterialAccountingStep& y) {
  MaterialAccountingStep step(x);
  step -= y;
  return step;
}

inline MaterialAccountingStep operator*(const MaterialAccountingStep& x, const MaterialAccountingStep& y) {
  MaterialAccountingStep step(x);
  step *= y;
  return step;
}

inline MaterialAccountingStep operator*(const MaterialAccountingStep& x, double y) {
  MaterialAccountingStep step(x);
  step *= y;
  return step;
}

inline MaterialAccountingStep operator*(double y, const MaterialAccountingStep& x) {
  MaterialAccountingStep step(x);
  step *= y;
  return step;
}

inline MaterialAccountingStep operator/(const MaterialAccountingStep& x, double y) {
  MaterialAccountingStep step(x);
  step /= y;
  return step;
}

#endif  // MaterialAccountingStep_h