kinematics

kinematics

Classes and functions for relativistic four-momentum kinematics.

class HelicityAdapter(transitions: Union[ReactionInfo, Iterable[Union[Topology, StateTransition]]])[source]

Bases: object

Converter for four-momenta to kinematic variable data.

The create_expressions method forms the bridge between four-momentum data for the decay you are studying and the kinematic variables that are in the HelicityModel. These are invariant mass (see get_invariant_mass_label()) and the \(\theta\) and \(\phi\) helicity angles (see get_helicity_angle_label()).

register_transition(transition: StateTransition) None[source]
register_topology(topology: Topology) None[source]
property registered_topologies: FrozenSet[Topology]
permutate_registered_topologies() None[source]

Register outgoing edge permutations of all registered_topologies.

See Extend kinematic variables.

create_expressions() Dict[str, Expr][source]
create_four_momentum_symbols(topology: Topology) FourMomenta[source]

Create a set of array-symbols for a Topology.

>>> from qrules.topology import create_isobar_topologies
>>> topologies = create_isobar_topologies(3)
>>> create_four_momentum_symbols(topologies[0])
{0: p0, 1: p1, 2: p2}
FourMomenta

A mapping of state IDs to their corresponding FourMomentumSymbol.

It’s best to create a dict of FourMomenta with create_four_momentum_symbols().

alias of Dict[int, FourMomentumSymbol]

FourMomentumSymbol

Array-Symbol that represents an array of four-momenta.

The array is assumed to be of shape \(n\times 4\) with \(n\) the number of events. The four-momenta are assumed to be in the order \(\left(E,\vec{p}\right)\). See also Energy, FourMomentumX, FourMomentumY, and FourMomentumZ.

class Energy(momentum: FourMomentumSymbol, **hints: Any)[source]

Bases: UnevaluatedExpression

Represents the energy-component of a FourMomentumSymbol.

\(E\left(p\right)=p\left[:, 0\right]\)

class FourMomentumX(momentum: FourMomentumSymbol, **hints: Any)[source]

Bases: UnevaluatedExpression

Component \(x\) of a FourMomentumSymbol.

\({p}_x=p\left[:, 1\right]\)

class FourMomentumY(momentum: FourMomentumSymbol, **hints: Any)[source]

Bases: UnevaluatedExpression

Component \(y\) of a FourMomentumSymbol.

\({p}_y=p\left[:, 2\right]\)

class FourMomentumZ(momentum: FourMomentumSymbol, **hints: Any)[source]

Bases: UnevaluatedExpression

Component \(z\) of a FourMomentumSymbol.

\({p}_z=p\left[:, 3\right]\)

class ThreeMomentumNorm(momentum: FourMomentumSymbol, **hints: Any)[source]

Bases: NumPyPrintable, UnevaluatedExpression

Norm of the three-momentum of a FourMomentumSymbol.

\(\left|\vec{p}\right|=\sqrt{\sum_{\mathrm{axis1}}{p\left[:, 1:\right]^{2}}}\)

from numpy import sqrt, sum

numpy.sqrt(sum(p[:, 1:]**2, axis=1))
class InvariantMass(momentum: FourMomentumSymbol, **hints: Any)[source]

Bases: UnevaluatedExpression

Invariant mass of a FourMomentumSymbol.

(1)\[ \begin{eqnarray} m_{p} & = & \sqrt[\mathrm{c}]{E\left(p\right)^{2} - \left|\vec{p}\right|^{2}} \end{eqnarray}\]
class Phi(momentum: FourMomentumSymbol, **hints: Any)[source]

Bases: UnevaluatedExpression

Azimuthal angle \(\phi\) of a FourMomentumSymbol.

(2)\[ \begin{eqnarray} \phi\left(p\right) & = & \operatorname{atan_{2}}{\left({p}_y,{p}_x \right)} \end{eqnarray}\]
class Theta(momentum: FourMomentumSymbol, **hints: Any)[source]

Bases: UnevaluatedExpression

Polar (elevation) angle \(\theta\) of a FourMomentumSymbol.

(3)\[ \begin{eqnarray} \theta\left(p\right) & = & \operatorname{acos}{\left(\frac{{p}_z}{\left|\vec{p}\right|} \right)} \end{eqnarray}\]
class BoostZMatrix(beta: Expr, n_events: Optional[Symbol] = None, **kwargs: Any)[source]

Bases: UnevaluatedExpression

Represents a Lorentz boost matrix in the \(z\)-direction.

Parameters
  • beta – Velocity in the \(z\)-direction, \(\beta=p_z/E\).

  • n_events – Number of events \(n\) for this matrix array of shape \(n\times4\times4\). Defaults to the len of beta.

This boost operates on a FourMomentumSymbol and looks like:

