Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-04-19 09:06:55

0001 // -*- C++ -*-
0002 #ifndef RIVET_Utils_HH
0003 #define RIVET_Utils_HH
0004 
0005 #include "Rivet/Tools/RivetSTL.hh"
0006 #include "Rivet/Tools/TypeTraits.hh"
0007 #include "Rivet/Tools/PrettyPrint.hh"
0008 #include "Rivet/Tools/Exceptions.hh"
0009 #include <ostream>
0010 #include <cctype>
0011 #include <cerrno>
0012 #include <stdexcept>
0013 #include <numeric>
0014 #include <limits>
0015 #include <climits>
0016 #include <cfloat>
0017 #include <cmath>
0018 #include <sstream>
0019 #include <functional>
0020 
0021 /// @defgroup utils Other utilities
0022 
0023 
0024 /// Macro to help mark code as deprecated to produce compiler warnings
0025 #ifndef DEPRECATED
0026 #if __GNUC__ && __cplusplus && RIVET_NO_DEPRECATION_WARNINGS == 0
0027   #define DEPRECATED(x) __attribute__((deprecated(x)))
0028 #else
0029   #define DEPRECATED(x)
0030 #endif
0031 #endif
0032 
0033 
0034 namespace Rivet {
0035 
0036 
0037   /// Convenient const for getting the double NaN value
0038   static constexpr double DBL_NAN = std::numeric_limits<double>::quiet_NaN();
0039 
0040 
0041   /// @defgroup strutils String utils
0042   /// @{
0043 
0044   /// @brief Exception class for throwing from lexical_cast when a parse goes wrong
0045   struct bad_lexical_cast : public std::runtime_error {
0046     bad_lexical_cast(const std::string& what) : std::runtime_error(what) {}
0047   };
0048 
0049   /// @brief Convert between any types via stringstream
0050   template<typename T, typename U>
0051   T lexical_cast(const U& in) {
0052     try {
0053       std::stringstream ss;
0054       ss << in;
0055       T out;
0056       ss >> out;
0057       return out;
0058     } catch (const std::exception& e) {
0059       throw bad_lexical_cast(e.what());
0060     }
0061   }
0062 
0063   /// @brief Convert any object to a string
0064   ///
0065   /// Just a convenience wrapper for the more general Boost lexical_cast
0066   template <typename T>
0067   inline string to_str(const T& x) {
0068     return lexical_cast<string>(x);
0069   }
0070 
0071   /// @brief Convert any object to a string
0072   ///
0073   /// An alias for to_str() with a more "Rivety" mixedCase name.
0074   template <typename T>
0075   inline string toString(const T& x) {
0076     return lexical_cast<string>(x);
0077   }
0078 
0079   /// Replace the first instance of patt with repl
0080   inline string& replace_first(string& str, const string& patt, const string& repl) {
0081     if (!contains(str, patt)) return str; //< contains from RivetSTL
0082     str.replace(str.find(patt), patt.size(), repl);
0083     return str;
0084   }
0085 
0086   /// @brief Replace all instances of patt with repl
0087   ///
0088   /// @note Finding is interleaved with replacement, so the second search happens after
0089   /// first replacement, etc. This could lead to infinite loops and other counterintuitive
0090   /// behaviours if not careful.
0091   inline string& replace_all(string& str, const string& patt, const string& repl) {
0092     if (!contains(str, patt)) return str; //< contains from RivetSTL
0093     while (true) {
0094       string::size_type it = str.find(patt);
0095       if (it == string::npos) break;
0096       str.replace(it, patt.size(), repl);
0097     }
0098     return str;
0099   }
0100 
0101 
0102   /// Case-insensitive string comparison function
0103   inline int nocase_cmp(const string& s1, const string& s2) {
0104     string::const_iterator it1 = s1.begin();
0105     string::const_iterator it2 = s2.begin();
0106     while ( (it1 != s1.end()) && (it2 != s2.end()) ) {
0107       if(::toupper(*it1) != ::toupper(*it2)) { // < Letters differ?
0108         // Return -1 to indicate smaller than, 1 otherwise
0109         return (::toupper(*it1) < ::toupper(*it2)) ? -1 : 1;
0110       }
0111       // Proceed to the next character in each string
0112       ++it1;
0113       ++it2;
0114     }
0115     size_t size1 = s1.size(), size2 = s2.size(); // Cache lengths
0116     // Return -1,0 or 1 according to strings' lengths
0117     if (size1 == size2) return 0;
0118     return (size1 < size2) ? -1 : 1;
0119   }
0120 
0121 
0122   /// Case-insensitive string equality function
0123   inline bool nocase_equals(const string& s1, const string& s2) {
0124     return nocase_cmp(s1, s2) == 0;
0125   }
0126 
0127 
0128   /// Convert a string to lower-case
0129   inline string toLower(const string& s) {
0130     string out = s;
0131     std::transform(out.begin(), out.end(), out.begin(), (int(*)(int)) std::tolower);
0132     return out;
0133   }
0134 
0135 
0136   /// Convert a string to upper-case
0137   inline string toUpper(const string& s) {
0138     string out = s;
0139     std::transform(out.begin(), out.end(), out.begin(), (int(*)(int)) std::toupper);
0140     return out;
0141   }
0142 
0143 
0144   /// Check whether a string @a start is found at the start of @a s
0145   inline bool startsWith(const string& s, const string& start) {
0146     if (s.length() < start.length()) return false;
0147     return s.substr(0, start.length()) == start;
0148   }
0149 
0150 
0151   /// Check whether a string @a end is found at the end of @a s
0152   inline bool endsWith(const string& s, const string& end) {
0153     if (s.length() < end.length()) return false;
0154     return s.substr(s.length() - end.length()) == end;
0155   }
0156 
0157 
0158   // Terminating version of strjoin, for empty fargs list
0159   inline string strcat() { return ""; }
0160   /// Make a string containing the concatenated string representations of each item in the variadic list
0161   template<typename T, typename... Ts>
0162   inline string strcat(T value, Ts... fargs) {
0163     const string strthis = lexical_cast<string>(value);
0164     const string strnext = strcat(fargs...);
0165     return strnext.empty() ? strthis : strthis + strnext;
0166   }
0167 
0168 
0169   /// Make a string containing the string representations of each item in v, separated by sep
0170   template <typename T>
0171   inline string join(const vector<T>& v, const string& sep=" ") {
0172     string rtn;
0173     for (size_t i = 0; i < v.size(); ++i) {
0174       if (i != 0) rtn += sep;
0175       rtn += to_str(v[i]);
0176     }
0177     return rtn;
0178   }
0179 
0180   /// Make a string containing the string representations of each item in v, separated by sep
0181   template <>
0182   inline string join(const vector<string>& v, const string& sep) {
0183     string rtn;
0184     for (size_t i = 0; i < v.size(); ++i) {
0185       if (i != 0) rtn += sep;
0186       rtn += v[i];
0187     }
0188     return rtn;
0189   }
0190 
0191   /// Make a string containing the string representations of each item in s, separated by sep
0192   template <typename T>
0193   inline string join(const set<T>& s, const string& sep=" ") {
0194     string rtn;
0195     for (const T& x : s) {
0196       if (rtn.size() > 0) rtn += sep;
0197       rtn += to_str(x);
0198     }
0199     return rtn;
0200   }
0201 
0202   /// Make a string containing the string representations of each item in s, separated by sep
0203   template <>
0204   inline string join(const set<string>& s, const string& sep) {
0205     string rtn;
0206     for (const string & x : s) {
0207       if (rtn.size() > 0) rtn += sep;
0208       rtn += x;
0209     }
0210     return rtn;
0211   }
0212 
0213   /// @brief Split a string on a specified separator string
0214   inline vector<string> split(const string& s, const string& sep) {
0215     vector<string> dirs;
0216     string tmp = s;
0217     while (true) {
0218       const size_t delim_pos = tmp.find(sep);
0219       if (delim_pos == string::npos) break;
0220       const string dir = tmp.substr(0, delim_pos);
0221       if (dir.length()) dirs.push_back(dir); // Don't insert "empties"
0222       tmp.replace(0, delim_pos+1, "");
0223     }
0224     if (tmp.length()) dirs.push_back(tmp); // Don't forget the trailing component!
0225     return dirs;
0226   }
0227 
0228   /// Left-pad the given string @a s to width @a width
0229   inline string lpad(const string& s, size_t width, const string& padchar=" ") {
0230     if (s.size() >= width) return s;
0231     return string(width - s.size(), padchar[0]) + s;
0232   }
0233 
0234   /// Right-pad the given string @a s to width @a width
0235   inline string rpad(const string& s, size_t width, const string& padchar=" ") {
0236     if (s.size() >= width) return s;
0237     return s + string(width - s.size(), padchar[0]);
0238   }
0239 
0240   /// @}
0241 
0242 
0243 
0244   /// @defgroup pathutils Path utils
0245   /// @{
0246 
0247   /// @brief Split a path string with colon delimiters
0248   ///
0249   /// Ignores zero-length substrings. Designed for getting elements of filesystem paths, naturally.
0250   inline vector<string> pathsplit(const string& path) {
0251     return split(path, ":");
0252   }
0253 
0254   /// @brief Join several filesystem paths together with the standard ':' delimiter
0255   ///
0256   /// Note that this does NOT join path elements together with a platform-portable
0257   /// directory delimiter, cf. the Python @c {os.path.join}!
0258   inline string pathjoin(const vector<string>& paths) {
0259     return join(paths, ":");
0260   }
0261 
0262   /// Operator for joining strings @a a and @a b with filesystem separators
0263   inline string operator / (const string& a, const string& b) {
0264     // Ensure that a doesn't end with a slash, and b doesn't start with one, to avoid "//"
0265     const string anorm = (a.find("/") != string::npos) ? a.substr(0, a.find_last_not_of("/")+1) : a;
0266     const string bnorm = (b.find("/") != string::npos) ? b.substr(b.find_first_not_of("/")) : b;
0267     return anorm + "/" + bnorm;
0268   }
0269 
0270   /// Get the basename (i.e. terminal file name) from a path @a p
0271   inline string basename(const string& p) {
0272     if (!contains(p, "/")) return p;
0273     return p.substr(p.rfind("/")+1);
0274   }
0275 
0276   /// Get the dirname (i.e. path to the penultimate directory) from a path @a p
0277   inline string dirname(const string& p) {
0278     if (!contains(p, "/")) return "";
0279     return p.substr(0, p.rfind("/"));
0280   }
0281 
0282   /// Get the stem (i.e. part without a file extension) from a filename @a f
0283   inline string file_stem(const string& f) {
0284     if (!contains(f, ".")) return f;
0285     return f.substr(0, f.rfind("."));
0286   }
0287 
0288   /// Get the file extension from a filename @a f
0289   inline string file_extn(const string& f) {
0290     if (!contains(f, ".")) return "";
0291     return f.substr(f.rfind(".")+1);
0292   }
0293 
0294   /// @}
0295 
0296 
0297 
0298   /// @defgroup contutils Container utils
0299   /// @{
0300 
0301   /// Return number of true elements in the container @a c .
0302   template <typename CONTAINER,
0303             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0304   inline unsigned int count(const CONTAINER& c) {
0305     // return std::count_if(std::begin(c), std::end(c), [](const typename CONTAINER::value_type& x){return bool(x);});
0306     unsigned int rtn = 0;
0307     for (const auto& x : c) if (bool(x)) rtn += 1;
0308     return rtn;
0309   }
0310 
0311   // /// Return number of elements in the container @a c for which @c f(x) is true.
0312   // template <typename CONTAINER>
0313   // inline unsigned int count(const CONTAINER& c, const std::function<bool(typename CONTAINER::value_type)>& f) {
0314   //   return std::count_if(std::begin(c), std::end(c), f);
0315   // }
0316 
0317   /// Return number of elements in the container @a c for which @c f(x) is true.
0318   template <typename CONTAINER, typename FN,
0319             //typename FN = bool(const typename std::decay_t<CONTAINER>::value_type&, const std::decay_t<T>&),
0320             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0321   inline unsigned int count(const CONTAINER& c, const FN& f) {
0322     return std::count_if(std::begin(c), std::end(c), f);
0323   }
0324 
0325 
0326 
0327   /// Return true if x is true for any x in container c, otherwise false.
0328   template <typename CONTAINER,
0329             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0330   inline bool any(const CONTAINER& c) {
0331     // return std::any_of(std::begin(c), std::end(c), [](const auto& x){return bool(x);});
0332     for (const auto& x : c) if (bool(x)) return true;
0333     return false;
0334   }
0335 
0336   // /// Return true if f(x) is true for any x in container c, otherwise false.
0337   // template <typename CONTAINER>
0338   // inline bool any(const CONTAINER& c, const std::function<bool(typename CONTAINER::value_type)>& f) {
0339   //   return std::any_of(std::begin(c), std::end(c), f);
0340   // }
0341 
0342   /// Return true if f(x) is true for any x in container c, otherwise false.
0343   template <typename CONTAINER, typename FN,
0344             //typename FN = double(const typename std::decay_t<CONTAINER>::value_type&, const std::decay_t<T>&),
0345             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0346   inline bool any(const CONTAINER& c, const FN& f) {
0347     return std::any_of(std::begin(c), std::end(c), f);
0348   }
0349 
0350 
0351 
0352   /// Return true if @a x is true for all @c x in container @a c, otherwise false.
0353   template <typename CONTAINER,
0354             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0355   inline bool all(const CONTAINER& c) {
0356     // return std::all_of(std::begin(c), std::end(c), [](const auto& x){return bool(x);});
0357     for (const auto& x : c) if (!bool(x)) return false;
0358     return true;
0359   }
0360 
0361   // /// Return true if @a f(x) is true for all @c x in container @a c, otherwise false.
0362   // template <typename CONTAINER>
0363   // inline bool all(const CONTAINER& c, const std::function<bool(typename CONTAINER::value_type)>& f) {
0364   //   return std::all_of(std::begin(c), std::end(c), f);
0365   // }
0366 
0367   /// Return true if @a f(x) is true for all @c x in container @a c, otherwise false.
0368   template <typename CONTAINER, typename FN,
0369             //typename FN = double(const typename std::decay_t<CONTAINER>::value_type&, const std::decay_t<T>&),
0370             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0371   inline bool all(const CONTAINER& c, const FN& f) {
0372     return std::all_of(std::begin(c), std::end(c), f);
0373   }
0374 
0375 
0376   /// Return true if @a x is false for all @c x in container @a c, otherwise false.
0377   template <typename CONTAINER,
0378             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0379   inline bool none(const CONTAINER& c) {
0380     // return std::none_of(std::begin(c), std::end(c), [](){});
0381     for (const auto& x : c) if (bool(x)) return false;
0382     return true;
0383   }
0384 
0385   // /// Return true if @a f(x) is false for all @c x in container @a c, otherwise false.
0386   // template <typename C>
0387   // inline bool none(const C& c, const std::function<bool(typename C::value_type)>& f) {
0388   //   return std::none_of(std::begin(c), std::end(c), f);
0389   // }
0390 
0391   /// Return true if @a f(x) is false for all @c x in container @a c, otherwise false.
0392   template <typename CONTAINER, typename FN,
0393             //typename FN = double(const typename std::decay_t<CONTAINER>::value_type&, const std::decay_t<T>&),
0394             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0395   inline bool none(const CONTAINER& c, const FN& f) {
0396     return std::none_of(std::begin(c), std::end(c), f);
0397   }
0398 
0399 
0400   // /// A single-container-arg version of std::transform, aka @c map
0401   // template <typename CONTAINER1, typename CONTAINER2>
0402   // inline const CONTAINER2& transform(const CONTAINER1& in, CONTAINER2& out,
0403   //                            const std::function<typename CONTAINER2::value_type(typename CONTAINER1::value_type)>& f) {
0404   //   out.clear(); out.resize(in.size());
0405   //   std::transform(in.begin(), in.end(), out.begin(), f);
0406   //   return out;
0407   // }
0408 
0409   /// A single-container-arg version of std::transform
0410   /// @todo Make the function template specific to ParticleBase
0411   /// or introduce C++20 concepts
0412   template <typename CONTAINER1, typename CONTAINER2,
0413             typename FN = typename std::decay_t<CONTAINER2>::value_type(
0414                           const typename std::decay_t<CONTAINER1>::value_type::ParticleBase&),
0415             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER1>> &&
0416                                          is_citerable_v<std::decay_t<CONTAINER2>> >>
0417   inline const CONTAINER2& transform(const CONTAINER1& in, CONTAINER2& out, FN&& f) {
0418     out.clear(); out.resize(in.size());
0419     std::transform(in.begin(), in.end(), out.begin(), std::forward<FN>(f));
0420     return out;
0421   }
0422 
0423   /// A single-container-arg, return-value version of std::transform, aka @c map
0424   /// @todo Make the function template polymorphic... or specific to ParticleBase
0425   template <typename CONTAINER1, typename RTN,
0426             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER1>> >>
0427   inline std::vector<RTN> transform(const CONTAINER1& in,
0428                                     const std::function<RTN(typename CONTAINER1::value_type::ParticleBase)>& f) {
0429     std::vector<RTN> out;
0430     transform(in, out, f);
0431     return out;
0432   }
0433 
0434   // /// A single-container-arg version of std::accumulate, aka @c reduce
0435   // template <typename CONTAINER1, typename T>
0436   // inline T accumulate(const CONTAINER1& in, const T& init, const std::function<T(typename CONTAINER1::value_type)>& f) {
0437   //   const T rtn = std::accumulate(in.begin(), in.end(), init, f);
0438   //   return rtn;
0439   // }
0440 
0441   /// A single-container-arg version of std::accumulate, aka @c reduce
0442   template <typename CONTAINER1, typename T, typename FN,
0443             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER1>> >>
0444   inline T accumulate(const CONTAINER1& in, const T& init, const FN& f) {
0445     const T rtn = std::accumulate(in.begin(), in.end(), init, f);
0446     return rtn;
0447   }
0448 
0449 
0450   /// @brief Generic sum function, adding @c x for all @c x in container @a c
0451   ///
0452   /// @note Default-constructs the return type -- not always possible! Supply an explicit start value if necessary.
0453   template <typename CONTAINER,
0454             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0455   inline typename CONTAINER::value_type sum(const CONTAINER& c) {
0456     typename CONTAINER::value_type rtn; //< default construct return type
0457     for (const auto& x : c) rtn += x;
0458     return rtn;
0459   }
0460 
0461   /// Generic sum function, adding @c x for all @c x in container @a c, starting with @a start
0462   ///
0463   /// @note It's more more flexible here to not use CONTAINER::value_type, allowing implicit casting to T.
0464   template <typename CONTAINER, typename T,
0465             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0466   inline T sum(const CONTAINER& c, const T& start) {
0467     T rtn = start;
0468     for (const auto& x : c) rtn += x;
0469     return rtn;
0470   }
0471 
0472   /// Generic sum function, adding @a fn(@c x) for all @c x in container @a c, starting with @a start
0473   template <typename CONTAINER, typename T,
0474             typename FN = T(const typename std::decay_t<CONTAINER>::value_type&),
0475             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0476   inline T sum(CONTAINER&& c, FN&& fn, const T& start=T()) {
0477     auto f = std::function(std::forward<FN>(fn));
0478     T rtn = start;
0479     for (const auto& x : c) rtn += fn(x);
0480     return rtn;
0481   }
0482 
0483 
0484   /// In-place generic sum function, adding @c x on to container @a out for all @c x in container @a c
0485   ///
0486   /// @note It's more more flexible here to not use CONTAINER::value_type, allowing implicit casting to T.
0487   template <typename CONTAINER, typename T,
0488             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0489   inline T& isum(const CONTAINER& c, T& out) {
0490     for (const auto& x : c) out += x;
0491     return out;
0492   }
0493 
0494   /// In-place generic sum function, adding @a fn(@c x) on to container @a out for all @c x in container @a c
0495   ///
0496   /// @note It's more more flexible here to not use CONTAINER::value_type, allowing implicit casting to T.
0497   template <typename CONTAINER, typename FN, typename T,
0498             //typename FN = double(const typename std::decay_t<CONTAINER>::value_type&, const std::decay_t<T>&),
0499             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0500   inline T& isum(const CONTAINER& c, const FN& f, T& out) {
0501     for (const auto& x : c) out += f(x);
0502     return out;
0503   }
0504 
0505 
0506 
0507   /// Filter a collection in-place, removing the subset that passes the supplied function
0508   ///
0509   /// @todo Use const std::function<bool(typename CONTAINER::value_type)>... but need polymorphism for ParticleBase
0510   template <typename CONTAINER, typename FN,
0511             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0512   inline CONTAINER& idiscard(CONTAINER& c, const FN& f) {
0513     const auto newend = std::remove_if(std::begin(c), std::end(c), f);
0514     c.erase(newend, c.end());
0515     return c;
0516   }
0517 
0518   /// Version with element-equality comparison in place of a function
0519   template <typename CONTAINER,
0520             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0521   inline CONTAINER& idiscard(CONTAINER& c, const typename CONTAINER::value_type& y) {
0522     return idiscard(c, [&](typename CONTAINER::value_type& x){ return x == y; });
0523   }
0524 
0525   /// Version with several element-equality comparisons in place of a function
0526   template <typename CONTAINER,
0527             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0528   inline CONTAINER& idiscard_if_any(CONTAINER& c, const CONTAINER& ys) {
0529     return idiscard(c, [&](typename CONTAINER::value_type& x){ return contains(ys, x); });
0530   }
0531 
0532 
0533   /// Filter a collection by copy, removing the subset that passes the supplied function
0534   ///
0535   /// @todo Use const std::function<bool(typename CONTAINER::value_type)>... but need polymorphism for ParticleBase
0536   template <typename CONTAINER, typename FN,
0537             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0538   inline CONTAINER discard(const CONTAINER& c, const FN& f) {
0539     CONTAINER rtn = c;
0540     return idiscard(rtn, f); ///< @todo More efficient would be copy_if with back_inserter...
0541   }
0542 
0543   /// Version with element-equality comparison in place of a function
0544   template <typename CONTAINER,
0545             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0546   inline CONTAINER discard(const CONTAINER& c, const typename CONTAINER::value_type& y) {
0547     return discard(c, [&](typename CONTAINER::value_type& x){ return x == y; });
0548   }
0549 
0550   /// Version with several element-equality comparisons in place of a function
0551   template <typename CONTAINER,
0552             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0553   inline CONTAINER discard_if_any(const CONTAINER& c, const CONTAINER& ys) {
0554     return discard(c, [&](typename CONTAINER::value_type& x){ return contains(ys, x); });
0555   }
0556 
0557 
0558   /// Filter a collection by copy into a supplied container, removing the subset that passes the supplied function
0559   ///
0560   /// @note New container will be replaced, not appended to
0561   ///
0562   /// @todo Use const std::function<bool(typename CONTAINER::value_type)>... but need polymorphism for ParticleBase
0563   template <typename CONTAINER, typename FN,
0564             //typename FN = bool(const typename std::decay_t<CONTAINER>::value_type&, const std::decay_t<T>&),
0565             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0566   inline CONTAINER& discard(const CONTAINER& c, const FN& f, CONTAINER& out) {
0567     out = discard(c, f);
0568     return out;
0569   }
0570 
0571   /// Version with element-equality comparison in place of a function
0572   template <typename CONTAINER,
0573             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0574   inline CONTAINER& discard(const CONTAINER& c, const typename CONTAINER::value_type& y, CONTAINER& out) {
0575     return discard(c, [&](typename CONTAINER::value_type& x){ return x == y; }, out);
0576   }
0577 
0578   /// Version with several element-equality comparisons in place of a function
0579   template <typename CONTAINER,
0580             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0581   inline CONTAINER& discard_if_any(const CONTAINER& c, const CONTAINER& ys, CONTAINER& out) {
0582     return discard(c, [&](typename CONTAINER::value_type& x){ return contains(ys, x); }, out);
0583   }
0584 
0585 
0586 
0587   /// Filter a collection in-place, keeping the subset that passes the supplied function
0588   ///
0589   /// @todo Use const std::function<bool(typename CONTAINER::value_type)>... but need polymorphism for ParticleBase
0590   template <typename CONTAINER, typename FN,
0591             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0592       inline CONTAINER& iselect(CONTAINER& c, const FN& f) {
0593     auto invf = [&](const typename CONTAINER::value_type& x){ return !f(x); };
0594     return idiscard(c, invf); //< yes, intentional!
0595   }
0596 
0597   // No single equality-comparison version for select, since that would be silly!
0598 
0599   /// Version with several element-equality comparisons in place of a function
0600   template <typename CONTAINER,
0601             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0602   inline CONTAINER& iselect_if_any(CONTAINER& c, const CONTAINER& ys) {
0603     return iselect(c, [&](typename CONTAINER::value_type& x){ return contains(ys, x); });
0604   }
0605 
0606 
0607   /// Filter a collection by copy, keeping the subset that passes the supplied function
0608   ///
0609   /// @todo Use const std::function<bool(typename CONTAINER::value_type)>... but need polymorphism for ParticleBase
0610   template <typename CONTAINER, typename FN,
0611             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0612   inline CONTAINER select(const CONTAINER& c, const FN& f) {
0613     CONTAINER rtn = c;
0614     return iselect(rtn, f); ///< @todo More efficient would be copy_if with back_inserter ... but is that equally container agnostic?
0615   }
0616 
0617   // No single equality-comparison version for select, since that would be silly!
0618 
0619   /// Version with several element-equality comparisons in place of a function
0620   template <typename CONTAINER,
0621             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0622   inline CONTAINER select_if_any(const CONTAINER& c, const CONTAINER& ys) {
0623     return select(c, [&](typename CONTAINER::value_type& x){ return contains(ys, x); });
0624   }
0625 
0626 
0627   /// Filter a collection by copy into a supplied container, keeping the subset that passes the supplied function
0628   ///
0629   /// @note New container will be replaced, not appended to
0630   ///
0631   /// @todo Use const std::function<bool(typename CONTAINER::value_type)>... but need polymorphism for ParticleBase
0632   template <typename CONTAINER, typename FN,
0633             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0634   inline CONTAINER& select(const CONTAINER& c, const FN& f, CONTAINER& out) {
0635     out = select(c, f);
0636     return out;
0637   }
0638 
0639   // No single equality-comparison version for select, since that would be silly!
0640 
0641   /// Version with several element-equality comparisons in place of a function
0642   template <typename CONTAINER,
0643             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0644   inline CONTAINER& select_if_any(const CONTAINER& c, const CONTAINER& ys, CONTAINER& out) {
0645     return select(c, [&](typename CONTAINER::value_type& x){ return contains(ys, x); }, out);
0646   }
0647 
0648 
0649 
0650   /// @brief Slice of the container elements cf. Python's [i:j] syntax
0651   ///
0652   /// The element at the @a j index is not included in the returned container.
0653   /// @a i and @a j can be negative, treated as backward offsets from the end of the container.
0654   template <typename CONTAINER,
0655             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0656   inline CONTAINER slice(const CONTAINER& c, int i, int j) {
0657     CONTAINER rtn;
0658     const size_t off1 = (i >= 0) ? i : c.size() + i;
0659     const size_t off2 = (j >= 0) ? j : c.size() + j;
0660     if (off1 > c.size() || off2 > c.size()) throw RangeError("Attempting to slice beyond requested offsets");
0661     if (off2 < off1) throw RangeError("Requested offsets in invalid order");
0662     rtn.resize(off2 - off1);
0663     std::copy(c.begin()+off1, c.begin()+off2, rtn.begin());
0664     return rtn;
0665   }
0666 
0667   /// @brief Tail slice of the container elements cf. Python's [i:] syntax
0668   ///
0669   /// Single-index specialisation of @c slice(c, i, j)
0670   template <typename CONTAINER,
0671             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0672   inline CONTAINER slice(const CONTAINER& c, int i) {
0673     return slice(c, i, c.size());
0674   }
0675 
0676   /// @brief Head slice of the @a n first container elements
0677   ///
0678   /// Negative @a n means to take the head excluding the @a n -element tail
0679   template <typename CONTAINER,
0680             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0681   inline CONTAINER head(const CONTAINER& c, int n) {
0682     // if (n > c.size()) throw RangeError("Requested head longer than container");
0683     if (n < 0) n = std::max(0, (int)c.size()+n);
0684     n = std::min(n, (int)c.size());
0685     return slice(c, 0, n);
0686   }
0687 
0688   /// @brief Tail slice of the @a n last container elements
0689   ///
0690   /// Negative @a n means to take the tail from after the @a n th element
0691   template <typename CONTAINER,
0692             typename = std::enable_if_t< is_citerable_v<std::decay_t<CONTAINER>> >>
0693   inline CONTAINER tail(const CONTAINER& c, int n) {
0694     // if (n > c.size()) throw RangeError("Requested tail longer than container");
0695     if (n < 0) n = std::max(0, (int)c.size()+n);
0696     n = std::min(n, (int)c.size());
0697     return slice(c, c.size()-n);
0698   }
0699 
0700 
0701   using std::min;
0702   using std::max;
0703 
0704   /// Find the minimum value in the vector
0705   inline double min(const vector<double>& in, double errval=DBL_NAN) {
0706     const auto e = std::min_element(in.begin(), in.end());
0707     return e != in.end() ? *e : errval;
0708   }
0709 
0710   /// Find the maximum value in the vector
0711   inline double max(const vector<double>& in, double errval=DBL_NAN) {
0712     const auto e = std::max_element(in.begin(), in.end());
0713     return e != in.end() ? *e : errval;
0714   }
0715 
0716   /// Find the minimum and maximum values in the vector
0717   inline pair<double,double> minmax(const vector<double>& in, double errval=DBL_NAN) {
0718     const auto e = std::minmax_element(in.begin(), in.end());
0719     const double rtnmin = e.first != in.end() ? *e.first : errval;
0720     const double rtnmax = e.second != in.end() ? *e.first : errval;
0721     return std::make_pair(rtnmin, rtnmax);
0722   }
0723 
0724 
0725   /// Find the minimum value in the vector
0726   inline int min(const vector<int>& in, int errval=-1) {
0727     const auto e = std::min_element(in.begin(), in.end());
0728     return e != in.end() ? *e : errval;
0729   }
0730 
0731   /// Find the maximum value in the vector
0732   inline int max(const vector<int>& in, int errval=-1) {
0733     const auto e = std::max_element(in.begin(), in.end());
0734     return e != in.end() ? *e : errval;
0735   }
0736 
0737   /// Find the minimum and maximum values in the vector
0738   inline pair<int,int> minmax(const vector<int>& in, int errval=-1) {
0739     const auto e = std::minmax_element(in.begin(), in.end());
0740     const double rtnmin = e.first != in.end() ? *e.first : errval;
0741     const double rtnmax = e.second != in.end() ? *e.first : errval;
0742     return std::make_pair(rtnmin, rtnmax);
0743   }
0744 
0745   /// @}
0746 
0747 
0748   /// @defgroup contcombutils Container-combinatorics utils
0749   /// @{
0750 
0751   /// @brief Return the index from a vector which best matches fn(c[i]) to the target value
0752   ///
0753   /// A NaN return from the function will be counted as a no-match, as will combinations giving
0754   /// a value outside the given @a minval .. @a maxval range. A -1 index is returned in the
0755   /// case that no valid match is found.
0756   template <typename CONTAINER,
0757             typename FN = double(const typename std::decay_t<CONTAINER>::value_type&),
0758             typename = isCIterable<CONTAINER>>
0759   inline int closestMatchIndex(CONTAINER&& c, FN&& fn,
0760                    double target, double minval=-DBL_MAX, double maxval=DBL_MAX) {
0761     auto f = std::function(std::forward<FN>(fn));
0762     int ibest = -1;
0763     double best = DBL_NAN;
0764     for (size_t i = 0; i < c.size(); ++i) {
0765       const double val = f(c[i]);
0766       if (isnan(val)) continue;
0767       if (val < minval || val > maxval) continue;
0768       if (isnan(best) || fabs(val-target) < fabs(best-target)) {
0769         best = val;
0770         ibest = i;
0771       }
0772     }
0773     return ibest;
0774   }
0775 
0776   /// @brief Return the indices from two vectors which best match fn(c1[i], c2[j]) to the target value
0777   ///
0778   /// A NaN return from the function will be counted as a no-match, as will combinations giving
0779   /// a value outside the given @a minval .. @a maxval range. A {-1,-1} pair is returned in the
0780   /// case that no valid match is found.
0781   ///
0782   /// @note No attempt is made to avoid duplication: the input lists are assumed independent and
0783   ///   the user should check the sanity of the result, e.g. c1[i] and c2[j] are not the same object.
0784   template <typename CONTAINER1, typename CONTAINER2,
0785             typename FN = double(const typename std::decay_t<CONTAINER1>::value_type&,
0786                                  const typename std::decay_t<CONTAINER2>::value_type&),
0787             typename = isCIterable<CONTAINER1, CONTAINER2>>
0788   inline pair<int,int> closestMatchIndices(CONTAINER1&& c1, CONTAINER2&& c2, FN&& fn,
0789                        double target, double minval=-DBL_MAX, double maxval=DBL_MAX) {
0790     auto f = std::function(std::forward<FN>(fn));
0791     pair<int,int> ijbest{-1,-1};
0792     double best = DBL_NAN;
0793     for (size_t i = 0; i < c1.size(); ++i) {
0794       for (size_t j = 0; j < c2.size(); ++j) {
0795         const double val = f(c1[i], c2[j]);
0796         if (isnan(val)) continue;
0797         if (val < minval || val > maxval) continue;
0798         if (isnan(best) || fabs(val-target) < fabs(best-target)) {
0799           best = val;
0800           ijbest = {i,j};
0801         }
0802       }
0803     }
0804     return ijbest;
0805   }
0806 
0807   /// @brief Return the index from a vector which best matches fn(c[i], x) to the target value
0808   ///
0809   /// A NaN return from the function will be counted as a no-match, as will combinations giving
0810   /// a value outside the given @a minval .. @a maxval range. A -1 index is returned in the
0811   /// case that no valid match is found.
0812   ///
0813   /// @note No attempt is made to avoid duplication: the inputs are assumed independent and
0814   ///   the user should check the sanity of the result, e.g. c[i] and x are not the same object.
0815   template <typename CONTAINER, typename T,
0816             typename FN = double(const typename std::decay_t<CONTAINER>::value_type&, const std::decay_t<T>&),
0817             typename = isCIterable<CONTAINER>>
0818   inline int closestMatchIndex(CONTAINER&& c, const T& x, FN&& fn,
0819                    double target, double minval=-DBL_MAX, double maxval=DBL_MAX) {
0820     pair<int,int> ijbest = closestMatchIndices(std::forward<CONTAINER>(c), vector<T>{x}, std::forward<FN>(fn), target, minval, maxval);
0821     return ijbest.first;
0822   }
0823 
0824   /// @brief Return the index from a vector which best matches fn(x, c[j]) to the target value
0825   ///
0826   /// A NaN return from the function will be counted as a no-match, as will combinations giving
0827   /// a value outside the given @a minval .. @a maxval range. A -1 index is returned in the
0828   /// case that no valid match is found.
0829   ///
0830   /// @note No attempt is made to avoid duplication: the inputs are assumed independent and
0831   ///   the user should check the sanity of the result, e.g. c[i] and x are not the same object.
0832   template <typename CONTAINER, typename T, typename FN, typename = isCIterable<CONTAINER>>
0833   inline int closestMatchIndex(T&& x, CONTAINER&& c, FN&& fn,
0834                    double target, double minval=-DBL_MAX, double maxval=DBL_MAX) {
0835     return closestMatchIndex(std::forward<CONTAINER>(c), std::forward<T>(x), std::forward<FN>(fn), target, minval, maxval);
0836   }
0837 
0838   /// @}
0839 
0840 
0841   /// @brief Get a parameter from a named environment variable, with automatic type conversion
0842   ///
0843   /// @note Return @a fallback if the variable is not defined, otherwise convert its string to the template type
0844   /// @todo Should the param name have to be specific to an analysis? Can specialise as an Analysis member fn.
0845   /// @ingroup utils
0846   template <typename T>
0847   T getEnvParam(const std::string name, const T& fallback) {
0848     char* env = getenv(name.c_str());
0849     return env ? lexical_cast<T>(env) : fallback;
0850   }
0851 
0852 
0853 
0854   /// @brief Get the contents of a vector between two indices
0855   ///
0856   /// @ingroup utils
0857   template<class T>
0858   vector<T> slice(const vector<T>& v, int startidx, int endidx) {
0859 
0860     if (startidx < 0 || endidx <= startidx || endidx >= v.size())
0861       throw RangeError("Requested start or end indices incompatible with given vector");
0862 
0863     auto start = v.begin() + startidx;
0864     auto end = v.begin() + endidx;
0865 
0866     vector<T> output(endidx - startidx);
0867 
0868     copy(start, end, output.begin());
0869 
0870     return output;
0871   }
0872 
0873 
0874 }
0875 
0876 #endif