Fork me on GitHub

Opened 12 years ago

Closed 12 years ago

#163 closed Bug (fixed)

no b-tagged jets // can't read hepmc file

Reported by: Benjamin Gerber Owned by:
Priority: major Milestone:
Component: Delphes miscellaneous Version: Delphes 3
Keywords: Cc:

Description

I use Herwig++-2.6.2 and Delphes 3.0.3 and I have two problems:

  1. I didn't change anything and it's the same for both detector cards, delphes_card_ATLAS.tcl and delphes_card_CMS.tcl, I don't get any b-tagged jets, no matter what input file I use (for Herwig). I also tried with setting the efficiencyformula to 1 and the two misidentification to 0, still no b-tagged jets.
  1. If I simulate e.g. 1'000 events with Herwig (combined with a slha-file) and feed the hepmc-file to Delphes, it works. But if I do exactly the same, but with 25'000 events, there are hundreds of lines (the whole terminal) with "Reading nameoffile.hepmc" untill the last one says: "ERROR: can't open nameoffile.hepmc".

Thx&Greets,
Beny

Attachments (1)

b.hepmc (1.4 MB ) - added by Benjamin Gerber 12 years ago.
created with LHC-Bottom.in from Herwig

Download all attachments as: .zip

Change History (28)

comment:1 by Pavel Demin, 12 years ago

Could you, please provide an example of your HepMC file with b-jets and the command(s) that you are using in 2?

by Benjamin Gerber, 12 years ago

Attachment: b.hepmc added

created with LHC-Bottom.in from Herwig

in reply to:  1 comment:2 by Benjamin Gerber, 12 years ago

Voila...

Commands in 2:

(for Herwig: "./Herwig++ read someinfile.in", then "Herwig++ run somefile.run")
for Delphes: "./DelphesHepMC delphes_card_ATLAS.tcl name.root name.hepmc"

comment:3 by Pavel Demin, 12 years ago

Thanks for the file. All gluons and quarks in this file have status 1:

P 10009 5 -8.9362277011480789e+00 2.0281696228266203e+01 1.1163856056512866e+02 1.1389473900262614e+02 4.2000000000000002e+00 1 0 0 0 1 1 501

Looks like the hadronization is switched off.

Normally, Delphes expects hadronized events as its input with only stable particles having status 1.

If you need to test b-tagging with events without hadronization, you can replace line 334 in readers/DelphesHepMC.cpp:

  else if(fStatus == 2)

with

  else if(fStatus == 1)

comment:4 by Pavel Demin, 12 years ago

Could you send me the HepMC file with 25k events? If it's too big I'm not sure that it would be possible to attach it to this ticket or send by e-mail. Could you put it on a web server or some file sharing service?

comment:5 by Benjamin Gerber, 12 years ago

Hm...I still don't get any b-tagged jets in the root file.

It's 6.1 GB. Can I send it per Skype (name: benygerber)?

comment:6 by Pavel Demin, 12 years ago

6.1 GB is too big.

Could you give me some instructions on how to generate such file (I never used Herwig++)?

The replacement line should be

if(fStatus == 1)

without 'else'.

Sorry for the confusion.

comment:7 by Pavel Demin, 12 years ago

I've just looked more carefully into your file and I don't understand why you use it as input for Delphes.

The file contains only quarks and gluons in the final state. No hadrons, no leptons.

Delphes does not do hadronization.

Could you, please explain what you trying to achieve?

comment:8 by Pavel Demin, 12 years ago

I've composed a large 6GB HepMC file sticking together many small ones and I can reproduce your problem. I'll try to fix it ASAP.

comment:9 by Benjamin Gerber, 12 years ago

Ok, I found now some b-tagged jets in the b.hepmc with this change.

The b.hepmc comes from the LHC-Bottom.in file which is the input for Herwig. It is a test file and I think it should more or less just produce b-quarks, but I don't understand what it exactly does.

