[Rivet] [Rivet-svn] r2394 - in trunk: include/Rivet src/Core src/Tools

Frank Siegert frank.siegert at durham.ac.uk
Thu Apr 8 23:31:12 BST 2010


Hi guys,

Here some explanations about this commit.
Getting frustrated that a simple (diphoton, no jets) analysis in Rivet is 
taking longer than the actual event generation in Sherpa, I started 
profiling our AnalysisHandler::analyze() phase today. What I found is 
quite simple, you can find kcachegrind input files for all of these in 
http://www.ippp.dur.ac.uk/~fsiegert/tmp/profile/ if you want to look at it 
yourself. The starting point (callgrid.out.0initial) was: 
D0_2010_S8570965::analyze() takes up only 11% of the runtime of 
AnalysisHandler::analyze()!!!

1. The getParticleNamesMap() function in ParticleName.hh is ridiculously 
expensive and called even more ridiculously often. It alone was 
responsible for >50% of analyze()!
I have changed that into a static map which is initialized upon 
AnalysisHandler construction. If somebody has a better idea for that I'm 
happy to hear it.
-> callgrind.out.1particlenames

2. Now the GenEvent copying in Rivet::Event takes up 30% of analyze().
I have changed that such that the Event is only copied if we have to 
rotate it or change the units. Otherwise the external const reference is 
used.
-> callgrind.out.2geneventref

3. Now Log::getLog takes up 20% of analyze().
This is caused by the map searching, not any actual output. I haven't 
solved that yet, as that would probably involve some rewrite. Instead I've 
used a simple hack (prof_getlog.patch) to disable the log searching for 
further profiling. If anybody has a good idea here, please shout.
-> callgrind.out.3getlog

4. The output of beam pairs to TRACE level in Event::_geNormAlignment 
takes up 15% now. I have removed it.
-> callgrind.out.4nooutputinalign

5. Now beam projection and sqrtS() are 20% of analyze().
I don't think we need to continuously monitor the event beams, but one 
initial check is enough. So I have removed that from the analyze() phase.
-> callgrind.out.5onlyinitialbeamcheck

Bottom line: D0_2010_S8570965::analyze() actually takes 94% of the total 
analyze() time now.

I have committed the proposed fixes for 1, 2, 4 and 5. Would somebody be 
willing to look for a proper solution to 3?

So far for today ...
Frank


