Experimental Typst web edition · Veto chapter pilot
3 Simulation and Analysis Framework
3.1 From incident particles to detector events
To understand the background measured by the detector, it is necessary to follow the particles from their source to the signals recorded by the electronics. The simulation therefore includes particle generation, transport through the shielding, detector response and event reconstruction (Figure 3.1). The reconstructed events can then be selected in the same way as experimental data. This chapter describes the tools used for this calculation and the developments made during this thesis. After following the event through signal formation and reconstruction, we define the observables and common time projection chamber (TPC) selection in Section 3.5. The comparison with calibration data then tests the reconstructed distributions and selection acceptance before these methods are applied to the veto and background studies.

Geant4 follows particles through the detector geometry, and the response model converts their interactions into channel waveforms. Simulated and measured waveforms are reconstructed with REST to obtain the TPC and veto observables used for event selection. The two TPC response models are shown in more detail in Figure 3.8.The calculation combines REST-for-Physics [86], Geant4 [87–89], ROOT [90], CRY [91] and the geometry tools developed within the IAXO collaboration. The main contributions are summarized in Table 3.1. They include a more efficient cosmic-ray generator, extensions to particle transport, a detailed Micromegas response and tools for large simulation campaigns and detector operation.
The measured signals described in Section 2 provide the experimental reference for the two TPC responses compared here. The developed Garfield++ response is introduced in Section 3.4, including the preliminary neutron extension in Section 3.4.3. Section 3.6 compares reconstructed calibration events with this model and the conventional energy-deposit response used in existing background productions.
| Stage | Contribution of this thesis |
|---|---|
| Particle injection | Published PDP method and reusable cosmic-source distributions. |
| Transport and geometry | Multithreading, track pruning, delayed subevents and integration of parameterized IAXO geometries. |
| Detector response | Garfield gas transport, resolved mesh field, individual electron collection and calibration with common REST reconstruction. |
| Production and inspection | Batch simulations, readout mapping and tools for detector operation and event inspection. |
3.2 REST-for-Physics
The comparison between simulation and data is simpler when both can be processed with the same tools. REST-for-Physics [86] (Rare Event Searches Toolkit for Physics) provides this common framework. Developed for rare-event experiments, it includes libraries for particle transport, detector response and data analysis, written in C++ and based on ROOT. It has already been used for the IAXO pathfinder detector at CAST, as described in the thesis of Cristina Margalejo Blasco [77].
Three programs cover the main tasks in this thesis. restG4 runs the Geant4 simulation, restManager applies a sequence of processing steps, and restRoot allows interactive inspection of the results. The geometry, source and processing steps are specified in RML configuration files, so different detector conditions can be studied without modifying the programs.
3.2.1 Data storage and analysis
Simulation and experimental data are stored in ROOT files [90], together with the detector and processing settings needed to interpret them. The tree structure of these files allows selected observables to be read efficiently, which is useful when comparing large samples or testing different selection criteria. Most figures in this thesis were produced in Python with Matplotlib or Plotly. The data were accessed through PyROOT, which also gives access to the REST classes, or Uproot [92], which reads ROOT files directly in Python.
3.2.2 Events and reconstruction

An event contains the information associated with one particle history or one acquisition trigger. As it passes through the processing chain, its description changes: particle interactions become electronic signals, and signals become reconstructed hits and tracks (Figure 3.2). REST provides corresponding event classes, such as TRestGeant4Event, TRestRawSignalEvent and TRestTrackEvent, derived from the common TRestEvent class.
The results are stored in two complementary trees. The EventTree contains the event at the chosen processing stage, while the AnalysisTree contains quantities such as energy, position, width and multiplicity. These observables can be used to draw control plots or change a selection without repeating the reconstruction. Earlier stages, such as waveforms or transport hits, are saved separately when they are needed for later reprocessing.
In the TPC, the reconstructed shape helps distinguish compact X-ray events from background. In the veto, pulse amplitudes, times and panel multiplicities identify activity associated with a TPC event. For both detectors, simulated waveforms are processed with the same algorithms and settings as measured waveforms. The original particle history remains available to study the interactions responsible for an event.
The TPC baseline selection is defined in Section 3.5.3 and illustrated in Figure 3.15. Section 3.6 tests the resulting observables and acceptance with calibration data. Section 5 and Section 4 then apply the TPC and veto selections to the background studies.
3.2.3 Readout geometry and channel mapping
A recorded signal must be associated with a position in the detector. In the Micromegas readout, connected pads form strips; in the veto, each channel corresponds to a scintillator panel. The readout description relates these elements to their positions, orientations and electronics identifiers, and is shared by simulation and data (Figure 3.3).
We developed readout-generation tools for the IAXO-D0 and IAXO-D1 Micromegas detectors and their veto systems. For the TPC, the description includes the shapes of individual pads, including those at the edges. This is needed to assign arriving electrons to pads in the detailed response model and to reconstruct positions from strip signals. For the veto, the mapping includes the orientation and sensitive depth of each panel.
We checked the mapping by comparing the sensitive volumes, acquisition channels and reconstructed positions. These checks help identify geometry or channel-assignment errors.
![Figure 3.3: Readout description of the Micromegas and veto detectors in REST [ 86 ]. The upper panel relates detector elements to the TRestDetectorReadout channel map. The lower panel shows how the same description is used to assign simulated charge to channels and reconstruct measured signals.](assets/cd4ab862ad157b6b462d7a55.png)
TRestDetectorReadout channel map. The lower panel shows how the same description is used to assign simulated charge to channels and reconstruct measured signals.3.3 Particle generation and radiation transport
3.3.1 Cosmic-ray source generators
Cosmic rays produce showers of secondary particles in the atmosphere. Their flux at the laboratory can be calculated with programs such as CRY [91] or CORSIKA [93]. The detector simulation then follows these particles through the shielding and active volumes. The source tools described here make this second calculation more efficient by reusing atmospheric distributions and directing particles toward the detector.
We developed a program to obtain cosmic-ray secondary distributions from CRY for the chosen site conditions [94]. For each species—muons, electrons, positrons, photons, protons and neutrons—it stores the joint kinetic-energy and zenith-angle distribution, , in a ROOT histogram. These histograms are then sampled by restG4. The particles are generated independently. This is a reasonable approximation for a compact detector: atmospheric muons and neutrons can be spread over kilometer scales, and coincidence measurements with meter-scale detectors found no significant neutron–muon correlation [95]. Coincidences between separate atmospheric particles are therefore omitted, while secondaries created by the same incident particle in the shielding are still followed together by Geant4.
Calculating the atmospheric distributions once makes it possible to reuse them for different shield thicknesses and veto geometries. They can also be compared with EXPACS calculations [96, 97] or HENSA measurements [98] before running the detector simulation. The background studies use the Guan distribution [99] for muons, CRY for photons, protons, electrons and positrons, and the measured HENSA spectrum [100] for the nominal outdoor-neutron energy distribution. CRY provides an alternative source for muons and neutrons. The HENSA energy spectrum is combined with an assumed angular distribution proportional to . The conversion from plane to spherical generation also affects the sampled angles; the distributions produced by the complete generator are examined in Section 5.4.1.

