Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Fix to prevent propagation of heavy neutrino from Geant4 #33299

Merged
merged 2 commits into from Apr 8, 2021
Merged

Fix to prevent propagation of heavy neutrino from Geant4 #33299

Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
c8d591c
proposed fix to avoid G4 propagation of HNL particle
mgratti Mar 24, 2021
86092ad
code checks
mgratti Mar 30, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
9 changes: 8 additions & 1 deletion SimG4Core/Generators/src/Generator.cc
View file
Open in desktop
Expand Up @@ -152,7 +152,8 @@ void Generator::HepMC2G4(const HepMC::GenEvent *evt_orig, G4Event *g4evt) {
for (pitr = (*vitr)->particles_begin(HepMC::children); pitr != (*vitr)->particles_end(HepMC::children); ++pitr) {
// For purposes of this function, the status is defined as follows:
// 1: particles are not decayed by generator
// 2: particles are decayed by generator but need to be propagated by
// 2: particles are decayed by generator but need to be propagated by GEANT
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@mgratti , to be more clear: status=1 means that Geant4 takes care on particle transport and decay particles according to the list of Geant4 decay channels; status=2 means that Geant4 takes care on particle transport but decay is predefined by generator - list of secondaries and their kinematics coming from PHYTHIA or other generator, Geant4 peaks up this list, transports particles to the decay point and generate secondaries according to the list. Status = 2 may be given not only to an exotic particles but to any normal as well, for example, for study of a rare decay.

This is a standard prescription of HepMC interface. Unfortunately, in the case of exotics generators use to provide status code whatever. It is why we have such complex logic inside this code. We do not fully control what values of the status code can be given and update this logic for each new case.

As I understand your needs, your long-lived particle do not visible in the detector but its decay products are visible, so you need displaced secondary vertex to be given to Geant4.

Let us discuss this case assuming that secondary vertex is defined and secondary particles are available. Here I have a difficulty with definition status=3, because some generators may provide such code. The worry is,g that your analysis will work fine but another analysis may be broken.

So, the question: if you claim a particle with abs(pdg) = 9900015 invisible, is it possible to provide Geant4 displaced vertex with final particles or you already tried and this is not working for you?

Copy link
Contributor Author

@mgratti mgratti Mar 31, 2021 •
edited

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@civanch

“status=2 means that Geant4 takes care on particle transport but decay is predefined by generator - list of secondaries and their kinematics coming from PHYTHIA or other generator, Geant4 picks up this list, transports particles to the decay point and generate secondaries according to the list. “

=> I see what you mean, indeed my heavy neutrino has status=2 and its displaced decay products have status=1 (output codes obtained from performing decay of the particle in EvtGen); my tests have shown that something must go wrong with the transportation of the particle when the long-lived particle is unknown to Geant4, like in my case (pdgId=9900015, g4code is null), because the simTracks associated to the decay products are not displaced.

“As I understand your needs, your long-lived particle do not visible in the detector but its decay products are visible, so you need displaced secondary vertex to be given to Geant4.”

=> yes

“Here I have a difficulty with definition status=3, because some generators may provide such code. The worry is,g that your analysis will work fine but another analysis may be broken.”

=> Iiuc, the pathologic situation that you describe could happen also for all other exotic particles; if they were generated by a generator different from Pythia, the rules for modifying their states that are currently implemented would not be valid anymore.

Also, I have searched for “9900015” in cmssw github repository and in the fragments repository without finding any occurrence.

“So, the question: if you claim a particle with abs(pdg) = 9900015 invisible, is it possible to provide Geant4 displaced vertex with final particles or you already tried and this is not working for you?”

=> If the status of the particle with abs(pdg) = 9900015 is equal to 3, then I am able to provide Geant4 the displaced vertex with the final particles. If instead the status of the particle is 2 (as originally set by EvtGen generator), the final state particles get attached to the production vertex of the particle with abs(pdg) = 9900015.

// 3: particles are decayed by generator and do not need to be propagated by GEANT
int status = (*pitr)->status();
int pdg = (*pitr)->pdg_id();
if (status > 3 && isExotic(pdg) && (!(isExoticNonDetectable(pdg)))) {
Expand All @@ -163,6 +164,9 @@ void Generator::HepMC2G4(const HepMC::GenEvent *evt_orig, G4Event *g4evt) {
// be propagated by GEANT, so do not change their status code.
status = 2;
}
if (status == 2 && abs(pdg) == 9900015) {
Copy link
Contributor

@civanch civanch Apr 3, 2021 •
edited

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@mgratti , would it be more correct to have following lines:

// heavy neutrino should not be propagated by Geant4
// but the vertex should be considered if is outside fiductial volume

if(std::abs(pdg) == 9900015) {
status = 2;
}

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

here if decay vertex of heavy neutrino will be inside fiductial volume (probability is low but still possible) then decay vertex of it will not be created but secondary particles will be created.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@civanch, Sorry for the delay in replying to your comments.

For my case, it does not make a difference because the production vertex of the HNL is anyway created because there is another status=1 particle that is a child of that vertex.

I thought of potential other use case: for example a signal with two HNLs produced at a vertex and no other status=1 particle (extremely rare). If both HNLs decay within the fiducial region, depending on the HNL status, the production vertex will (status=3) or will not (status=2) be saved as a G4 vertex. However, this will have no impact in the end; in fact the HNL will have status=3 in the second part (to avoid G4 propagation) and, as a result, no particle will be associated to this G4 vertex.

status = 3;
}

// Particles which are not decayed by generator
if (status == 1) {
Expand Down Expand Up @@ -244,6 +248,9 @@ void Generator::HepMC2G4(const HepMC::GenEvent *evt_orig, G4Event *g4evt) {
if (status > 3 && isExotic(pdg) && (!(isExoticNonDetectable(pdg)))) {
status = hasDecayVertex ? 2 : 1;
}
if (status == 2 && abs(pdg) == 9900015) {
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@mgratti , why not to add following lines:

// heavy neutrino should not be propagated by Geant4
if(std::abs(pdg) == 9900015) status = 3;

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@civanch, sure, this clearly preserves my case.

status = 3;
}

// this particle has predefined decay but has no vertex
if (2 == status && !hasDecayVertex) {
Expand Down