Document Geant4 scoring and its FLUKA counterpart

docs/scoring/README.md

@@ -0,0 +1,190 @@
+---
+sort: 9
+title: Geant4 scoring
+---
+
+# Geant4 scoring
+
+Geant4 can accumulate radiation maps — fluence, dose, 1 MeV neutron equivalent
+fluence, hadron fluence above 20 MeV — on a mesh laid over the ALICE geometry,
+without any change to detector code. The mesh is configured with Geant4 macro
+commands, the same way FLUKA users configure `USRBIN` cards, and the result is
+written as a text file at the end of the run.
+
+This page explains how to switch it on in `o2-sim`, how to write the macro, and
+how to read the output. The [FLUKA counterpart](flukacomparison.md) page gives
+the card-by-card translation and the differences you have to know about before
+comparing numbers.
+
+## Quick start
+
+Write a macro file `scoring.in`:
+
+```
+/score/create/cylinderMesh fine
+/score/mesh/cylinderSize 10. 30. cm
+/score/mesh/nBin 100 60 1
+/score/quantity/cellFlux flux percm2
+/score/quantity/doseDeposit dose Gy
+/score/close
+```
+
+Append it to the standard ALICE Geant4 configuration and run:
+
+```bash
+cat $O2_ROOT/share/Detectors/gconfig/g4config.in scoring.in > g4scoring.in
+
+o2-sim-serial -e TGeant4 -g pythia8pp -n 100 \
+  --configKeyValues "G4.g4scoring=true;G4.configMacroFile=$PWD/g4scoring.in"
+```
+
+The run writes `fine.worker<pid>.txt` into the working directory, one file per
+mesh.
+
+`G4.configMacroFile` **replaces** the standard `g4config.in`, it does not add to
+it. Always concatenate, otherwise the run loses the ALICE physics settings
+(range cuts, PAI in the TRD, biasing, field steppers).
+
+## Meshes
+
+Each mesh is a block between `/score/create/...` and `/score/close`. Several
+meshes can coexist; they are parallel worlds and do not disturb the transport.
+
+| Command | Meaning |
+|---|---|
+| `/score/create/cylinderMesh <name>` | cylindrical mesh, axis along z |
+| `/score/create/boxMesh <name>` | box mesh |
+| `/score/mesh/cylinderSize <R> <dz> <unit>` | outer radius and **half** length |
+| `/score/mesh/boxSize <dx> <dy> <dz> <unit>` | **half** lengths |
+| `/score/mesh/nBin <nR> <nZ> <nPhi>` | bins of a cylinder mesh |
+| `/score/mesh/nBin <nX> <nY> <nZ>` | bins of a box mesh |
+| `/score/mesh/translate/xyz <x> <y> <z> <unit>` | move the mesh off the origin |
+| `/score/close` | end the mesh definition |
+
+A mesh is centred on the origin unless it is translated. Bins are equidistant.
+A single phi bin gives the phi-averaged map, which is what a radiation study
+normally wants.
+
+## Quantities and filters
+
+Inside a mesh, every `/score/quantity/...` line adds one scorer, and the
+`/score/filter/...` lines that follow apply to the scorer above them.
+
+| Command | Quantity |
+|---|---|
+| `/score/quantity/cellFlux <name> percm2` | track length per cell volume, i.e. fluence |
+| `/score/quantity/doseDeposit <name> Gy` | absorbed dose |
+| `/score/quantity/energyDeposit <name> MeV` | deposited energy |
+| `/score/quantity/nOfStep <name>` | number of steps, useful for debugging |
+
+| Filter | Selects |
+|---|---|
+| `/score/filter/particle <fname> <p1> <p2> ...` | those Geant4 particle names |
+| `/score/filter/particleWithKineticEnergy <fname> <emin> <emax> <unit> <p1> ...` | those particles in that kinetic energy range |
+| `/score/filter/charged <fname>` | all charged particles |
+| `/score/filter/neutral <fname>` | all neutral particles |
+
+Particle names are the Geant4 ones: `neutron`, `proton`, `anti_proton`, `pi+`,
+`pi-`, `kaon0L`, `gamma`, `e-`, `opticalphoton`, and so on.
+
+A scorer without a filter counts everything, **including optical photons**.
+Those dominate the total fluence around the FT0 and are not transported by
+FLUKA at all, so an unfiltered `cellFlux` is not a quantity you can compare
+between the two engines. Filter, or read the filtered scorers instead.
+
+## Output
+
+Each mesh is dumped at the end of the run to `<mesh>.worker<pid>.txt`:
+
+```
+# mesh name: fine
+# primitive scorer name: flux
+# iZ, iPHI, iR, total(value) [percm2], total(val^2), entry
+0,0,0,1856.089717,405107.423708,10
+```
+
+The columns are the bin indices, the sum over the whole run, the sum of
+squares, and the number of contributing steps. Two things to keep in mind:
+
+- The value is summed over **all events of the run**, not per event. Divide by
+  the number of events yourself.
+- `total(val^2)` is the sum of squares of the per-step contributions, so it is
+  not the statistical error of the bin. Use several runs with different seeds
+  if you need an uncertainty.
+
+Parallel `o2-sim` writes one file per worker and merges them into `<mesh>.txt`
+at the end of the run. If a merge has to be redone by hand, for example after
+collecting dumps from several jobs:
+
+```bash
+o2-sim-merge-g4scoring <directory> <expected number of workers>
+```
+
+The second argument is optional; when given, the merger fails loudly if a
+worker file is missing, which is the common way to lose a few percent of the
+statistics without noticing.
+
+## 1 MeV neutron equivalent fluence
+
+Geant4 has no built-in NIEL scorer. The ALICE Geant4 build takes a per-step
+weight from a user function, and `o2-sim` feeds it with RD50 damage weights, so
+that a `cellFlux` scorer becomes a damage-weighted fluence:
+
+```
+/score/quantity/cellFlux neq percm2 true
+```
+
+The trailing `true` is the `scoreweighted` flag. It only exists for `cellFlux`.
+
+Switch the weights on with:
+
+```bash
+--configKeyValues "G4.g4scoring=true;G4.g4fluenceweight=true;G4.fluenceWeightFile=$PWD/rd50_niel.csv;G4.configMacroFile=$PWD/g4scoring.in"
+```
+
+The weight file is a CSV of PDG code, kinetic energy in MeV and damage weight:
+
+```
+# pdg,ekin[MeV],weight
+2112,1.025000e-10,1.575000e-02
+```
+
+The same curves ship with O2 as ROOT graphs in
+`$O2_ROOT/share/Detectors/gconfig/data/rd50_niel.root`, under the names
+`neutronDW`, `protonDW`, `pionDW` and `electronDW`. If you do not have the text
+file, write it out from there — the loader only accepts CSV.
+
+Only four tables are read — neutron (2112), proton (2212), pion (211) and
+electron (11) — and they are applied like this:
+
+| Particle | Weight used |
+|---|---|
+| neutron | neutron table |
+| e+, e- | electron table, zero if the table is absent |
+| all other baryons and antibaryons | proton table |
+| all other mesons | pion table |
+| everything else | zero |
+
+Outside the tabulated energy range the weight is clamped to the first or last
+tabulated value.
+
+The consequence for a comparison: **do not filter the NIEL scorer down to
+neutrons, protons and pions.** FLUKA's `SI1MEVNE` weights every hadron and
+electron, so an ALICE inner-tracker map scored only on n, p and pi comes out
+about 18% low. Kaons alone account for most of it.
+
+## A ready example
+
+[`g4scoring_alice.in`](g4scoring_alice.in) is the macro used for the ALICE
+FLUKA/Geant4 radiation comparison: a fine mesh over the inner tracker and a
+coarse one over the cavern, each with fluence, 1 MeV n eq, hadrons above
+20 MeV, charged fluence, energy deposit and dose.
+
+```bash
+cat $O2_ROOT/share/Detectors/gconfig/g4config.in g4scoring_alice.in > g4full.in
+
+o2-sim -e TGeant4 -g pythia8pp -n 1000 -j 8 \
+  --configKeyValues "G4.g4scoring=true;G4.g4fluenceweight=true;G4.fluenceWeightFile=$PWD/rd50_niel.csv;G4.configMacroFile=$PWD/g4full.in"
+```
+
+{% include list.liquid all=true %}