CRY supplies energy–angle histograms for each particle species; the HENSA neutron spectrum and Guan muon distribution provide alternative inputs. These sources are combined with the detector geometry before transport and response simulation. Arrows connect source preparation, configuration and detector simulation.3.3.2 Probability Distribution Projection cosmic-ray generator
In conventional cosmic-ray simulations, particles are generated over a large plane above the detector. Most trajectories miss the geometry, and their transport contributes little to the background calculation. For a compact detector such as IAXO-D0, this can account for a large fraction of the computing time.
The Probability Distribution Projection (PDP) method was developed and published during this thesis to reduce this cost [101]. It generates only trajectories that cross a sphere enclosing the geometry. For a given direction, their starting points lie within the ellipse obtained by projecting the sphere onto a horizontal plane (Figure 3.5). If is the distribution of horizontal-plane crossings per unit zenith angle, the distribution for the sphere is
(3.1)
The factor accounts for the different projected areas. To obtain the simulated exposure, the generation rate must be calculated using the same convention. For intensity per unit energy, solid angle, time and area perpendicular to the direction of incidence, the rates through a horizontal area and a sphere of radius are
(3.2)
The integrals cover the chosen incident hemisphere and energy range; for azimuthal symmetry, . An input already defined per unit area perpendicular to the ray does not require the plane-to-sphere correction again.

Geant4 transport at the blue entry points. (c) The corresponding starting positions on the horizontal plane form an ellipse. In practice, positions are sampled on an equivalent disk perpendicular to the incident direction.![Figure 3.6: Summary of the IAXO-D0 benchmark reported in [ 101 ]. Left: saved events per second of computation for three conventional generation-disk radii and PDP. Center: the corresponding physical event rate. Right: relative rate difference from PDP, Δ PDP = 100 ( 𝑁 ̇ / 𝑁 ̇ PDP − 1 ) , in percent. The first two panels use the statistical errors quoted in Table 2; the third propagates the MC and PDP rate errors assuming independent samples. PDP defines the zero reference, and its shared uncertainty correlates the three MC differences.](assets/4250ebfa947d2dad1110990b.png)
The sampled zenith angles and distances from the detector axis agree with conventional generation within statistical uncertainties when the generation plane is sufficiently large. The IAXO-D0 benchmark also shows why this convergence matters (Figure 3.6). Small conventional generation disks are fast, but underestimate the physical event rate. Relative to PDP, the deficits are about 6.4%, 1.0% and 0.5% as the disk radius increases. The analytical projection and the trajectory-distribution checks support using PDP as the reference; the zero in the difference panel follows from this choice. PDP gives about three times more saved events per unit computing time than the fastest of these configurations, and about 37 times more than conventional generation at comparable physical accuracy. This improvement made the large cosmic-ray samples used in the following chapters more practical.
3.3.3 Detector geometry
Comparing shielding and veto designs requires geometries that can be changed easily while keeping the detector dimensions and materials consistent. The IAXO-D0 and IAXO-D1 geometries were therefore built from parameterized descriptions and exported to the Geometry Description Markup Language (GDML) [102], which can be read by Geant4 and ROOT.
Some veto geometries contain more than 400 components. The collaboration’s Kotlin tool, gdml.kt, constructs these geometries from reusable elements rather than requiring each volume to be written by hand. During this thesis, we contributed to its development and testing, and integrated the generated geometries with particle transport, readout mapping and event visualization.
The geometry files are versioned, so the configuration used in each simulation can be recovered. The restG4 interface was also extended to resolve volume names changed during GDML assembly parsing and to accept logical-volume names. This allows transport settings to be applied consistently to repeated components. The same detector description is used for the channel map and event viewer: a panel with a simulated energy deposit can be highlighted and checked against its reconstructed signal (Figure 3.7).

3.3.4 Particle transport with restG4
Geant4 [87–89] calculates particle interactions and transport through matter. REST supplies the application, restG4, and the Geant4Lib library needed to configure the calculation and store the resulting histories. Several extensions were developed during this thesis to reduce computing and storage costs and to study delayed activity.
3.3.4.1 Particle histories and delayed activity
Each simulated event begins with one or more primary particles. Geant4 follows their trajectories in steps and also tracks the secondary particles produced by their interactions. The stored history contains the positions, momenta, times, energies and interaction processes needed to understand how energy reaches the active detector volumes.
Some secondaries produce signals long after the initial interaction. A stacking action was added to identify new tracks and separate long-lived activity into delayed subevents. Their times and parent identifiers are preserved, so prompt and delayed signals can be studied separately while remaining associated with the same primary history.
3.3.4.2 Multithreading and track pruning
Independent particle histories can be processed in parallel. Multithreading was added to restG4 to use several processor cores, together with synchronized output and clean interruption of a running simulation.
Track pruning reduces the amount of stored information. Only tracks related to the chosen sensitive volumes are kept, including the trajectories through passive material that led to those interactions. This reduced file sizes by factors of 1–10 in the tested configurations while preserving the checked TPC and veto observables. The retained information depends on the response to be calculated: scintillator quenching may need step lengths and recoil identities, while delayed coincidences need times and parent associations. These quantities must therefore be included when preparing a sample for such studies.
| Transport extension | Purpose |
|---|---|
| Multithreading | Process independent particle histories in parallel. |
| Track pruning | Reduce file sizes while keeping trajectories related to the sensitive volumes. |
| Delayed subevents | Separate late activity while keeping its relation to the parent history. |
| Geometry-name resolution | Apply volume settings consistently to repeated detector components. |
3.4 TPC detector response and event reconstruction
We developed a Garfield++ detector response that follows primary electrons through drift, amplification and individual collection before forming channel waveforms [72]. This is the response adopted for subsequent simulations; Section 3.4.3 describes the neutron conversion intended for the final response and marks its current development status. The REST detector response used in the existing background productions provides a comparison. Both responses enter the same reconstruction as measured waveforms (Figure 3.8); the physical processes underlying the signals are described in Section 2.

