Forward declaration to allow less code verbosity

Hector/Beamline.h

@@ -1,7 +1,7 @@
 #ifndef Hector_Beamline_h
 #define Hector_Beamline_h
 
-#include "Hector/Elements/ElementBase.h"
+#include "Hector/Elements/ElementBaseFwd.h"
 #include "Hector/Elements/Marker.h"
 
 #include <map>
@@ -21,7 +21,7 @@ namespace hector {
     /// Build a beamline from a longitudinal size and a interaction point position
     /// \param[in] length Longitudinal length of the beamline
     /// \param[in] ip Position of the interaction point
-    Beamline(double length, const std::shared_ptr<element::ElementBase>& ip = nullptr);
+    Beamline(double length, const element::ElementPtr& ip = nullptr);
     ~Beamline();
 
     /// Compute all drifts between each element in the beamline
@@ -35,42 +35,42 @@ namespace hector {
     void dump(std::ostream& os, bool show_drifts = true) const;
 
     /// Set the position of the interaction point
-    void setInteractionPoint(std::shared_ptr<element::ElementBase> ip) { ip_ = ip; }
+    void setInteractionPoint(const element::ElementPtr& ip) { ip_ = ip; }
     /// Retrieve the position of the interaction point
-    const std::shared_ptr<element::ElementBase>& interactionPoint() const { return ip_; }
+    const element::ElementPtr& interactionPoint() const { return ip_; }
 
     /// Add a new element in the beamline
     /// \param[in] elem Element to be copied and added to the beamline
-    void add(const std::shared_ptr<element::ElementBase> elem);
+    void add(const element::ElementPtr& elem);
     /// Get the full beamline content (vector of elements)
-    const Elements& elements() const { return elements_; }
+    const element::Elements& elements() const { return elements_; }
     /// Get the full beamline content (vector of elements)
-    Elements& elements() { return elements_; }
+    element::Elements& elements() { return elements_; }
     /// Retrieve a beamline element given its name
     /// \param[in] name Name of the element to be retrieved
-    std::shared_ptr<element::ElementBase>& get(const std::string& name);
+    element::ElementPtr& get(const std::string& name);
     /// Retrieve a beamline element given its name
     /// \param[in] name Name of the element to be retrieved
-    const std::shared_ptr<element::ElementBase> get(const std::string& name) const;
+    const element::ElementPtr& get(const std::string& name) const;
     /// Retrieve a beamline element given its s-position
     /// \param[in] s s-position of the element (computed wrt the interaction point)
-    std::shared_ptr<element::ElementBase>& get(double s);
+    element::ElementPtr& get(double s);
     /// Retrieve a beamline element given its s-position
     /// \param[in] s s-position of the element (computed wrt the interaction point)
-    const std::shared_ptr<element::ElementBase> get(double s) const;
+    const element::ElementPtr& get(double s) const;
     /// Find an element by name
-    std::vector<std::shared_ptr<element::ElementBase> > find(const std::string&);
+    element::Elements find(const std::string&);
     /// Number of elements in the beamline
     unsigned short numElements() const { return elements_.size(); }
 
     /// Iterator to the first element in the beamline
-    Elements::iterator begin() { return elements_.begin(); }
+    element::Elements::iterator begin() { return elements_.begin(); }
     /// Iterator to the last element in the beamline
-    Elements::iterator end() { return elements_.end(); }
+    element::Elements::iterator end() { return elements_.end(); }
     /// Constant iterator to the first element in the beamline
-    const Elements::const_iterator begin() const { return elements_.begin(); }
+    const element::Elements::const_iterator begin() const { return elements_.begin(); }
     /// Constant iterator to the last element in the beamline
-    const Elements::const_iterator end() const { return elements_.end(); }
+    const element::Elements::const_iterator end() const { return elements_.end(); }
 
     /// Add a new marker element to the collection
     void addMarker(const element::Marker& marker);
@@ -102,10 +102,10 @@ namespace hector {
     /// Beamline maximal length (in m)
     double max_length_;
     /// Pointer to the interaction point
-    std::shared_ptr<element::ElementBase> ip_;
+    element::ElementPtr ip_;
 
     /// List of elements defining the beamline
-    Elements elements_;
+    element::Elements elements_;
     /// List of markers in the beamline
     MarkersMap markers_;
   };

Hector/Elements/ElementBase.h

@@ -181,8 +181,6 @@ namespace hector {
       }
     };
   }  // namespace element