(4)\[\begin{split}\boldsymbol{B_z}\left(\beta\right) = \left[\begin{matrix}\frac{1}{\sqrt{1 - \beta^{2}}} & 0 & 0 & - \frac{\beta}{\sqrt{1 - \beta^{2}}}\\0 & 1 & 0 & 0\\0 & 0 & 1 & 0\\- \frac{\beta}{\sqrt{1 - \beta^{2}}} & 0 & 0 & \frac{1}{\sqrt{1 - \beta^{2}}}\end{matrix}\right]\end{split}\]

In TensorWaves, this class is expressed in a computational backend and it should operate on four-momentum arrays of rank-2. As such, this boost matrix becomes a rank-3 matrix. When using NumPy as backend, the computation looks as follows:

from numpy import array, ones, sqrt, zeros

def _lambdifygenerated(b):
    x0 = 1/numpy.sqrt(1 - b**2)
    x1 = len(b)
    return (array(
            [
                [x0, zeros(x1), zeros(x1), -b*x0],
                [zeros(x1), ones(x1), zeros(x1), zeros(x1)],
                [zeros(x1), zeros(x1), ones(x1), zeros(x1)],
                [-b*x0, zeros(x1), zeros(x1), x0],
            ]
        ).transpose((2, 0, 1)))

Note that this code was generated with sympy.lambdify with cse=True. The repetition of numpy.ones() is still bothersome, but these sub-nodes is also extracted by sympy.cse if the expression is nested further down in an expression tree, for instance when boosting a FourMomentumSymbol \(p\) in the \(z\)-direction:

(5)\[\boldsymbol{B_z}\left(\beta\right) \boldsymbol{R_y}\left(\theta\right) \boldsymbol{R_z}\left(\phi\right) p\]

which in numpy code becomes:

from numpy import array, cos, einsum, ones, sin, sqrt, zeros

def _lambdifygenerated(beta, p, phi, theta):
    x0 = 1/numpy.sqrt(1 - beta**2)
    x1 = len(p)
    x2 = ones(x1)
    x3 = zeros(x1)
    return (einsum("...ij,...jk,...kl,...l->...i", array(
            [
                [x0, x3, x3, -beta*x0],
                [x3, x2, x3, x3],
                [x3, x3, x2, x3],
                [-beta*x0, x3, x3, x0],
            ]
        ).transpose((2, 0, 1)), array(
            [
                [x2, x3, x3, x3],
                [x3, numpy.cos(theta), x3, numpy.sin(theta)],
                [x3, x3, x2, x3],
                [x3, -numpy.sin(theta), x3, numpy.cos(theta)],
            ]
        ).transpose((2, 0, 1)), array(
            [
                [x2, x3, x3, x3],
                [x3, numpy.cos(phi), -numpy.sin(phi), x3],
                [x3, numpy.sin(phi), numpy.cos(phi), x3],
                [x3, x3, x3, x2],
            ]
        ).transpose((2, 0, 1)), p))
class RotationYMatrix(angle: Expr, n_events: Optional[Symbol] = None, **hints: Any)[source]

Bases: UnevaluatedExpression

Rotation matrix around the \(y\)-axis for a FourMomentumSymbol.

Parameters
  • angle – Angle with which to rotate, see e.g. Phi and Theta.

  • n_events – Number of events \(n\) for this matrix array of shape \(n\times4\times4\). Defaults to the len of angle.

The matrix for a rotation over angle \(\alpha\) around the \(y\)-axis operating on FourMomentumSymbol looks like:

(6)\[\begin{split}\boldsymbol{R_y}\left(\alpha\right) = \left[\begin{matrix}1 & 0 & 0 & 0\\0 & \cos{\left(\alpha \right)} & 0 & \sin{\left(\alpha \right)}\\0 & 0 & 1 & 0\\0 & - \sin{\left(\alpha \right)} & 0 & \cos{\left(\alpha \right)}\end{matrix}\right]\end{split}\]

See RotationZMatrix for the computational code.

class RotationZMatrix(angle: Expr, n_events: Optional[Symbol] = None, **hints: Any)[source]

Bases: UnevaluatedExpression

Rotation matrix around the \(z\)-axis for a FourMomentumSymbol.

Parameters
  • angle – Angle with which to rotate, see e.g. Phi and Theta.

  • n_events – Number of events \(n\) for this matrix array of shape \(n\times4\times4\). Defaults to the len of angle.

The matrix for a rotation over angle \(\alpha\) around the \(z\)-axis operating on FourMomentumSymbol looks like:

(7)\[\begin{split}\boldsymbol{R_z}\left(\alpha\right) = \left[\begin{matrix}1 & 0 & 0 & 0\\0 & \cos{\left(\alpha \right)} & - \sin{\left(\alpha \right)} & 0\\0 & \sin{\left(\alpha \right)} & \cos{\left(\alpha \right)} & 0\\0 & 0 & 0 & 1\end{matrix}\right]\end{split}\]

In TensorWaves, this class is expressed in a computational backend and it should operate on four-momentum arrays of rank-2. As such, this boost matrix becomes a rank-3 matrix. When using NumPy as backend, the computation looks as follows:

from numpy import array, cos, ones, sin, zeros