The other 6.1 GB hepmc file is from LHC-MSSM.in file which includes a SLHA file to simulate a SUSY model. My goal is just to redo the m(gluino,neutralino)=(800,750) and (500,150) plot from: http://inspirehep.net/record/1198427/files/MET_MC_vs_Data_inc_signal.png ,i.e. from http://arxiv.org/abs/ARXIV:1211.1167 , to check if Herwig&Delphes work correctly. I have both SLHA files and for the (800,750) I made, as they did, 10'000 events. This works. But for the (500,150) I need 25'000 events but there suddenly I have the "can't read" error.

Ah ok...so it is a bug or what do you think?

anyway, thx already!

comment:10 by Pavel Demin, 12 years ago

Yes, the problem with 6GB files is a bug.

It appears on 32-bit systems with files that are greater than 2GB.

It's a known problem with a known solution:

  • replace ftell with ftello
  • add compilation option -D_FILE_OFFSET_BITS=64

I'll try make a release for the fixed version by tomorrow.

comment:11 by Pavel Demin, 12 years ago

My point is that the events generated with your LHC-Bottom.in card contain only partons and no stable particles.

This type of events should not be used as Delphes input.

Last edited 12 years ago by Pavel Demin (previous) (diff)

comment:13 by cbpark, 12 years ago

I'm not sure whether this is relevant in this ticket but want to say my experience. I have tested Delphes 3 (both 3.0.3 and 3.0.5) to see the effect of the b-tagging efficiency.

For a simple example, a HepMC event file containing the SM top-pair process has been produced using Pythia 8, and then I set the b-tagging efficiency in the delphes_card_ATLAS.tcl to 1.0 and the mis-tagging efficiency was set to 0.0. But eventually, I found that there was no b-tagged jet no matter how large the number of events.

So, I looked closely into the BTagging module. Simply saying, it is testing the existence of partonArray at first, and then uses the b-tagging efficiency formula when the parton in the array has particle ID 5. The partonArray is filled when the particle status is 2 in the code, DelphesHepMC.cpp. Although the convention for the particle status code in Pythia 8 is dissimilar, Pythia 8 converts the status code into the convention of HepMC when producing the hepmc output file (See here).

The problem is that the notion of the status code 2 in Pythia 8 seems not to be matched that in HepMC. The b-quark parton in the hepmc file produced in Pythia 8 has the status code 23 instead of 2. So, I changed the lines in DelphesHepMC.cpp

  else if(fStatus == 2)
  {
    fPartonOutputArray->Add(candidate);
  }

into

  else if(fStatus == 23)
  {
    fPartonOutputArray->Add(candidate);
  }

then found that the b-tagging efficiency formula is applying well (if the efficiency were set to be 1.0, I can see 1 or 2 b-jets in every event).

One bad thing is that Herwig++ seems to be using another different convention for that. I'm not sure whether this issue should go to Pythia 8 or Herwig++ (ThePEG).

FYI, I used Pythia 8175, which is the latest public version.

comment:14 by Pavel Demin, 12 years ago

Thank you very much for the detailed analysis of the problem.

One solution for this problem could be to not rely on the status codes but on the PDG ID codes (1-5, 21 and 15 for tau-tagging).

Here is the code:

Int_t pdgCode = TMath::Abs(candidate->PID);

...

  else if(pdgCode < 5 || pdgCode == 21 || pdgCode == 15)
  {
    fPartonOutputArray->Add(candidate);
  }

What do you think?

comment:15 by Benjamin Gerber, 12 years ago

@pavel: Thank you very much! Works fine now.

@b-tag-topic: I'm very new in this area (I'll do simulation of SUSY models as Masters thesis), so I don't really know how all this codes work atm.

in reply to:  14 ; comment:16 by cbpark, 12 years ago

Thanks a lot! I have just tested your suggestion (by the way, it should be like else if (pdgCode < 6 || pdgCode == 21 || pdgCode == 15), right?) with the same hepmc data. It seems to be working as well. Still, I wonder whether this part affects the lepton, i.e., electron or muon, and photon isolation part (I'm asking because I don't know the whole code flow).

Another issue is that the jets have Charge, which has a value 1 or -1. Is this an electromagnetic charge or something else? If it's the EM charge, it looks quite strange as the jet does not possess it in reality.

