Weights are lost form lhe to hepmc

Dear Manuel,

as far as the HepMC file is concerned, the multiple weights
are not propagated to the HepMC file with our shower interface.
Weights are propagated only using the 'HwU' interface, namely
by specifying in shower_card.dat the following line

ANALYSE = HwU.o your_analysis.o

where 'your_analisis.f' is a fortran analysis file that has to be
located in MCatNLO/PY8Analyzer. There are many examples
in that directory, including some analyses for the process you
generated.

As for the weights not appearing in the summary, what do you
mean exactly? Your summary.txt file does not contain a statement
like
 --------------------------------------------------------------
      Scale variation (computed from LHE events):
          Dynamical_scale_choice -1 (envelope of 9 values):
              1.002e+05 pb +6.0% -9.8%
   --------------------------------------------------------------
?

I tried with the default run_card and it does for me.
Please, let me know.

Best.
Paolo

Dear Paolo,

Thank you very much for your answer. It has been very useful.

I have set my run_card to default and I get the right summary file, I must have made some mistake earlier.

As for propagating the weights, I am afraid some new issues have arosen. When I try using ANALYSE = HwU.o py8an_HwU_pp_V.o, I get the error below. Would you know what's causing it?

As an aside, is there an equivalent to HwU for mcatnlo_pyan_pp_V_hepmc.f? And do you happen to know how to include an HwU in a Delphes analysis?

Thank you again for your help

Best regards,

Manuel

g++ -O -I/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/HEPTools/pythia8/include \
  -I/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/wnlo/MCatNLO/include Pythia82.cc -o Pythia8.exe HwU.o py8an_HwU_pp_V.o \
  -L/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/HEPTools/pythia8/lib -lpythia8 \
  -lHepMCfio -lgfortran -I/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/HEPTools/hepmc/include -L/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/HEPTools/hepmc/lib -lHepMC -L../lib -L/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/HEPTools/pythia8//lib -L/home/manuel/Softwarenoroot/MG5_aMC_v2.6.1/MG5_aMC_v2_6_1/HEPTools/boost/lib -L/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/HEPTools/zlib/lib -lpythia8 -lboost_iostreams -lz -ldl -lstdc++ \

In file included from /home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/wnlo/MCatNLO/include/LHEFRead.h:1:0,
                 from Pythia82.cc:13:
/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/wnlo/MCatNLO/include/LHEF.h: In member function ‘bool LHEF::Reader::getline()’:
/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/wnlo/MCatNLO/include/LHEF.h:1974:46: error: cannot convert ‘std::basic_istream<char>’ to ‘bool’ in return
     return ( std::getline(file, currentLine) );
                                              ^
Makefile:49: recipe for target 'Pythia82' failed
make: *** [Pythia82] Error 1
Pythia8 compilation did not succeed, exiting

Dear Manuel,
I tried with that analysis, and in my setup it is working.

A possible reason is that I’ve installed Pythia (v8223)
manually, i.e. not from MadGraph. Could you try and do
the same?

Recall to install HepMC before installing Pythia, and then
specify its location while installing Pythia through the
option --with-hepmc2.
Once you’ve done that, you should also update the
input/mg5_configuration.txt, and the
<your_process>/Cards/amcatnlo_configuration.txt
files to point to the new Pythia location.

As for your latter question, ‘mcatnlo_pyan_pp_V_hepmc.f’
may be a misleading name, since it is not actually related
to an HepMC output, but only to the fact that it uses some
HepMC internal information to manage momenta and particle
identities. That analysis (and the ones similar to that)
outputs a histogram file in top-drawer format, which of
course can be used, but being a very old tool, is not
something one commonly uses these days.

As for Delphes, I don’t this there is an easy way to do that
with HwU. Maybe with other tools like MadAnalysis5 it is
possible. I’ll anyway try and see if the original issue,
namely including multiple weights in our interface’s HepMC
files, can be easily fixed.

Cheers.
Paolo