docs/scoring/flukacomparison.md

@@ -0,0 +1,114 @@
+---
+sort: 1
+title: The same scoring in FLUKA
+---
+
+# The same scoring in FLUKA
+
+FLUKA scores the same maps with `USRBIN` cards. In `o2-sim` the cards come from
+a file:
+
+```bash
+o2-sim-serial -e TFluka -g pythia8pp -n 100 \
+  --configKeyValues "FlukaParam.scoringFile=$PWD/scoring.inp"
+```
+
+The file is appended to the FLUKA input that `o2-sim` generates, so it contains
+only the cards, in the fixed FLUKA column format:
+
+```
+*...+....1....+....2....+....3....+....4....+....5....+....6....+....7...+...8
+USRBIN           11.      201.       21.      10.0        0.      30.0fflux
+USRBIN            0.        0.     -30.0      100.        1.       60.&
+USRBIN           11.      236.       22.      10.0        0.      30.0fneq
+USRBIN            0.        0.     -30.0      100.        1.       60.&
+```
+
+Each pair of cards is one binning: type, quantity, output unit, `Rmax`, x of the
+axis, `zmax`, name; then `Rmin`, y of the axis, `zmin`, `nR`, `nPhi`, `nZ`.
+
+Two details worth spelling out:
+
+- Type `11` is an r-phi-z mesh scored with the track-length apportioning
+  algorithm. Type `1` uses the old point-wise algorithm and is **not allowed**
+  for `SI1MEVNE` and `HADGT20M`.
+- A positive output unit writes formatted text to `fort.<unit>`, a negative one
+  writes binary that has to go through `usbsuw`/`usbrea` first. Units below 21
+  are reserved. Text output is easier to read back and costs nothing at ALICE
+  mesh sizes.
+
+## Translation table
+
+| Quantity | FLUKA | Geant4 |
+|---|---|---|
+| all-particle fluence | `201` (ALL-PART) | `/score/quantity/cellFlux f percm2` |
+| charged fluence | `202` (ALL-CHAR) | `cellFlux` + `/score/filter/charged` |
+| neutral fluence | `203` (ALL-NEUT) | `cellFlux` + `/score/filter/neutral` |
+| neutron fluence | `8` (NEUTRON) | `cellFlux` + `/score/filter/particle f neutron` |
+| photon fluence | `7` (PHOTON) | `cellFlux` + `/score/filter/particle f gamma` |
+| deposited energy | `208` (ENERGY), GeV/cm3 | `/score/quantity/energyDeposit e MeV` |
+| dose | `228` (DOSE), GeV/g | `/score/quantity/doseDeposit d Gy` |
+| 1 MeV n eq (Si) | `236` (SI1MEVNE) | weighted `cellFlux`, see [Geant4 scoring](README.md) |
+| hadrons above 20 MeV | `237` (HADGT20M) | `cellFlux` + `particleWithKineticEnergy` filter |
+
+`HADGT20M` has no Geant4 equivalent either, so the particle list has to be
+written out:
+
+```
+/score/quantity/cellFlux had20 percm2
+/score/filter/particleWithKineticEnergy had20filter 20. 10000000. MeV proton anti_proton neutron anti_neutron pi+ pi- kaon+ kaon- kaon0L kaon0S lambda anti_lambda sigma+ sigma- xi- xi0 omega-
+```
+
+Filtering works the other way round in the two codes. In FLUKA a binning scores
+one generalised particle and `AUXSCORE` can restrict it further, but not by
+energy. In Geant4 a scorer starts from everything and is narrowed by filters,
+including energy filters.
+
+## Normalisation
+
+This is the most common source of a wrong ratio.
+
+- FLUKA normalises a `USRBIN` to **one unit of primary weight**. The text
+  output carries the total weight of the primaries in its header; multiply by
+  it to get the sum over the run.
+- Geant4 sums over the **whole run** and divides by nothing. Divide by the
+  number of events yourself.
+
+So per event: `fluka_value * total_weight / nevents` against
+`g4_value / nevents`.
+
+## Settings to match before comparing
+
+- **Same kinematics.** Generate once with `--noGeant`, then replay the same
+  file in both engines with `-g extkinO2 --extKinFile o2sim_Kine.root`.
+  Otherwise the comparison also contains the generator.
+- **Low-energy neutrons.** FLUKA transports neutrons down to thermal energies
+  by default. Geant4 does not, unless an HP physics list is used and the
+  neutron cut is lowered:
+  `G4.physicsmode=kFTFP_BERT_HP_optical;SimCutParams.lowneut=true;GlobalSimProcs.CUTNEU=5.e-12`.
+- **Optical photons.** Geant4 produces and transports Cherenkov photons in the
+  FT0; FLUKA in the ALICE setup does not. They show up in any unfiltered
+  fluence scorer and swamp it locally. Dose, 1 MeV n eq and hadron fluence are
+  unaffected.
+- **Same geometry.** `--detectorList` and `--skipModules` have to match, and
+  hit creation should be switched off or on in both.
+
+## Pitfall: black-hole regions in the FLUKA geometry
+
+FLUKA gets one region per ROOT volume number. Two volumes that share a *name*
+share a number, and only one of them gets its material assigned. If the winner
+is a `TGeoVolumeAssembly`, its placeholder material index `-1` becomes
+`BLCKHOLE`, and a real piece of detector silently deletes every particle that
+enters it.
+
+Check this in every FLUKA run directory before you trust a map:
+
+```bash
+grep -c "ASSIGNMAT       -1.0" flukaMat.inp   # 0 when healthy
+```
+
+This affected the ten MFT PEEK support disks in ALICE 2 and pulled the FLUKA
+fluence in the barrel down by up to a factor two — for a long time this was read
+as a physics disagreement between the two codes. The geometry is fixed by
+giving the support volume a name of its own; make sure your O2 version has that
+fix.

