From 287f265b503a9f73841f2cf2bcb5cb7a112924e3 Mon Sep 17 00:00:00 2001 From: Juraj Smiesko Date: Wed, 19 Aug 2026 14:28:51 +0200 Subject: [PATCH 1/6] Flesh out Gen/ee tutorial: WHIZARD mumuH card, Pythia8 k4Gen step Adapts the FCC-config golden card (winter2023 branch) for e+e- -> mu mu H, stripping its internal legacy-PYTHIA6 shower/decay to keep the two-step WHIZARD (LHE) -> Pythia8 (k4Gen Gaudi, H -> b b) chain explicit for students, and adds the runnable mumuH.sin, mumuH_Hbb.cmd, and pythia_mumuH.py files alongside the README. --- Gen/ee/README.md | 193 ++++++++++++++++++++++++++++++++++++++++- Gen/ee/mumuH.sin | 29 +++++++ Gen/ee/mumuH_Hbb.cmd | 7 ++ Gen/ee/pythia_mumuH.py | 28 ++++++ 4 files changed, 254 insertions(+), 3 deletions(-) create mode 100644 Gen/ee/mumuH.sin create mode 100644 Gen/ee/mumuH_Hbb.cmd create mode 100644 Gen/ee/pythia_mumuH.py diff --git a/Gen/ee/README.md b/Gen/ee/README.md index 20ebf86..830c8c4 100644 --- a/Gen/ee/README.md +++ b/Gen/ee/README.md @@ -1,5 +1,192 @@ # Gen - FCC-ee -WHIZARD generation of e+e- -> mu mu H (output HEPMC), followed by Pythia -showering and H -> b b decay. Golden run cards live in FCCConfig; the actual -card used here comes from FCCDirac. +WHIZARD generates the hard process e+e- -> mu+ mu- H, then Pythia8 showers, +hadronizes, and decays H -> b b. See the callout below for why the process +is e+e- -> mu+ mu- H and not e+e- -> Z H. + +## Environment setup + +Everything needed (WHIZARD, Pythia8, the Key4hep/Gaudi tools) comes from the +Key4hep stack: + +``` +source /cvmfs/sw.hsf.org/key4hep/setup.sh +``` + +(`/cvmfs/fcc.cern.ch/sw/latest/setup.sh` provides the same Key4hep stack +plus FCC-specific tools on top.) + +## Step 1: WHIZARD - hard process + +Check WHIZARD is available: + +``` +which whizard +``` + +The golden card for this process, +[`wzp6_ee_mumuH_Hbb_ecm240.sin`](https://github.com/HEP-FCC/FCC-config/blob/winter2023/FCCee/Generator/Whizard/v3.0.3/wzp6_ee_mumuH_Hbb_ecm240.sin), +lives on the `winter2023` campaign branch of `FCC-config` (each production +campaign gets its own branch off `main` — the golden card wasn't on `main` +itself). That card does everything in one go: it showers, hadronizes, and +forces H -> b b internally through WHIZARD's built-in interface to the +*legacy* PYTHIA6, and writes `stdhep`. That's how production samples are +actually made, but it skips the explicit, separate Pythia8 step this +tutorial wants students to see and configure themselves. + +> WHIZARD *can* be linked against PYTHIA8 at build time +> (`--enable-pythia8`), but as of the current WHIZARD manual (covering up +> to v3.4.3) that integration isn't finished: the "Parton shower and +> hadronization from PYTHIA8" manual section is an empty stub, and the +> documented values for `$shower_method` are only `"WHIZARD"` (in-house) +> and `"PYTHIA6"` — `"PYTHIA8"` isn't one of them. So showering with +> Pythia8 as a separate step (Step 2 below) isn't just the pedagogically +> cleaner choice, it's the only well-supported way to bring Pythia8 into +> this chain at all. + +So `mumuH.sin` below keeps the golden card's physics (process, energy, beam +spread and ISR settings) but strips the internal shower/hadronization/decay +lines and switches the output to LHEf, the same bare-generator pattern used +by the existing tutorial's own +[`Z_mumu.sin`](https://fccsw.web.cern.ch/fccsw/share/gen/whizard/Zpole/Z_mumu.sin) +example. The card is [`mumuH.sin`](mumuH.sin) in this directory: + +``` +model = SM + +# Center of mass energy +sqrts = 240 GeV + +process mumuH = e1, E1 => e2, E2, H + +beams = e1, E1 => gaussian => isr +?keep_beams = true +?keep_remnants = true + +gaussian_spread1 = 0.185% +gaussian_spread2 = 0.185% + +?isr_handler = true +$isr_handler_mode = "recoil" +isr_alpha = 0.0072993 +isr_mass = 0.000511 + +# Golden-card value (production-quality precision, slow to integrate live): +# integrate (mumuH) { iterations = 10:100000:"gw", 5:200000:"" } +# Classroom default (faster, lower precision): +integrate (mumuH) { iterations = 3:2000:"gw" } + +n_events = 1000 + +$lhef_version = "3.0" +sample_format = lhef +simulate (mumuH) { $sample = "mumuH" } +``` + +Run it in its own directory: + +``` +mkdir -p test_whizard/mumuH && cd test_whizard/mumuH +cp ../../mumuH.sin . +whizard mumuH.sin +``` + +This produces `mumuH.lhe`. + +### Why mu mu H, not Z H? + +WHIZARD is asked to generate e+e- -> mu+ mu- H directly, rather than +e+e- -> Z H with Z -> mu mu. The final state looks the same, but the +generator-level parent history doesn't: in the mu mu H process the two +muons' parent particles are the incoming e+/e-, not a Z boson. This matters +downstream once students truth-match particles in the Analysis stage — +don't be surprised the muons have no Z parent in the event record. + +## Step 2: Pythia8 - shower, hadronize, decay H -> b b + +FCC-ee tooling doesn't call Pythia8 as a bare standalone binary — it goes +through the Gaudi components in +[`key4hep/k4Gen`](https://github.com/key4hep/k4Gen) (`PythiaInterface` + +`GenAlg`), run with `k4run`, the same framework used later for Delphes and +FCCAnalyses. `PythiaInterface` reads a `.cmd` card, which can point at an +external LHE file exactly like the working example shipped in `k4Gen` itself +([`data/Pythia_LHEinput.cmd`](https://github.com/key4hep/k4Gen/blob/main/k4Gen/data/Pythia_LHEinput.cmd)). + +[`mumuH_Hbb.cmd`](mumuH_Hbb.cmd) in this directory combines that LHE-reading +pattern with the H -> b b decay-forcing settings from +[`FCC-config`'s `p8_ee_H_Hbb_ecm125.cmd`](https://github.com/HEP-FCC/FCC-config/blob/main/FCCee/Generator/Pythia8/p8_ee_H_Hbb_ecm125.cmd): + +``` +! Read in the WHIZARD LHEf file +Beams:frameType = 4 +Beams:LHEF = mumuH.lhe + +! Force H -> b b +25:onMode = off +25:onIfAny = 5 +``` + +Steering script [`pythia_mumuH.py`](pythia_mumuH.py), adapted from `k4Gen`'s +own +[`options/pythia.py`](https://github.com/key4hep/k4Gen/blob/main/k4Gen/options/pythia.py) +example — swap in the LHE-reading card above and write out EDM4hep rather +than plain HepMC (`HepMCFileWriter`'s own docstring says it's for debugging, +not for actual event storage): + +```python +from Gaudi.Configuration import * +from edm4hep import labels as e4_labels + +from Configurables import EventDataSvc +from k4FWCore import ApplicationMgr, IOSvc +ApplicationMgr().EvtSel = 'NONE' +ApplicationMgr().EvtMax = 1000 +ApplicationMgr().ExtSvc += ["RndmGenSvc", EventDataSvc("EventDataSvc")] + +from Configurables import PythiaInterface +pythia8gentool = PythiaInterface() +pythia8gentool.pythiacard = "mumuH_Hbb.cmd" + +from Configurables import GenAlg +pythia8gen = GenAlg("Pythia8") +pythia8gen.SignalProvider = pythia8gentool +pythia8gen.hepmc.Path = "hepmc" +ApplicationMgr().TopAlg += [pythia8gen] + +from Configurables import HepMCToEDMConverter +hepmc_converter = HepMCToEDMConverter() +hepmc_converter.hepmc.Path = "hepmc" +hepmc_converter.hepmcStatusList = [] +hepmc_converter.GenParticles.Path = e4_labels.MCParticles +ApplicationMgr().TopAlg += [hepmc_converter] + +iosvc = IOSvc() +iosvc.Output = "mumuH_Hbb.root" +``` + +Copy both files next to the LHE file produced in Step 1 and run: + +``` +cp ../../mumuH_Hbb.cmd ../../pythia_mumuH.py . +k4run pythia_mumuH.py +``` + +This produces `mumuH_Hbb.root`, an EDM4hep file with the showered, +hadronized, H -> b b decayed event record. + +> **Confirmed handoff to `Sim/ee`:** this stage hands off an EDM4hep ROOT +> file with an `MCParticles` collection, not a HepMC file and not a combined +> Pythia8+Delphes run. Checked `key4hep/k4SimDelphes`: its Gaudi component +> `k4SimDelphesAlg` reads a generic `edm4hep::MCParticleCollection` (data +> path `"GenParticles"`) — its own example steering script feeds it from a +> particle gun, not Pythia8, proving it doesn't care how that collection was +> produced. `DelphesPythia8_EDM4HEP` (which runs Pythia8 internally) is a +> separate, alternative standalone entry point in the same repo, not the +> only way in. So `Sim/ee` should use `k4SimDelphesAlg` via `k4run`, reading +> the `MCParticles` collection from `mumuH_Hbb.root` produced here — the +> Gen/Sim split as designed is correct. + +## What's next + +The showered, H -> b b decayed sample (`mumuH_Hbb.root`) is the input to +Delphes fast simulation with the FCC-ee IDEA card — see `Sim/ee`. diff --git a/Gen/ee/mumuH.sin b/Gen/ee/mumuH.sin new file mode 100644 index 0000000..af102f9 --- /dev/null +++ b/Gen/ee/mumuH.sin @@ -0,0 +1,29 @@ +model = SM + +# Center of mass energy +sqrts = 240 GeV + +process mumuH = e1, E1 => e2, E2, H + +beams = e1, E1 => gaussian => isr +?keep_beams = true +?keep_remnants = true + +gaussian_spread1 = 0.185% +gaussian_spread2 = 0.185% + +?isr_handler = true +$isr_handler_mode = "recoil" +isr_alpha = 0.0072993 +isr_mass = 0.000511 + +# Golden-card value (production-quality precision, slow to integrate live): +# integrate (mumuH) { iterations = 10:100000:"gw", 5:200000:"" } +# Classroom default (faster, lower precision): +integrate (mumuH) { iterations = 3:2000:"gw" } + +n_events = 1000 + +$lhef_version = "3.0" +sample_format = lhef +simulate (mumuH) { $sample = "mumuH" } diff --git a/Gen/ee/mumuH_Hbb.cmd b/Gen/ee/mumuH_Hbb.cmd new file mode 100644 index 0000000..b83cfa5 --- /dev/null +++ b/Gen/ee/mumuH_Hbb.cmd @@ -0,0 +1,7 @@ +! Read in the WHIZARD LHEf file +Beams:frameType = 4 +Beams:LHEF = mumuH.lhe + +! Force H -> b b +25:onMode = off +25:onIfAny = 5 diff --git a/Gen/ee/pythia_mumuH.py b/Gen/ee/pythia_mumuH.py new file mode 100644 index 0000000..40377e2 --- /dev/null +++ b/Gen/ee/pythia_mumuH.py @@ -0,0 +1,28 @@ +from Gaudi.Configuration import * +from edm4hep import labels as e4_labels + +from Configurables import EventDataSvc +from k4FWCore import ApplicationMgr, IOSvc +ApplicationMgr().EvtSel = 'NONE' +ApplicationMgr().EvtMax = 1000 +ApplicationMgr().ExtSvc += ["RndmGenSvc", EventDataSvc("EventDataSvc")] + +from Configurables import PythiaInterface +pythia8gentool = PythiaInterface() +pythia8gentool.pythiacard = "mumuH_Hbb.cmd" + +from Configurables import GenAlg +pythia8gen = GenAlg("Pythia8") +pythia8gen.SignalProvider = pythia8gentool +pythia8gen.hepmc.Path = "hepmc" +ApplicationMgr().TopAlg += [pythia8gen] + +from Configurables import HepMCToEDMConverter +hepmc_converter = HepMCToEDMConverter() +hepmc_converter.hepmc.Path = "hepmc" +hepmc_converter.hepmcStatusList = [] +hepmc_converter.GenParticles.Path = e4_labels.MCParticles +ApplicationMgr().TopAlg += [hepmc_converter] + +iosvc = IOSvc() +iosvc.Output = "mumuH_Hbb.root" From 27e93f0427144c36ebdfd52d8b5f39faf949e765 Mon Sep 17 00:00:00 2001 From: Juraj Smiesko Date: Wed, 19 Aug 2026 14:50:44 +0200 Subject: [PATCH 2/6] Fix ISR/FSR double-counting and vertex smearing in Gen/ee Pythia8 step mumuH_Hbb.cmd now bases its Beams:/PartonLevel settings on FCC-config's actual LHE-reading card (p8_ee_default.cmd) instead of an unrelated production card, turning off Pythia8's own ISR/FSR since WHIZARD already applies it. Vertex/time smearing moves to the Gaudi GaussSmearVertex tool in pythia_mumuH.py to avoid double-smearing against Pythia8's own Beams:allowVertexSpread. Also documents the golden card's PYTHIA6 settings (Higgs mass/width, fragmentation tune, Bose-Einstein correlations, long-lived particle stability) that still need porting, as open TODOs. --- Gen/ee/README.md | 66 +++++++++++++++++++++++++++++++++++++++--- Gen/ee/mumuH_Hbb.cmd | 15 ++++++++++ Gen/ee/pythia_mumuH.py | 13 +++++++++ 3 files changed, 90 insertions(+), 4 deletions(-) diff --git a/Gen/ee/README.md b/Gen/ee/README.md index 830c8c4..d213286 100644 --- a/Gen/ee/README.md +++ b/Gen/ee/README.md @@ -112,14 +112,28 @@ FCCAnalyses. `PythiaInterface` reads a `.cmd` card, which can point at an external LHE file exactly like the working example shipped in `k4Gen` itself ([`data/Pythia_LHEinput.cmd`](https://github.com/key4hep/k4Gen/blob/main/k4Gen/data/Pythia_LHEinput.cmd)). -[`mumuH_Hbb.cmd`](mumuH_Hbb.cmd) in this directory combines that LHE-reading -pattern with the H -> b b decay-forcing settings from -[`FCC-config`'s `p8_ee_H_Hbb_ecm125.cmd`](https://github.com/HEP-FCC/FCC-config/blob/main/FCCee/Generator/Pythia8/p8_ee_H_Hbb_ecm125.cmd): +[`mumuH_Hbb.cmd`](mumuH_Hbb.cmd) in this directory is based on +[`FCC-config`'s `p8_ee_default.cmd`](https://github.com/HEP-FCC/FCC-config/blob/winter2023/FCCee/Generator/Pythia8/p8_ee_default.cmd) +(`winter2023` branch) — the actual card FCC-config uses to read WHIZARD LHE +output into Pythia8 — with H -> b b decay-forcing added on top. Note +`PartonLevel:ISR/FSR = off`: WHIZARD's `isr_handler` in Step 1 already +generated the initial-state radiation, so leaving Pythia8's own parton-level +ISR/FSR on would double-count it: ``` ! Read in the WHIZARD LHEf file Beams:frameType = 4 Beams:LHEF = mumuH.lhe +Beams:setProductionScalesFromLHEF = off +Beams:allowMomentumSpread = off + +! ISR/FSR already handled by WHIZARD's isr_handler in Step 1 - keep off here +! to avoid double-counting radiation in Pythia8's parton-level shower +PartonLevel:ISR = off +PartonLevel:FSR = off + +Check:epTolErr = 1e-1 +LesHouches:matchInOut = off ! Force H -> b b 25:onMode = off @@ -131,10 +145,15 @@ own [`options/pythia.py`](https://github.com/key4hep/k4Gen/blob/main/k4Gen/options/pythia.py) example — swap in the LHE-reading card above and write out EDM4hep rather than plain HepMC (`HepMCFileWriter`'s own docstring says it's for debugging, -not for actual event storage): +not for actual event storage). Beamspot vertex/time smearing is done here, +via the Gaudi `GaussSmearVertex` tool wired into `GenAlg`, using the same +FCC-ee IDEA beamspot values as `p8_ee_default.cmd`'s `Beams:sigmaVertex*` +settings — not via Pythia8's own `Beams:allowVertexSpread`, which would +apply it a second time on top of this: ```python from Gaudi.Configuration import * +from GaudiKernel import SystemOfUnits as units from edm4hep import labels as e4_labels from Configurables import EventDataSvc @@ -143,6 +162,17 @@ ApplicationMgr().EvtSel = 'NONE' ApplicationMgr().EvtMax = 1000 ApplicationMgr().ExtSvc += ["RndmGenSvc", EventDataSvc("EventDataSvc")] +# Beamspot vertex/time smearing (FCC-ee IDEA values, from FCC-config's +# p8_ee_default.cmd Beams:sigmaVertex{X,Y,Z}/sigmaTime). Done here via the +# Gaudi VertexSmearingTool rather than Pythia8's own Beams:allowVertexSpread, +# so it isn't applied twice. +from Configurables import GaussSmearVertex +smeartool = GaussSmearVertex() +smeartool.xVertexSigma = 5.96e-3 * units.mm +smeartool.yVertexSigma = 23.8e-6 * units.mm +smeartool.zVertexSigma = 0.397 * units.mm +smeartool.tVertexSigma = 10.89 * units.mm + from Configurables import PythiaInterface pythia8gentool = PythiaInterface() pythia8gentool.pythiacard = "mumuH_Hbb.cmd" @@ -150,6 +180,7 @@ pythia8gentool.pythiacard = "mumuH_Hbb.cmd" from Configurables import GenAlg pythia8gen = GenAlg("Pythia8") pythia8gen.SignalProvider = pythia8gentool +pythia8gen.VertexSmearingTool = smeartool pythia8gen.hepmc.Path = "hepmc" ApplicationMgr().TopAlg += [pythia8gen] @@ -186,6 +217,33 @@ hadronized, H -> b b decayed event record. > the `MCParticles` collection from `mumuH_Hbb.root` produced here — the > Gen/Sim split as designed is correct. +## Open TODOs + +The golden card's internal PYTHIA6 shower/hadronization step +(`$ps_PYTHIA_PYGIVE` in +[`wzp6_ee_mumuH_Hbb_ecm240.sin`](https://github.com/HEP-FCC/FCC-config/blob/winter2023/FCCee/Generator/Whizard/v3.0.3/wzp6_ee_mumuH_Hbb_ecm240.sin)) +sets several physics parameters that have **not** been ported to +`mumuH_Hbb.cmd`, because PYTHIA6 parameter names don't map mechanically onto +Pythia8 settings and doing this properly needs someone to work out the +correct Pythia8-native equivalents (or confirm Pythia8 defaults are close +enough for a teaching sample): + +- **Higgs mass and width** — the golden card sets `PMAS(25,1)=125.` and + `PMAS(25,2)=0.4143E-02` (4.143 MeV) explicitly. `mumuH_Hbb.cmd` doesn't + set `25:m0` / `25:mWidth`, so Pythia8's own defaults apply instead. +- **Hadronization / fragmentation tune** — the golden card carries a full + set of Lund string parameters (`PARJ(1,2,3,4,11-17,21,41,42,54,55)`, + `MSTJ(11)`, `MSTP(3)`). None of this has been translated into a Pythia8 + tune; Pythia8 defaults are used instead. +- **Bose-Einstein correlations** — the golden card turns these on with a + specific tune (`MSTP(151)=1`, `PARP(151-154)`). Not enabled in + `mumuH_Hbb.cmd` (off by Pythia8 default). +- **Long-lived particle stability treatment** — the golden card sets + `MSTJ(22)=4` with `PARJ(73)=2250`, `PARJ(74)=2500`, controlling which + particles get left stable (for the detector to handle) based on decay + length. Not addressed in `mumuH_Hbb.cmd`; this could matter for how + Delphes sees long-lived particles like K_S/Lambda downstream. + ## What's next The showered, H -> b b decayed sample (`mumuH_Hbb.root`) is the input to diff --git a/Gen/ee/mumuH_Hbb.cmd b/Gen/ee/mumuH_Hbb.cmd index b83cfa5..ca8ba35 100644 --- a/Gen/ee/mumuH_Hbb.cmd +++ b/Gen/ee/mumuH_Hbb.cmd @@ -1,6 +1,21 @@ +! Based on FCC-config's p8_ee_default.cmd (winter2023 campaign), the actual +! card used to read WHIZARD LHE output into Pythia8, plus H -> b b +! decay-forcing on top. Vertex/time smearing is done separately in the +! Gaudi steering script (pythia_mumuH.py), not here. + ! Read in the WHIZARD LHEf file Beams:frameType = 4 Beams:LHEF = mumuH.lhe +Beams:setProductionScalesFromLHEF = off +Beams:allowMomentumSpread = off + +! ISR/FSR already handled by WHIZARD's isr_handler in Step 1 - keep off here +! to avoid double-counting radiation in Pythia8's parton-level shower +PartonLevel:ISR = off +PartonLevel:FSR = off + +Check:epTolErr = 1e-1 +LesHouches:matchInOut = off ! Force H -> b b 25:onMode = off diff --git a/Gen/ee/pythia_mumuH.py b/Gen/ee/pythia_mumuH.py index 40377e2..d6cf205 100644 --- a/Gen/ee/pythia_mumuH.py +++ b/Gen/ee/pythia_mumuH.py @@ -1,4 +1,5 @@ from Gaudi.Configuration import * +from GaudiKernel import SystemOfUnits as units from edm4hep import labels as e4_labels from Configurables import EventDataSvc @@ -7,6 +8,17 @@ ApplicationMgr().EvtMax = 1000 ApplicationMgr().ExtSvc += ["RndmGenSvc", EventDataSvc("EventDataSvc")] +# Beamspot vertex/time smearing (FCC-ee IDEA values, from FCC-config's +# p8_ee_default.cmd Beams:sigmaVertex{X,Y,Z}/sigmaTime). Done here via the +# Gaudi VertexSmearingTool rather than Pythia8's own Beams:allowVertexSpread, +# so it isn't applied twice. +from Configurables import GaussSmearVertex +smeartool = GaussSmearVertex() +smeartool.xVertexSigma = 5.96e-3 * units.mm +smeartool.yVertexSigma = 23.8e-6 * units.mm +smeartool.zVertexSigma = 0.397 * units.mm +smeartool.tVertexSigma = 10.89 * units.mm + from Configurables import PythiaInterface pythia8gentool = PythiaInterface() pythia8gentool.pythiacard = "mumuH_Hbb.cmd" @@ -14,6 +26,7 @@ from Configurables import GenAlg pythia8gen = GenAlg("Pythia8") pythia8gen.SignalProvider = pythia8gentool +pythia8gen.VertexSmearingTool = smeartool pythia8gen.hepmc.Path = "hepmc" ApplicationMgr().TopAlg += [pythia8gen] From 6de0fe51208a701356faa3cb30149a34ad1a4107 Mon Sep 17 00:00:00 2001 From: Juraj Smiesko Date: Wed, 19 Aug 2026 17:07:56 +0200 Subject: [PATCH 3/6] Test Gen/ee pipeline end-to-end, recommend VSCodium, add tooling TODOs Ran WHIZARD + Pythia8 (k4Gen) against the real Key4hep stack: WHIZARD step works and gives a cross section matching the known sigma(ee->ZH)*BR(Z->mumu) value, but the Pythia8 step crashes inside Pythia8's own resonanceShowers/iBotCopyId machinery when decaying and showering the bare Higgs resonance from the LHE file. Documents the crash, backtrace, and ruled-out fixes as a blocking Open TODO. Also recommends VSCodium over VSCode in the top-level README, and adds Open TODOs for background-sample generation, a jupytext-based Jupyter export, and a custom TextMate grammar for WHIZARD/Pythia8 syntax highlighting. --- Gen/ee/README.md | 57 ++++++++++++++++++++++++++++++++++++++++++++++++ README.md | 6 ++--- 2 files changed, 60 insertions(+), 3 deletions(-) diff --git a/Gen/ee/README.md b/Gen/ee/README.md index d213286..165e683 100644 --- a/Gen/ee/README.md +++ b/Gen/ee/README.md @@ -205,6 +205,35 @@ k4run pythia_mumuH.py This produces `mumuH_Hbb.root`, an EDM4hep file with the showered, hadronized, H -> b b decayed event record. +> **Status: Step 2 currently crashes — do not treat as working yet.** +> Actually ran the full chain against the real Key4hep stack +> (`/cvmfs/sw.hsf.org/key4hep/setup.sh`, WHIZARD 3.1.5, Pythia8 8.315). +> Step 1 works: `whizard mumuH.sin` runs cleanly and gives +> sigma ~6.8 fb, matching the known sigma(ee->ZH)*BR(Z->mumu) ~ 200 fb * +> 3.37% ~ 6.7 fb at 240 GeV. Step 2 (`k4run pythia_mumuH.py`) crashes on +> the first event with: +> ``` +> Pythia8 FATAL Standard std::exception is caught in sysExecute +> Pythia8 ERROR vector::_M_range_check: __n (which is 9) >= this->size() (which is 9) +> ``` +> A `gdb` backtrace (`catch throw`) traces this to +> `Pythia8::PartonLevel::resonanceShowers()` -> +> `Pythia8::Particle::iBotCopyId()`, i.e. inside Pythia8's own machinery for +> decaying+showering a resonance (the Higgs) found mid-event, thrown while +> walking a daughter-index chain that runs past the end of the +> not-yet-fully-built event record. +> +> Ruled out (still crashes identically with each of these changed, alone or +> combined): `PartonLevel:FSR` on vs off, `LesHouches:matchInOut` on vs off, +> `PartonLevel:earlyResDec = on`, WHIZARD's `?keep_beams` true vs false, +> LHEF version 2.0 vs 3.0, and stripping WHIZARD's `` +> `sqme_prc` block from the LHE file. The crash appears to be specifically +> about Pythia8 decaying+showering a *bare, undecayed resonance* (H) that +> arrives alone in an LHE event with no other colored partons — needs +> either a working combination of Pythia8 settings we haven't found yet, a +> different Pythia8 version, or input from someone with deeper Pythia8 +> internals experience. + > **Confirmed handoff to `Sim/ee`:** this stage hands off an EDM4hep ROOT > file with an `MCParticles` collection, not a HepMC file and not a combined > Pythia8+Delphes run. Checked `key4hep/k4SimDelphes`: its Gaudi component @@ -219,6 +248,34 @@ hadronized, H -> b b decayed event record. ## Open TODOs +- **Blocking: Step 2 crashes when actually run** — see the "Status" callout + above. `k4run pythia_mumuH.py` throws inside Pythia8's own + `resonanceShowers`/`iBotCopyId` machinery on the first event; several + plausible fixes were tried and ruled out. Needs resolving before this + tutorial step can be called working. +- **Background samples** — undecided whether this tutorial includes any + background processes alongside the mumuH signal, and if so how they'd be + generated. One tentative idea floated: generate backgrounds with pure + Pythia8 (no WHIZARD step), since Pythia8 alone can produce e.g. generic + qqbar/WW/ZZ final states without needing WHIZARD's matrix-element + machinery. Not decided or attempted. +- **Jupyter notebook export** — Jupyter itself was ruled out as the primary + authoring format (too complicated with FCCAnalyses, per planning notes), + but a one-way export of this markdown material to `.ipynb` (e.g. via + [`jupytext`](https://jupytext.readthedocs.io/)) could still be useful for + students who'd rather work in a notebook. Not attempted — would need some + markup convention to mark which fenced code blocks are meant to be + executable Python cells versus illustrative shell/Sindarin/Pythia8-card + snippets, since jupytext doesn't know the difference on its own. +- **Syntax highlighting for WHIZARD/Pythia8 files** — neither WHIZARD's + Sindarin (`.sin`) format nor Pythia8's `.cmd` cards have a grammar in + GitHub's Linguist (so no fenced-code-block language tag lights them up on + GitHub) or an existing VSCodium/VSCode extension. Writing a small custom + TextMate grammar for one or both (packaged as a minimal VSCodium + extension, or bundled in this repo) would fix local editing at least; + GitHub rendering would still fall back to a closest-fit generic tag (e.g. + `ini`-ish for the Pythia8 cards) or plain text. Not started. + The golden card's internal PYTHIA6 shower/hadronization step (`$ps_PYTHIA_PYGIVE` in [`wzp6_ee_mumuH_Hbb_ecm240.sin`](https://github.com/HEP-FCC/FCC-config/blob/winter2023/FCCee/Generator/Whizard/v3.0.3/wzp6_ee_mumuH_Hbb_ecm240.sin)) diff --git a/README.md b/README.md index 1e65e59..63d2d24 100644 --- a/README.md +++ b/README.md @@ -36,9 +36,9 @@ Analysis/ ## Format Each section is introduced with brief slides, then students work through the -markdown material directly in this repo (VSCode + extensions recommended for an -all-in-one setup). Students present their solutions; reference solutions are -also provided here for offline use. +markdown material directly in this repo (VSCodium + extensions recommended +for an all-in-one setup). Students present their solutions; reference +solutions are also provided here for offline use. ## Branches From 61877ade8ef480ee06765391ed1c2e4cc599f471 Mon Sep 17 00:00:00 2001 From: Juraj Smiesko Date: Thu, 20 Aug 2026 16:25:48 +0200 Subject: [PATCH 4/6] Fix ?keep_beams crash and correct ISR/FSR handling in Gen/ee The Step 2 (Pythia8) crash was caused by WHIZARD's ?keep_beams = true, which FCC-config's own (unheeded) maintainer comment on a sibling card already warns about. Fixed by no longer setting ?keep_beams at all in mumuH.sin, relying on its default (false). Verified with a standalone, zero-Gaudi Pythia8 reproducer and a full 1000-event k4run pipeline run. Also corrects PartonLevel:FSR in mumuH_Hbb.cmd: it must stay on (Pythia8 default) since it's what showers the H -> b b decay products before hadronization, not just a source of double-counted radiation as previously assumed - verified by dumping the Pythia8 event record with it on and off. PartonLevel:ISR stays off, confirmed correct the same way. Also removes explicit defaults (?keep_beams, ?keep_remnants) from mumuH.sin and documents mH, ISR/FSR, and the ISR-photon-visibility caveat in Gen/ee/README.md. --- Gen/ee/README.md | 148 +++++++++++++++++++++++++++++++++---------- Gen/ee/mumuH.sin | 2 - Gen/ee/mumuH_Hbb.cmd | 7 +- 3 files changed, 117 insertions(+), 40 deletions(-) diff --git a/Gen/ee/README.md b/Gen/ee/README.md index 165e683..7cc8eec 100644 --- a/Gen/ee/README.md +++ b/Gen/ee/README.md @@ -60,8 +60,6 @@ sqrts = 240 GeV process mumuH = e1, E1 => e2, E2, H beams = e1, E1 => gaussian => isr -?keep_beams = true -?keep_remnants = true gaussian_spread1 = 0.185% gaussian_spread2 = 0.185% @@ -93,6 +91,26 @@ whizard mumuH.sin This produces `mumuH.lhe`. +Two settings are deliberately left out rather than pinned explicitly, +relying on their defaults: + +- **`mH` (Higgs mass)** — not set, because WHIZARD's `SM` model already + defaults it to 125 GeV (`parameter mH = 125` in `SM.mdl`), matching the + golden card's own `PMAS(25,1)=125.` and confirmed empirically: generated + LHE events carry the Higgs at exactly 125.0 GeV. Note the sibling + `wzp6_ee_mumuH_ecm240.sin` card instead sets `mH = 125.1 GeV` explicitly + — a minor inconsistency between FCC-config's own cards, not something + reconciled here (see the Higgs mass/width item in Open TODOs). +- **`?keep_beams` / `?keep_remnants`** — not set, because their defaults + (`false` and `true` respectively) are exactly what's needed here. + `?keep_beams = true` is what caused the Step 2 crash (see the Status + callout below) — it writes the original beam particles into the LHE + record as extra entries, which WHIZARD's own manual explicitly warns + against for reading into PYTHIA. Leaving it at its default `false` avoids + that. `?keep_remnants` only has any effect when `?keep_beams = true` + (per the manual), so with `?keep_beams` at its default it's inert either + way and isn't worth setting. + ### Why mu mu H, not Z H? WHIZARD is asked to generate e+e- -> mu+ mu- H directly, rather than @@ -115,10 +133,7 @@ external LHE file exactly like the working example shipped in `k4Gen` itself [`mumuH_Hbb.cmd`](mumuH_Hbb.cmd) in this directory is based on [`FCC-config`'s `p8_ee_default.cmd`](https://github.com/HEP-FCC/FCC-config/blob/winter2023/FCCee/Generator/Pythia8/p8_ee_default.cmd) (`winter2023` branch) — the actual card FCC-config uses to read WHIZARD LHE -output into Pythia8 — with H -> b b decay-forcing added on top. Note -`PartonLevel:ISR/FSR = off`: WHIZARD's `isr_handler` in Step 1 already -generated the initial-state radiation, so leaving Pythia8's own parton-level -ISR/FSR on would double-count it: +output into Pythia8 — with H -> b b decay-forcing added on top: ``` ! Read in the WHIZARD LHEf file @@ -127,10 +142,11 @@ Beams:LHEF = mumuH.lhe Beams:setProductionScalesFromLHEF = off Beams:allowMomentumSpread = off -! ISR/FSR already handled by WHIZARD's isr_handler in Step 1 - keep off here -! to avoid double-counting radiation in Pythia8's parton-level shower +! ISR already handled by WHIZARD's isr_handler in Step 1 - keep off here to +! avoid double-counting radiation. FSR stays on (Pythia8 default): it's what +! actually showers the H -> b b decay products before hadronization - with +! it off, b/bbar go straight into string fragmentation with zero shower. PartonLevel:ISR = off -PartonLevel:FSR = off Check:epTolErr = 1e-1 LesHouches:matchInOut = off @@ -140,6 +156,50 @@ LesHouches:matchInOut = off 25:onIfAny = 5 ``` +`p8_ee_default.cmd` itself sets `PartonLevel:FSR = off` too — appropriate +*there* because its use case is WHIZARD LHE files with the full final state, +including any colored partons, already present at the matrix-element level +(nothing left for Pythia8 to decay+shower). That's not our case: our LHE +only has e+e- -> mu+ mu- H with the Higgs undecayed, and we rely on Pythia8 +itself to decay H -> b b and shower the result. Checked this directly by +dumping the Pythia8 event record: with `FSR = off`, the b/bbar from the +Higgs decay went straight into hadronization with no shower emissions at +all (unrealistic jets); with `FSR` left on (Pythia8's default), the event +record shows proper gluon-emission branching off the b/bbar before +hadronization, as expected. + +`ISR = off` was checked the same way, and unlike FSR this one didn't reveal +a hidden dependency: with `ISR = on`, Pythia8 adds a genuine extra +initial-state radiation photon (visible in the event record, status `-43`, +radiated directly off one of the incoming leptons by Pythia8's own spacelike +shower) on top of the energy redistribution WHIZARD's `isr_handler` already +applied to the beam momenta — real double-counted radiation. With +`ISR = off`, that extra photon doesn't appear and nothing else changes; ISR +only concerns the initial state, so there's no FSR-style gating issue here. + +> **Caveat: the ISR photon itself is invisible downstream, as a side +> effect of the `?keep_beams` crash fix, not a deliberate physics choice.** +> WHIZARD's `?isr_handler` (`$isr_handler_mode = "recoil"`) still generates +> a real ISR photon per beam internally and correctly recoils the visible +> system against it — checked directly on generated events: the mu+ mu- H +> system's total energy varies event-by-event (~235.5-239.7 GeV instead of +> a fixed 240) and sometimes carries non-zero net transverse momentum +> (e.g. 0.14 GeV), exactly the signature of a real recoil against a missing +> photon. But the WHIZARD manual classifies these radiated ISR photons as +> "beam remnants" (`?keep_remnants` docs: *"for ISR and/or beamstrahlung +> spectra, the radiated photons are considered as beam remnants"*), and +> remnants are gated by the same `?keep_beams` flag that had to be set to +> `false` to fix the Step 2 crash. So the photon itself is never written +> into the event record — confirmed structurally too: the old +> `?keep_beams = true` events had explicit outgoing photon lines; the +> fixed `?keep_beams = false` events don't. For this tutorial's actual +> measurements (mu mu recoil mass, H -> b b dijet mass) this doesn't +> matter, since both only depend on the visible mu/mu/b/bbar kinematics, +> which already correctly reflect the ISR-induced energy loss and recoil. +> It does mean there's no possibility of ever seeing a reconstructed ISR +> photon downstream in `Sim/ee`/`Analysis` — that option is gone as an +> unavoidable consequence of the crash fix, not a deliberate simplification. + Steering script [`pythia_mumuH.py`](pythia_mumuH.py), adapted from `k4Gen`'s own [`options/pythia.py`](https://github.com/key4hep/k4Gen/blob/main/k4Gen/options/pythia.py) @@ -205,34 +265,57 @@ k4run pythia_mumuH.py This produces `mumuH_Hbb.root`, an EDM4hep file with the showered, hadronized, H -> b b decayed event record. -> **Status: Step 2 currently crashes — do not treat as working yet.** -> Actually ran the full chain against the real Key4hep stack -> (`/cvmfs/sw.hsf.org/key4hep/setup.sh`, WHIZARD 3.1.5, Pythia8 8.315). -> Step 1 works: `whizard mumuH.sin` runs cleanly and gives +> **Status: fixed — Step 2 now runs cleanly.** Actually ran the full chain +> against the real Key4hep stack (`/cvmfs/sw.hsf.org/key4hep/setup.sh`, +> WHIZARD 3.1.5, Pythia8 8.315). Step 1: `whizard mumuH.sin` gives > sigma ~6.8 fb, matching the known sigma(ee->ZH)*BR(Z->mumu) ~ 200 fb * -> 3.37% ~ 6.7 fb at 240 GeV. Step 2 (`k4run pythia_mumuH.py`) crashes on -> the first event with: +> 3.37% ~ 6.7 fb at 240 GeV. Step 2 (`k4run pythia_mumuH.py`) initially +> crashed on the first event with: > ``` > Pythia8 FATAL Standard std::exception is caught in sysExecute > Pythia8 ERROR vector::_M_range_check: __n (which is 9) >= this->size() (which is 9) > ``` -> A `gdb` backtrace (`catch throw`) traces this to +> A `gdb` backtrace (`catch throw`) traced this to > `Pythia8::PartonLevel::resonanceShowers()` -> -> `Pythia8::Particle::iBotCopyId()`, i.e. inside Pythia8's own machinery for -> decaying+showering a resonance (the Higgs) found mid-event, thrown while -> walking a daughter-index chain that runs past the end of the -> not-yet-fully-built event record. +> `Pythia8::Particle::iBotCopyId()`, inside Pythia8's own machinery for +> decaying+showering a resonance (the Higgs) found mid-event. A zero-Gaudi +> standalone C++ reproducer (`Pythia pythia; pythia.readFile(...); +> pythia.init(); pythia.next();`) hit the identical crash, ruling out the +> 2026-08-19 meeting's "needs to be run as Gaudi functionals" lead — this +> was never a Gaudi/k4Gen issue. +> +> **Root cause, found via a warning comment on a sibling FCC-config card** +> ([`wzp6_ee_mumuH_ecm240.sin`](https://github.com/HEP-FCC/FCC-config/blob/winter2023/FCCee/Generator/Whizard/v3.0.3/wzp6_ee_mumuH_ecm240.sin) +> has `?keep_beams = true # do not use this option, makes Pythia crash`, +> left in place but unheeded in that card since it only ever exercises the +> internal PYTHIA6 path): our `mumuH.sin` had `?keep_beams = true` too. A +> clean retest (the very first attempt to test this flag turned out to be +> unreliable, from a contaminated debugging script, and wrongly seemed to +> rule it out) confirms `?keep_beams = false` fixes it completely — a full +> 1000-event run now completes with a single non-fatal +> "energy-momentum not quite conserved" warning and produces a valid +> `mumuH_Hbb.root` (`MCParticles` collection, 102 particles in the first +> event, fully hadronized). `mumuH.sin` no longer sets `?keep_beams` at +> all, relying on its default of `false` (see the Environment/Step 1 +> section above). +> +> For the record, before finding the real cause, these were also tried and +> made no difference on their own (the actual fix was always `?keep_beams`): +> `PartonLevel:FSR` on vs off, `LesHouches:matchInOut` on vs off, +> `PartonLevel:earlyResDec = on`, LHEF version 2.0 vs 3.0, stripping +> WHIZARD's `` `sqme_prc` block, forcing H -> b b vs leaving the +> decay unforced, and forcing H completely stable (`25:mayDecay = off`). > -> Ruled out (still crashes identically with each of these changed, alone or -> combined): `PartonLevel:FSR` on vs off, `LesHouches:matchInOut` on vs off, -> `PartonLevel:earlyResDec = on`, WHIZARD's `?keep_beams` true vs false, -> LHEF version 2.0 vs 3.0, and stripping WHIZARD's `` -> `sqme_prc` block from the LHE file. The crash appears to be specifically -> about Pythia8 decaying+showering a *bare, undecayed resonance* (H) that -> arrives alone in an LHE event with no other colored partons — needs -> either a working combination of Pythia8 settings we haven't found yet, a -> different Pythia8 version, or input from someone with deeper Pythia8 -> internals experience. +> One earlier claim in this callout was wrong and is corrected here: an +> initial test seemed to show the crash was specific to particle ID 25 +> (Higgs) as the LHE resonance, based on a WHIZARD LHE with a Z boson in +> H's place reading in cleanly — but that Z-boson test card never had +> `?keep_beams = true` set in the first place, so it wasn't a like-for-like +> comparison. Redone properly (identical card, only H swapped for Z, +> `?keep_beams = true` kept in both): the Z version crashes identically to +> the H version. So the crash is generic to `?keep_beams = true` plus any +> downstream resonance decay, exactly matching the (unheeded) FCC-config +> warning comment, not something Higgs-specific. > **Confirmed handoff to `Sim/ee`:** this stage hands off an EDM4hep ROOT > file with an `MCParticles` collection, not a HepMC file and not a combined @@ -248,11 +331,6 @@ hadronized, H -> b b decayed event record. ## Open TODOs -- **Blocking: Step 2 crashes when actually run** — see the "Status" callout - above. `k4run pythia_mumuH.py` throws inside Pythia8's own - `resonanceShowers`/`iBotCopyId` machinery on the first event; several - plausible fixes were tried and ruled out. Needs resolving before this - tutorial step can be called working. - **Background samples** — undecided whether this tutorial includes any background processes alongside the mumuH signal, and if so how they'd be generated. One tentative idea floated: generate backgrounds with pure diff --git a/Gen/ee/mumuH.sin b/Gen/ee/mumuH.sin index af102f9..0580d83 100644 --- a/Gen/ee/mumuH.sin +++ b/Gen/ee/mumuH.sin @@ -6,8 +6,6 @@ sqrts = 240 GeV process mumuH = e1, E1 => e2, E2, H beams = e1, E1 => gaussian => isr -?keep_beams = true -?keep_remnants = true gaussian_spread1 = 0.185% gaussian_spread2 = 0.185% diff --git a/Gen/ee/mumuH_Hbb.cmd b/Gen/ee/mumuH_Hbb.cmd index ca8ba35..1ecc7bd 100644 --- a/Gen/ee/mumuH_Hbb.cmd +++ b/Gen/ee/mumuH_Hbb.cmd @@ -9,10 +9,11 @@ Beams:LHEF = mumuH.lhe Beams:setProductionScalesFromLHEF = off Beams:allowMomentumSpread = off -! ISR/FSR already handled by WHIZARD's isr_handler in Step 1 - keep off here -! to avoid double-counting radiation in Pythia8's parton-level shower +! ISR already handled by WHIZARD's isr_handler in Step 1 - keep off here to +! avoid double-counting radiation. FSR stays on (Pythia8 default): it's what +! actually showers the H -> b b decay products before hadronization - with +! it off, b/bbar go straight into string fragmentation with zero shower. PartonLevel:ISR = off -PartonLevel:FSR = off Check:epTolErr = 1e-1 LesHouches:matchInOut = off From 7d28eecee1fd0ae5b31e10840091327cb7b83e94 Mon Sep 17 00:00:00 2001 From: Juraj Smiesko Date: Fri, 21 Aug 2026 08:34:31 +0200 Subject: [PATCH 5/6] Add ErrorMax to Gen/ee Pythia8 step for large-scale runs A small fraction of events (~1 in a few thousand) hit an unrecoverable Pythia8-level failure that retrying doesn't fix. Without raising Gaudi's default per-algorithm ErrorMax = 1, a single such event aborts the entire run - only shows up at scale (a 1000-event test never hits it), confirmed by generating and processing 10,000 events. ErrorMax = 20 lets a handful be skipped instead, verified with a clean 10,000/10,000-event run. --- Gen/ee/README.md | 19 +++++++++++++++++++ Gen/ee/pythia_mumuH.py | 7 +++++++ 2 files changed, 26 insertions(+) diff --git a/Gen/ee/README.md b/Gen/ee/README.md index 7cc8eec..4a1a7d1 100644 --- a/Gen/ee/README.md +++ b/Gen/ee/README.md @@ -242,6 +242,13 @@ pythia8gen = GenAlg("Pythia8") pythia8gen.SignalProvider = pythia8gentool pythia8gen.VertexSmearingTool = smeartool pythia8gen.hepmc.Path = "hepmc" +# A small fraction of events (~1 in a few thousand) hit an unrecoverable +# Pythia8-level failure (e.g. energy-momentum conservation check) that +# retrying doesn't fix. Without raising this, Gaudi's default ErrorMax = 1 +# means a single such event aborts the whole run - only shows up at scale +# (not in a 1000-event test). ErrorMax = 20 lets a handful be skipped +# instead. Deprecated on some Gaudi versions but still functional. +pythia8gen.ErrorMax = 20 ApplicationMgr().TopAlg += [pythia8gen] from Configurables import HepMCToEDMConverter @@ -317,6 +324,18 @@ hadronized, H -> b b decayed event record. > downstream resonance decay, exactly matching the (unheeded) FCC-config > warning comment, not something Higgs-specific. +> **Tested at 10,000-event scale, found and fixed a second, scale-dependent +> issue.** A 1,000-event run only ever showed one non-fatal warning, but at +> 10,000 events a handful (~3 in 10,000) hit a Pythia8-level failure +> (energy-momentum check) that retrying doesn't recover from — it's the +> same event failing deterministically every time, not a transient issue. +> Without any fix, Gaudi's default per-algorithm `ErrorMax = 1` means even +> one such event aborts the *entire* run (job exits with an error, most +> events lost) — this never shows up at 1,000 events, only at scale. Fixed +> by setting `pythia8gen.ErrorMax = 20` in `pythia_mumuH.py`, so a handful +> of individually-unrecoverable events get skipped rather than ending the +> job. Verified: a clean 10,000/10,000-event run with this set. + > **Confirmed handoff to `Sim/ee`:** this stage hands off an EDM4hep ROOT > file with an `MCParticles` collection, not a HepMC file and not a combined > Pythia8+Delphes run. Checked `key4hep/k4SimDelphes`: its Gaudi component diff --git a/Gen/ee/pythia_mumuH.py b/Gen/ee/pythia_mumuH.py index d6cf205..58e29cc 100644 --- a/Gen/ee/pythia_mumuH.py +++ b/Gen/ee/pythia_mumuH.py @@ -28,6 +28,13 @@ pythia8gen.SignalProvider = pythia8gentool pythia8gen.VertexSmearingTool = smeartool pythia8gen.hepmc.Path = "hepmc" +# A small fraction of events (~1 in a few thousand) hit an unrecoverable +# Pythia8-level failure (e.g. energy-momentum conservation check) that +# retrying doesn't fix. Without raising this, Gaudi's default ErrorMax = 1 +# means a single such event aborts the whole run - only shows up at scale +# (not in a 1000-event test). ErrorMax = 20 lets a handful be skipped +# instead. Deprecated on some Gaudi versions but still functional. +pythia8gen.ErrorMax = 20 ApplicationMgr().TopAlg += [pythia8gen] from Configurables import HepMCToEDMConverter From 211714f99a59ed2d1a6f050d2c68be4029b4b21b Mon Sep 17 00:00:00 2001 From: Juraj Smiesko Date: Fri, 21 Aug 2026 15:36:05 +0200 Subject: [PATCH 6/6] Resolve unported PYTHIA6 tune flags: match features, not numbers For each PYTHIA6 setting the golden card had enabled, checked whether Pythia8's own default already covers it rather than porting the tuned numeric values (risky - the two codes don't always define equivalent- sounding parameters identically). Most items needed no change: Bowler fragmentation for b/c quarks is already Pythia8's default, and Higgs mass/width already closely match. Bose-Einstein correlations were checked and deliberately left off - they don't affect the mu mu recoil mass and only marginally affect the H -> b b dijet mass, not worth the added complexity here. The long-lived-particle setting went through a wrong attempt (mapped to a proper-lifetime cutoff) before checking the actual PYTHIA 6.4 manual, which shows MSTJ(22)=4 is a geometric decay-vertex-position cutoff instead - Pythia8's default already matches its practical outcome, confirmed on real generated events (K_S0/Lambda decay behavior), so no setting is needed there either. --- Gen/ee/README.md | 105 ++++++++++++++++++++++++++++++++++--------- Gen/ee/mumuH_Hbb.cmd | 8 ++++ 2 files changed, 92 insertions(+), 21 deletions(-) diff --git a/Gen/ee/README.md b/Gen/ee/README.md index 4a1a7d1..a4f3e5f 100644 --- a/Gen/ee/README.md +++ b/Gen/ee/README.md @@ -154,6 +154,9 @@ LesHouches:matchInOut = off ! Force H -> b b 25:onMode = off 25:onIfAny = 5 + +! No extra long-lived-particle or Bose-Einstein settings here - see the +! Open TODOs section below for why. ``` `p8_ee_default.cmd` itself sets `PartonLevel:FSR = off` too — appropriate @@ -376,27 +379,87 @@ hadronized, H -> b b decayed event record. The golden card's internal PYTHIA6 shower/hadronization step (`$ps_PYTHIA_PYGIVE` in [`wzp6_ee_mumuH_Hbb_ecm240.sin`](https://github.com/HEP-FCC/FCC-config/blob/winter2023/FCCee/Generator/Whizard/v3.0.3/wzp6_ee_mumuH_Hbb_ecm240.sin)) -sets several physics parameters that have **not** been ported to -`mumuH_Hbb.cmd`, because PYTHIA6 parameter names don't map mechanically onto -Pythia8 settings and doing this properly needs someone to work out the -correct Pythia8-native equivalents (or confirm Pythia8 defaults are close -enough for a teaching sample): - -- **Higgs mass and width** — the golden card sets `PMAS(25,1)=125.` and - `PMAS(25,2)=0.4143E-02` (4.143 MeV) explicitly. `mumuH_Hbb.cmd` doesn't - set `25:m0` / `25:mWidth`, so Pythia8's own defaults apply instead. -- **Hadronization / fragmentation tune** — the golden card carries a full - set of Lund string parameters (`PARJ(1,2,3,4,11-17,21,41,42,54,55)`, - `MSTJ(11)`, `MSTP(3)`). None of this has been translated into a Pythia8 - tune; Pythia8 defaults are used instead. -- **Bose-Einstein correlations** — the golden card turns these on with a - specific tune (`MSTP(151)=1`, `PARP(151-154)`). Not enabled in - `mumuH_Hbb.cmd` (off by Pythia8 default). -- **Long-lived particle stability treatment** — the golden card sets - `MSTJ(22)=4` with `PARJ(73)=2250`, `PARJ(74)=2500`, controlling which - particles get left stable (for the detector to handle) based on decay - length. Not addressed in `mumuH_Hbb.cmd`; this could matter for how - Delphes sees long-lived particles like K_S/Lambda downstream. +sets several physics parameters. Rather than trying to port PYTHIA6's exact +tuned numbers onto Pythia8 (risky — PYTHIA6 and Pythia8 don't always define +equivalent-sounding parameters the same way, and getting this wrong would be +worse than just using Pythia8's own tune), the resolution here was to match +which *features* the golden card had enabled, using Pythia8's own default +values for how each behaves. Checked against Pythia8's manual (version- +matched to the actual stack) parameter by parameter: + +- **Bose-Einstein correlations** (`MSTP(151)=1`, `PARP(151-154)`) — + deliberately left off (Pythia8 default), after checking whether it + actually matters for this tutorial's two measurements: it doesn't touch + muons at all (only identical-boson pairs like pi/K), so the mu mu + recoil mass is completely unaffected; for the H -> b b dijet mass it + would only be a small within-jet momentum redistribution, since the + algorithm is explicitly designed to conserve overall jet 4-momentum + (per Pythia8's own manual). Checked it doesn't introduce a higher event + failure rate either (1000/1000 events succeed with it on, same as off). + Not worth the added complexity for what it buys here — if this + changes (e.g. this tutorial ever adds a jet-substructure exercise), it + can be turned on via `HadronLevel:BoseEinstein = on`, leaving + `BoseEinstein:lambda`/`QRef` at Pythia8 defaults (PYTHIA6's + `PARP(151-154)` is a 4-parameter tune and Pythia8's is only 3, so the + golden card's specific numbers were never going to be portable anyway). +- **Long-lived particle stability** (`MSTJ(22)=4`, `PARJ(73)=2250`, + `PARJ(74)=2500`) — **left off (Pythia8 default), corrected after an + earlier wrong attempt at this.** First pass mapped this to + `ParticleDecays:limitTau0 = on` (a pure proper-lifetime/ctau cutoff), + reasoning it was a lifetime threshold like PYTHIA6's `PARJ(71)`-based + `MSTJ(22)=2` option. That was the wrong PYTHIA6 option: checked the + actual PYTHIA 6.4 manual, and `MSTJ(22)=4` is a **geometric** cutoff — + "a particle is decayed only if the decay vertex is within a cylindrical + volume with radius `PARJ(73)` ... and extent to `±PARJ(74)`" — i.e. + decay-vertex *position* (momentum-dependent), not proper lifetime at + all. The correct Pythia8 equivalent is `ParticleDecays:limitCylinder` + (`xyMax`/`zMax`), not `limitTau0`. Checked the practical difference on + real generated events (500-event sample): with the wrong `limitTau0 = + on` setting, 100% of K_S0/Lambda came out *stable*; with it left off + (Pythia8 default), 100% come out *decayed* — matching what the golden + card's actual cylinder cutoff does in practice, since at these energies + essentially every K_S0/Lambda decays well inside the golden card's + ~2.25m/2.5m volume anyway. So the fix was to remove the setting + entirely, not port a different one — Pythia8's own default already + matches the golden card's practical behavior for this tutorial's + particle content. +- **Fragmentation function for b/c quarks** (`MSTJ(11)=3`, PYTHIA6's + "Bowler" option) — **no change needed**: Pythia8's manual states + outright that "for massive quarks..., the Bowler modification to the + Lund FF is the default choice." So this golden-card setting is already + Pythia8's default behavior with zero configuration. +- **Higgs mass and width** (`PMAS(25,1)=125.`, `PMAS(25,2)=0.4143E-02`, + i.e. 4.143 MeV) — no change: Pythia8's own default (125.0 GeV, + 4.08 MeV width) already matches closely. +- **Lund `a`/`b`/`sigma` and diquark/meson-multiplet tune** + (`PARJ(1,2,3,4,11-17,21,41,42)`, `MSTP(3)`) — no change: these are + baseline numeric tuning knobs Pythia8 always applies some value for + (`StringZ:aLund`/`bLund`, `StringPT:sigma`, `StringFlav:...`), not + on/off features to toggle, so under this approach they're left at + Pythia8's own defaults rather than guessing at cross-code equivalence. + As a sanity check on how different the actual numbers are: the golden + card's `PARJ(41)/(42)` (Lund a/b = 0.11/0.52) versus Pythia8's defaults + (`aLund`/`bLund` = 0.68/0.98) differ substantially — a real tune + difference, not a rounding one, which is exactly why these weren't + guessed at. +- **`MSTJ(28)=0`** (disables PYTHIA6's own tau decay, deferring to an + external tool like TAUOLA for correct tau spin/polarization + correlations) — checked, and it's not as inapplicable as first assumed: + taus do appear in this chain, not from the H -> b b decay itself, but + from semitauonic B-hadron decays (B -> D(*) tau nu, a real ~2-3%-per-B + branching fraction) after the b/bbar hadronize - confirmed empirically, + ~8.5% of a 200-event sample had at least one tau. TAUOLA's value is + correct tau polarization for analyses where tau decay kinematics *are* + the observable (e.g. Z/H -> tau tau); here the taus are secondary, + buried inside b-jets, and nothing in this tutorial's analysis (mu mu + recoil mass, H -> b b dijet mass) is sensitive to tau polarization, so + Pythia8's own native tau decay treatment should be adequate without + TAUOLA. If more precise B-hadron decay modeling ever mattered, there's + an available upgrade path already in this pipeline: + `PythiaInterface`'s `doEvtGenDecays` option (currently unset/off, + see `pythia_mumuH.py`) would route B-hadron decays through EvtGen + instead of Pythia8's simpler built-in table, improving tau kinematics + as a side effect too. Not enabled here. ## What's next diff --git a/Gen/ee/mumuH_Hbb.cmd b/Gen/ee/mumuH_Hbb.cmd index 1ecc7bd..0b5bc8d 100644 --- a/Gen/ee/mumuH_Hbb.cmd +++ b/Gen/ee/mumuH_Hbb.cmd @@ -21,3 +21,11 @@ LesHouches:matchInOut = off ! Force H -> b b 25:onMode = off 25:onIfAny = 5 + +! No extra long-lived-particle or Bose-Einstein settings here - see +! Gen/ee/README.md's Open TODOs for why: the golden card's MSTJ(22)=4 is +! a geometric (detector-cylinder) decay-vertex cutoff, not a proper- +! lifetime one, and Pythia8's default (off) already matches its practical +! outcome for K_S0/Lambda at these energies. Bose-Einstein correlations +! (golden card's MSTP(151)=1) deliberately left off too - see the same +! section for why.