tomopt.inference package

Submodules

tomopt.inference.scattering module

class tomopt.inference.scattering.GenScatterBatch(mu, volume)

Bases: ScatterBatch

Class for computing scattering information from the true hits via incoming/outgoing trajectory fitting.

Warning

This class is intended for diagnostic purposes only. The tracks and scatter variables carry no gradient w.r.t. detector parameters (except z position).

Linear fits are performed separately to all hits associated with layer groups, as indicated by the pos attribute of the layers which recorded hits. Currently, the inference methods expect detectors above the passive layer to have pos==’above’, and those below the passive volume to have pos==’below’. Trajectory fitting is performed using an analytic likelihood minimisation, but no uncertainties on the hits are considered.

Important

The current separation of hits into above and below groups does not allow for e.g. a third set of detectors, since this split is based on the value of the n_hits_above attribute.

One instance of this class should created for each MuonBatch. As part of the initialisation, muons will be filtered using _filter_scatters() in order to avoid NaN/Inf values. This results in direct, in-place changes to the MuonBatch.

Since many variables of the scattering can be inferred, but not all are required for further inference downstream, variables, and their uncertainties, are computed on a lazy basis, with memoisation: the values are only computed on the first request (if at all) and then stored in case of further requests.

The dtheta, dphi, and total scattering variables are computed under the assumption of small angular scatterings. An assumption is necessary here, since there is a loss of information in the when the muons undergo scattering in theta and phi: since theta is [0,pi] a negative scattering in theta will always results in a positive theta, but phi can become phi+pi. When inferring the angular scattering, one cannot precisely tell whether instead a large scattering in phi occurred. The total scattering (total_scatter) is the quadrature sum of dtheta and dphi, and all three are computed under both hypotheses. The final values of these are chosen using the hypothesis which minimises the total amount of scattering. This assumption has been tested and found to be good.

Parameters:

• mu (MuonBatch) – muons with hits to infer on

• volume (Volume) – volume through which the muons travelled

class tomopt.inference.scattering.ScatterBatch(mu, volume)

Bases: object

Class for computing scattering information from the hits via incoming/outgoing trajectory fitting.

Linear fits are performed separately to all hits associated with layer groups, as indicated by the pos attribute of the layers which recorded hits. Currently, the inference methods expect detectors above the passive layer to have pos==’above’, and those below the passive volume to have pos==’below’. Trajectory fitting is performed using an analytic likelihood minimisation, which considers uncertainties and efficiencies on the hits in x and y.

Important

The current separation of hits into above and below groups does not allow for e.g. a third set of detectors, since this split is based on the value of the n_hits_above attribute.

One instance of this class should created for each MuonBatch. As part of the initialisation, muons will be filtered using _filter_scatters() in order to avoid NaN/Inf gradients or values. This results in direct, in-place changes to the MuonBatch.

Since many variables of the scattering can be inferred, but not all are required for further inference downstream, variables, and their uncertainties, are computed on a lazy basis, with memoisation: the values are only computed on the first request (if at all) and then stored in case of further requests.

The dtheta, dphi, and total scattering variables are computed under the assumption of small angular scatterings. An assumption is necessary here, since there is a loss of information in the when the muons undergo scattering in theta and phi: since theta is [0,pi] a negative scattering in theta will always results in a positive theta, but phi can become phi+pi. When inferring the angular scattering, one cannot precisely tell whether instead a large scattering in phi occurred. The total scattering (total_scatter) is the quadrature sum of dtheta and dphi, and all three are computed under both hypotheses. The final values of these are chosen using the hypothesis which minimises the total amount of scattering. This assumption has been tested and found to be good.

Parameters:

• mu (MuonBatch) – muons with hits to infer on

• volume (Volume) – volume through which the muons travelled

property above_gen_hits: Tensor | None

Returns: (muons,hits,xyz) tensor of true hits in the “above” detectors

property above_hit_effs: Tensor | None

Returns: (muons,hits,effs) tensor of hit efficiencies in the “above” detectors

property above_hit_uncs: Tensor | None

Returns: (muons,hits,xyz) tensor of uncertainties on hits in the “above” detectors

property above_hits: Tensor | None

Returns: (muons,hits,xyz) tensor of recorded hits in the “above” detectors

property below_gen_hits: Tensor | None

Returns: (muons,hits,xyz) tensor of true hits in the “below” detectors

property below_hit_effs: Tensor | None

Returns: (muons,hits,eff) tensor of hit efficiencies in the “below” detectors

property below_hit_uncs: Tensor | None