Magboltz supplies the gas-transport coefficients. Measured acquisitions calibrate the effective gain and noise and provide the pulse shape. Simulated and measured waveforms enter the same REST reconstruction. The REST detector response is included for comparison.3.4.1 Gas-transport tables and response models
The gas mixture, pressure, temperature and electric field determine the transport coefficients used by both responses. We developed gas-cli to prepare Magboltz tables for these conditions [103]. Figure 3.9 compares the drift velocity and diffusion coefficients as the isobutane fraction changes.
![Figure 3.9: Dependence of Micromegas gas transport on the isobutane concentration, calculated with Garfield++ / Magboltz gas tables at 20 ∘ C [ 72 ]. The upper row shows argon–isobutane mixtures, and the lower row shows equal-partial-pressure xenon–neon mixtures with isobutane. The horizontal axis is the reduced drift field 𝐸 𝑑 / 𝑝 , with pressure in bar. Dotted lines mark representative detector-response conditions: the midpoint of the approximately 133 – 153 V cm − 1 argon drift-field range at 𝑝 = 1.4 bar , and 𝐸 𝑑 ≃ 111 V cm − 1 for xenon–neon at 𝑝 = 1.05 bar .](assets/7c287742b07b06e3c565d80a.png)
Garfield++/Magboltz gas tables at [72]. The upper row shows argon–isobutane mixtures, and the lower row shows equal-partial-pressure xenon–neon mixtures with isobutane. The horizontal axis is the reduced drift field , with pressure in bar. Dotted lines mark representative detector-response conditions: the midpoint of the approximately argon drift-field range at , and for xenon–neon at .Increasing the isobutane fraction generally reduces transverse diffusion in these scans, giving a narrower cloud after the same drift distance. Changes in drift velocity and longitudinal diffusion also affect pulse timing. The dotted lines identify reference fields, rather than an optimized operating range or the conditions of every calibration run. The argon–isobutane curve illustrates the transport study; the conservative IAXO-D1 background reference in Section 5 uses isobutane. The xenon–neon–isobutane mixture belongs to the separate BabyIAXO response studies.
The Garfield response uses Heed to generate primary ionization and then calculates drift, multiplication and collection of individual electrons. In the REST response, Geant4 supplies gas energy deposits, which REST converts into primary charge with statistical fluctuations. Drift is described through transport coefficients and amplification through effective gain parameters before the charge is assigned to channels and converted into sampled pulses. Its effective response is calibrated against measured signals, making this calculation practical for large background samples.
3.4.2 Ionization and electron drift
Geant4 transports the incident particle through the chamber materials and entrance window. For photons and the charged particles supported by the present Heed interface, the position, direction, energy and time are recorded at first gas entry, and Geant4 transport is stopped before any gas interaction. Heed then calculates the interaction and ionization [104], avoiding a second calculation of the same gas interaction. Each ionization site supplies primary electrons with their positions and creation times. Neutrons and other particles outside this interface require the treatment described in Section 3.4.3.
Bulk drift samples diffusion and attachment through macroscopic transport coefficients, without following every electron collision. The microscopic mesh calculation starts from each electron’s drift endpoint and accumulated time. These calculations describe successive parts of the electron’s path.
3.4.3 Neutrons and primary-electron generation
The present Heed interface supplies primary electrons for photons, electrons, positrons and muons directly to Garfield [104]. Their number is not sampled again from the deposited energy. Neutrons and other particles outside this interface must remain in Geant4 inside the gas, so their collisions, recoils, captures and secondary particles are calculated there. An additional conversion from gas deposits to primary electrons connects these interactions to the same Garfield drift, amplification and waveform calculation (Figure 3.10).
3.4.3.0.1 TODO / PRELIMINARY.
Complete and validate the neutron deposit-to-electron interface. The description below specifies the method intended for the final version of this thesis; the current response comparisons do not test this neutron extension.

Heed calculation. The dashed branch identifies the neutron deposit-to-electron interface, marked preliminary in this draft. Both routes supply electron positions and creation times to the same drift, amplification and waveform calculation.The deposit-based interface applies the primary-charge yield and fluctuation relations in Equation (3.8) and Equation (3.9) to the electron-equivalent energy . For a recoil stopped in the gas,
(3.3)
where the ionization quenching factor is defined in Section 2.1.3 and depends on recoil species, energy and gas composition [69]. The interface samples an integer electron number according to the charge-yield and fluctuation model and places the electrons along the recorded gas steps, retaining their creation times. As an existing example, the local REST conversion used in the timing study starts from a Poisson draw (the approximation), replacing a zero count by one electron for a positive deposit. Other REST configurations specify Fano fluctuations; the sampling prescription must therefore be identified for each response. The mean pair-creation energy determines the charge yield, rather than the initial kinetic energy assigned to each electron.
3.4.4 Focusing and amplification at the Micromegas
The local field describes the nominal pattern of circular mesh holes above a diamond readout pad. Figure 3.11 shows the field around the holes and through the amplification gap. For the illustrated Ar/isobutane 99/1 mixture at , the drift field is , while the amplification gap has a nominal field of , over 400 times stronger.

