Fork me on GitHub

Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed
  • readers/DelphesPythia8.cpp

    r3e276d4 ra446115  
    2020#include <iostream>
    2121#include <sstream>
     22#include <string>
    2223
    2324#include <signal.h>
     
    3839#include "classes/DelphesClasses.h"
    3940#include "classes/DelphesFactory.h"
     41#include "classes/DelphesLHEFReader.h"
    4042
    4143#include "ExRootAnalysis/ExRootTreeWriter.h"
     
    149151  char appName[] = "DelphesPythia8";
    150152  stringstream message;
     153  FILE *inputFile = 0;
    151154  TFile *outputFile = 0;
    152155  TStopwatch readStopWatch, procStopWatch;
    153156  ExRootTreeWriter *treeWriter = 0;
    154157  ExRootTreeBranch *branchEvent = 0;
     158  ExRootTreeBranch *branchEventLHEF = 0, *branchWeightLHEF = 0;
    155159  ExRootConfReader *confReader = 0;
    156160  Delphes *modularDelphes = 0;
    157161  DelphesFactory *factory = 0;
    158162  TObjArray *stableParticleOutputArray = 0, *allParticleOutputArray = 0, *partonOutputArray = 0;
     163  TObjArray *stableParticleOutputArrayLHEF = 0, *allParticleOutputArrayLHEF = 0, *partonOutputArrayLHEF = 0;
     164  DelphesLHEFReader *reader = 0;
    159165  Long64_t eventCounter, errorCounter;
    160166  Long64_t numberOfEvents, timesAllowErrors;
     
    205211    partonOutputArray = modularDelphes->ExportArray("partons");
    206212
    207     modularDelphes->InitTask();
    208 
    209     // Initialize pythia
     213    // Initialize Pythia
    210214    pythia = new Pythia8::Pythia;
    211215
     
    221225    numberOfEvents = pythia->mode("Main:numberOfEvents");
    222226    timesAllowErrors = pythia->mode("Main:timesAllowErrors");
     227
     228    inputFile = fopen(pythia->word("Beams:LHEF").c_str(), "r");
     229    if(inputFile)
     230    {
     231      reader = new DelphesLHEFReader;
     232      reader->SetInputFile(inputFile);
     233
     234      branchEventLHEF = treeWriter->NewBranch("EventLHEF", LHEFEvent::Class());
     235      branchWeightLHEF = treeWriter->NewBranch("WeightLHEF", LHEFWeight::Class());
     236
     237      allParticleOutputArrayLHEF = modularDelphes->ExportArray("allParticlesLHEF");
     238      stableParticleOutputArrayLHEF = modularDelphes->ExportArray("stableParticlesLHEF");
     239      partonOutputArrayLHEF = modularDelphes->ExportArray("partonsLHEF");
     240    }
     241
     242    modularDelphes->InitTask();
    223243
    224244    pythia->init();
     
    234254    for(eventCounter = 0; eventCounter < numberOfEvents && !interrupted; ++eventCounter)
    235255    {
     256      while(reader && reader->ReadBlock(factory, allParticleOutputArrayLHEF,
     257        stableParticleOutputArrayLHEF, partonOutputArrayLHEF) && !reader->EventReady());
     258
    236259      if(!pythia->next())
    237260      {
    238261        // If failure because reached end of file then exit event loop
    239         if (pythia->info.atEndOfFile())
     262        if(pythia->info.atEndOfFile())
    240263        {
    241264          cerr << "Aborted since reached end of Les Houches Event File" << endl;
     
    244267
    245268        // First few failures write off as "acceptable" errors, then quit
    246         if (++errorCounter < timesAllowErrors) continue;
    247         cerr << "Event generation aborted prematurely, owing to error!" << endl;
    248         break;
     269        if(++errorCounter > timesAllowErrors)
     270        {
     271          cerr << "Event generation aborted prematurely, owing to error!" << endl;
     272          break;
     273        }
     274
     275        modularDelphes->Clear();
     276        reader->Clear();
     277        continue;
    249278      }
    250279
     
    258287      procStopWatch.Stop();
    259288
     289      if(reader)
     290      {
     291        reader->AnalyzeEvent(branchEventLHEF, eventCounter, &readStopWatch, &procStopWatch);
     292        reader->AnalyzeWeight(branchWeightLHEF);
     293      }
     294
    260295      treeWriter->Fill();
    261296
    262297      treeWriter->Clear();
    263298      modularDelphes->Clear();
     299      if(reader) reader->Clear();
    264300
    265301      readStopWatch.Start();
     
    277313    cout << "** Exiting..." << endl;
    278314
     315    delete reader;
    279316    delete pythia;
    280317    delete modularDelphes;
Note: See TracChangeset for help on using the changeset viewer.