-  /// List of elements
-  typedef std::vector<std::shared_ptr<element::ElementBase> > Elements;
   /// Human-readable printout of the properties of an element
   std::ostream& operator<<(std::ostream&, const element::ElementBase&);
   /// Human-readable printout of the properties of an element

Hector/Elements/ElementBaseFwd.h

@@ -0,0 +1,18 @@
+#ifndef Hector_Elements_ElementBaseFwd_h
+#define Hector_Elements_ElementBaseFwd_h
+
+#include <memory>
+#include <vector>
+
+namespace hector {
+  /// Collection of beamline elements
+  namespace element {
+    class ElementBase;
+    /// A smart pointer to a beamline element
+    typedef std::shared_ptr<element::ElementBase> ElementPtr;
+    /// List of elements
+    typedef std::vector<ElementPtr> Elements;
+  }  // namespace element
+}  // namespace hector
+
+#endif

Hector/ParticleStoppedException.h

@@ -3,6 +3,7 @@
 
 #include "Hector/Exception.h"
 #include "Hector/Elements/ElementBase.h"
+#include "Hector/Elements/ElementBaseFwd.h"
 
 namespace hector {
   /// Specific exception for particles stopped in the course of a beamline propagation
@@ -15,19 +16,19 @@ namespace hector {
     /// \param[in] additional_info Any other useful information to help the debugging
     ParticleStoppedException(const std::string& from,
                              ExceptionType type = ExceptionType::undefined,
-                             const element::ElementBase* elem = nullptr)
+                             const element::ElementPtr& elem = nullptr)
         : Exception(from, type, 10000), elem_(elem) {
       message_ << "Particle stopped";
       if (elem)
         message_ << " at " << elem->name() << " (" << elem->type() << ")";
       message_ << ".\n";
     }
     /// Retrieve the beamline element that stopped the particle
-    const element::ElementBase* stoppingElement() const { return elem_; }
+    const element::ElementPtr& stoppingElement() const { return elem_; }
 
   private:
     /// Beamline element that stopped the particle
-    const element::ElementBase* elem_;
+    const element::ElementPtr& elem_;
   };
 }  // namespace hector
 

src/Beamline.cc

@@ -21,7 +21,7 @@ namespace hector {
       setElements(rhs);
   }
 
-  Beamline::Beamline(double length, const std::shared_ptr<element::ElementBase>& ip)
+  Beamline::Beamline(double length, const element::ElementPtr& ip)
       : max_length_(length + 5.),  // artificially increase the size to include next elements
         ip_(ip) {}
 
@@ -36,7 +36,7 @@ namespace hector {
     markers_.insert(std::pair<double, element::Marker>(marker.s(), marker));
   }
 
