Fix to prevent propagation of heavy neutrino from Geant4

[mgratti]

PR description:

This PR implements a minor change in order to prevent Geant4 from performing the propagation of a heavy neutrino (pdg=9900015) with status 2 outside of the beam pipe.
Heavy neutrinos of status 2 are the result of the event generation by Pythia+EvtGen, where EvtGen handles the decay of the heavy neutrino.

Additional information documenting the need for this fix can be found in the hypernews thread

This is meant to be used for an MC production using the b-parking conditions, so a backport to 10_2_X will be submitted too, once this PR is merged.

PR validation:

Followed the instructions here:
http://cms-sw.github.io/PRWorkflow.html

The report of runTheMatrix.py -l limited -i all --ibeos is attached:
runTheMatrix_short.log

Code checks have been performed and included as a separate commit.

[cmsbuild]

+code-checks

Logs: https://cmssdt.cern.ch/SDT/code-checks/cms-sw-PR-33299/21823

• This PR adds an extra 20KB to repository

[cmsbuild]

A new Pull Request was created by @mgratti (Maria Giulia Ratti) for master.

It involves the following packages:

SimG4Core/Generators

@cmsbuild, @civanch, @mdhildreth can you please review it and eventually sign? Thanks.
@makortel, @cvuosalo, @rovere, @fabiocos, @slomeo this is something you requested to watch as well.
@silviodonato, @dpiparo, @qliphy you are the release manager for this.

cms-bot commands are listed here

[civanch]

@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?

[mgratti]

@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.

[civanch]

please test

[cmsbuild]

+1

Summary: https://cmssdt.cern.ch/SDT/jenkins-artifacts/pull-request-integration/PR-c06ed5/13984/summary.html
COMMIT: 86092ad
CMSSW: CMSSW_11_3_X_2021-04-02-1800/slc7_amd64_gcc900
User test area: For local testing, you can use /cvmfs/cms-ci.cern.ch/week0/cms-sw/cmssw/33299/13984/install.sh to create a dev area with all the needed externals and cmssw changes.

Comparison Summary

Summary:

• No significant changes to the logs found
• Reco comparison results: 4 differences found in the comparisons
• DQMHistoTests: Total files compared: 37
• DQMHistoTests: Total histograms compared: 2640868
• DQMHistoTests: Total failures: 7
• DQMHistoTests: Total nulls: 0
• DQMHistoTests: Total successes: 2640839
• DQMHistoTests: Total skipped: 22
• DQMHistoTests: Total Missing objects: 0
• DQMHistoSizes: Histogram memory added: 0.0 KiB( 36 files compared)
• Checked 155 log files, 37 edm output root files, 37 DQM output files
• TriggerResults: no differences found

[civanch]

@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;
}

[civanch]

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.

[mgratti]

@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.

[civanch]

@mgratti , why not to add following lines:

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

[mgratti]

@civanch, sure, this clearly preserves my case.

[civanch]

@mgratti ,

please, try out if my proposal is working for you. I was able to understand what is going on slowly, in order to be more effective next time it is important to add comments to your modifications. Please, check if this variant work for long distance of heavy neutrino and also for short one like 1 mm.

If you agree with these modifications, please update this PR and after propose backport PRs to Run-2 simulation production: CMSSW_10_6_X, CMSSW_9_3_X, and _CMSSW_7_1_X.

Note, that Generator.cc may be slightly different in previous releases, at least, lines with printouts may be different, so use cut/paste to add your fixes.

[civanch]

+1

let us merge this for the time being

[cmsbuild]

This pull request is fully signed and it will be integrated in one of the next master IBs (tests are also fine). This pull request will now be reviewed by the release team before it's merged. @silviodonato, @dpiparo, @qliphy (and backports should be raised in the release meeting by the corresponding L2)

[silviodonato]

+1

[mgratti]

@civanch, In relation to your comment, I can say that your proposed modifications still preserve my use case and slightly improve the logic, at the cost of having a different status of the HNL in different parts of the code.
Please, let me know if you would like me to open a PR that implements your suggestions.
Otherwise, I will proceed with back-porting in 10_2_X.

[civanch]

@mgratti , please check if one or both solutions work fine for short decay length (let say below 2 cm) as for long decay length.

[mgratti]

@civanch, I checked that both solutions (mine and yours) work fine for short and long decay length.

[civanch]

@mgratti , then proceed with backport of your solution may be pick up my comments.