Compute JES/JER-varied MELA probabilities - #400
Conversation
namapane
left a comment
namapane
left a comment
There was a problem hiding this comment.
Thanks @GlennnLiu for this development! I have put some comments/suggestions.
@pianonick411 : you may want to have a look to changes in MELAProbHelper. As far as I can see, there is nothing affecting the computation of LHE probabilities.
| ADD_ALLEVENTS = getConf("ADD_ALLEVENTS", IsSIGNAL) # if true, add a separate tree with gen-level variables for all events (not just those passing the candidate selection); by default, this is done for signal samples | ||
| ADD_LHE_PROB = getConf("ADD_LHE_PROB", ADD_ALLEVENTS) # Add LHE angles and probabilities. This is in general the case whenever ADD_ALLEVENTS is true (ie for signals) | ||
| JES_SPLITTING = getConf("JES_SPLITTING", True) # Whether to split JES variations into 11 components (if false, only up/down variations are produced, by summing all components in quadrature) | ||
| COMPUTE_JET_VARIATIONS_MELA = getConf("COMPUTE_JET_VARIATIONS_MELA", True) # Compute MELA probabilities for JES/JER shifted jets |
There was a problem hiding this comment.
Can we leave the default to False to avoid changing the current behaviour for other users?
For example, VBS people use some jet probabilities (pyFragments/VBS_probs.py) but may not want to activate this, especially given the large overhead and the fact that they implemented variations in their preprocessor-equivalent code.
We could add COMPUTE_JET_VARIATIONS_MELA=True it in the STXS pyFragment like we do for bestCandByMELA to have it on automatically for our processing without having to change the csv files.
(It is now possible to add default variables in the format line of csv fles, so we could simply modify 1 line in each csv file, but I prefer the PyFragment way as it concentrates analysis-specific customizations in a single place).
| if recoProbFiller is not None and COMPUTE_JET_VARIATIONS_MELA: | ||
| recoProbFiller.setJetVariations(jetCorrector) | ||
| elif recoProbFiller is not None: | ||
| print("***RecoProbFiller: jet-varied MELA computation disabled", flush=True) |
There was a problem hiding this comment.
I would actually skip this message - the logic of initialization messages is that everything that is set up issues a message (and the message in this case is already issued by RecoProbFiller), so no message -> not done
| py = pt * math.sin(phi) | ||
| pz = pt * math.sinh(eta) | ||
| energy = math.sqrt(max(mass * mass + px * px + py * py + pz * pz, 0.)) | ||
| return Mela.SimpleParticle_t(pdgId, px, py, pz, energy) |
There was a problem hiding this comment.
The rest of the code goes through p4() which internally uses TLorentzVector [here]. The math is obviously the same but the different implementation in C++ and python may give small numerical effects and in principle it's not so nice to do nominal in one way and variations in another. I would have a small preference to adopt an as close as possible implementation in the two cases, although any numerical effect is in principle already covered by your validation tests.
| selected = [] | ||
| for idx, jet in enumerate(jets): | ||
| pt = getattr(jet, variation + "_pt") | ||
| mass = getattr(jet, variation + "_mass") |
There was a problem hiding this comment.
Since jets is the full list of jets with no selection, it would make sense to:
- check
jet.jetIdfirst so that no computation is done for jets that are discarded anyhow - move the creation of the variable names (
variation + "_pt", variation + "_mass") outside the loop so that it is done only once per variation and not for every jet for every variation - move the extraction of mass below the current L137 so that it is not done for jets that fail other cuts
| mass = getattr(jet, variation + "_mass") | ||
| leptonPt = jet.ZZLepEF * jet.pt | ||
| overlapsLeptons = pt <= 0. or leptonPt / pt > 0.5 | ||
| if overlapsLeptons or jet.jetId != 6 or pt <= jet.ptThreshold: |
There was a problem hiding this comment.
One potentially dangerous thing here is that if for some reason somebody (other than you) decides to change ID or cleaning criteria in jetFiller.py, it is practically certain (s)he will not realize a matching change is required here. There is no quick workaround, but at the very least we should put a big reminder in jetFiller.py where ID cut and EFthreshold are defined.
| if overlapsLeptons or jet.jetId != 6 or pt <= jet.ptThreshold: | ||
| continue | ||
| selected.append((pt, idx, mass)) | ||
| selected.sort(key=lambda item: item[0], reverse=True) |
There was a problem hiding this comment.
selected.sort(key=itemgetter(0), reverse=True) (with from operator import itemgetter) would save the overhead of defining and calling the lambda.
[Note: the code in jetFiller.py picks the 2 elements manually avoiding sort on purpose, to avoid the cost of sorting all elements when we only need the first two. Of course that does not make sense here]
| 0, pt, jet.eta, jet.phi, mass | ||
| )) | ||
| for idx in (aCand.extraLep1Idx, aCand.extraLep2Idx): | ||
| if idx < 0: |
There was a problem hiding this comment.
extra leptons don't change across jet variations: the extra leptons' SimpleParticle_t could be built once for each candidate and cached before the loop on variations, so that they do not need to be remade 22 times.
|
A general comment on your very nice comparison plots: |
2 participants
Dear Experts,
I (and codex) managed to address the computation of JES/JER-varied MELA probabilities. Only 24 * 11 MELA probabilities are additionally selected; no intermediate variables like varied leading / subleading jet indices are saved.
Tests have been performed with 2022EE ggH125, VBFH125 and partially ZZTo4l. For ggH125, VBFH125, I performed with or without JES/JER-varied probabilities. The runtime and file sizes are reported below:
We can see that the runtime is significantly slower, but in any case this is the price we have to pay, no matter computing them here or in the preprocessor.
Besides, comparison plots of the up/down variations with the nominal ones are made. Plots showing the distributions are shown in [1], and ratios of variations divided by the nominal ones are shown in [2], for the three processes.
Please let me know if the results are acceptable.
Besides, in
nanoZZ4lAnalysis.py, jetVetoMap is put after JERC.Thanks!
Yours Sincerely,
Geliang
[1] https://geliu.web.cern.ch/HZZ4l_Run3/JetVariedMELA/jet_varied_mela/
[2] https://geliu.web.cern.ch/HZZ4l_Run3/JetVariedMELA/jet_varied_mela_ratio/