def _lambdifygenerated(a):
    x0 = len(a)
    return (array(
            [
                [ones(x0), zeros(x0), zeros(x0), zeros(x0)],
                [zeros(x0), numpy.cos(a), -numpy.sin(a), zeros(x0)],
                [zeros(x0), numpy.sin(a), numpy.cos(a), zeros(x0)],
                [zeros(x0), zeros(x0), zeros(x0), ones(x0)],
            ]
        ).transpose((2, 0, 1)))

See also the note that comes with Equation (5).

compute_helicity_angles(four_momenta: FourMomenta, topology: Topology) Dict[str, Expr][source]

Formulate expressions for all helicity angles in a topology.

Formulate expressions (Expr) for all helicity angles appearing in a given Topology. The expressions are given in terms of FourMomenta The expressions returned as values in a dict, where the keys are defined by get_helicity_angle_label().

Example

>>> from qrules.topology import create_isobar_topologies
>>> topologies = create_isobar_topologies(3)
>>> topology = topologies[0]
>>> four_momenta = create_four_momentum_symbols(topology)
>>> angles = compute_helicity_angles(four_momenta, topology)
>>> angles["theta_1+2"]
Theta(p1 + p2)
compute_invariant_masses(four_momenta: FourMomenta, topology: Topology) Dict[str, Expr][source]

Compute the invariant masses for all final state combinations.

get_helicity_angle_label(topology: Topology, state_id: int) Tuple[str, str][source]

Generate labels that can be used to identify helicity angles.

The generated subscripts describe the decay sequence from the right to the left, separated by commas. Resonance edge IDs are expressed as a sum of the final state IDs that lie below them (see determine_attached_final_state()). The generated label does not state the top-most edge (the initial state).

Example

The following two allowed isobar topologies for a 1-to-5-body decay illustrates how the naming scheme results in a unique label for each of the eight edges in the decay topology. Note that label only uses final state IDs, but still reflects the internal decay topology.

>>> from qrules.topology import create_isobar_topologies
>>> topologies = create_isobar_topologies(5)
>>> topology = topologies[0]
>>> for i in topology.intermediate_edge_ids | topology.outgoing_edge_ids:
...     phi_label, theta_label = get_helicity_angle_label(topology, i)
...     print(f"{i}: '{phi_label}'")
0: 'phi_0,0+3+4'
1: 'phi_1,1+2'
2: 'phi_2,1+2'
3: 'phi_3,3+4,0+3+4'
4: 'phi_4,3+4,0+3+4'
5: 'phi_0+3+4'
6: 'phi_1+2'
7: 'phi_3+4,0+3+4'
>>> topology = topologies[1]
>>> for i in topology.intermediate_edge_ids | topology.outgoing_edge_ids:
...     phi_label, theta_label = get_helicity_angle_label(topology, i)
...     print(f"{i}: '{phi_label}'")
0: 'phi_0,0+1'
1: 'phi_1,0+1'
2: 'phi_2,2+3+4'
3: 'phi_3,3+4,2+3+4'
4: 'phi_4,3+4,2+3+4'
5: 'phi_0+1'
6: 'phi_2+3+4'
7: 'phi_3+4,2+3+4'

Some labels explained:

  • phi_1+2: edge 6 on the left topology, because for this topology, we have \(p_6=p_1+p_2\).

  • phi_2+3+4: edge 6 right, because for this topology, \(p_6=p_2+p_3+p_4\).

  • phi_1,1+2: edge 1 left, because 1 decays from \(p_6=p_1+p_2\).

  • phi_1,0+1: edge 1 right, because it decays from \(p_5=p_0+p_1\).

  • phi_4,3+4,2+3+4: edge 4 right, because it decays from edge 7 (\(p_7=p_3+p_4\)), which comes from edge 6 (\(p_7=p_2+p_3+p_4\)).

As noted, the top-most parent (initial state) is not listed in the label.

../_images/graphviz_2.svg

topologies[0]

../_images/graphviz_3.svg

topologies[1]

get_invariant_mass_label(topology: Topology, state_id: int) str[source]

Generate an invariant mass label for a state (edge on a topology).

Example

In the case shown in Figure topologies[0], the invariant mass of state \(5\) is \(m_{034}\), because \(p_5=p_0+p_3+p_4\):

>>> from qrules.topology import create_isobar_topologies
>>> topologies = create_isobar_topologies(5)
>>> get_invariant_mass_label(topologies[0], state_id=5)
'm_034'

Naturally, the ‘invariant’ mass label for a final state is just the mass of the state itself:

>>> get_invariant_mass_label(topologies[0], state_id=1)
'm_1'
determine_attached_final_state(topology: Topology, state_id: int) List[int][source]

Determine all final state particles of a transition.

These are attached downward (forward in time) for a given edge (resembling the root).

Example

For edge 5 in Figure topologies[0], we get:

>>> from qrules.topology import create_isobar_topologies
>>> topologies = create_isobar_topologies(5)
>>> determine_attached_final_state(topologies[0], state_id=5)
[0, 3, 4]