The microscopic calculation follows multiplication and losses by attachment or to the conductors, retaining each electron’s position relative to the holes. The mean collected gain includes injections that produce no collected charge. Figure 3.12 follows one event from chamber-scale drift to microscopic amplification. The calculation uses an ideal periodic geometry; the complete chamber field and variations in manufactured dimensions remain to be studied.

3.4.5 Charge collection and channel waveforms
Each electron reaching the anode is assigned to the pad containing its landing position, and connected pads contribute to the same channel. This retains sharing between strips when an avalanche spans several pads, instead of assigning all its charge to the centroid. Electrons lost before collection or landing outside the active readout do not contribute to channel charge.
The Python waveform calculation adds pulses delayed by the collected electrons’ arrival times, samples the result and converts it into integer analog-to-digital converter (ADC) values. Measured acquisitions provide the shared pulse template, channel delays and residual noise; the development energy peak fixes the conversion from collected charge to ADC amplitude. The effective template may include gas and illumination time spread, while induced currents and the complete electronic circuit are not calculated separately. Figure 3.13 illustrates the resulting waveforms and reconstruction.

3.4.6 Common event reconstruction
For the measured calibration run used in the response comparison, REST removes slow baseline variation with a 75-sample moving average, limiting the influence of large excursions. It then estimates the common fluctuation at each time sample from the median of the channels and subtracts its variation about the baseline. TPC and veto channels are treated separately, and the median reduces the influence of pulses confined to a few channels.
The saved measured waveforms have already undergone these corrections. We generate simulated waveforms at the same stage, using the collected-electron response, empirical pulse and residual noise obtained from calibration data. Both samples enter the same REST reconstruction with the same readout mapping and settings, without repeating the common-mode correction. The remaining baseline and noise are estimated from samples 20–199, approximately at the sampling period, before detecting pulses above noise and fitting the accepted signals. The channel fixes each strip hit’s transverse position, the fitted peak time its drift coordinate and the amplitude its weight. Summed hit amplitudes supply track energies; the two main tracks provide the calibrated event energy used for selection (Section 3.5.3), rather than waveform integrals or deposited truth energy. Weak channels can be lost through thresholding or unsuccessful fits.
Hits are clustered independently in the XZ and YZ projections. All reconstructed tracks contribute to multiplicity, while the track with the largest summed amplitude in each view supplies the main-track observables: the weighted center and widths, energy, and charge and hit-count balance. These projections do not give a unique three-dimensional track. The longitudinal width reflects pulse times after shaping and fitting as well as gas diffusion.
Section 3.5 defines these observables and the fixed baseline, and Section 3.6 compares their distributions and selection acceptance with calibration data. No correction is applied to reconstructed observables after waveform reconstruction; the additional gain calibration in Section 3.6.3 acts on the signal before this stage.
The IAXO-shaped cathode aperture provides a simple check of the reconstructed transverse positions. Figure 3.14 compares measured calibration events with a REST detector response simulation using a similar source position. The IAXO lettering is clearly visible in both samples.

3.5 Event observables and multivariate selection
The reconstruction provides the quantities needed to identify X-ray candidates and compare measured events with simulations. Energy, position and track multiplicity define the initial candidate population; the reconstructed charge distribution then provides further discrimination. This section introduces these observables, fixes the common TPC baseline and describes how a multivariate selection combines the shape information. The following section tests whether the simulated response reproduces both the observables and the fraction of calibration events passing that baseline.
3.5.1 Energy, position and track multiplicity
The most energetic track in each strip view supplies a charge estimate and a position in that projection. The two charge estimates are combined and converted to energy using the accepted calibration for the run or simulated sample. Their amplitude-weighted positions give the transverse track center. These quantities determine whether the candidate lies in the energy window and central detector region used for the background comparison.
The number of reconstructed tracks describes whether the hits form one connected track or several in each view. A localized X-ray conversion usually forms one track per view, while extended interactions and separated secondary deposits can form additional tracks. Track multiplicity is therefore useful for selection, although it does not directly count the physical particles that interacted. Section 3.5.3 specifies the exact requirements on these reconstructed quantities.
3.5.2 TPC observables used by the BDT
The TPC boosted decision tree (BDT) uses fifteen quantities derived from the main track in each view, listed in Table 3.3. Pulse amplitudes weight the hit positions when calculating the track center and widths. The transverse widths describe the spread of charge over strips, while the longitudinal width is obtained from the fitted pulse times and drift velocity. It therefore includes the effects of signal formation and reconstruction as well as electron diffusion. Hit counts describe how many reconstructed strip hits contribute to a track, and the balances compare the charge and shape reconstructed in the two views. Skewness measures departures from a symmetric charge distribution.
| Property | BDT inputs |
|---|---|
| Track widths | , , , and |
| Width balance | |
| Charge balance | Main-track fraction Energy balance |
| Track asymmetry | Absolute transverse and longitudinal skewness, and |
| Strip hits | , , their minimum and maximum, and |
Calibrated energy, radius and number of tracks enter the baseline selection rather than the BDT. Particle identity and simulated interaction information are also excluded, so the classifier can be applied to reconstructed data. The main-track energy fraction is retained in the existing classifier, although the one-track requirement makes it essentially unity at this stage.
3.5.3 Universal TPC baseline selection
Throughout this thesis, the TPC baseline denotes the conjunction
(3.4)
The four criteria in Table 3.4 apply to laboratory data and reconstructed simulations for every source, including , muons and neutrons. They use measured or emulated detector observables rather than simulation truth. The criteria remain fixed when the gas, voltage, available readout channels or calibrated response change between runs. The term fiducial refers specifically to .
| Stage | Required condition |
|---|---|
| Valid, | The native maximum-track validity flags for XZ, YZ and their combination are finite and each . The two projection energies are finite and ; their combined energy , the calibrated energy and the coordinates are finite. |
| Energy, | , including both endpoints; uses the previously accepted calibration for that run or sample. |
| Multiplicity, | Exactly one reconstructed track in each projection: NTracksX and NTracksY. |
| Fiducial, | about the fixed readout origin . This is a radius of , not a diameter, with geometric area . The boundary itself is excluded. |
The maximum-energy XZ and YZ tracks supply , and their respective amplitude-weighted centers , . The energy scale is fixed before evaluating the selection. These quantities describe the main tracks, rather than all event hits or the true gas deposit.
3.5.3.0.1 Reconstruct first; select afterwards.
All channels admitted by the documented readout and waveform processing are reconstructed before applying the baseline. The native direct hits-to-track chain clusters XZ and YZ hits separately: an edge joins two hits if their distance in that projection is strictly less than , and connected components form tracks. Connectivity is transitive, so is not a maximum track diameter. An isolated hit is a track; there is no exemption for a low-energy secondary track. No radial hit mask is applied before counting them. A second track outside the central circle therefore still fails . One track in each view does not establish that only one physical particle interacted.
The validity checks establish usable observables, not fit quality or particle identity. The spatial requirement concerns only the maximum-track center: some hits may lie outside the circle, and no additional selection on is imposed. Widths, energy or hit-count balance, pulse shape, learned classifiers and veto activity are not baseline requirements. Undefined unused topology features do not remove an otherwise accepted event. Waveform thresholds, channel masks and drift conversion are reconstruction settings applied before this event selection and documented separately.
3.5.3.0.2 Selection sequence and scope of efficiencies.
Report the stages in the order input, valid reconstruction, energy, multiplicity and fiducial, retaining the failures at each stage in the accounting. The arrows in Figure 3.15 follow this reporting order; an event must satisfy all four criteria. Every quoted acceptance states its numerator and denominator. An acceptance among saved or valid reconstructed events is not an absolute efficiency per generated photon; the latter also requires primaries with no saved or reconstructed event. Diagnostic plots may deliberately precede the baseline, provided their population is named. Additional veto or topology selections are separate stages and cannot be absorbed into the word “baseline.”