> On 18 Feb 2019, at 14:12, Manuel Szewc <email address hidden> wrote:
>
> Question #678584 on MadGraph5_aMC@NLO changed:
> https://answers.launchpad.net/mg5amcnlo/+question/678584
>
> Status: Answered => Open
>
> Manuel Szewc is still having a problem:
> Dear Paolo,
>
> Thank you very much for your answer. It has been very useful.
>
> I have set my run_card to default and I get the right summary file, I
> must have made some mistake earlier.
>
> As for propagating the weights, I am afraid some new issues have arosen.
> When I try using ANALYSE = HwU.o py8an_HwU_pp_V.o, I get the error
> below. Would you know what's causing it?
>
> As an aside, is there an equivalent to HwU for
> mcatnlo_pyan_pp_V_hepmc.f? And do you happen to know how to include an
> HwU in a Delphes analysis?
>
> Thank you again for your help
>
> Best regards,
>
> Manuel
>
> g++ -O -I/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/HEPTools/pythia8/include \
> -I/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/wnlo/MCatNLO/include Pythia82.cc -o Pythia8.exe HwU.o py8an_HwU_pp_V.o \
> -L/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/HEPTools/pythia8/lib -lpythia8 \
> -lHepMCfio -lgfortran -I/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/HEPTools/hepmc/include -L/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/HEPTools/hepmc/lib -lHepMC -L../lib -L/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/HEPTools/pythia8//lib -L/home/manuel/Softwarenoroot/MG5_aMC_v2.6.1/MG5_aMC_v2_6_1/HEPTools/boost/lib -L/home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/HEPTools/zlib/lib -lpythia8 -lboost_iostreams -lz -ldl -lstdc++ \
>
> In file included from /home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/wnlo/MCatNLO/include/LHEFRead.h:1:0,
> from Pythia82.cc:13:
> /home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/wnlo/MCatNLO/include/LHEF.h: In member function ‘bool LHEF::Reader::getline()’:
> /home/manuel/Softwarenoroot/MG5_aMC_v2_6_5/wnlo/MCatNLO/include/LHEF.h:1974:46: error: cannot convert ‘std::basic_istream<char>’ to ‘bool’ in return
> return ( std::getline(file, currentLine) );
> ^
> Makefile:49: recipe for target 'Pythia82' failed
> make: *** [Pythia82] Error 1
> Pythia8 compilation did not succeed, exiting
>
> --
> You received this question notification because you are subscribed to
> the question.

Dear Paolo

I am very sorry for the late response.

I tried installing manually Pythia8 with hepmc but it returned the same problem as before. Could there be a problem with hepmc2? When installing it I noticed that the testHepMC.sh test failed.

I will start looking at MadAnalysis then!

Best regards

Manuel

Dear Paolo

Thank you very much for all your help! I will try with different Pythia installations.

Best regards,

Manuel

Thanks Paolo Torrielli, that solved my question.

Hi Paolo,
concering the question of weight propagation from LHE to HEPMC in NLO runs: Can the additional weights from the LHE file be transferred to the HEPMC file 1:1 or do they have to be changed because of PS matching/merging?
Thanks,
Peter

Dear Peter,
yes, the weights can be transferred 1:1 to the HEPMC file.
Cheers.
Paolo

> On 25 Oct 2019, at 14:17, Peter Galler <email address hidden> wrote:
>
> Question #678584 on MadGraph5_aMC@NLO changed:
> https://answers.launchpad.net/mg5amcnlo/+question/678584
>
> Peter Galler posted a new comment:
> Hi Paolo,
> concering the question of weight propagation from LHE to HEPMC in NLO runs: Can the additional weights from the LHE file be transferred to the HEPMC file 1:1 or do they have to be changed because of PS matching/merging?
> Thanks,
> Peter
>
> --
> You received this question notification because you are subscribed to
> the question.