Replying to pavel:

Thank you very much for the detailed analysis of the problem.

One solution for this problem could be to not rely on the status codes but on the PDG ID codes (1-5, 21 and 15 for tau-tagging).

Here is the code:

Int_t pdgCode = TMath::Abs(candidate->PID);

...

  else if(pdgCode < 5 || pdgCode == 21 || pdgCode == 15)
  {
    fPartonOutputArray->Add(candidate);
  }

What do you think?

comment:17 by Pavel Demin, 12 years ago

Thanks for the test and for the correction. You are right it should be pdgCode < 6.

The jet's charge is set in the TauTagging module as gRandom->Uniform() > 0.5 ? 1 : -1.

Normally, this charge is only relevant for tau-jets. May be we should set it only for jets with TauTag == 1.

Last edited 12 years ago by Pavel Demin (previous) (diff)

comment:18 by Benjamin Gerber, 12 years ago

Maybe a stupid question, but why do you take quarks "pgdCode<6" and gluons "pdgCode==21" for tau-tagging? I mean, why not just "pdgCode==15"?

comment:19 by Pavel Demin, 12 years ago

'pdgCode < 6' and 'pdgCode == 21' are for b-tagging

'pdgCode == 15' is for tau-tagging

comment:20 by Benjamin Gerber, 12 years ago

So do I get this right: fPartonOutuputArray checks for a jet and then somewhere else is checked if pgdCode==5 --> b-tagging and pdgCode==15 --> tau-tagging?

comment:21 by Pavel Demin, 12 years ago

fPartonOutuputArray contains all quarks, gluons and tau leptons in the event.

Then the !BTagging and TauTagging modules use the particles from fPartonOutuputArray to match them with jets to define the jet's origin.

Last edited 12 years ago by Pavel Demin (previous) (diff)

in reply to:  16 comment:22 by Pavel Demin, 12 years ago

Replying to cbpark:

Still, I wonder whether this part affects the lepton, i.e., electron or muon, and photon isolation part (I'm asking because I don't know the whole code flow).

This modification should not affect anything except the BTagging and TauTagging modules.

There is a data flow diagram at WorkBook/DataFlowDiagram

comment:23 by Archana, 12 years ago

Has this issue been resolved? I find the same problem with Delphes-3.0.5 and the version downloaded from svn on April 4. I find no b-tagged jets in data generated with Pythia in hepmc format. At the parton level, I do find bs. The fix suggested by cbpark seems to work-

else if(fStatus == 23)

{

fPartonOutputArray->Add(candidate);

}

But I am curious why this hasnt been implemented in updates from svn. Would you suggest changing the code manually to fix this issue?

comment:24 by Pavel Demin, 12 years ago

The problem with the (fStatus == 23) fix is that it is generator dependent. For example, it won't work with older Pythia versions.

Normally, the (pdgCode <= 5 || pdgCode == 21 || pdgCode == 15) fix is in the svn repository and we're planning to release it in the coming days.

Could you, please check if the following prerelease works better for you?

http://cp3.irmp.ucl.ac.be/downloads/Delphes-3.0.6.pre2.tar.gz

comment:25 by Archana, 12 years ago

The pre-release 3.0.6 version seems to working fine. I find that the number of b-jets roughly matches with the "manual" fix to 3.0.5

comment:26 by Pavel Demin, 12 years ago

Thanks for the test.

Currently, the repository for the Delphes 3.x code is in 'branches/ModularDelphes'. We'll try to move this repository to the trunk when the 3.0.6 release is ready.

We've just discovered another problem with the 3.0.6.pre2 prerelease and fixed it in the new prerelease:

http://cp3.irmp.ucl.ac.be/downloads/Delphes-3.0.6.pre3.tar.gz

comment:27 by Pavel Demin, 12 years ago

Resolution: fixed
Status: newclosed

Delphes 3.0.6 is released.

Subversion trunk contains now the Delphes 3.x code.

Thanks a lot for your suggestions and tests.

Note: See TracTickets for help on using tickets.