3.5.4 Calibration and background distributions
Figure 3.16 shows six of these observables for the simulated calibration and cosmic-background candidates used to check the conventional-response classifier. The calibration sample contains 11 499 events from the original test partition of one merged simulation file. This retrospective partition separates events by their identifiers; it is not an independent simulation production or a blinded validation sample. The background sample contains 172 candidates from 129 cosmic-production groups excluded from classifier development: 141 neutron, 30 proton and one electron candidate. Both populations pass the common baseline and have finite values for all 15 classifier inputs, before applying the BDT or veto selection. Each candidate is counted equally; these shapes do not represent source-normalized background rates. These samples illustrate the populations entering the conventional-response selector. The measured-data comparison in Section 3.6 checks how well the response models reproduce the reconstructed shapes.

The remaining cosmic candidates often occupy fewer strips and have smaller widths than the calibration events: the median minimum hit count is two, compared with nine for . Their transverse widths are also less balanced between views. The energy balance overlaps substantially, showing why no single quantity fully separates the two populations. These are backgrounds that already passed the compact-event baseline; extended tracks rejected earlier do not appear in this comparison. The BDT combines the differences and their correlations into one selection.
3.5.5 Multivariate event selection
Combining reconstructed observables can improve the separation of signal and background compared with individual selection criteria. A multivariate classifier assigns a score to the reconstructed observables , and events are selected by placing a threshold on this score. Two methods are used here: binned log-likelihood ratios and boosted decision trees. The TPC classifier uses the shape quantities introduced above, while the veto classifier uses its own pulse and panel observables.
3.5.5.0.1 Binned log-likelihood ratio (“log-odds”).
Normalized histograms of the training samples estimate how each observable is distributed in signal () and background (). Their contributions are combined as
(3.5)
Values more common in signal contribute positively to the score. Small smoothing terms avoid infinite ratios in empty bins. The method is easy to interpret, although it neglects correlations between observables within each class, as in the naive-Bayes approximation [105].
3.5.5.0.2 Boosted decision trees (BDTs).
A decision tree separates events through successive selection criteria on their observables. Boosting combines many trees. In gradient boosting, successive trees are fitted to improve the combined prediction [106]. This allows correlations and more complex boundaries to be used, but also makes overfitting possible. Applications include particle-identification studies for MiniBooNE [107], rejection of misreconstructed pulses in XENON1T [108], and selection of rare decays in CMS and LHCb [109]. The classifiers used in this thesis are implemented in scikit-learn [105].