docs/scoring/g4scoring_alice.in

@@ -0,0 +1,33 @@
+# ALICE radiation scoring: two cylindrical meshes, phi averaged.
+# Append this to $O2_ROOT/share/Detectors/gconfig/g4config.in and pass the
+# result as G4.configMacroFile. The scorers mirror the FLUKA USRBIN cards
+# 201 (fluence), 236 (1 MeV n eq), 237 (hadrons > 20 MeV), 202 (charged)
+# and 208 (energy deposit).
+
+# inner tracker region: R < 10 cm, |z| < 30 cm, 1 mm x 1 cm bins
+/score/create/cylinderMesh fine
+/score/mesh/cylinderSize 10. 30. cm
+/score/mesh/nBin 100 60 1
+/score/quantity/cellFlux flux percm2
+/score/quantity/cellFlux neq percm2 true
+/score/quantity/cellFlux had20 percm2
+/score/filter/particleWithKineticEnergy had20filter 20. 10000000. MeV proton anti_proton neutron anti_neutron pi+ pi- kaon+ kaon- kaon0L kaon0S lambda anti_lambda sigma+ sigma- xi- xi0 omega-
+/score/quantity/cellFlux charged percm2
+/score/filter/charged chargedfilter
+/score/quantity/energyDeposit edep MeV
+/score/quantity/doseDeposit dose Gy
+/score/close
+
+# whole cavern: R < 5 m, |z| < 5 m, 2.5 cm x 5 cm bins
+/score/create/cylinderMesh coarse
+/score/mesh/cylinderSize 500. 500. cm
+/score/mesh/nBin 200 200 1
+/score/quantity/cellFlux flux percm2
+/score/quantity/cellFlux neq percm2 true
+/score/quantity/cellFlux had20 percm2
+/score/filter/particleWithKineticEnergy had20filter 20. 10000000. MeV proton anti_proton neutron anti_neutron pi+ pi- kaon+ kaon- kaon0L kaon0S lambda anti_lambda sigma+ sigma- xi- xi0 omega-
+/score/quantity/cellFlux charged percm2
+/score/filter/charged chargedfilter
+/score/quantity/energyDeposit edep MeV
+/score/quantity/doseDeposit dose Gy
+/score/close

docs/transport/engines.md

@@ -7,6 +7,8 @@ title: Transport engines
 
 ## GEANT4
 
+Radiation maps (fluence, dose, 1 MeV n eq) can be scored on a mesh, see [Geant4 scoring](../scoring/).
+
 ### Scaling hadronic cross sections
 
 Hadronic cross sections can be scaled by passing a specific Geant4 configuration to the transport simulation.
@@ -16,6 +18,8 @@ Hadronic cross sections can be scaled by passing a specific Geant4 configuration
 
 Informtation about FLUKA.
 
+For `USRBIN` scoring in `o2-sim` and how it maps to Geant4, see [the same scoring in FLUKA](../scoring/flukacomparison.md).
+
 ## GEANT3
 
 Information about GEANT3.