blackhole at projects.hepforge.org, Thursday 08 April 2010:
> Author: fsiegert
> Date: Thu Apr  8 23:26:49 2010
> New Revision: 2394
> 
> Log:
> Performance improvements that will e.g. bring the actual analysis from taking up 11% of the AnalysisHandler::analyze runtime up to 94%.
> 
> Deleted:
>    trunk/src/Tools/ParticleName.cc
> Modified:
>    trunk/include/Rivet/AnalysisHandler.hh
>    trunk/include/Rivet/Event.hh
>    trunk/include/Rivet/ParticleName.hh
>    trunk/src/Core/AnalysisHandler.cc
>    trunk/src/Core/Event.cc
> 
> Modified: trunk/include/Rivet/AnalysisHandler.hh
> ==============================================================================
> --- trunk/include/Rivet/AnalysisHandler.hh	Thu Apr  8 00:17:33 2010	(r2393)
> +++ trunk/include/Rivet/AnalysisHandler.hh	Thu Apr  8 23:26:49 2010	(r2394)
> @@ -269,6 +269,8 @@
>      /// fact, it should not even be implemented.
>      AnalysisHandler(const AnalysisHandler&);
>  
> +    static void initializeParticleNames();
> +
>    };
>  
>  
> 
> Modified: trunk/include/Rivet/Event.hh
> ==============================================================================
> --- trunk/include/Rivet/Event.hh	Thu Apr  8 00:17:33 2010	(r2393)
> +++ trunk/include/Rivet/Event.hh	Thu Apr  8 23:26:49 2010	(r2394)
> @@ -30,8 +30,8 @@
>      /// The copy constructor.
>      Event(const Event& e);
>  
> -    /// Copy assignment operator
> -    Event& operator=(const Event& e);
> +    /// The destructor
> +    ~Event();
>  
>      //@}
>  
> @@ -39,9 +39,7 @@
>    public:
>  
>      /// Return the generated event obtained from an external event generator.
> -    const GenEvent& genEvent() const {
> -      return _genEvent;
> -    }
> +    const GenEvent& genEvent() const;
>  
>      /// The weight associated with the event.
>      double weight() const {
> @@ -82,10 +80,14 @@
>  
>  
>    private:
> +    void _geNormAlignment();
> +
>  
>      /// @brief The generated event, obtained from an external generator.
>      /// Note that it is only externally accessible via a const reference.
> -    GenEvent _genEvent;
> +    GenEvent const& _genEvent;
> +
> +    GenEvent* _modGenEvent;
>  
>      /// The set of Projection objects applied so far.
>      mutable std::set<ConstProjectionPtr> _projections;
> 
> Modified: trunk/include/Rivet/ParticleName.hh
> ==============================================================================
> --- trunk/include/Rivet/ParticleName.hh	Thu Apr  8 00:17:33 2010	(r2393)
> +++ trunk/include/Rivet/ParticleName.hh	Thu Apr  8 23:26:49 2010	(r2394)
> @@ -87,48 +87,13 @@
>    typedef std::map<std::string, PdgId> ParticleNameMapR;
>  
>  
> -  /// Function which returns a map from beam particle enums to the corresponding name strings.
> -  inline ParticleNameMap getParticleNamesMap() {
> -    ParticleNameMap bpmap;
> -    bpmap[ELECTRON] = "ELECTRON";
> -    bpmap[POSITRON] = "POSITRON";
> -    bpmap[PROTON] = "PROTON";
> -    bpmap[ANTIPROTON] = "ANTIPROTON";
> -    bpmap[PHOTON] = "PHOTON";
> -    bpmap[NEUTRON] = "NEUTRON";
> -    bpmap[ANTINEUTRON] = "ANTINEUTRON";
> -    bpmap[MUON] = "MUON";
> -    bpmap[ANTIMUON] = "ANTIMUON";
> -    bpmap[NU_E] = "NU_E";
> -    bpmap[NU_EBAR] = "NU_EBAR";
> -    bpmap[NU_MU] = "NU_MU";
> -    bpmap[NU_MUBAR] = "NU_MUBAR";
> -    bpmap[NU_TAU] = "NU_TAU";
> -    bpmap[NU_TAUBAR] = "NU_TAUBAR";
> -    bpmap[PIPLUS] = "PIPLUS";
> -    bpmap[PIMINUS] = "PIMINUS";
> -    bpmap[TAU] = "TAU";
> -    bpmap[WPLUSBOSON] = "WPLUSBOSON";
> -    bpmap[WMINUSBOSON] = "WMINUSBOSON";
> -    bpmap[ZBOSON] = "ZBOSON";
> -    bpmap[HIGGS] = "HIGGS";
> -    bpmap[ANTITAU] = "ANTITAU";
> -    bpmap[PHOTOELECTRON] = "PHOTOELECTRON";
> -    bpmap[PHOTOPOSITRON] = "PHOTOPOSITRON";
> -    bpmap[PHOTOMUON] = "PHOTOMUON";
> -    bpmap[PHOTOANTIMUON] = "PHOTOANTIMUON";
> -    bpmap[PHOTOTAU] = "PHOTOTAU";
> -    bpmap[PHOTOANTITAU] = "PHOTOANTITAU";
> -    bpmap[ANY] = "*";
> -    return bpmap;
> -  }
> +  static ParticleNameMap s_pnames;
>  
>  
>    /// Function which returns a map from beam particle name strings to the corresponding enums.
>    inline ParticleNameMapR getParticleNamesRMap() {
> -    ParticleNameMap bpmap = getParticleNamesMap();
>      ParticleNameMapR bpmapr;
> -    for (ParticleNameMap::const_iterator bp = bpmap.begin(); bp != bpmap.end(); ++bp) {
> +    for (ParticleNameMap::const_iterator bp = s_pnames.begin(); bp != s_pnames.end(); ++bp) {
>        bpmapr[bp->second] = bp->first;
>      }
>      return bpmapr;
> @@ -143,8 +108,7 @@
>    /// the ParticleName enum.
>    inline ParticleNameList getParticleNameEnums() {
>      ParticleNameList names;
> -    ParticleNameMap bpmap = getParticleNamesMap();
> -    for (ParticleNameMap::const_iterator bp = bpmap.begin(); bp != bpmap.end(); ++bp) {
> +    for (ParticleNameMap::const_iterator bp = s_pnames.begin(); bp != s_pnames.end(); ++bp) {
>        names.push_back(bp->first);
>      }
>      return names;
> @@ -161,8 +125,7 @@
>    /// Function which returns a vector of all the beam particle name strings.
>    inline std::vector<std::string> getParticleNames() {
>      vector<string> names;
> -    ParticleNameMap bpmap = getParticleNamesMap();
> -    for (ParticleNameMap::const_iterator bp = bpmap.begin(); bp != bpmap.end(); ++bp) {
> +    for (ParticleNameMap::const_iterator bp = s_pnames.begin(); bp != s_pnames.end(); ++bp) {
>        names.push_back(bp->second);
>      }
>      return names;
> @@ -171,14 +134,14 @@
>  
>    /// Print a ParticleName as a string.
>    inline std::string toString(const ParticleName& p) {
> -    return getParticleNamesMap()[p];
> +    return s_pnames[p];
>    }
>  
>  
>    /// Print a PdgId as a named string.
>    inline std::string toParticleName(PdgId p) {
> -    if (getParticleNamesMap().find(p) != getParticleNamesMap().end()) {
> -      return getParticleNamesMap()[p];
> +    if (s_pnames.find(p) != s_pnames.end()) {
> +      return s_pnames[p];
>      }
>      ostringstream ss;
>      ss << p;
> 
> Modified: trunk/src/Core/AnalysisHandler.cc
> ==============================================================================
> --- trunk/src/Core/AnalysisHandler.cc	Thu Apr  8 00:17:33 2010	(r2393)
> +++ trunk/src/Core/AnalysisHandler.cc	Thu Apr  8 23:26:49 2010	(r2394)
> @@ -1,6 +1,7 @@
>  // -*- C++ -*-
>  #include "Rivet/Rivet.hh"
>  #include "Rivet/Tools/Logging.hh"
> +#include "Rivet/ParticleName.hh"
>  #include "Rivet/AnalysisHandler.hh"
>  #include "Rivet/RivetAIDA.hh"
>  #include "Rivet/Analysis.hh"
> @@ -20,6 +21,7 @@
>    {
>      _theAnalysisFactory = createAnalysisFactory();
>      _setupFactories(basefilename, storetype);
> +    initializeParticleNames();
>    }
>  
>  
> @@ -30,6 +32,7 @@
>        _theAnalysisFactory(&afac) 
>    {
>      _setupFactories(basefilename, storetype);
> +    initializeParticleNames();
>    }
>  
>  
> @@ -86,15 +89,15 @@
>      }
>      // Proceed with event analysis
>      assert(_initialised);
> -    // Ensure that beam details match those from first event
> -    const BeamPair beams = Rivet::beamIds(ge);
> -    const double sqrts = Rivet::sqrtS(ge);
> -    if (!compatible(beams, _beams) || !fuzzyEquals(sqrts, sqrtS())) {
> -      getLog() << Log::ERROR << "Event beams mismatch: "
> -               << beams << " @ " << sqrts/GeV << " GeV" << " vs. first beams "
> -               << this->beams() << " @ " << this->sqrtS()/GeV << " GeV" << endl;
> -      exit(1);
> -    }
> +//     // Ensure that beam details match those from first event @todo This is too expensive
> +//     const BeamPair beams = Rivet::beamIds(ge);
> +//     const double sqrts = Rivet::sqrtS(ge);
> +//     if (!compatible(beams, _beams) || !fuzzyEquals(sqrts, sqrtS())) {
> +//       getLog() << Log::ERROR << "Event beams mismatch: "
> +//                << beams << " @ " << sqrts/GeV << " GeV" << " vs. first beams "
> +//                << this->beams() << " @ " << this->sqrtS()/GeV << " GeV" << endl;
> +//       exit(1);
> +//     }
>  
>      
>      Event event(ge);
> @@ -365,4 +368,39 @@
>    }
>  
>  
> +  void AnalysisHandler::initializeParticleNames() {
> +    Rivet::s_pnames[ELECTRON] = "ELECTRON";
> +    Rivet::s_pnames[POSITRON] = "POSITRON";
> +    Rivet::s_pnames[PROTON] = "PROTON";
> +    Rivet::s_pnames[ANTIPROTON] = "ANTIPROTON";
> +    Rivet::s_pnames[PHOTON] = "PHOTON";
> +    Rivet::s_pnames[NEUTRON] = "NEUTRON";
> +    Rivet::s_pnames[ANTINEUTRON] = "ANTINEUTRON";
> +    Rivet::s_pnames[MUON] = "MUON";
> +    Rivet::s_pnames[ANTIMUON] = "ANTIMUON";
> +    Rivet::s_pnames[NU_E] = "NU_E";
> +    Rivet::s_pnames[NU_EBAR] = "NU_EBAR";
> +    Rivet::s_pnames[NU_MU] = "NU_MU";
> +    Rivet::s_pnames[NU_MUBAR] = "NU_MUBAR";
> +    Rivet::s_pnames[NU_TAU] = "NU_TAU";
> +    Rivet::s_pnames[NU_TAUBAR] = "NU_TAUBAR";
> +    Rivet::s_pnames[PIPLUS] = "PIPLUS";
> +    Rivet::s_pnames[PIMINUS] = "PIMINUS";
> +    Rivet::s_pnames[TAU] = "TAU";
> +    Rivet::s_pnames[WPLUSBOSON] = "WPLUSBOSON";
> +    Rivet::s_pnames[WMINUSBOSON] = "WMINUSBOSON";
> +    Rivet::s_pnames[ZBOSON] = "ZBOSON";
> +    Rivet::s_pnames[HIGGS] = "HIGGS";
> +    Rivet::s_pnames[ANTITAU] = "ANTITAU";
> +    Rivet::s_pnames[PHOTOELECTRON] = "PHOTOELECTRON";
> +    Rivet::s_pnames[PHOTOPOSITRON] = "PHOTOPOSITRON";
> +    Rivet::s_pnames[PHOTOMUON] = "PHOTOMUON";
> +    Rivet::s_pnames[PHOTOANTIMUON] = "PHOTOANTIMUON";
> +    Rivet::s_pnames[PHOTOTAU] = "PHOTOTAU";
> +    Rivet::s_pnames[PHOTOANTITAU] = "PHOTOANTITAU";
> +    Rivet::s_pnames[ANY] = "*";
> +  }
> +
> +
> +
>  }
> 
> Modified: trunk/src/Core/Event.cc
> ==============================================================================
> --- trunk/src/Core/Event.cc	Thu Apr  8 00:17:33 2010	(r2393)
> +++ trunk/src/Core/Event.cc	Thu Apr  8 23:26:49 2010	(r2394)
> @@ -7,15 +7,6 @@
>  namespace Rivet {
>  
>  
> -  // Convert the GenEvent to use conventional units (GeV, mm)
> -  void _geNormUnits(GenEvent& ge) {
> -    // Specify units if supported
> -    #ifdef HEPMC_HAS_UNITS
> -    ge.use_units(HepMC::Units::GEV, HepMC::Units::MM);
> -    #endif
> -  }
> -
> -
>    void _geRot180x(GenEvent& ge) {
>      for (HepMC::GenEvent::particle_iterator ip = ge.particles_begin(); ip != ge.particles_end(); ++ip) {
>        const HepMC::FourVector& mom = (*ip)->momentum();
> @@ -32,12 +23,12 @@
>    // (proton or electron on +ve z-axis?)
>    // For example, FHerwig only produces DIS events in the
>    // unconventional orientation and has to be corrected
> -  void _geNormAlignment(GenEvent& ge) {
> -    if (!ge.valid_beam_particles()) return;
> +  void Event::_geNormAlignment() {
> +    if (!_genEvent.valid_beam_particles()) return;
>      typedef pair<HepMC::GenParticle*, HepMC::GenParticle*> GPPair;
> -    GPPair bps = ge.beam_particles();
> +    GPPair bps = _genEvent.beam_particles();
>      const BeamPair beamids = make_pdgid_pair(bps.first->pdg_id(), bps.second->pdg_id());
> -    Log::getLog("Rivet.Event") << Log::TRACE << "Beam IDs: " << beamids << endl;
> +    //Log::getLog("Rivet.Event") << Log::TRACE << "Beam IDs: " << beamids << endl;
>      const HepMC::GenParticle* plusgp = 0;
>      bool rot = false;
>  
> @@ -66,13 +57,14 @@
>                                     << bps.first->pdg_id() << "@pz=" << bps.first->momentum().pz()/GeV << ", "
>                                     << bps.second->pdg_id() << "@pz=" << bps.second->momentum().pz()/GeV << endl;
>        }
> -      _geRot180x(ge);
> +      if (!_modGenEvent) _modGenEvent = new GenEvent(_genEvent);
> +      _geRot180x(*_modGenEvent);
>      }
>    }
>  
>  
>    Event::Event(const GenEvent& ge)
> -    : _genEvent(ge), _weight(1.0)
> +    : _genEvent(ge), _modGenEvent(NULL), _weight(1.0)
>    {
>      // Set the weight if there is one, otherwise default to 1.0
>      if (!_genEvent.weights().empty()) {
> @@ -80,10 +72,16 @@
>      }
>  
>      // Use Rivet's preferred units if possible
> -    _geNormUnits(_genEvent);
> +    #ifdef HEPMC_HAS_UNITS
> +    if (_genEvent.momentum_unit()!=HepMC::Units::GEV ||
> +        _genEvent.length_unit()!=HepMC::Units::MM) {
> +      if (!_modGenEvent) _modGenEvent = new GenEvent(ge);
> +      _modGenEvent->use_units(HepMC::Units::GEV, HepMC::Units::MM);
> +    }
> +    #endif
>   
>      // Use the conventional alignment
> -    _geNormAlignment(_genEvent);
> +    _geNormAlignment();
>  
>      // Debug printout to check that copying/magling has worked
>      //_genEvent.print();
> @@ -91,17 +89,20 @@
>  
>  
>    Event::Event(const Event& e)
> -    : _genEvent(e._genEvent),
> +    : _genEvent(e._genEvent), _modGenEvent(e._modGenEvent),
>        _weight(e._weight)
>    {
>      //
>    }
>  
>  
> -  Event& Event::operator=(const Event& e) {
> -    _genEvent = e._genEvent;
> -    _weight = e._weight;
> -    return *this;
> +  Event::~Event() {
> +    if (_modGenEvent) delete _modGenEvent;
> +  }
> +
> +  const GenEvent& Event::genEvent() const {
> +    if (_modGenEvent) return *_modGenEvent;
> +    return _genEvent;
>    }
>  
>  
> _______________________________________________
> Rivet-svn mailing list
> Rivet-svn at projects.hepforge.org
> http://www.hepforge.org/lists/listinfo/rivet-svn
> 


More information about the Rivet mailing list