3.5.5.0.3 Selection and performance.
Training, threshold choice and final evaluation use separate samples. Splitting by complete runs or simulated histories helps avoid using closely related events for both training and testing [105]. The partitions used for the TPC and veto studies, and the limitations of the archived calibration samples, are described in their respective chapters. For a threshold , events with are accepted. The signal and background acceptances are and , and the background rejection is . A tighter threshold generally rejects more background but also loses signal. For example, a working point with and keeps 90% of the signal sample and rejects 99% of the background sample. These acceptances refer to the samples entering the classifier, rather than the full detector efficiency. Applying the result to measured data also requires an adequate description of the detector response.
3.5.5.0.4 Comparing classifiers.
Scanning the threshold shows the background rejection achievable at each signal acceptance. This efficiency curve is called the receiver operating characteristic (ROC), a term also used in HEP tools such as TMVA [110]. Its standard form plots signal acceptance against background acceptance ; HEP plots often show the equivalent background rejection against signal acceptance instead. The latter directly answers a useful experimental question: how much background can be rejected while retaining a chosen fraction of signal?
The area under the ROC curve (AUC) summarizes the separation over all thresholds [110]. With on the horizontal axis and on the vertical axis, an uninformative score gives a diagonal curve with AUC , while perfect separation gives AUC . A larger AUC indicates better overall separation, but it does not determine the best operating point for the experiment. For rare-event searches, the background rejection at the required signal acceptance remains the more direct measure of selection performance.
Section 4 applies these methods to distinguish neutron-associated veto activity from calibration-triggered activity; the control acceptance measures the fraction not tagged by the veto. Section 5 sets the TPC working point and evaluates the background surviving it, using the observables described in Section 3.5.2.
3.6 Response calibration and comparison with data
The refined response was compared with IAXO-D1 acquisitions in Ar/isobutane 99/1 to test reconstructed event properties and baseline acceptance. We first calibrated the mean amplitude and compared the reconstructed distributions, then examined an additional gain correction at the reference condition. These comparisons test the calibrated response as a whole, including waveform inputs derived from measured data.
3.6.1 Calibration across operating conditions
The gas, geometry and fields determine electron transport for each operating condition. Two electronics groups use shared effective pulse shapes and channel delays learned from 28 whole development runs. Each acquisition provides its residual-noise sample, trigger placement and reconstruction settings. In this initial comparison, only the conversion from collected electrons to ADC amplitude is fitted within each run, using the energy-peak position in its development sample. Widths, multiplicity, spatial distributions and selection acceptance are not used in that fit. The per-run amplitude calibration therefore accompanies a shared empirical waveform response, as summarized in Table 3.5.
| Level | Response information | Treatment |
|---|---|---|
| Physical condition | Gas, geometry and fields; electron transport | Fixed inputs |
| Electronics group | Effective pulse and channel-delay map | Shared; development data |
| Acquisition | Residual noise, trigger placement and reconstruction settings | Measured and fixed |
| Run amplitude | ADC per collected electron | Fitted to development peak |
| Additional corrections | Pulse stretch, added jitter and event-gain spread | Fixed at unity or zero |
The comparison covers 40 acquisitions, including twelve whole checking runs excluded from the shared-pulse calibration. Measured and simulated checking blocks are separate from amplitude fitting; previously inspected acquisitions and shared simulation histories make this a consistency check. Before the additional R02729 gain variation discussed below, the median simulation/data width ratios in the twelve checking runs are 0.977 transversely and 0.943 longitudinally. Nevertheless, simulation has greater baseline acceptance in 38 of 40 acquisitions, with a median absolute difference of 16.1 percentage points in the twelve checking runs (Figure 3.18). Agreement of central widths therefore does not establish agreement of the selection efficiency.

The larger R02729 sample below examines the response distributions and their effect on selection in more detail.
3.6.2 A representative calibration run
To compare the response distributions with larger samples, we selected run R02729, recorded in nominal Ar/isobutane 99/1 at . It contains over 62 000 saved events, with a calibration-peak variation of about 1% across twelve time blocks. We used 10 000 measured waveforms spread across its checking blocks and 15 000 fresh gas-entry histories for the individual-electron Garfield model. A new conventional REST production supplies 12 000 gas-interacting histories, of which 200 are reserved for the amplitude calibration and the remaining 11 800 for the comparison. The two simulations use the reference chamber geometry and source position, with the same nominal gas, pressure and drift field.
REST starts from gas energy deposits, while Garfield starts from particles entering the gas. The comparison is therefore conditional on these populations and does not measure the efficiency for an emitted photon to produce a recorded event.
The measured and simulated waveforms undergo the same REST reconstruction, with the same readout, energy conversion and settings. Both simulations use residual noise from the measured development sample. The conventional calculation retains its analytic pulse and additional charge smearing common to the whole event, while Garfield uses microscopic avalanches and the shared empirical pulse described above. The avalanches fluctuate, but no additional event-by-event gain variation is included in this Garfield comparison. The amplitude calibration and separation of development and checking events follow the procedure above; widths, positions and acceptance are not fitted. The previously inspected acquisition provides a consistency check of the complete calibrated chains.
Figure 3.19 compares the reconstructed energy, position, widths and balances between strip views before the baseline selection criteria. Energy and radius provide response checks alongside the shape variables used by the BDT.

Garfield better reproduces the longitudinal-width distribution and the charge balance between the two strip views. For example, the median is in data, with conventional REST and with Garfield. The transverse-width distribution changes more modestly, and both models describe the hit-count balance reasonably well. The agreement is not uniform: the Garfield energy peak remains too narrow, conventional REST describes the energy distribution better, and both models place events somewhat closer to the detector center than the measurement. Table 3.6 tests how these differences affect acceptance when the common selection criteria are applied together.
Section 5 develops the BDT selection using the separate calibration and background productions.
3.6.3 Energy response and additional gain variations
The remaining difference in the energy spectrum motivates an additional gain calibration. Primary ionization and avalanche multiplication both contribute to the energy width [73]. Garfield calculates these processes for the reference gas and mesh, and the waveform calculation adds the measured residual noise. The REST detector response also applies a common gain fluctuation to all hits in an event. We introduced a similar additional term in the Garfield waveforms to describe the remaining energy spread. Unlike fluctuations of separate electron avalanches, this term changes the signal from the whole charge cloud together.
The conversion from collected electrons to signal amplitude is multiplied by one positive factor for the whole event:
(3.6)
Here sets the spread and the subtraction in the exponent gives . We therefore keep the mean amplitude calibration while adding a gain variation between events. The reconstructed peak need not remain at exactly the same position after digitization, pulse detection and fitting. The same factor applies to the signal in both strip views before ADC conversion and reconstruction; the noise is added separately. This term describes an effective variation of the measured response, rather than changing the simulated ionization, electron trajectories or avalanches. The spectrum alone does not identify whether the missing spread comes from gain variations, the electronics response or other unresolved effects.
We fitted to the width of the main energy peak in separate measured and simulated development events, keeping the existing mean amplitude conversion, pulse and noise fixed. The measured sample contains valid reconstructions from 10 000 events. The simulation uses 2000 fresh development histories from the reference source production, separate from the 15 000 checking histories. We fit a Gaussian in the same window in the two samples, before applying the baseline or classifier selection, and compare its width in . The peak position is checked separately. The fitted spread corresponds to an additional gain RMS of about 9%. Separate checking events test the resulting spectrum after the parameter is fixed. This is a refinement motivated by the observed discrepancy, rather than a blind validation of the original response.
Figure 3.20 compares the measured and REST spectra with Garfield before and after the additional gain calibration. The two Garfield samples use the same checking histories. The measured, REST and Garfield checking samples contain, respectively, , and valid events; the Garfield count is unchanged after the gain correction.