Dear Paolo,
great! Then until this is implemented in MadGraph the following procedure seems to work.
I copied the the relevant code piece ("HEPMC2HACK") from
MadGraph/HEPTools/MG5aMC_PY8_interface/MG5aMC_PY8_interface.cc
(which is the Pythia driver for LO runs where multiple weights in HEPMC work)
into the file
MadGraph/Template/NLO/MCatNLO/srcPythia8/Pythia82_hep.cc
and this seems to do the trick.
I include the listing for the entire modified file Pythia82_hep.cc below.
If one doesn't want to mess up the MG installation one can also modify
the same file in the local run directory of the process under study, i.e.
<your_favourite_process>/MCatNLO/srcPythia8/Pythia82_hep.cc

*********** Pythia82_hep.cc (modified) *********************************

// Driver for Pythia 8. Reads an input file dynamically created on
// the basis of the inputs specified in MCatNLO_MadFKS_PY8.Script
#include "Pythia8/Pythia.h"
#include "Pythia8Plugins/HepMC2.h"
#include "Pythia8Plugins/aMCatNLOHooks.h"
#include "Pythia8Plugins/CombineMatchingInput.h"
#include "HepMC/GenEvent.h"
#include "HepMC/IO_GenEvent.h"

using namespace Pythia8;

string convertFromAMCATNLO( string input) {
  string output="";
  // Count number of blanks in weight name.
  int appearances = 0;
  for(int n = input.find(" ", 0); n != int(string::npos);
          n = input.find(" ", n)) {
    appearances++;
    n++;
  }
  // Cut string by position of blanks.
  vector <string> pieces;
  for(int i =0; i < appearances;++i) {
    int n = input.find(" ", 0);
    pieces.push_back(input.substr(0,n));
    input = input.substr(n+1,input.size());
  }
  // Now pieces contains the details of the weight name.
  for(int i = 4; i < int(pieces.size()); ++i) {
    // Make upper case.
    transform(pieces[i].begin(),pieces[i].end(),pieces[i].begin(), ::toupper);
    // Add piece to weight name, with "_" as separator.
    output+=pieces[i] + "_";
  }
  return output;
}