Returns: (muons,hits,xyz) tensor of uncertainties on hits in the “below” detectors

property below_hits: Tensor | None

Returns: (muons,hits,xyz) tensor of recorded hits in the “below” detectors

property dphi: Tensor

Returns: (muons,1) delta phi between incoming & outgoing muons

property dphi_unc: Tensor

Returns: (muons,1) uncertainty on dphi

property dtheta: Tensor

Returns: (muons,1) delta theta between incoming & outgoing muons

property dtheta_unc: Tensor

Returns: (muons,1) uncertainty on dtheta

property dtheta_xy: Tensor

Returns: (muons,xy) delta theta_xy between incoming & outgoing muons in the zx and zy planes

property dtheta_xy_unc: Tensor

Returns: (muons,xy) uncertainty on dtheta_xy

property dxy: Tensor

Returns: (muons,xy) distances in x & y from PoCA to incoming|outgoing muons

property dxy_unc: Tensor

Returns: (muons,xy) uncertainty on dxy

property gen_hits: Tensor | None

Returns: (muons,hits,xyz) tensor of true hits

static get_muon_trajectory(hits, uncs, lw)

Fits a linear trajectory to a group of hits, whilst considering their uncertainties on their xy positions. No uncertainty is considered for z positions of hits. The fit is performed via an analytical likelihood-maximisation.

Important

Muons with <2 hits have NaN trajectory

Parameters:

• hits (Tensor) – (muons,hits,xyz) tensor of hit positions

• uncs (Tensor) – (muons,hits,(unc x,unc y,0)) tensor of hit uncertainties

• lw (Tensor) – length and width of the passive layers of the volume

Returns:

(muons,xyz) fitted-vector directions start: (muons,xyz) initial point of fitted-vector

Return type:

vec

get_scatter_mask()

Return type:

Tensor

Returns:

(muons) Boolean tensor where True indicates that the PoCA of the muon is located within the passive volume

property hit_effs: Tensor | None

Returns: (muons,hits,eff) tensor of hit efficiencies

property hit_uncs: Tensor | None

Returns: (muons,hits,xyz) tensor of uncertainties on hits

property hits: Dict[str, Dict[str, Tensor]]

Returns: Dictionary of hits, as returned by get_hits()

property n_hits_above: int | None

Returns: Number of hits per muon in the “above” detectors

property n_hits_below: int | None

Returns: Number of hits per muon in the “below” detectors

property phi_in: Tensor

Returns: (muons,1) phi of incoming muons

property phi_in_unc: Tensor

Returns: (muons,1) uncertainty on phi_in

property phi_out: Tensor

Returns: (muons,1) phi of outgoing muons

property phi_out_unc: Tensor

Returns: (muons,1) uncertainty on phi_out

plot_scatter(idx, savename=None)

Plots representation of hits and fitted trajectories for a single muon.

Parameters:

• idx (int) – index of muon to plot

• savename (Optional[Path]) – optional path to save figure to

Return type:

None

property poca_xyz: Tensor

Returns: (muons,xyz) xyz location of PoCA

property poca_xyz_unc: Tensor

Returns: (muons,xyz) uncertainty on poca_xyz

property reco_hits: Tensor | None

Returns: (muons,hits,xyz) tensor of recorded hits

property theta_in: Tensor

Returns: (muons,1) theta of incoming muons

property theta_in_unc: Tensor

Returns: (muons,1) uncertainty on theta_in

property theta_msc: Tensor

Returns: (muons,1) theta_msc; the total amount of angular scattering

property theta_msc_unc: Tensor

Returns: (muons,1) uncertainty on total_scatter

property theta_out: Tensor

Returns: (muons,1) theta of outgoing muons

property theta_out_unc: Tensor

Returns: (muons,1) uncertainty on theta_out

property theta_xy_in: Tensor

Returns: (muons,xy) decomposed theta and phi of incoming muons in the zx and zy planes

property theta_xy_in_unc: Tensor

Returns: (muons,xy) uncertainty on theta_xy_in

property theta_xy_out: Tensor

Returns: (muons,xy) decomposed theta and phi of outgoing muons in the zx and zy planes

property theta_xy_out_unc: Tensor

Returns: (muons,xy) uncertainty on theta_xy_out

property total_scatter: Tensor

Returns: (muons,1) theta_msc; the total amount of angular scattering

property total_scatter_unc: Tensor

Returns: (muons,1) uncertainty on total_scatter

property track_in: Tensor | None

Returns: (muons,xyz) incoming xyz vector