Before the additional gain variation, the Garfield response already reproduces the main calibration and escape-peak structures, alongside the improvements in drift width and charge balance shown in Figure 3.19. Its main energy peak is nevertheless too narrow. In the separate checking sample, the additional gain spread changes the fitted Gaussian width from to , close to the measured . The fitted center and the high-energy tail still differ: 0.76% of corrected Garfield events lie above , compared with 9.91% of measured events. With normalization to all valid events, the smaller simulated tail also contributes to the larger fraction in the main peak. The microscopic calculation therefore provides a useful description of the charge cloud response, while the additional empirical term improves the central energy spread. Section 3.6.4 evaluates the baseline acceptance after this correction, including the effect of the remaining spectral differences.
3.6.4 Impact of the detector response on event selection
The response distributions determine how many reconstructed events survive the baseline selection. A broader energy peak can move events across the energy limits, while secondary tracks and positions near the edge of the fiducial region affect the multiplicity and spatial selection criteria. We therefore apply the unchanged requirements of Section 3.5.3 to the same R02729 checking samples used above. Table 3.6 compares the selected fractions of valid reconstructed events; no response parameter or selection criterion is adjusted using these checking results.
| Sample | Baseline / valid events | Retention [%] |
|---|---|---|
| Measured calibration | ||
| REST detector response | ||
| Garfield detector response | ||
| Garfield + gain variation |
The REST and original Garfield acceptances exceed the measured value by 6.7 and 10.4 percentage points, respectively. The additional gain variation reduces the Garfield discrepancy to 7.2 percentage points, but does not remove it. The improved energy width therefore accompanies a remaining difference in the selected fraction, for which the energy tails, positions and track multiplicities also matter. These selected fractions must be checked before transferring the response to the background studies in Section 5.
3.6.5 Limitations and further checks
Separate tests of electron assignment and noiseless waveform reconstruction check charge accounting, timing and pulse addition. These establish that the calculation is implemented consistently, while its physical accuracy requires comparison with measurements. The empirical pulse can include gas or illumination time spread, and the amplitude calibration can absorb both avalanche gain and electronics conversion. Consequently, the agreement above does not independently determine transport coefficients or absolute gain. A change in the field grid altered the collected charge per injected primary by about 5.8%, so gain convergence also requires further study.
Further checks must address the discrepancy in selection acceptance as well as the gain convergence.
3.6.6 Computational cost
The Garfield detector response follows individual electrons through drift and avalanche multiplication, while the REST detector response uses parameters for diffusion and effective gain. We compared the calculation times for calibration photons, muons crossing the gas, and a muon background sample transported through the shielding and veto geometry. The comparison quantifies the cost of the additional detail used in the response validation above. Both calculations use Ar/isobutane 99/1 at and the same central region, with fixed response parameters. The Garfield calculation omits ionization before entry into this region and gas-to-surroundings transport feedback; equal admitted areas do not ensure identical microscopic paths. They save channel waveforms after pulse formation, noise and ADC conversion; common event reconstruction is excluded. Each calculation uses a long continuous batch, spreading the initialization cost over many events.
| Simulation | Response | Events | Simulation time [min] | Time per event [s] |
|---|---|---|---|---|
| Calibration photons | REST | |||
| Garfield | ||||
| Gas-crossing muons | REST | |||
| Garfield | ||||
| Muon background | REST | |||
| Garfield |
Table 3.7 gives the generated particle counts and total wall time, including initialization and output. Figure 3.21 shows the fractions of CPU time spent in each stage and the CPU time per generated particle. Particles producing no gas signal remain included, so the background comparison also accounts for the transport of muons that do not produce a TPC event.
The microscopic gas calculation accounts for about 99% of CPU time for calibration photons and 96% for gas-crossing muons. The REST calibration calculation instead spends most of its time forming and saving waveforms. For the muon background, only 151 of the incident particles have a primary or secondary reaching the gas; 44 events produce waveforms with REST and 46 with Garfield. Transport outside the gas takes 57% of the complete CPU time with REST and 48% with Garfield. Including the remaining Geant4 event processing, these fractions rise to 95% and 80%, respectively. In this background sample, the complete Garfield calculation takes 18% more wall time than REST, although its individual gas responses are much more expensive. The cost of the microscopic response therefore depends strongly on how often the simulated background produces a gas signal.
A single worker runs on the same processor core in both cases. The neutron guard rejects unsupported gas interactions; its success in these muon samples does not implement a general neutron charge response. The timing comparison is conditional on these source populations and response boundaries.

3.7 Simulation production
The background studies require large samples for several particle sources, shielding configurations and veto designs. These samples are generated in independent jobs following the chain in Figure 3.1. The simulation settings are maintained in the collaboration’s repository [111]. To obtain background rates, the exposure is calculated from the source distribution and the number of generated primary particles, including those that leave no recorded detector event.
3.7.1 Batch production
HTCondor [112] distributes the jobs across the DESY NAF-IAXO computing environment, with files stored in dCache [113]. The restG4ToCondor.py tool developed in this thesis prepares and submits the jobs and collects their logs and outputs. The workflow is shown in Figure 3.22.