int main() {
  Pythia pythia;

  string inputname="Pythia8.cmd",outputname="Pythia8.hep";

  pythia.readFile(inputname.c_str());

  //Create UserHooks pointer for the FxFX matching. Stop if it failed. Pass pointer to Pythia.
  CombineMatchingInput* combined = NULL;
  UserHooks* matching = NULL;

  int nAbort=10;
  int nPrintLHA=1;
  int iAbort=0;
  int iPrintLHA=0;
  int iEventtot=pythia.mode("Main:numberOfEvents");
  int iEventshower=pythia.mode("Main:spareMode1");
  string evt_norm=pythia.word("Main:spareWord1");

  //FxFx merging
  bool isFxFx=pythia.flag("JetMatching:doFxFx");
  if (isFxFx) {
    matching = combined->getHook(pythia);
    if (!matching) {
      std::cout << " Failed to initialise jet matching structures.\n"
                << " Program stopped.";
      return 1;
    }
    pythia.setUserHooksPtr(matching);
    int nJmax=pythia.mode("JetMatching:nJetMax");
    double Qcut=pythia.parm("JetMatching:qCut");
    double PTcut=pythia.parm("JetMatching:qCutME");
    if (Qcut <= PTcut || Qcut <= 0.) {
      std::cout << " \n";
      std::cout << "Merging scale (shower_card.dat) smaller than pTcut (run_card.dat)"
                << Qcut << " " << PTcut << "\n";
      return 0;
    }
  }

  pythia.init();

  HepMC::Pythia8ToHepMC ToHepMC;
  HepMC::IO_GenEvent ascii_io(outputname.c_str(), std::ios::out);
  // Do not store cross section information, as this will be done manually.
  ToHepMC.set_store_pdf(false);
  ToHepMC.set_store_proc(false);
  ToHepMC.set_store_xsec(false);

  // Cross section an error.
  double sigmaTotal = 0.;
  double errorTotal = 0.;

  for (int iEvent = 0; ; ++iEvent) {
    if (!pythia.next()) {
      if (++iAbort < nAbort) continue;
      break;
    }
    // the number of events read by Pythia so far
    int nSelected=pythia.info.nSelected();

    if (nSelected > iEventshower) break;
    if (pythia.info.isLHA() && iPrintLHA < nPrintLHA) {
      pythia.LHAeventList();
      pythia.info.list();
      pythia.process.list();
      pythia.event.list();
      ++iPrintLHA;
    }

    HepMC::GenEvent* hepmcevt = new HepMC::GenEvent();
    double evtweight = pythia.info.weight();
    double normhepmc;
    // ALWAYS NORMALISE HEPMC WEIGHTS TO SUM TO THE CROSS SECTION
    if (evt_norm != "sum") {
      normhepmc = 1. / double(iEventshower);
    } else {
      normhepmc = double(iEventtot) / double(iEventshower);
    }
    sigmaTotal += evtweight*normhepmc;
    //hepmcevt->weights().push_back(evtweight*normhepmc);

    // Loop through weightgroups.
    for ( std::map<std::string,LHAweightgroup>::const_iterator
      it_wg = pythia.info.weightgroups->begin();
      it_wg != pythia.info.weightgroups->end(); ++it_wg ) {
      // Loop through weights in the weightgroup.
      for ( std::map<std::string,LHAweight>::const_iterator
        it_w = it_wg->second.weights.begin();
        it_w != it_wg->second.weights.end(); ++it_w ) {
        // Get value of weight indexed by the present key.
        std::map<std::string,double>::const_iterator it_value
          = pythia.info.weights_detailed->find(it_w->first);
        double w = it_value->second; //*merging_weight*vetoWeights[iCut];
        // Get name.
        string name;
        for ( std::map<string,string>::const_iterator
          it_att = it_w->second.attributes.begin();
          it_att != it_w->second.attributes.end(); ++it_att)
          name += it_att->first + "=" + it_att->second + "_";

        // In aMC@NLO, the weight names should be extracted from the
        // weight tag contents, not its attributes.
        if (it_w->second.attributes.size()==0)
          name = convertFromAMCATNLO(it_w->second.contents);

        //if (pythia.info.lhaStrategy() == 3)
        // w *= norm_event_wgt / pythia.info.eventWeightLHEF;
        //// Weighted events.
        //else w *= 1 / double(1e9*nEvent);

        //w *=hepmcWeightRescaling;

        // Add event weight.
        hepmcevt->weights().push_back(w*normhepmc, name);

      } // Done looping through weights.
    } // Done looping through weight groups.

    ToHepMC.fill_next_event( pythia, hepmcevt );
    // Add the weight of the current event to the cross section.
    // Report cross section to hepmc
    HepMC::GenCrossSection xsec;
    xsec.set_cross_section( sigmaTotal, pythia.info.sigmaErr() );
    hepmcevt->set_cross_section( xsec );
    // Write the HepMC event to file. Done with it.
    ascii_io << hepmcevt;
    delete hepmcevt;
  }

  pythia.stat();
  if (isFxFx){
    std::cout << " \n";
    std::cout << "*********************************************************************** \n";
    std::cout << "*********************************************************************** \n";
    std::cout << "Cross section, including FxFx merging is: "
              << sigmaTotal << "\n";
    std::cout << "*********************************************************************** \n";
    std::cout << "*********************************************************************** \n";
  }

  return 0;
}

*********** End of Pythia82_hep.cc (modified) *********************************

Maybe that can help someone who needs a quick fix. No guarantees, though!

Cheers,
Peter

Dear Peter,
many thanks for this implementation.

I’ve talked to one of the Pythia authors, who warned
me about potential weight-ordering/convention issues
that may arise by following straight the way you propose.

I’ll discuss with him in detail about this in the following
days, and we will get back to you as soon as possible.

Thanks,
Cheers.
Paolo

