Generate transitions

A Partial Wave Analysis starts by defining an amplitude model that describes the reaction process that is to be analyzed. Such a model is generally very complex and requires a fair amount of effort by the analyst (you). This gives a lot of room for mistakes.

The ‘expert system’ is responsible to give you advice on the form of an amplitude model, based on the problem set you define (initial state, final state, allowed interactions, intermediate states, etc.). Internally, the system propagates the quantum numbers through the reaction graph while satisfying the specified conservation rules. How to control this procedure is explained in more detail below.

Afterwards, the amplitude model of the expert system can be exported into TensorWaves. The model can for instance be used to generate a data set (toy Monte Carlo) for this reaction and to optimize its parameters to resemble an actual data set as good as possible. For more info on that see Formulate amplitude model.

Note

Simple channels can be treated with the generate_transitions() façade function. This notebook shows how to treat more complicated cases with the StateTransitionManager.

1. Define the problem set

We first define the boundary conditions of our physics problem, such as initial state, final state, formalism type, etc. and pass all of that information to the StateTransitionManager. This is the main user interface class of qrules.

By default, the StateTransitionManager loads all particles from the PDG. The qrules would take a long time to check the quantum numbers of all these particles, so in this notebook, we use a smaller subset of relatively common particles.

from qrules import InteractionType, StateTransitionManager

stm = StateTransitionManager(
    initial_state=["J/psi(1S)"],
    final_state=["gamma", "pi0", "pi0"],
    formalism_type="helicity",
)

2. Prepare Problem Sets

Create all ProblemSet’s using the boundary conditions of the StateTransitionManager instance. By default it uses the isobar model (tree of two-body decays) to build Topology’s. Various InitialFacts are created for each topology based on the initial and final state. Lastly some reasonable default settings for the solving process are chosen. Remember that each interaction node defines its own set of conservation laws.

The StateTransitionManager (STM) defines three interaction types:

Interaction

Strength

strong

\(60\)

electromagnetic (EM)

\(1\)

weak

\(10^{-4}\)

By default, all three are used in the preparation stage. The create_problem_sets() method of the STM generates graphs with all possible combinations of interaction nodes. An overall interaction strength is assigned to each graph and they are grouped according to this strength.

problem_sets = stm.create_problem_sets()

3. Find solutions

If you are happy with the default settings generated by the StateTransitionManager, just start with solving directly!

This step takes about 23 sec on an Intel(R) Core(TM) i7-6820HQ CPU of 2.70GHz running, multi-threaded.

result = stm.find_solutions(problem_sets)

The find_solutions() method returns a Result object from which you can extract the transitions. Now, you can use get_intermediate_particles() to print the names of the intermediate states that the StateTransitionManager found:

print("found", len(result.transitions), "solutions!")
result.get_intermediate_particles().names
found 270 solutions!
{'a(0)(1450)0',
 'a(0)(980)0',
 'a(1)(1260)0',
 'a(1)(1640)0',
 'b(1)(1235)0',
 'f(0)(1370)',
 'f(0)(1500)',
 'f(0)(1710)',
 'f(0)(500)',
 'f(0)(980)',
 'f(1)(1285)',
 'f(1)(1420)',
 'h(1)(1170)',
 'omega(1420)',
 'omega(1650)',
 'omega(782)',
 'phi(1020)',
 'phi(1680)',
 'rho(1450)0',
 'rho(1700)0',
 'rho(770)0'}

Now we have a lot of solutions that are actually heavily suppressed (they involve two weak decays).

In general, you can modify the ProblemSets returned by create_problem_sets() directly, but the STM also comes with functionality to globally choose the allowed interaction types. So, go ahead and disable the EM and InteractionType.WEAK interactions:

stm.set_allowed_interaction_types([InteractionType.STRONG])
problem_sets = stm.create_problem_sets()
result = stm.find_solutions(problem_sets)