An initial validation campaign comprised approximately 600 jobs and 40 000 saved cosmic-neutron and cosmic-muon events. Larger campaigns followed, as described in Section 5. Before combining the results, the files, simulation settings and exposure are checked (Table 3.8).
| Check | Purpose |
|---|---|
| Output files | Verify that the job finished and the events can be read. |
| Simulation settings | Identify the source, geometry, gas and response used. |
| Exposure | Retain generated-primary counts and the corresponding source integral. |
| Event grouping | Keep repeated responses and delayed subevents linked to their original particle. |
3.7.2 Reprocessing and saved information
Saving intermediate stages allows part of the calculation to be repeated when a model or analysis setting changes (Table 3.9). The starting point depends on the change: the Garfield response requires particle information at gas entry, while saved charge arrivals allow new waveforms to be generated. Each change also requires the following stages to be repeated.
| Change | Saved information | Repeat from |
|---|---|---|
| Source, chamber materials or geometry | Source definition and geometry | Particle transport through the detector |
| Gas or bulk-drift conditions (Garfield response) | Particles at first gas entry | Primary ionization, drift and amplification |
| Local mesh field or microscopic transport | Bulk-drift boundary arrivals | Amplification and pad collection |
| Pulse response, noise or ADC conversion | Individual pad arrivals, with channel and time | Channel waveform formation |
| Pulse fitting, clustering or energy calibration | Corrected channel waveforms | Shared REST reconstruction |
| Event selection | Reconstructed observables and normalization records | Selection criteria and physical rate calculation |
Heed. Reusing bulk-drift arrivals assumes unchanged bulk conditions and handover surface; a mesh-field change that alters bulk transport requires repeating drift. Changing the gas in the REST response also requires recalculating its gas interactions in Geant4. Gas deposits cannot replace gas-entry records, and reconstructed tracks cannot recover earlier waveforms.Repeated responses and delayed subevents remain linked to their original particle for statistical analysis. Cosmic-production grouping also applies to classifier training and testing, so that related events do not appear in both samples. The retrospective partition of the archived calibration sample is described in Section 3.5.2. Threshold calibration and final evaluation are recorded separately.
Software versions, input files, geometry, seeds and analysis settings are saved to make the calculation reproducible. Automated software tests were migrated to GitHub Actions during this work. Physical validation relies on the comparisons with detector data and the source-normalization studies presented in the following chapters.
3.8 Supporting software and detector operation
We developed several tools to inspect simulated events, acquire waveforms and monitor the detector during operation. Figure 3.23 summarizes their roles.

Uproot and fsspec allow the analysis tools to read local and remote ROOT files.3.8.1 Browser-based event inspection
We developed a browser-based viewer to display simulated particle tracks in the IAXO geometry. It uses three.js to render event information exported from REST, allowing events to be inspected without a ROOT installation. The design was inspired by the Phoenix event-display project [114].

The display helps identify the interactions responsible for a detector event and check the geometry and readout mapping. Time playback makes it possible to follow prompt cascades and delayed neutron captures. The viewer was also used to prepare the event illustrations in Section 5.
3.8.2 Acquisition and online waveform inspection
We developed feminos-daq from the legacy mclient software for the AGET/Feminos electronics [115]. It adds ROOT output, monitoring and a waveform viewer, with a CMake build. Incoming UDP frames are buffered while events are decoded and written, reducing the effect of short processing delays (Figure 3.25).

feminos-daq. A growing queue indicates that writing cannot keep up with reception. Periodic flushing reduces the amount of data left unwritten after an unexpected stop.Event information and waveform samples are stored without custom ROOT dictionaries, allowing them to be read with ROOT or Uproot and converted to REST events.

feminos-viewer application. Left: waveforms from a single event, with channels colored by their readout position. Right: the corresponding hit map in the detector plane.The Python application feminos-viewer displays waveforms, channel activity and readout positions from local or remote files, including files still being written. During detector operation it was used to identify malformed events, inactive or noisy channels and incorrect readout mappings.
3.8.3 High-voltage control and operating conditions
We developed the hvps library [116] to provide a common Python interface to CAEN and iseg serial high-voltage supplies. It allows voltages, currents and channel status to be read, channels to be switched and voltage ramps to be configured. Command-line tools, Node.js bindings and a Node-RED node support scripts and slow-control dashboards. We also developed a graphical prototype for the CAEN N1471H supply.
These quantities help monitor the detector and identify periods affected by trips or unstable conditions (Table 3.10). Modeling gas transport also requires the electrode geometry and gas conditions. The library was distributed through PyPI and npm.
| Operating quantity | Use in detector operation or response interpretation |
|---|---|
| Set and monitored voltage | Check electrode operating points and determine fields with the gap geometry. |
| Current and channel status | Identify trips or abnormal operating periods. |
| Ramp and interlock settings | Control voltage changes and channel protection. |
3.8.4 Uproot and fsspec integration
During a three-month IRIS-HEP fellowship, we contributed to the integration of fsspec into Uproot, the Python library for reading ROOT files. Previously, local files and remote protocols such as HTTP and XRootD required separate implementations in Uproot. The integration, released in Uproot v5.2.0, provides a common file-access interface and adds support for protocols handled by fsspec, such as SSH, cloud object stores and WebDAV.
This development allowed feminos-viewer to open remote DAQ files over SSH during acquisition and the background-analysis scripts to read simulation files directly from NAF-IAXO dCache storage. The experience with the ROOT format also informed the design of the DAQ output.
3.8.5 Python prototyping with Geant4
We developed geant4-python-application for small transport studies in notebooks and scripts. The package runs Geant4 in separate worker processes and was demonstrated with a photon-transmission calculation. The main simulation campaigns used restG4.
3.9 Summary
The framework described in this chapter follows particles from their generation to the reconstructed detector event. We extended restG4 and Geant4Lib with multithreading, track pruning, delayed-subevent handling and geometry-name resolution. The PDP method improved cosmic-particle injection, and the HTCondor tools allowed large simulation samples to be produced.
The Garfield response follows gas ionization, electron drift, amplification and individual charge collection before the common REST reconstruction. Measured pulses and residual noise complete the waveform response. The planned deposit-to-electron conversion extends this calculation to neutron gas interactions; its implementation and validation are marked preliminary in this draft. The reconstructed energy, position and multiplicity define the common TPC baseline. Widths, hit counts, charge balance and track asymmetry provide the inputs to the TPC BDT. Calibration and background distributions illustrate their role in event selection. After an empirical waveform calibration, typical reconstructed event widths are close to measurements, but the simulation overestimates baseline selection acceptance (Section 3.6). Section 5 fixes the additional BDT working point and evaluates signal acceptance, background rejection and source-normalized rates. Event viewers and readout mapping tools supported the simulation studies, while acquisition and high-voltage software supported detector operation. The Uproot/fsspec contribution and Python interfaces made local and remote data accessible for analysis and prototyping.
The following chapters apply these tools to the shielding and veto studies (Section 4) and the background model (Section 5).