-  void Beamline::add(const std::shared_ptr<element::ElementBase> elem) {
+  void Beamline::add(const element::ElementPtr& elem) {
     const double new_size = elem->s() + elem->length();
     if (new_size > max_length_ && max_length_ < 0.)
       throw H_ERROR << "Element " << elem->name() << " is too far away for this beamline!\n"
@@ -73,7 +73,7 @@ namespace hector {
               << "Hector will fix the overlap by splitting the earlier.";
       const double prev_length = prev_elem->length();
 
-      std::shared_ptr<element::ElementBase> next_elem = nullptr;
+      element::ElementPtr next_elem = nullptr;
       // check if one needs to add an extra piece to the previous element
       if (elem->s() + elem->length() < prev_elem->s() + prev_elem->length()) {
         const std::string prev_name = prev_elem->name();
@@ -106,35 +106,34 @@ namespace hector {
     std::sort(elements_.begin(), elements_.end(), element::ElementsSorter());
   }
 
-  const std::shared_ptr<element::ElementBase> Beamline::get(const std::string& name) const {
-    for (size_t i = 0; i < elements_.size(); ++i) {
-      const auto elem = elements_.at(i);
+  const element::ElementPtr& Beamline::get(const std::string& name) const {
+    for (const auto& elem : elements_)
       if (elem->name().find(name) != std::string::npos)
         return elem;
-    }
-    return 0;
+    H_WARNING << "Beamline element \"" << name << "\" not found.";
+    return *elements_.end();
   }
 
-  std::shared_ptr<element::ElementBase>& Beamline::get(const std::string& name) {
-    for (auto& elem : elements_) {
+  element::ElementPtr& Beamline::get(const std::string& name) {
+    for (auto& elem : elements_)
       if (elem->name().find(name) != std::string::npos)
         return elem;
-    }
+    H_WARNING << "Beamline element \"" << name << "\" not found.";
     return *elements_.end();
   }
 
-  const std::shared_ptr<element::ElementBase> Beamline::get(double s) const {
-    for (size_t i = 0; i < elements_.size(); ++i) {
-      const auto elem = elements_.at(i);
+  const element::ElementPtr& Beamline::get(double s) const {
+    for (const auto& elem : elements_) {
       if (elem->s() > s)
         continue;
       if (elem->s() <= s && elem->s() + elem->length() >= s)
         return elem;
     }
-    return 0;
+    H_WARNING << "Beamline has no element at s=" << s << ".";
+    return *elements_.end();
   }
 
-  std::shared_ptr<element::ElementBase>& Beamline::get(double s) {
+  element::ElementPtr& Beamline::get(double s) {
     for (auto& elem : elements_) {
       if (elem->s() > s)
         continue;
@@ -144,11 +143,11 @@ namespace hector {
     return *elements_.end();
   }
 
-  std::vector<std::shared_ptr<element::ElementBase> > Beamline::find(const std::string& regex) {
+  element::Elements Beamline::find(const std::string& regex) {
     try {
       std::regex rgx_search(regex);
       std::cmatch m;
-      std::vector<std::shared_ptr<element::ElementBase> > out;
+      element::Elements out;
       for (auto& elem : elements_)
         if (std::regex_search(elem->name().c_str(), m, rgx_search))
           out.emplace_back(elem);
@@ -162,7 +161,7 @@ namespace hector {
     Matrix out = DiagonalMatrix(6, 1);
 
     for (const auto& elem : elements_) {
-      const auto mat = elem->matrix(eloss, mp, qp);
+      const auto& mat = elem->matrix(eloss, mp, qp);
       H_INFO << "Multiplication by transfer matrix of element \"" << elem->name() << "\".\n"
              << " value: " << mat;
       out = out * mat;
@@ -174,14 +173,15 @@ namespace hector {
   double Beamline::length() const {
     if (elements_.empty())
       return 0.;
-    const auto elem = *elements_.rbegin();
-    return (elem->s() + elem->length());
+    const auto& elem = *elements_.rbegin();
+    return elem->s() + elem->length();
   }
 
   void Beamline::dump(std::ostream& os, bool show_drifts) const {
-    os << "=============== Beamline dump ===============\n"
-       << " interaction point: " << ip_->name() << " at " << ip_->position() << "\n"
-       << " length: " << length() << " m\n"
+    os << "=============== Beamline dump ===============\n";
+    if (ip_)
+      os << " interaction point: " << ip_->name() << " at " << ip_->position() << "\n";
+    os << " length: " << length() << " m\n"
        << " elements list: ";
     for (const auto& elemPtr : elements_) {
       if (!show_drifts && elemPtr->type() == element::aDrift)

src/Particle.cc

@@ -48,8 +48,8 @@ namespace hector {
 
     //PrintInfo( Form( "Interpolating for s = %.2f between %.2f and %.2f", s, lower_it->first, upper_it->first ) );
 
-    const StateVector sv_before = lower_it->second, sv_after = upper_it->second;
-    const TwoVector in = sv_before.position(), out = sv_after.position();
+    const StateVector &sv_before = lower_it->second, &sv_after = upper_it->second;
+    const TwoVector &in = sv_before.position(), &out = sv_after.position();
 
     const double drift_length = upper_it->first - lower_it->first;
     if (drift_length == 0)

src/Propagator.cc

@@ -23,7 +23,7 @@ namespace hector {
       return;
     }
     try {
-      for (Elements::const_iterator it = beamline_->begin() + 1; it != beamline_->end(); ++it) {
+      for (auto it = beamline_->begin() + 1; it != beamline_->end(); ++it) {
         // extract the previous and the current element in the beamline
         const auto prev_elem = *(it - 1), elem = *it;
         if (elem->s() > s_max)
@@ -69,14 +69,14 @@ namespace hector {
 
         const TwoVector pos_prev_elem(part.stateVectorAt(prev_elem->s()).position());
         if (!aper->contains(pos_prev_elem))
-          throw ParticleStoppedException(__PRETTY_FUNCTION__, ExceptionType::warning, prev_elem.get())
+          throw ParticleStoppedException(__PRETTY_FUNCTION__, ExceptionType::warning, prev_elem)
               << "Entering at " << pos_prev_elem << ", s = " << prev_elem->s() << " m\n\t"
               << "Aperture centre at " << aper->position() << "\n\t"
               << "Distance to aperture centre: " << (aper->position() - pos_prev_elem).mag() * 1.e2 << " cm.";
         // has passed through the element?
         //std::cout << prev_elem->s()+prev_elem->length() << "\t" << part.stateVectorAt( prev_elem->s()+prev_elem->length() ).position() << std::endl;
         if (!aper->contains(part.stateVectorAt(prev_elem->s() + prev_elem->length()).position()))
-          throw ParticleStoppedException(__PRETTY_FUNCTION__, ExceptionType::warning, prev_elem.get())
+          throw ParticleStoppedException(__PRETTY_FUNCTION__, ExceptionType::warning, prev_elem)
               << "Did not pass aperture " << aper->type() << ".";
       }
     } catch (const ParticleStoppedException&) {
@@ -90,7 +90,7 @@ namespace hector {
   }
 
   bool Propagator::stopped(Particle& part, double s_max) const {
-    for (Elements::const_iterator it = beamline_->begin() + 1; it != beamline_->end(); ++it) {
+    for (auto it = beamline_->begin() + 1; it != beamline_->end(); ++it) {
       // extract the previous and the current element in the beamline
       const auto prev_elem = *(it - 1), elem = *it;
       if (s_max > 0 && elem->s() > s_max)
@@ -110,7 +110,7 @@ namespace hector {
   }
 
   Particle::Position Propagator::propagateThrough(const Particle::Position& ini_pos,
-                                                  const std::shared_ptr<element::ElementBase> elem,
+                                                  const element::ElementPtr elem,
                                                   double eloss,
                                                   int qp) const {
     try {

test/hitmap.cxx

@@ -69,7 +69,7 @@ int main(int argc, char* argv[]) {
   //hector::BeamProducer::Xscanner gun( num_particles, hector::Parameters::get()->beamEnergy(), 0., 0.01 );
 
   unsigned short num_stopped = 0;
-  map<const hector::element::ElementBase*, unsigned short> stopped_at;
+  map<const hector::element::ElementPtr, unsigned short> stopped_at;
   for (size_t i = 0; i < num_particles; ++i) {
     if ((int)(i * (double)num_particles / num_particles) % 1000 == 0)
       H_INFO << ">>> Generating particle " << i << " / " << num_particles;