phasespace package

phasespace.phasespace module

Implementation of the Raubold and Lynch method to generate n-body events.

The code is based on the GENBOD function (W515 from CERNLIB), documented in:

  1. James, Monte Carlo Phase Space, CERN 68-15 (1968)

class phasespace.phasespace.GenParticle(name: str, mass: Callable | int | float)[source]

Bases: object

Representation of a particle.

Instances of this class can be combined with each other to build decay chains, which can then be used to generate phase space events through the generate or generate_tensor method.

A GenParticle must have
  • a name, which is ensured not to clash with any others in

    the decay chain.

  • a mass, which can be either a number or a function to generate it according to

    a certain distribution. The returned ~`tf.Tensor` needs to have shape (nevents,). In this case, the particle is not considered as having a fixed mass and the has_fixed_mass method will return False.

It may also have:

  • Children, ie, decay products, which are also GenParticle instances.

Parameters:
  • name (str) – Name of the particle.

  • mass (float, Tensor`, callable) – Mass of the particle. If it’s a float, it get converted to a tf.constant.

generate(n_events: int | Tensor | Variable, boost_to: Tensor | None = None, normalize_weights: bool = True, seed: int | Generator | None = None) tuple[Tensor, dict[str, Tensor]][source]

Generate normalized n-body phase space as tensorflow tensors.

Any TensorFlow tensor can always be converted to a numpy array with the method numpy().

Events are generated in the rest frame of the particle, unless boost_to is given.

Note

In this method, the event weights are returned normalized to their maximum.

Parameters:
  • n_events (int) – Number of events to generate.

  • boost_to (optional) – Momentum vector of shape (x, 4), where x is optional, to where the resulting events will be boosted. If not specified, events are generated in the rest frame of the particle.

  • normalize_weights (bool, optional) – Normalize the event weight to its max?

  • seed (SeedLike) – The seed can be a number or a tf.random.Generator that are used as a seed to create a random number generator inside the function or directly as the random number generator instance, respectively.

Returns:

Result of the generation, which varies with the value of normalize_weights:

  • If True, the tuple elements are the normalized event weights as a tensor of shape (n_events, ), and the momenta generated particles as a dictionary of tensors of shape (4, n_events) with particle names as keys.

  • If False, the tuple weights are the unnormalized event weights as a tensor of shape (n_events, ), the maximum per-event weights as a tensor of shape (n_events, ) and the momenta generated particles as a dictionary of tensors of shape (4, n_events) with particle names as keys.

Return type:

tuple

Raises:
  • tf.errors.InvalidArgumentError – If the the decay is kinematically forbidden.

  • ValueError – If n_events and the size of boost_to don’t match.

  • See GenParticle.generate_unnormalized.

generate_tensor(n_events: int, boost_to=None, normalize_weights: bool = True)[source]

Generate normalized n-body phase space as numpy arrays.

Events are generated in the rest frame of the particle, unless boost_to is given.

Note

In this method, the event weights are returned normalized to their maximum.

Parameters:
  • n_events (int) – Number of events to generate.

  • boost_to (optional) – Momentum vector of shape (x, 4), where x is optional, to where the resulting events will be boosted. If not specified, events are generated in the rest frame of the particle.

  • normalize_weights (bool, optional) – Normalize the event weight to its max

Returns:

Result of the generation, which varies with the value of normalize_weights:

  • If True, the tuple elements are the normalized event weights as an array of shape (n_events, ), and the momenta generated particles as a dictionary of arrays of shape (4, n_events) with particle names as keys.

  • If False, the tuple weights are the unnormalized event weights as an array of shape (n_events, ), the maximum per-event weights as an array of shape (n_events, ) and the momenta generated particles as a dictionary of arrays of shape (4, n_events) with particle names as keys.

Return type:

tuple

Raises:
  • tf.errors.InvalidArgumentError – If the the decay is kinematically forbidden.

  • ValueError – If n_events and the size of boost_to don’t match.

  • See GenParticle.generate_unnormalized.

get_mass(min_mass: tf.Tensor = None, max_mass: tf.Tensor = None, n_events: tf.Tensor | tf.Variable = None, seed: SeedLike = None) tf.Tensor[source]

Get the particle mass.

If the particle is resonant, the mass function will be called with the min_mass, max_mass, n_events and optionally a seed parameter.

Parameters:
  • min_mass (tensor) – Lower mass range. Defaults to None, which is only valid in the case of fixed mass.

  • max_mass (tensor) – Upper mass range. Defaults to None, which is only valid in the case of fixed mass.

  • () (n_events) – Number of events to produce. Has to be specified if the particle is resonant.

  • seed (SeedLike) – The seed can be a number or a tf.random.Generator that are used as a seed to create a random number generator inside the function or directly as the random number generator instance, respectively.

Returns:

Mass of the particles, either a scalar or shape (nevents,)

Return type:

Tensor`

Raises:

ValueError – If the mass is requested and has not been set.

property has_children

Does the particle have children?

Type:

bool

property has_fixed_mass

Is the mass a callable function?

Type:

bool

property has_grandchildren

Does the particle have grandchildren?

Type:

bool

set_children(*children)[source]

Assign children.

Parameters:

children (GenParticle) – Two or more children to assign to the current particle.

Returns:

self

Raises:
  • ValueError – If there is an inconsistency in the parent/children relationship, ie,

  • if children were already set, if their parent was or if less than two children were given.

  • KeyError – If there is a particle name clash.

  • RuntimeError – If generate was already called before.

class phasespace.phasespace.Particle[source]

Bases: object

Deprecated Particle class.

Renamed to GenParticle.

phasespace.phasespace.generate_decay(*args, **kwargs)[source]

Deprecated.

phasespace.phasespace.nbody_decay(mass_top: float, masses: list, top_name: str = '', names: list = None)[source]

Shortcut to build an n-body decay of a GenParticle.

If the particle names are not given, the top particle is called ‘top’ and the children ‘p_{i}’, where i corresponds to their position in the masses sequence.

Parameters:
  • mass_top (tensor, list) – Mass of the top particle. Can be a list of 4-vectors.

  • masses (list) – Masses of the child particles.

  • name_top (str, optional) – Name of the top particle. If not given, the top particle is named top.

  • names (list, optional) – Names of the child particles. If not given, they are build as ‘p_{i}’, where i is given by their ordering in the masses list.

Returns:

Particle decay.

Return type:

GenParticle

Raises:

ValueError – If the length of masses and names doesn’t match.

phasespace.phasespace.pdk(a, b, c)[source]

Calculate the PDK (2-body phase space) function.

Based on Eq. (9.17) in CERN 68-15 (1968).

Parameters:
  • a (Tensor`) – \(M_{i+1}\) in Eq. (9.17).

  • b (Tensor`) – \(M_{i}\) in Eq. (9.17).

  • c (Tensor`) – \(m_{i+1}\) in Eq. (9.17).

Returns:

~`tf.Tensor`

phasespace.phasespace.process_list_to_tensor(lst)[source]

Convert a list to a tensor.

The list is converted to a tensor and transposed to get the proper shape.

Note

If lst is a tensor, nothing is done to it other than convert it to tf.float64.

Parameters:

lst (list) – List to convert.

Returns:

~`tf.Tensor`

phasespace.kinematics module

Basic kinematics.

phasespace.kinematics.beta(vector)[source]

Calculate beta of a given 4-vector.

Parameters:

vector – Input Lorentz momentum vector.

phasespace.kinematics.boost_components(vector)[source]

Get the boost components of a given 4-vector.

Parameters:

vector – Input Lorentz momentum vector.

phasespace.kinematics.lorentz_boost(vector, boostvector)[source]

Perform Lorentz boost.

Parameters:
  • vector – 4-vector to be boosted

  • boostvector – Boost vector. Can be either 3-vector or 4-vector, since only spatial components are used.

phasespace.kinematics.lorentz_vector(space, time)[source]

Make a Lorentz vector from spatial and time components.

Parameters:
  • space – 3-vector of spatial components.

  • time – Time component.

phasespace.kinematics.mass(vector)[source]

Calculate mass scalar for Lorentz 4-momentum.

Parameters:

vector – Input Lorentz momentum vector.

phasespace.kinematics.metric_tensor()[source]

Metric tensor for Lorentz space (constant).

phasespace.kinematics.scalar_product(vec1, vec2)[source]

Calculate scalar product of two 3-vectors.

Parameters:
  • vec1 – First vector.

  • vec2 – Second vector.

phasespace.kinematics.spatial_component(vector)[source]

Extract spatial components of the input Lorentz vector.

Parameters:

vector – Input Lorentz vector (where indexes 0-2 are space, index 3 is time).

phasespace.kinematics.time_component(vector)[source]

Extract time component of the input Lorentz vector.

Parameters:

vector – Input Lorentz vector (where indexes 0-2 are space, index 3 is time).

phasespace.kinematics.x_component(vector)[source]

Extract spatial X component of the input Lorentz or 3-vector.

Parameters:

vector – Input vector.

phasespace.kinematics.y_component(vector)[source]

Extract spatial Y component of the input Lorentz or 3-vector.

Parameters:

vector – Input vector.

phasespace.kinematics.z_component(vector)[source]

Extract spatial Z component of the input Lorentz or 3-vector.

Parameters:

vector – Input vector.

phasespace.fromdecay module

class phasespace.fromdecay.genmultidecay.GenMultiDecay(gen_particles: list[tuple[float, GenParticle]])[source]

Bases: object

A GenParticle-type container that can handle multiple decays.

Parameters:

gen_particles – All the GenParticles and their corresponding probabilities. The list must be of the format [[probability, GenParticle instance], [probability, …

DEFAULT_MASS_FUNC = 'relbw'
MASS_WIDTH_TOLERANCE = 0.01
classmethod from_dict(dec_dict: dict, mass_converter: dict[str, Callable] = None, tolerance: float = None, particle_model_map: dict[str, str] = None)[source]

Create a GenMultiDecay instance from a dict in the DecayLanguage package format.

Parameters:
  • dec_dict – The input dict from which the GenMultiDecay object will be created from. A dec_dict has the same structure as the dicts used in DecayLanguage, see the examples below.

  • mass_converter – A dict with mass function names and their corresponding mass functions. These functions should take the particle mass and the mass width as inputs and return a mass function that phasespace can understand. This dict will be combined with the predefined mass functions in this package. See the Example below or the tutorial for how to use this parameter.

  • tolerance – Minimum mass width of the particle to use a mass function instead of assuming the mass to be constant. If None, the default value, defined by the class variable MASS_WIDTH_TOLERANCE, will be used. This value can be customized if desired.

  • particle_model_map – A dict where the key is a particle name and the value is a mass function name. If a particle is specified in the particle_model_map, then all appearances of this particle in dec_dict will get the same mass function. This way, one does not have to manually add the zfit parameter to every place where this particle appears in dec_dict. If the zfit parameter is specified for a particle which is also included in particle_model_map, the zfit parameter mass function name will be prioritized.

Returns:

The created GenMultiDecay object.

Examples

DecayLanguage is usually used to create a dict that can be understood by GenMultiDecay. >>> from decaylanguage import DecFileParser >>> from phasespace.fromdecay import GenMultiDecay

Parse a .dec file to create a DecayLanguage dict describing a D*+ particle that can decay in 2 different ways: D*+ -> D0(->K- pi+) pi+ or D*+ -> D+ gamma.

>>> parser = DecFileParser('example_decays.dec')
>>> parser.parse()
>>> dst_chain = parser.build_decay_chains("D*+")
>>> dst_chain
{'D*+': [{'bf': 0.984,
    'fs': [{'D0': [{'bf': 1.0,
            'fs': ['K-', 'pi+']}]},
        'pi+']},
    {'bf': 0.016,
     'fs': ['D+', 'gamma']}]}

If the D0 particle should have a mass distribution of a gaussian when it decays, one can pass the particle_model_map parameter to from_dict: >>> dst_gen = GenMultiDecay.from_dict(dst_chain, particle_model_map={“D0”: “gauss”}) This will then set the mass function of D0 to a gaussian for all its decays.

If more custom control is required, e.g., if D0 can decay in multiple ways and one of the decays should have a specific mass function, a zfit parameter can be added to its decay dict: >>> dst_chain[“D*+”][0][“fs”][0][“D0”][0][“zfit”] = “gauss” >>> dst_chain {‘D*+’: [{‘bf’: 0.984,

‘fs’: [{‘D0’: [{‘bf’: 1.0,

‘fs’: [‘K-’, ‘pi+’], ‘zfit’: ‘gauss’}]},

‘pi+’]},

{‘bf’: 0.016, ‘fs’: [‘D+’, ‘gamma’]}]}

This dict can then be passed to GenMultiDecay.from_dict: >>> dst_gen = GenMultiDecay.from_dict(dst_chain) This will now convert make the D0 particle have a gaussian mass function, only when it decays into K- and pi+. In this case, there are no other ways that D0 can decay, so using particle_model_map is a cleaner and easier option.

If the decay of the D0 particle should be modelled by a mass distribution that does not come with the package, a custom distribution can be created: >>> def custom_gauss(mass, width): >>> particle_mass = tf.cast(mass, tf.float64) >>> particle_width = tf.cast(width, tf.float64) >>> def mass_func(min_mass, max_mass, n_events): >>> min_mass = tf.cast(min_mass, tf.float64) >>> max_mass = tf.cast(max_mass, tf.float64) >>> # Use a zfit PDF >>> pdf = zfit.pdf.Gauss(mu=particle_mass, sigma=particle_width, obs=””) >>> iterator = tf.stack([min_mass, max_mass], axis=-1) >>> return tf.vectorized_map( >>> lambda lim: pdf.sample(1, limits=(lim[0], lim[1])), iterator >>> ) >>> return mass_func

Once again change the distribution in the dst_chain dict. Here, it is changed to “custom_gauss” but any name can be used. >>> dst_chain[“D*+”][0][“fs”][0][“D0”][0][“zfit”] = “custom_gauss”

One can then pass the custom_gauss function and its name (in this case “custom_gauss”) as a dict`to `from_dict as the mass_converter parameter: >>> dst_gen = GenMultiDecay.from_dict(dst_chain, mass_converter={“custom_gauss”: custom_gauss})

Notes

For a more in-depth tutorial, see the tutorial on GenMultiDecay in the [documentation](https://phasespace.readthedocs.io/en/stable/GenMultiDecay_Tutorial.html).

generate(n_events: int, normalize_weights: bool = True, **kwargs) tuple[list[Tensor], list[Tensor]] | tuple[list[Tensor], list[Tensor], list[Tensor]][source]

Generate four-momentum vectors from the decay(s).

Parameters:
  • n_events – Total number of events combined, for all the decays.

  • normalize_weights – Normalize weights according to all events generated. This also changes the return values. See the phasespace documentation for more details.

  • kwargs – Additional parameters passed to all calls of GenParticle.generate

Returns:

The arguments returned by GenParticle.generate are returned. See the phasespace documentation for details. However, instead of being 2 or 3 tensors, it is 2 or 3 lists of tensors, each entry in the lists corresponding to the return arguments from the corresponding GenParticle instances in self.gen_particles. Note that when normalize_weights is True, the weights are normalized to the maximum of all returned events.

phasespace.fromdecay.mass_functions.breitwigner_factory(mass, width)[source]
phasespace.fromdecay.mass_functions.gauss_factory(mass, width)[source]
phasespace.fromdecay.mass_functions.relativistic_breitwigner_factory(mass, width)[source]