ThePEG/Utilities/Selector.tcc

0001 // -*- C++ -*-
0002 //
0003 // Selector.tcc is a part of ThePEG - Toolkit for HEP Event Generation
0004 // Copyright (C) 1999-2019 Leif Lonnblad
0005 //
0006 // ThePEG is licenced under version 3 of the GPL, see COPYING for details.
0007 // Please respect the MCnet academic guidelines, see GUIDELINES for details.
0008 //
0009
0010 namespace ThePEG {
0011
0012 template <typename T, typename WeightType>
0013 WeightType Selector<T,WeightType>::erase(const T & t) {
0014   Selector<T,WeightType> newSelector;
0015   WeightType oldsum = WeightType();
0016   for ( iterator it = theMap.begin();
0017         it != theMap.end(); ++it ) {
0018     WeightType r = it->first - oldsum;
0019     oldsum = it->first;
0020     if ( it->second != t ) newSelector.insert(r, it->second);
0021   }
0022   theMap.swap(newSelector.theMap);
0023   return theSum = newSelector.theSum;
0024 }
0025
0026 template <typename T, typename WeightType>
0027 const T & Selector<T,WeightType>::
0028 select(double rnd, double * remainder) const {
0029   if ( rnd <= 0 )
0030     throw range_error("Random number out of range in Selector::select.");
0031   const_iterator it = theMap.upper_bound(rnd*theSum);
0032   if ( it == theMap.end() )
0033     throw range_error("Empty Selector, or random number out of range "
0034                       "in Selector::select");
0035   if ( remainder ) {
0036     if ( it == theMap.begin() )
0037       *remainder = rnd*theSum/(it->first);
0038     else {
0039       const_iterator prit = it;
0040       --prit;
0041       *remainder = (rnd*theSum - prit->first)/(it->first - prit->first);
0042     }
0043   }
0044   return it->second;
0045 }
0046
0047 template <typename T, typename WeightType>
0048 T & Selector<T,WeightType>::
0049 select(double rnd, double * remainder) {
0050   if ( rnd <= 0 )
0051     throw range_error("Random number out of range in Selector::select.");
0052   iterator it = theMap.upper_bound(rnd*theSum);
0053   if ( it == theMap.end() )
0054     throw range_error("Empty Selector, or random number out of range "
0055                       "in Selector::select");
0056   if ( remainder ) {
0057     if ( it == theMap.begin() )
0058       *remainder = rnd*theSum/(it->first);
0059     else {
0060       const_iterator prit = it;
0061       --prit;
0062       *remainder = (rnd*theSum - prit->first)/(it->first - prit->first);
0063     }
0064   }
0065   return it->second;
0066 }
0067
0068 template <typename T, typename WeightType>
0069 template <typename OStream>
0070 void Selector<T,WeightType>::output(OStream & os, DimensionT) const {
0071   os << ounit(theSum,WeightType::baseunit()) << size();
0072   for ( const_iterator it = theMap.begin(); it != theMap.end(); ++it )
0073     os << ounit(it->first,WeightType::baseunit()) << it->second;
0074 }
0075
0076 template <typename T, typename WeightType>
0077 template <typename IStream>
0078 void Selector<T,WeightType>::input(IStream & is, DimensionT) {
0079   typedef typename MapType::value_type value_type;
0080   clear();
0081   T t;
0082   WeightType weightsum;
0083   long n;
0084   is >> iunit(theSum,WeightType::baseunit()) >> n;
0085   while ( n-- ) {
0086     is >> iunit(weightsum,WeightType::baseunit()) >> t;
0087     theMap.insert(theMap.end(), value_type(weightsum, t));
0088   }
0089 }
0090
0091 template <typename T, typename WeightType>
0092 template <typename OStream>
0093 void Selector<T,WeightType>::output(OStream & os, StandardT) const {
0094   os << theSum << size();
0095   for ( const_iterator it = theMap.begin(); it != theMap.end(); ++it )
0096     os << it->first << it->second;
0097 }
0098
0099 template <typename T, typename WeightType>
0100 template <typename IStream>
0101 void Selector<T,WeightType>::input(IStream & is, StandardT) {
0102   typedef typename MapType::value_type value_type;
0103   clear();
0104   T t;
0105   WeightType weightsum;
0106   long n;
0107   is >> theSum >> n;
0108   while ( n-- ) {
0109     is >> weightsum >> t;
0110     theMap.insert(theMap.end(), value_type(weightsum, t));
0111   }
0112 }
0113
0114 }