> On 25 Oct 2019, at 23:53, Peter Galler <email address hidden> wrote:
>
> Question #678584 on MadGraph5_aMC@NLO changed:
> https://answers.launchpad.net/mg5amcnlo/+question/678584
>
> Peter Galler posted a new comment:
> Dear Paolo,
> great! Then until this is implemented in MadGraph the following procedure seems to work.
> I copied the the relevant code piece ("HEPMC2HACK") from
> MadGraph/HEPTools/MG5aMC_PY8_interface/MG5aMC_PY8_interface.cc
> (which is the Pythia driver for LO runs where multiple weights in HEPMC work)
> into the file
> MadGraph/Template/NLO/MCatNLO/srcPythia8/Pythia82_hep.cc
> and this seems to do the trick.
> I include the listing for the entire modified file Pythia82_hep.cc below.
> If one doesn't want to mess up the MG installation one can also modify
> the same file in the local run directory of the process under study, i.e.
> <your_favourite_process>/MCatNLO/srcPythia8/Pythia82_hep.cc
>
> *********** Pythia82_hep.cc (modified) *********************************
>
> // Driver for Pythia 8. Reads an input file dynamically created on
> // the basis of the inputs specified in MCatNLO_MadFKS_PY8.Script
> #include "Pythia8/Pythia.h"
> #include "Pythia8Plugins/HepMC2.h"
> #include "Pythia8Plugins/aMCatNLOHooks.h"
> #include "Pythia8Plugins/CombineMatchingInput.h"
> #include "HepMC/GenEvent.h"
> #include "HepMC/IO_GenEvent.h"
>
> using namespace Pythia8;
>
> string convertFromAMCATNLO( string input) {
> string output="";
> // Count number of blanks in weight name.
> int appearances = 0;
> for(int n = input.find(" ", 0); n != int(string::npos);
> n = input.find(" ", n)) {
> appearances++;
> n++;
> }
> // Cut string by position of blanks.
> vector <string> pieces;
> for(int i =0; i < appearances;++i) {
> int n = input.find(" ", 0);
> pieces.push_back(input.substr(0,n));
> input = input.substr(n+1,input.size());
> }
> // Now pieces contains the details of the weight name.
> for(int i = 4; i < int(pieces.size()); ++i) {
> // Make upper case.
> transform(pieces[i].begin(),pieces[i].end(),pieces[i].begin(), ::toupper);
> // Add piece to weight name, with "_" as separator.
> output+=pieces[i] + "_";
> }
> return output;
> }
>
>
> int main() {
> Pythia pythia;
>
> string inputname="Pythia8.cmd",outputname="Pythia8.hep";
>
> pythia.readFile(inputname.c_str());
>
> //Create UserHooks pointer for the FxFX matching. Stop if it failed. Pass pointer to Pythia.
> CombineMatchingInput* combined = NULL;
> UserHooks* matching = NULL;
>
> int nAbort=10;
> int nPrintLHA=1;
> int iAbort=0;
> int iPrintLHA=0;
> int iEventtot=pythia.mode("Main:numberOfEvents");
> int iEventshower=pythia.mode("Main:spareMode1");
> string evt_norm=pythia.word("Main:spareWord1");
>
> //FxFx merging
> bool isFxFx=pythia.flag("JetMatching:doFxFx");
> if (isFxFx) {
> matching = combined->getHook(pythia);
> if (!matching) {
> std::cout << " Failed to initialise jet matching structures.\n"
> << " Program stopped.";
> return 1;
> }
> pythia.setUserHooksPtr(matching);
> int nJmax=pythia.mode("JetMatching:nJetMax");
> double Qcut=pythia.parm("JetMatching:qCut");
> double PTcut=pythia.parm("JetMatching:qCutME");
> if (Qcut <= PTcut || Qcut <= 0.) {
> std::cout << " \n";
> std::cout << "Merging scale (shower_card.dat) smaller than pTcut (run_card.dat)"
> << Qcut << " " << PTcut << "\n";
> return 0;
> }
> }
>
> pythia.init();
>
> HepMC::Pythia8ToHepMC ToHepMC;
> HepMC::IO_GenEvent ascii_io(outputname.c_str(), std::ios::out);
> // Do not store cross section information, as this will be done manually.
> ToHepMC.set_store_pdf(false);
> ToHepMC.set_store_proc(false);
> ToHepMC.set_store_xsec(false);
>
> // Cross section an error.
> double sigmaTotal = 0.;
> double errorTotal = 0.;
>
> for (int iEvent = 0; ; ++iEvent) {
> if (!pythia.next()) {
> if (++iAbort < nAbort) continue;
> break;
> }
> // the number of events read by Pythia so far
> int nSelected=pythia.info.nSelected();
>
> if (nSelected > iEventshower) break;
> if (pythia.info.isLHA() && iPrintLHA < nPrintLHA) {
> pythia.LHAeventList();
> pythia.info.list();
> pythia.process.list();
> pythia.event.list();
> ++iPrintLHA;
> }
>
> HepMC::GenEvent* hepmcevt = new HepMC::GenEvent();
> double evtweight = pythia.info.weight();
> double normhepmc;
> // ALWAYS NORMALISE HEPMC WEIGHTS TO SUM TO THE CROSS SECTION
> if (evt_norm != "sum") {
> normhepmc = 1. / double(iEventshower);
> } else {
> normhepmc = double(iEventtot) / double(iEventshower);
> }
> sigmaTotal += evtweight*normhepmc;
> //hepmcevt->weights().push_back(evtweight*normhepmc);
>
> // Loop through weightgroups.
> for ( std::map<std::string,LHAweightgroup>::const_iterator
> it_wg = pythia.info.weightgroups->begin();
> it_wg != pythia.info.weightgroups->end(); ++it_wg ) {
> // Loop through weights in the weightgroup.
> for ( std::map<std::string,LHAweight>::const_iterator
> it_w = it_wg->second.weights.begin();
> it_w != it_wg->second.weights.end(); ++it_w ) {
> // Get value of weight indexed by the present key.
> std::map<std::string,double>::const_iterator it_value
> = pythia.info.weights_detailed->find(it_w->first);
> double w = it_value->second; //*merging_weight*vetoWeights[iCut];
> // Get name.
> string name;
> for ( std::map<string,string>::const_iterator
> it_att = it_w->second.attributes.begin();
> it_att != it_w->second.attributes.end(); ++it_att)
> name += it_att->first + "=" + it_att->second + "_";
>
> // In aMC@NLO, the weight names should be extracted from the
> // weight tag contents, not its attributes.
> if (it_w->second.attributes.size()==0)
> name = convertFromAMCATNLO(it_w->second.contents);
>
> //if (pythia.info.lhaStrategy() == 3)
> // w *= norm_event_wgt / pythia.info.eventWeightLHEF;
> //// Weighted events.
> //else w *= 1 / double(1e9*nEvent);
>
> //w *=hepmcWeightRescaling;
>
> // Add event weight.
> hepmcevt->weights().push_back(w*normhepmc, name);
>
> } // Done looping through weights.
> } // Done looping through weight groups.
>
> ToHepMC.fill_next_event( pythia, hepmcevt );
> // Add the weight of the current event to the cross section.
> // Report cross section to hepmc
> HepMC::GenCrossSection xsec;
> xsec.set_cross_section( sigmaTotal, pythia.info.sigmaErr() );
> hepmcevt->set_cross_section( xsec );
> // Write the HepMC event to file. Done with it.
> ascii_io << hepmcevt;
> delete hepmcevt;
> }
>
> pythia.stat();
> if (isFxFx){
> std::cout << " \n";
> std::cout << "*********************************************************************** \n";
> std::cout << "*********************************************************************** \n";
> std::cout << "Cross section, including FxFx merging is: "
> << sigmaTotal << "\n";
> std::cout << "*********************************************************************** \n";
> std::cout << "*********************************************************************** \n";
> }
>
> return 0;
> }
>
>
> *********** End of Pythia82_hep.cc (modified)
> *********************************
>
> Maybe that can help someone who needs a quick fix. No guarantees,
> though!
>
> Cheers,
> Peter
>
> --
> You received this question notification because you are subscribed to
> the question.

Dear Paolo,
thank you very much for following this up.
Cheers,
Peter