property track_out: Tensor | None

Returns: (muons,xyz) outgoing xyz vector

property track_start_in: Tensor | None

Returns: (muons,xyz) initial point of incoming xyz vector

property track_start_out: Tensor | None

Returns: (muons,xyz) initial point of outgoing xyz vector

property xyz_in: Tensor

Returns: (muons,xyz) inferred xy position of muon at the z-level of the top of the passive volume

property xyz_in_unc: Tensor

Returns: (muons,xyz) uncertainty on xyz_in

property xyz_out: Tensor

Returns: (muons,xyz) inferred xy position of muon at the z-level of the bottom of the passive volume

property xyz_out_unc: Tensor

Returns: (muons,xyz) uncertainty on xyz_out

tomopt.inference.volume module

class tomopt.inference.volume.AbsIntClassifierFromX0(partial_x0_inferrer, volume, output_probs=True, class2float=None)

Bases: AbsVolumeInferrer

Abstract base class for inferring integer targets through multiclass classification from voxelwise X0 predictions. Inheriting classes must provide a way to convert voxelwise X0s into class probabilities of the required dimension. Requires a basic inferrer for providing the voxelwise X0 predictions. Optionally, the predictions can be returns as the raw class predictions, or the most probable class. In case of the latter, this class can be optionally be converted to a float value via a user-provided processing function.

Parameters:

• partial_x0_inferrer (Type[AbsX0Inferrer]) – (partial) class to instatiate to provide the voxelwise X0 predictions

• volume (Volume) – volume through which the muons will be passed

• output_probs (bool) – if True, will return the per-class probabilites, otherwise will return the argmax of the probabilities, over the last dimension

• class2float (Optional[Callable[[Tensor, Volume], Tensor]]) – optional function to convert class indices to a floating value

add_scatters(scatters)

Appends a new set of muon scatter vairables. When get_prediction() is called, the prediction will be based on all ScatterBatch s added up to that point

Return type:

None

compute_efficiency(scatters)

Compuates the per-muon efficiency according to the method implemented by the X0 inferrer.

Parameters:

scatters (ScatterBatch) – scatter batch containing muons whose efficiency should be computed

Return type:

Tensor

Returns:

(muons) tensor of muon efficiencies

get_prediction()

Computes the predicions for the volume. If class probabilities were requested during initialisation, then these will be returned. Otherwise the most probable class will be returned, and this will be converted to a float value if class2float is not None.

Returns:

(*) volume prediction

Return type:

pred

abstract x02probs(vox_preds)

Inheriting classes must override this method to convert voxelwise X0 predictions to class probabilities

Parameters:

vox_preds (Tensor) – (z,x,y) tensor of voxelwise X0 predictions

Return type:

Tensor

Returns:

(*) tensor of class probabilities

class tomopt.inference.volume.AbsVolumeInferrer(volume)

Bases: object

Abstract base class for volume inference.

Inheriting classes are expected to be fed multiple ScatterBatch s, via add_scatters(), for a single Volume and return a volume prediction based on all of the muon batches when get_prediction() is called.

Parameters:

volume (Volume) – volume through which the muons will be passed

add_scatters(scatters)

Appends a new set of muon scatter variables. When get_prediction() is called, the prediction will be based on all ScatterBatch s added up to that point

Return type:

None

abstract compute_efficiency(scatters)

Inheriting classes must override this method to provide a computation of the per-muon efficiency, given the individual muon hit efficiencies.

Return type:

Tensor

abstract get_prediction()

Inheriting classes must override this method to provide a prediction computed using the added scatter batches. E.g. the sum of muon efficiencies.

Return type:

Optional[Tensor]

class tomopt.inference.volume.AbsX0Inferrer(volume)

Bases: AbsVolumeInferrer

Abstract base class for inferring the X0 of every voxel in the passive volume.

The inference is based on the PoCA approach of assigning the entirety of the muon scattering to a single point, and the X0 computation is based on inversion of the PDG scattering model described in https://pdg.lbl.gov/2019/reviews/rpp2018-rev-passage-particles-matter.pdf.

Once all scatter batches have been added, the inference proceeds thusly:

• For each muon i, a probability p_ij, is computed according to the probability that the PoCA was located in voxel j.

• These probabilities are computed by integrating over the voxel the PDF of 3 uncorrelated Gaussians centred on the PoCA, with scales equal the uncertainty on the PoCA position in x,y,z.

• p_ij is multiplied by muon efficiency e_i to compute a muon/voxel weight w_ij.