print("found", len(result.transitions), "solutions!")
result.get_intermediate_particles().names
found 84 solutions!
{'b(1)(1235)0',
 'f(0)(1370)',
 'f(0)(1500)',
 'f(0)(1710)',
 'f(0)(500)',
 'f(0)(980)',
 'rho(1450)0',
 'rho(1700)0',
 'rho(770)0'}

Now note that, since a \(\gamma\) particle appears in one of the interaction nodes, the expert system knows that this node must involve EM interactions! Because the node can be an effective interaction, the weak interaction cannot be excluded, as it contains only a subset of conservation laws.

Since only the strong interaction was supposed to be used, this results in a warning and the STM automatically corrects the mistake.

Once the EM interaction is included, this warning disappears. Be aware, however, that the EM interaction is now available globally. Hence, there now might be solutions in which both nodes are electromagnetic.

stm.set_allowed_interaction_types([InteractionType.STRONG, InteractionType.EM])
problem_sets = stm.create_problem_sets()
result = stm.find_solutions(problem_sets)

print("found", len(result.transitions), "solutions!")
result.get_intermediate_particles().names
found 174 solutions!
{'a(0)(1450)0',
 'a(0)(980)0',
 'b(1)(1235)0',
 'f(0)(1370)',
 'f(0)(1500)',
 'f(0)(1710)',
 'f(0)(500)',
 'f(0)(980)',
 'h(1)(1170)',
 'omega(1420)',
 'omega(1650)',
 'omega(782)',
 'phi(1020)',
 'phi(1680)',
 'rho(1450)0',
 'rho(1700)0',
 'rho(770)0'}

Great! Now we selected only the strongest contributions. Be aware, though, that there are more effects that can suppress certain decays, like small branching ratios. In this example, the initial state \(J/\Psi\) can decay into \(\pi^0 + \rho^0\) or \(\pi^0 + \omega\).

decay

branching ratio

\(\omega\to\gamma+\pi^0\)

0.0828

\(\rho^0\to\gamma+\pi^0\)

0.0006

Unfortunately, the \(\rho^0\) mainly decays into \(\pi^0+\pi^0\), not \(\gamma+\pi^0\) and is therefore suppressed. This information is currently not known to the expert system, but it is possible to hand the expert system a list of allowed intermediate states.

# particles are found by name comparison,
# i.e. f2 will find all f2's and f all f's independent of their spin
stm.set_allowed_intermediate_particles(["f(0)", "f(2)"])

result = stm.find_solutions(problem_sets)

print("found", len(result.transitions), "solutions!")
result.get_intermediate_particles().names
found 30 solutions!
{'f(0)(1370)', 'f(0)(1500)', 'f(0)(1710)', 'f(0)(500)', 'f(0)(980)'}

Now we have selected all amplitudes that involve f states. The warnings appear only to notify the user that the list of solutions is not exhaustive: for certain edges in the graph, no suitable particle was found (since only f states were allowed).

import graphviz

from qrules import io

dot = io.asdot(result, collapse_graphs=True, render_node=False)
graphviz.Source(dot)
../_images/reaction_33_0.svg

4. Export generated transitions

The Result, StateTransitionGraph, and Topology can be serialized to and from a dict with io.asdict() and io.fromdict():

from qrules import io

graph = result.transitions[0]
io.asdict(graph.topology)
{'nodes': [0, 1],
 'edges': {-1: {'ending_node_id': 0},
  0: {'originating_node_id': 0},
  1: {'originating_node_id': 1},
  2: {'originating_node_id': 1},
  3: {'originating_node_id': 0, 'ending_node_id': 1}}}

This also means that the Result can be written to JSON or YAML format with io.write() and loaded again with io.load():

io.write(result, "transitions.json")
imported_result = io.load("transitions.json")
assert imported_result == result

Handy if it takes a lot of computation time to re-generate the transitions!

Warning

It’s not possible to pickle a Result, because StateTransitionGraph makes use of Generic.

Tip

The ampform package can formulate amplitude models based on the state transitions created by qrules. See Formulate amplitude model.