• Inversion of the PDG model gives: \(X_0 = \left(\frac{0.0136}{p^{\mathrm{rms}}}\right)^2\frac{\delta z}{\cos\left(\bar{\theta}^{\mathrm{rms}}\right)}\frac{2}{\theta^{\mathrm{rms}}_{\mathrm{tot.}}}\)

• In order to account for the muon weights and compute different X0s for the voxels whilst using the whole muon population:

  • Weighted RMSs are computed for each of the scattering terms in the right-hand side of the equation.

  • In addition to the muon weight w_ij, the variances of the squared values of the scattering variables is used to divide w_ij.

• The result is a set of X0 predictions X0_j.

Important

Inversion of the PDG model does NOT account for the natural log term.

Important

To simplify the computation code, this class relies heavily on lazy computation and memoisation; be careful if calling private methods manually.

Parameters:

volume (Volume) – volume through which the muons will be passed

get_prediction()

Computes the predicted X0 per voxel as a (z,x,y) tensor via PDG scatter-model inversion for the provided scatter batches.

Returns:

(z,x,y) voxelwise X0 predictions

Return type:

pred

property muon_efficiency: Tensor

Returns: (muons,1) tensor of the efficiencies of the muons

property muon_mom: Tensor

Returns: (muons,1) tensor of the momenta of the muons

property muon_mom_unc: Tensor

Returns: (muons,1) tensor of the uncertainty on the momenta of the muons

property muon_poca_xyz: Tensor

Returns: (muons,xyz) tensor of PoCA locations

property muon_poca_xyz_unc: Tensor

Returns: (muons,xyz) tensor of PoCA location uncertainties

property muon_probs_per_voxel_zxy: Tensor

Warning

Integration tested only

TODO: Don’t assume that poca_xyz uncertainties are uncorrelated TODO: Improve efficiency: currently CDFs are computed multiple times at the same points; could precompute x,y,z probs once, and combine in triples :returns: (muons,z,x,y) tensor of probabilities that the muons’ PoCAs were located in the given voxels.

property muon_theta_in: Tensor

Returns: (muons,1) tensor of the thetas of the incoming muons

property muon_theta_in_unc: Tensor

Returns: (muons,1) tensor of the uncertainty on the theta of the incoming muons

property muon_theta_out: Tensor

Returns: (muons,1) tensor of the thetas of the outgoing muons

property muon_theta_out_unc: Tensor

Returns: (muons,1) tensor of the uncertainty on the theta of the outgoing muons

property muon_total_scatter: Tensor

Returns: (muons,1) tensor of total angular scatterings

property muon_total_scatter_unc: Tensor

Returns: (muons,1) tensor of uncertainties on the total angular scatterings

property n_mu: int

Returns: Total number muons included in the inference

property vox_zxy_x0_pred_uncs: Tensor

Warning

Not recommended for use: long calculation; not unit-tested

Returns:

(z,x,y) tensor of uncertainties on voxelwise X0s

property vox_zxy_x0_preds: Tensor

Returns: (z,x,y) tensor of voxelwise X0 predictions

static x0_from_scatters(deltaz, total_scatter, theta_in, theta_out, mom)

Computes the X0 of a voxel, by inverting the PDG scattering model in terms of the scattering variables

Important

Inversion of the PDG model does NOT account for the natural log term.

Parameters:

• deltaz (float) – height of the voxels

• total_scatter (Tensor) – (voxels,1) tensor of the (RMS of the) total angular scattering of the muon(s)

• theta_in (Tensor) – (voxels,1) tensor of the (RMS of the) theta of the muon(s), as inferred using the incoming trajectory/ies

• theta_out (Tensor) – (voxels,1) tensor of the (RMS of the) theta of the muon(s), as inferred using the outgoing trajectory/ies

• mom (Tensor) – (voxels,1) tensor of the (RMS of the) momentum/a of the muon(s)

Return type:

Tensor

Returns:

(voxels,1) estimated X0 in metres

class tomopt.inference.volume.DenseBlockClassifierFromX0s(n_block_voxels, partial_x0_inferrer, volume, use_avgpool=True, cut_coef=10000.0, ratio_offset=-1.0, ratio_coef=1.0)

Bases: AbsVolumeInferrer

Class for inferreing the presence of a small amount of denser material in the passive volume.

Transforms voxel-wise X0 preds into binary classification statistic under the hypothesis of a small, dense block against a light-weight background. This test statistic, s is computed as:

\[r = 2 \frac{\bar{X0}_{0,\mathrm{bkg}} - \bar{X0}_{0,\mathrm{blk}}}{\bar{X0}_{0,\mathrm{bkg}} + \bar{X0}_{0,\mathrm{blk}}} s = \sigma\!(a(r+b))\]

where \(\bar{X0}_{0,\mathrm{blk}}\) is the mean X0 of the N lowest X0 voxels, and \(\bar{X0}_{0,\mathrm{bkg}}\) is the mean X0 of the remaining voxels. a and b are rescaling coefficients and offsets.

This results in a differentiable value constrained beween 0 and 1, with values near 0 indicating that no relatively dense material is present, and values nearer 1 indicating that it is present. In case it is expected that the dense material forms a contiguous block, the voxelwise X0s can be blurred via a stride-1 kernel-size-3 average pooling.

In actuality, the “cut” on X0s into background and block is implemented as a sigmoid weight, centred at the necessary kth value of the X0. This means that the test statisitc is also differentiable w.r.t. the cut.

Parameters:

• n_block_voxels (int) – number of voxels expected to be occupied by the dense material, if present

• partial_x0_inferrer (Type[AbsX0Inferrer]) – (partial) class to instatiate to provide the voxelwise X0 predictions

• volume (Volume) – volume through which the muons will be passed

• use_avgpool (bool) – wether to blur voxelwise X0 predicitons with a stride-1 kernel-size-3 average pooling useful when the dense material is expected to form a contiguous block

• cut_coef (float) – the “sharpness” of the sigmoid weight that splits voxels into block and background. Higher values results in a sharper cut.

• ratio_offset (float) – additive constant for the X0 ratio

• ratio_coef (float) – multiplicative coefficient for the offset X0 ratio

add_scatters(scatters)

Appends a new set of muon scatter vairables. When get_prediction() is called, the prediction will be based on all ScatterBatch s added up to that point

Return type:

None

compute_efficiency(scatters)

Compuates the per-muon efficiency according to the method implemented by the X0 inferrer.

Parameters:

scatters (ScatterBatch) – scatter batch containing muons whose efficiency should be computed

Return type:

Tensor

Returns:

(muons) tensor of muon efficiencies

get_prediction()

Computes the test statistic for the volume, with values near 0 indicating that no relatively dense material is present, and values nearer 1 indicating that it is present.

Returns:

(1,1,1) volume prediction

Return type:

pred

class tomopt.inference.volume.PanelX0Inferrer(volume)

Bases: AbsX0Inferrer

Class for inferring the X0 of every voxel in the passive volume using hits recorded by PanelDetectorLayer s.

The inference is based on the PoCA approach of assigning the entirety of the muon scattering to a single point, and the X0 computation is based on inversion of the PDG scattering model described in https://pdg.lbl.gov/2019/reviews/rpp2018-rev-passage-particles-matter.pdf.

Once all scatter batches have been added, the inference proceeds thusly:

• For each muon i, a probability p_ij, is computed according to the probability that the PoCA was located in voxel j.

• These probabilities are computed by integrating over the voxel the PDF of 3 uncorrelated Gaussians centred on the PoCA, with scales equal the uncertainty on the PoCA position in x,y,z.

• p_ij is multiplied by muon efficiency e_i to compute a muon/voxel weight w_ij.

• Inversion of the PDG model gives: \(X_0 = \left(\frac{0.0136}{p^{\mathrm{rms}}}\right)^2\frac{\delta z}{\cos\left(\bar{\theta}^{\mathrm{rms}}\right)}\frac{2}{\theta^{\mathrm{rms}}_{\mathrm{tot.}}}\)

• In order to account for the muon weights and compute different X0s for the voxels whilst using the whole muon population:

  • Weighted RMSs are computed for each of the scattering terms in the right-hand side of the equation.

  • In addition to the muon weight w_ij, the variances of the squared values of the scattering variables is used to divide w_ij.

• The result is a set of X0 predictions X0_j.

Important

Inversion of the PDG model does NOT account for the natural log term.

Important

To simplify the computation code, this class relies heavily on lazy computation and memoisation; be careful if calling private methods manually.

Parameters:

volume (Volume) – volume through which the muons will be passed

# TODO: refactor this to be provided to volume inference as a callable

compute_efficiency(scatters)

Computes the per-muon efficiency, given the individual muon hit efficiencies, as the probability of at least two hits above and below the passive volume.

Parameters:

scatters (ScatterBatch) – scatter batch containing muons whose efficiency should be computed

Return type:

Tensor

Returns:

(muons) tensor of muon efficiencies