pymule reference guide

This section describes all functions and classes in pymule. Most users will not have to view this.

Working with errors

pymule.errortools.addplots(a: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, b: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, sa: float = 1.0, sb: float = 1.0) → <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>

adds or subtracts two plots

Parameters:

  • a – Nx3 numpy matrix; the first plot

  • b – Nx3 numpy matrix; the second plot

  • sa – float, optional; the coefficient of the first plot

  • sb – float, optional; the coefficient of the second plot

Returns:

a Nx3 numpy matrix with \(s_a\cdot a + s_b\cdot b\)

Note

a and b must share x values, otherwise entries are dropped

Example:

subtract two plots a and b

>>> addplots(a, b, sb=-1)
Example:

Given the LO plots thetaLO and the NLO corrections thetadNLO, we calculate the \(K\) factor as either

>>> thetaNLO = addplots(thetaLO, thetadNLO)
pymule.errortools.chisq(values: ~typing.List[~typing.List[float] | <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>] | <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → float

calculates the \(\chi^2/\textrm{d.o.f.}\) of numbers

Parameters:

value – Nx2 numpy matrix or list of lists; the values as [[y1, dy1], [y2, dy2], ...]

Returns:

float; the \(\chi^2/\textrm{d.o.f.} = \frac{1}{n} \sum_{n=1}^n(\frac{y_i-\bar y}{\delta y_i})^2\) with the average value \(\bar y\)

Example:

a good example

>>> chisq([[20.0, 0.8],
...        [21.6, 0.9],
...        [18.7, 1.2]])
1.3348808062205872

and a bad example

>>> chisq([[16.2, 0.8],
...        [22.9, 0.9],
...        [8.81, 1.2]])
30.173852184366673
pymule.errortools.combineNplots(func: ~typing.Callable[[<MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>], <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>], plots: ~typing.List[<MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>] | <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>

combines a list of plots using a function

Parameters:

  • func – callable with two arguments; the function to combine the plots

  • plots – list of Nx3 numpy matrices

Returns:

a Nx3 numpy matrix \(f(p_0, f(p_1, f(p_2, \cdots)))\)

pymule.errortools.combineplots(a: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, b: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, yfunc: ~typing.Callable[[<MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>], <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>], efunc: ~typing.Callable[[<MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>], <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>], tol: float = 1e-08) → <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>

combines two plots using functions for the value and the error

Parameters:

  • a – Nx3 numpy matrix; the first plot

  • b – Nx3 numpy matrix; the second plot

  • yfunc – callable; a function to calculate the value yfunc(a, b)

  • efunc – callable; a function to calculate the squared error efunc(a, da**2 b, db**2)

  • tol – float, optional; the difference at which values are considered equal

Returns:

Nx3 numpy matrix; the combined plot np.array([[x1, yfunc(..), sqrt(efunc(..))], ..])

Note

a and b must share x values, up to the tolerance tol, otherwise values may be dropped

Example:

Add two plots A and B

>>> combineplots(A, B,
...              lambda a, b: a+b,
...              lambda a, da2, b, db2: da2 + db2))

Calculate a K factor

>>> combineplots(dnlo, lo,
...              lambda a, b: 1 + a/b,
...              lambda a, da2, b, db2: db2 * a**2 / b**4 + da2 / b**2)
pymule.errortools.dividenumbers(a: ~typing.List[float] | <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, b: ~typing.List[float] | <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → __getitem__()' id='140397293428560'>

divides numbers

Parameters:

  • a – list of floats; the numerator with error [a, da]

  • b – list of floats; the denominator with error [b, db]

Returns:

the result of the division a/b [y, dy]

Example:

Divide \((2.3\pm0.1) / (45\pm0.01)\)

>>> dividenumbers([2.3, 0.1], [45., 0.01])
array([0.05111111, 0.00222225])
pymule.errortools.divideplots(a: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, b: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, offset: float = 0.0) → <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>

divides two plots

Parameters:

  • a – Nx3 numpy matrix; the numerator plot

  • b – Nx3 numpy matrix; the denominator plot

  • offset – float, optional; shifts the result

Returns:

a Nx3 numpy matrix with \(a/b + offset\)

Note

a and b must share x values, otherwise entries are dropped

Example:

Given the LO plots thetaLO and the NLO corrections thetadNLO, we calculate the \(K\) factor as either

>>> thetaNLO = addplots(thetaLO, thetadNLO)
>>> thetaK = divideplots(thetaNLO, thetaLO)
>>> thetaK = divideplots(thetadNLO, thetaLO, offset=+1.)
pymule.errortools.integratehistogram(hist: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → float

integrates a histogram

Parameters:

hist – Nx3 numpy matrix; the histogram to integrate \(d\sigma/dx\) as np.array([[x1, y1, e1], [x2, y2, e2], ...])

Returns:

float; the integrated histogram \(\int d\sigma/dx dx\) without error estimate

Example:

Integrate a histogram

>>> hist
array([[          -inf, 0.00000000e+00, 0.00000000e+00],
       [5.00000000e-02, 4.77330751e+01, 2.26798977e-01],
       [1.50000000e-01, 7.40641192e+01, 2.36498021e-01],
       ...,
       [8.85000000e+00, 1.67513948e+00, 1.16218116e-01],
       [8.95000000e+00, 0.00000000e+00, 0.00000000e+00],
       [           inf, 0.00000000e+00, 0.00000000e+00]])
>>> integratehistogram(hist)
4188.519369660588
pymule.errortools.mergebins(p: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, n: int) → <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>

merges n adjacent bins into one larger bin, reducing the uncertainty.

Parameters:

  • p – Nx3 numpy matrix; the plot

  • n – int; how many bins to merge

Returns:

a (N/n)x3 numpy matrix

Note

This process loses len(p)%n bins at the end of the histogram

Example:

merge five bins

>>> len(p)
200
>>> len(mergebins(p, 5))
40

Bins may be lost

>>> len(p)
203
>>> len(mergebins(p, 5))
40
pymule.errortools.mergenumbers(values: ~typing.List[~typing.List[float] | <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>] | <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → __getitem__()' id='140397293428560'>

statistically combines values with uncertainties

Parameters:

values – Nx2 numpy matrix or list of lists; the values as [[y1, dy1], [y2, dy2], ...]

Returns:

either answer as numpy array [y, dy] or tuple of \(chi^2\) and answer

Example:

If quiet is not specified this will print the \(chi^2\)

>>> mergenumbers([[20.0, 0.8],
...               [21.6, 0.9],
...               [18.7, 1.2]])
1.3348808062205872
array([20.30718232,  0.53517179])

Otherwise, it will return it

>>> mergenumbers([[20.0, 0.8],
...               [21.6, 0.9],
...               [18.7, 1.2]], quiet=True)
(1.3348808062205872, array([20.30718232,  0.53517179]))
pymule.errortools.mergeplots(ps: ~typing.List[<MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>] | <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>

statistically combines a list of plots

Parameters:

ps – list of Nx3 numpy matrices; the plots to combine as [np.array([[x1, y1, e1], [x2, y2, e2], ...]), ...]

Returns:

a Nx3 numpy matrix

Example:

Load a number of vegas files and merge them

>>> data = [
...     importvegas(i)['thetae']
...     for i in glob.glob('out/em2em0*')
... ]
>>> mergeplots(data)
pymule.errortools.plusnumbers(*args: ~typing.List[float] | <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → __getitem__()' id='140397293428560'>

adds numbers and errors

Parameters:

yi – list of floats; a number with error [yi, dyi]

Returns:

the result of the addition [y, dy]

Example:

Adding \((10\pm1)+(20\pm0.5)+(-5\pm2)\)

>>> plusnumbers([10, 1], [20, 0.5], [-5, 2])
array([25.        ,  2.29128785])
pymule.errortools.printnumber(x: ~typing.List[float] | <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, prec: int = 0) → str

returns a string representation of a number with uncertainties to one significant digit

Parameters:

  • x – a list with two floats; the number as [x, dx]

  • prec – int, otpional; number of extra signficant figures

Returns:

str; the formatted string

Example:

printing \(53.2\pm0.1\) to one significant figure

>>> printnumber([53.2, 0.1])
"53.2(1)"
pymule.errortools.scaleplot(a: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, sx: float, sy: float | None = None) → <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>

rescales a plot such that the integrated plot remains unchanged, i.e. rescale \(x\to x/s\) and \(y\to y\cdot s\). This is useful to, for example, change units.

Parameters:

  • a – Nx3 numpy matrix; the plot

  • sx – float; the inverse scale factor for the x direction

  • sy – float, optional; if present, sy will be used for the y direction instead of sx

Returns:

a Nx3 numpy matrix

Example:

rescaling units from rad to mrad

>>> scaleplot(data, 1e-3)
pymule.errortools.timesnumbers(a: ~typing.List[float] | <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, b: ~typing.List[float] | <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → __getitem__()' id='140397293428560'>

multiplies numbers

Parameters:

  • a – list of floats; the first factor with error [a, da]

  • b – list of floats; the second factor with error [b, db]

Returns:

the result of the multiplication a*b [y, dy]

Example:

Divide \((0.5\pm0.02) * (45\pm0.01)\)

>>> timesnumbers([0.5,0.02], [45, 0.1])
array([22.5       ,  0.90138782])

Working with abstract records

class pymule.record.Chi2

An object to represent a \(\chi^2\) value.

This can be created using Chi2.from_single() or Chi2.from_merge(). This class supports addition and float.

Attribute type:

an enum to indicate the type of \(\chi^2\) this is

Attribute value:

the value of \(\chi^2\) for single and mergers

Attribute children:

the children Chi2 objects for sums and mergers

class Type(value, names=None, *, module=None, qualname=None, type=None, start=1, boundary=None)

The \(\chi^2\) type enum

Attribute NONE:

empty entry

Attribute SINGLE:

a single entry

Attribute SUM:

a sum of entries

Attribute MERGE:

a statistical merge of entries

classmethod from_merge(chi: float, *s: Chi2) → Chi2

combines multiple \(\chi^2\) into a single one

Parameters:

  • chi – the value of \(\chi^2\) of the merge, must be positive

  • s – a sequence of Chi2 objects

Returns:

a new Chi2 object

Return type:

Chi2

Raises:

  • ValueError – the parameter chi is not positive or no Chi2 is be provided in s

  • TypeError – the parameter chi is not a number or the objects provided are not Chi2

Example:

Provide a merged set

>>> c1 = Chi2.from_single(1.02)
>>> c2 = Chi2.from_single(1.41)
>>> c3 = Chi2.from_single(2.02)
>>> tot = Chi2.from_merge(1.3, c1, c2, c3)
classmethod from_single(v: float) → Chi2

creates a single Chi2 object

Parameters:

v – the value of \(\chi^2\), must be positive

Returns:

a new Chi2 object

Return type:

Chi2

Raises:

  • ValueError – the parameter v is not positive

  • TypeError – the parameter v is not a number

Example:

Create a \(\chi^2 = 1.02\)

>>> Chi2.from_single(1.02)
<Chi^2: single 1.02>
class pymule.record.Record

A McMule record, either corresponding to a single file or multiple.

This should usually not be instantiated directly but using either VegasRecord or XiRecord This supports addition, subtraction, and scalar multiplication.

Properties can be accessed as attributes or as keys. Histograms are accessed as keys.

Example:

Assuming that r0, rV, and rR contains records of leading order, virtual, and real corrections (as obtained using pymule.vegas.FileRecord.read() or mergefks())

>>> total = r0 + alpha * rV + alpha * rR

Add the VP contributions from rA to this by

>>> total += alpha * rA

Print the cross section and values

>>> print(total.value)
[1.00646123e-04 1.06589201e-08]

>>> total.histograms
dict_keys(['thetae', 'thetam'])

>>> total.chi2a
<Chi^2: single 1.0530713790136834>

>>> total['thetae']
array([[           -inf,  0.00000000e+00,  0.00000000e+00],
       [ 1.75000000e-03, -1.92701245e-03,  1.27052491e-05],
       [ 5.25000000e-03,  1.57686769e-03,  2.08601107e-05],
       [ 8.75000000e-03,  3.10660221e-03,  2.45339189e-05],
       [ 1.22500000e-02,  3.90154737e-03,  2.68166764e-05],
       [ 1.57500000e-02,  5.04369602e-03,  2.64601744e-05],
       [ 1.92500000e-02,  7.87371009e-03,  2.37634926e-05],
       [ 2.27500000e-02,  8.54989271e-03,  1.14554010e-05],
       [ 2.62500000e-02,  4.47774416e-04,  1.48451350e-06],
       [            inf,  7.29295875e-08,  1.24592720e-09]])
property chi2a: Chi2

the \(\chi^2\) estimate of the integrator or combination using Chi2

property histograms: list[str]

a list of histograms included in this record

classmethod merge(*s: Record) → Record

statistically merges a set of runs, combining cross sections, histograms, and run-time information.

Parameters:

s – a sequence of Record objects (or subclasses like VegasRecord)

Returns:

a new record containing cross sections, distributions, and run-time information.

Return type:

Record

Raises:

  • ValueError – no record is provided

  • TypeError – not all recrods are Record objects

Note

The result only contains histograms that are in all records. Entries that appear less often are dropped

Example:

Assuming that r1, r2, and r3 are statistically independent records that should be merged

>>> r = Record.merge(r1, r2, r3)

alternatively

>>> records = [r1, r2, r3]
>>> r = Record.merge(*records)
property time: float

the run time required by the record in seconds, either for a single job or combinations

property value: __getitem__()' id='140397293428560'>

the best estimate for the cross section and its error as np.array([y, e])

Working with vegas records

class pymule.vegas.BareRecordV4

A subclass of FileRecord for v4 bare files.

class pymule.vegas.FileRecord

A subclass of Record to track McMule output files.

Additionally to the attributes of Record this contains all relevant properties of a McMule run

Attribute version:

the file version

Attribute inttype:

the integer type used in McMule in the format for struct

Attribute sha:

the sha hash provided by the file

Attribute sampler:

the sampler used, currently either vegas or foam

Attribute iteration:

the current iteration of McMule

Attribute randy:

the current random seed

Attribute msg:

the integration message

Example:

Load a record stored at /path/to/record.mcmule and print its value

>>> rec = FileRecord.read("/path/to/record.mcmule")
>>> print(rec.value)
[1.00646123e-04 1.06589201e-08]
classmethod read(filename: str = '', fp: IO[bytes] | None = None) → FileRecord

opens a McMule file, infers the correct file version, and creates a subclass of FileRecord

Parameters:

  • filename – file name to open, optional

  • fp – file pointer to read, optional

Returns:

a recording corresponding to the file

Return type:

subclass of FileRecord

Note

Either filename xor fp need to be specified

Note

If less than two iterations have been completed, no histograms will be returned

Example:

Load a file for the muon decay

>>> FileRecord.read('m2ennRR_mu-e_S0000068031X0.50000D0.50000_ITMX080x150M_O12.vegas')
<Record m2ennRR_mu-e_S0000068031X0.50000D0.50000_ITMX080x150M_O12.vegas xsec = 0.00010065(1), 12 histograms>
class pymule.vegas.FoamRecordV4

A subclass of FileRecord to track Foam files. Additionally to the attributes of Record and FileRecord this contains relevant information for foam

Attribute active:

mask of active cells

Attribute parents:

list of parents of cells or -1 for the root cell

Attribute offset:

offset of the bottom right corner of each cell

Attribute length:

size of each cell

Attribute probability:

probability that a cell is selected for sampling

Attribute nbins:

number of bins in used in cell division

write(filename: str = '', fp: IO[bytes] | None = None) → None

writes this record to a file or file pointer

Parameters:

  • filename – file name to open, optional

  • fp – file pointer to write to, optional

Note

Either filename xor fp need to be specified

Example:

save a random run to disk

>>> rec.write("out.vegas")
class pymule.vegas.VegasRecord

A subclass of FileRecord to track Vegas files. Additionally to the attributes of Record and FileRecord this contains relevant information for vegas

Attribute ndo:

number of subdivisions on an axis

Attribute xi:

the vegas grid, location of the \(i\)-th division on the \(j\)-th axis, normalised to lie between 0 and 1.

class pymule.vegas.VegasRecordV1

A subclass of VegasRecord for v1 Vegas files

class pymule.vegas.VegasRecordV2

A subclass of VegasRecord for v2 Vegas files

class pymule.vegas.VegasRecordV3

A subclass of VegasRecord for v3 Vegas files

class pymule.vegas.VegasRecordV4

A subclass of VegasRecord for v4 Vegas files. Since this is the most recent version it supports writing as well as reading

write(filename: str = '', fp: IO[bytes] | None = None) → None

writes this record to a file or file pointer

Parameters:

  • filename – file name to open, optional

  • fp – file pointer to write to, optional

Note

Either filename xor fp need to be specified

Example:

save a random run to disk

>>> rec.write("out.vegas")
pymule.vegas.guess_version(fp: IO[bytes], inttype: str = 'i') → tuple[int, str]

infers version of the vegas file using either the version string (since v3) or the file length (v1 and v2).

Parameters:

  • fp – file pointer

  • inttype – either 'i' or 'q', optional; the integer type to use for v1 or v2, inferred otherwise

Returns:

tuple of version number and integer type

pymule.vegas.read_record(fp: IO[bytes], typ: str) → int | float | bytes | List[int | float]

reads a single FORTRAN record

Parameters:

  • fp – file pointer

  • typ –

    the type to read, everything that struct understands. Examples are

    • i: 32 bit signed integer (standard integer in Fortran)

    • I: 32 bit unsigned integer

    • q: 64 bit signed integer (integer*8 in Fortran)

    • c: 8 bit charater (character in Fortran)

    • d: 64 bit double precision (real(kind=prec) in Fortran,

      with default prec)

Records are data structures that are build as follows:

  • 4 byte header: length of the record as a 32 bit unsigned integer, called l1

  • body of length l1

  • 4 byte footer: a repetition of l1 to make sure the record is properly closed.

Records can contain multiple variables.

pymule.vegas.write_record(fp: IO[bytes], typ: str, content: int | float | bytes | List[int | float] | str) → None

writes a single FORTRAN record

Parameters:

  • fp – file pointer

  • typ –

    the type to read, everything that struct understands. Examples are

    • i: 32 bit signed integer (standard integer in Fortran)

    • I: 32 bit unsigned integer

    • q: 64 bit signed integer (integer*8 in Fortran)

    • c: 8 bit charater (character in Fortran)

    • d: 64 bit double precision (real(kind=prec) in Fortran,

      with default prec)

  • content – scalar, list, str or bytes; the data to write

Records are data structures that are build as follows:

  • 4 byte header: length of the record as a 32 bit unsigned integer, called l1

  • body of length l1

  • 4 byte footer: a repetition of l1 to make sure the record is properly closed.

Records can contain multiple variables.

Working with records of data

class pymule.loader.FolderLoader(path: str)

A Loader than can act on folders using os.listdir

classmethod can_open(path: str) → bool

checks if this loader can open a given path

Parameters:

path – the path to check

Returns:

a bool to indicate whether this loader can open the path

list() → list[str]

lists all files known to this loader

Returns:

an iterable of all files known

open1(file: str) → FileRecord

opens a file

Parameters:

file – the file path to be opened

Returns:

the loaded record

Return type:

FileRecord

class pymule.loader.GitLabLoader(path: str)

A Loader than can read tar balls from GitLab

classmethod can_open(path: str) → bool

checks if this loader can open a given path

Parameters:

path – the path to check

Returns:

a bool to indicate whether this loader can open the path

list() → list[str]

lists all files known to this loader

Returns:

an iterable of all files known

open1(file: str) → FileRecord

opens a file

Parameters:

file – the file path to be opened

Returns:

the loaded record

Return type:

FileRecord

class pymule.loader.Loader(path: str)

An abstract class to load vegas files from various sources.

Implementations need to provide Loader.can_open(), Loader.list(), and Loader.open1()

abstract classmethod can_open(path: str) → bool

checks if this loader can open a given path

Parameters:

path – the path to check

Returns:

a bool to indicate whether this loader can open the path

static chain(loaders: list[Loader], pattern: Pattern[str]) → Generator[tuple[Any, FileRecord], None, None]

chains a list of loaders together using their Loader.open()

Parameters:

  • loaders – a list of loaders to chain

  • pattern – the pattern to match in each loader

Returns:

a generator of FileRecord objects that match the pattern

Example:
>>> for rec in Loader.chain([l1, l2, l3], pattern):
...     print(rec)
abstract list() → Generator[str, None, None] | List[str]

lists all files known to this loader

Returns:

an iterable of all files known

open(pattern: Pattern[str]) → Generator[tuple[Any, FileRecord], None, None]

checks all files known to this loader and opens all that match the provided pattern

Parameters:

pattern – a regular expression to match the filename

Returns:

a generator of FileRecord objects that match the pattern

Example:
>>> for rec in loader.open(pattern):
...     print(rec)
abstract open1(file: str) → FileRecord

opens a file

Parameters:

file – the file path to be opened

Returns:

the loaded record

Return type:

FileRecord

class pymule.loader.Options

Options used for loading and merging records

Attribute merge:

a dict of histograms {'name': n} to merge n bins in the histogram name. defaults to no merging

Attribute sanitycheck:

callable; a function that, given a FileRecord, whether to include the file in the output (return True) or to skip (return False).

Attribute flavour:

str; the flavour to load, defaults to everything Initialised to everything, i.e. .*.

Attribute obs:

str; the observable to load (the bit after the O), defaults to everything Initialised to everything, i.e. ''.

Attribute folderp:

str; a regular expression to match directory structures of a tar file, defaults to everything Initialised to everything, i.e. .*.

Attribute KNOWN_LOADERS:

a list of Loader subclasses to try

add_folder(folder: str) → None

adds the first Loader that can open a folder or file

Parameters:

folder – the path to folder or file try add

Raises:

NotImplemented – no loader can open this object

reset() → None

resets this loader to the default state

class pymule.loader.TarLoader(path: str)

A Loader than can read tar balls, potentially compressed using BZip2 or GZip

classmethod can_open(path: str) → bool

checks if this loader can open a given path

Parameters:

path – the path to check

Returns:

a bool to indicate whether this loader can open the path

list() → list[str]

lists all files known to this loader

Returns:

an iterable of all files known

open1(file: str) → FileRecord

opens a file

Parameters:

file – the file path to be opened

Returns:

the loaded record

Return type:

FileRecord

class pymule.loader.XiRecord

A collection of McMule records, keyed by seed or \(\xi_c\)

This should be created using XiRecord.load() or sigma() which loads all matching records into this object. Default options can be changed using the options attribute. Records can then be merged or modified. This supports addition, subtraction, and scalar multiplication.

Specific values of \(\xi_c\) can be accessed as keys.

Attribute options:

An instance of Options with the defaults to use. This can be modified

Example:

Load all records for piece em2emFEE in folder out

>>> XiRecord.options.folder = "out/"
>>> rf = XiRecord.load('em2emFMM')
>>> rf
<XiRecord set em2emFMM, xic = [0.3, 0.5, 1.0, 0.1], 12 histograms>
>>> rr = XiRecord.load('em2emRMM', obs="0")
>>> tot = rf + rr
>>> tot
<XiRecord set em2emFMM + em2emRMM, xic = [0.3, 0.5, 1.0, 0.1], 12 histograms>
>>> tot.xic
{(0.1, 0.1), (0.3, 0.3), (0.5, 0.5), (1.0, 1.0)}
>>> tot[1.,1.]
<Record em2emFMM xsec = -0.000232279(3), 12 histograms>
as_record() → Record

returns the record if there is only one \(\xi_c\) record

returns:

the record

rtype:

Record

raises IndexError:

data is not FKS merged

Example:

Open a leading order calculation

>>> XiRecord.load('em2em0')
<XiRecord set em2em0, xic = [1.0], 12 histograms>
>>> XiRecord.load('em2em0').as_record()
<Record em2em0 xsec = 0.0058559441(2), 12 histograms>
property histograms: set[str]

the histograms present in in all records. Entries that appear less often are dropped

classmethod load(piece: str, folder: str | None = None, **kwargs: Any) → XiRecord

imports all matching records files

Parameters:

  • piece – str; which_piece to load

  • folder – str, optional; file name, optional; folder or tarball to search for vegas files

  • merge – a dict of histograms {'name': n} to merge n bins in the histogram name. defaults to no merging

  • sanitycheck – callable; a function that, given a FileRecord, whether to include the file in the output (return True) or to skip (return False).

  • flavour – str; the flavour to load, defaults to everything Initialised to everything, i.e. .*.

  • obs – str; the observable to load (the bit after the O), defaults to everything Initialised to everything, i.e. ''.

  • folderp – str; a regular expression to match directory structures of a tar file, defaults to everything Initialised to everything, i.e. .*.

Returns:

an in-memory representation of the records

Return type:

XiRecord

Note

Use setup() or modify XiRecord.options to set the defaults. Arguments provided here override the defaults

Example:
>>> XiRecord.load('em2emFMM')
<XiRecord set em2emFMM, xic = [0.3, 0.5, 1.0, 0.1], 12 histograms>
>>> XiRecord.load('em2emRMM', obs="0")
<XiRecord set em2emRMM, xic = [0.3, 0.5, 1.0, 0.1], 12 histograms>
property merged_seeds: dict[tuple[float, float], Record]

the random-seed-merged Record for each \(\xi_c\) combination

Example:
>>> rr
<XiRecord set em2emRMM, xic = [0.3, 0.5, 1.0, 0.1], 12 histograms>
>>> rr.merged_seeds
{(0.3, 0.3): <Record em2emRMM xsec = 0.00037124(3), 12 histograms>,
 (0.5, 0.5): <Record em2emRMM xsec = 0.00002654(2), 12 histograms>,
 (1.0, 1.0): <Record em2emRMM xsec = -0.0004411(3), 12 histograms>,
 (0.1, 0.1): <Record em2emRMM xsec = 0.00111258(3), 12 histograms>}
property seeds: set[int]

a set of seeds in this record

Example:
>>> rr.seeds
{96012, 96018, 96019}
property time: float

the run time required by the record in seconds, either for a single job or combinations

property xic: set[tuple[float, float]]

a set of pairs of \(\xi_c\) in this record

Example:
>>> rr.xic
{(0.1, 0.1), (0.3, 0.3), (0.5, 0.5), (1.0, 1.0)}
pymule.loader.mergefks(*sets: XiRecord, **kwargs: XiRecord) → Record

performs the FKS merge

Parameters:

sets – random-seed-merged results of type XiRecord (usually from sigma())

Returns:

the FKS-merged final set containing cross sections, distributions, and run-time information.

Return type:

Record

Note

Optional argument anyxi (or anything starting with anyxi): Sometimes it is necessary to merge \(\xi_c\)-dependent runs (such as a counter term) and \(\xi_c\)-independent runs (such as the one-loop term).

Example:

Load the LO results for the muon decay using sigma()

>>> mergefks(sigma("m2enn0"))
>>> <Record m2enn0 xsec = 2211500(2), 1 histograms>

Load the NLO results

>>> mergefks(sigma("m2ennV"), sigma("m2ennR"))
<Record xsec = -4143205.8581450013, 1 histograms>

Load the NNLO results where m2ennNF does not depend on \(\xi_c\)

>>> mergefks(
...     sigma("m2ennFF"), sigma("m2ennRF"), sigma("m2ennRR"),
...     anyxi=sigma("m2ennNF")
... )
<Record xsec = -131.13(9), 1 histograms>
pymule.loader.multiintersect(lists: list[list[T] | set[T]]) → set[T]

finds elements that are common to all lists. This is used to find a list of FKS parameters of a given run.

Parameters:

list – list of lists \(l_1\), \(l_2\), …, \(l_n\)

Returns:

the list \(l_1 \cap l_2 \cap \cdots \cap l_n\)

pymule.loader.setup(folder: str | None = None, **kwargs: Any) → None

sets the default arguments for sigma().

Parameters:

  • folder – str, optional; file name, optional; folder or tarball to search for vegas files Initialised to current directory (.).

  • flavour – str, optional; the flavour to load, defaults to everything Initialised to everything, i.e. .*.

  • obs – str, optional; the observable to load (the bit after the O), defaults to everything Initialised to everything, i.e. ''.

  • folderp – str, optional; a regular expression to match directory structures of a tar file, defaults to everything Initialised to everything, i.e. .*.

  • merge – dict, optional: a dict of histograms {'name': n} to merge n bins in the histogram name. Initialised to to no merging, i.e. {}

  • sanitycheck – callable, optional; a function that, given a FileRecord, whether to include the file in the output (return True) or to skip (return False). Initialised to lambda x : True, i.e. include everything.

Example:

Setup some folders, ensure that /tmp/mcmule exists

>>> setup(folder="path/to/data.tar.bz2", cachefolder="/tmp/mcmule")
Example:

Restrict observable

>>> setup(obs="3")
Example:

Drop runs with a \(\chi^2 > 10\)

>>> setup(sanitycheck=lambda x : x['chi2a'] < 10)
pymule.loader.sigma(piece: str, **kwargs: Any) → XiRecord

loads a which_piece and returns a XiRecord

Parameters:

  • piece – str; which_piece to load

  • folder – str, optional; file name, optional; folder or tarball to search for vegas files Initialised to current directory (.).

  • flavour – str, optional; the flavour to load, defaults to everything Initialised to everything, i.e. .*.

  • obs – str, optional; the observable to load (the bit after the O), defaults to everything. Initialised to everything, i.e. ''.

  • folderp – str, optional; a regular expression to match directory structures of a tar file, defaults to everything. Initialised to everything, i.e. .*.

  • merge – dict, optional: a dict of histograms {'name': n} to merge n bins in the histogram name. Initialised to to no merging, i.e. {}

  • sanitycheck – callable, optional; a function that, given a FileRecord, whether to include the file in the output (return True) or to skip (return False). Initialised to lambda x : True, i.e. include everything.

Returns:

a dict with the tuples of FKS parameters as keys and vegas datasets as values.

Note

Use setup() to set the defaults. Arguments provided here override the defaults

Example:

Load the leading order muon decay

>>> sigma("m2enn0")
>>> <XiRecord set m2enn0, xic = [1.0], 1 histograms>

Load only observable O3

>>> sigma("m2enn0", obs="3")
>>> <XiRecord set m2enn0, xic = [1.0], 1 histograms>

Working with \(\xi_c\) data

pymule.xicut.get_errorbands(x: <MagicMock name='mock.typing.NDArray.__getitem__().__or__()' id='140397273413136'>, coeff: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, covar: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, ndata: int, cf: float = 0.9) → <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>

evaluates the errorbands of the fit obtained by get_val()

Parameters:

  • x – iterable; values of \(\xi_c\) to evaluate

  • coeff – list; coefficient list [a_0, a_1, ..., a_n]

  • covar – list; the covariance matrix

  • ndata – int; the number of data points used in the fit, required for the \(t\) value estimation

  • cl – float, optional; the confidence level used

Returns:

list; the values and errors of the fit at the presented values as [x, y, y-dy, y+dy]

pymule.xicut.get_val(xs: ~typing.Iterable[float], coeff: <MagicMock name='mock.typing.NDArray.__getitem__().__ror__()' id='140397261699856'>) → list[float]

evaluates the fit obtained by get_val()

Parameters:

  • xs – iterable; values of \(\xi_c\) to evaluate

  • coeff – list; coefficient list [a_0, a_1, ..., a_n]

Returns:

list; the values of the fit at the presented values.

pymule.xicut.mergefkswithplot(*psets: XiRecord, showfit: tuple[bool, bool] = (True, True), xlim: tuple[float, float] = (-7.0, 0.0)) → tuple[Any, Record]

performs and FKS merge like mergefks() but it also produces a \(\xi_c\) independence plot.

Note

In contrast mergefks(), here phase-space partioned results need to be summed first.

Parameters:

  • sets – a sequence of XiRecord (or sums thereof, usually from sigma()), starting with the lowest particle number and going up

  • showfit – [bool, bool], optional; whether to show the fit lines in the overview plot (first element) and the zoomed in plot (second element)

  • xlim – tuple of floats, optional; upper and lower bounds for \(\log\xi_c\)

Returns:

a figure and the FKS-merged Record containing cross sections, distributions, and run-time information.

Example:

In the partioned muon-electron scattering case

>>> fig, res = mergefkswithplot(
...     alpha*(sigma('em2emFEE')+sigma('em2emFMM')+sigma('em2emFEM')),
...     alpha*(
...         sigma('em2emREE15') + sigma('em2emREE35') +
...         sigma('em2emRMM') + sigma('em2emREM')
...     )
... )
pymule.xicut.myfit(data: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, n: int) → tuple[<MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>]

performs a log-polynomial \(\sum_{i=0}^n a_i \log(\xi_c)^i\) fit

Parameters:

  • data – numpy array; The different \(\xi_c\) values in the format np.array([[xi1, y1, e1], [xi2, y2, e2], ...]).

  • n – the degree of the polynomial

Result:

the coefficients and covariant matrix

pymule.xicut.xiresidue(pset: XiRecord, n: int, xlim: tuple[float, float] = (-7.0, 0.0)) → tuple[~typing.Any, <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>]

creates a residue plot for a \(\xi_c\) fit

Parameters:

  • sets – a XiRecord

  • n – int; order of the fit, 1 at NLO, 2 at NNLO

  • xlim – tuple of floats, optional; upper and lower bounds for \(\log\xi_c\)

Returns:

a figure and the fit coefficients as a matrix

Working with plots

pymule.plot.errorband(p, ax=None, col='default', underflow=False, overflow=False, linestyle='solid')

plots an errorband of a compatible histogram

Parameters:

  • p – Nx3 numpy matrix; the histogram to plot as np.array([[x1, y1, e1], [x2, y2, e2], ...])

  • ax – axes, optional: the axes object to use, defaults to gca() which may create a new axes.

  • col – the colour to be used for the plot. Per default matplotlib decides using the order specified in colours

  • underflow – bool, optional; whether to plot the underflow bin. Either logical or number indicating the how much bigger it shall be

  • overflow – bool, optional; whether to plot the overflow bin. Either logical or number indicating the how much bigger it shall be

  • linestyle – str, optional; which line style to use

Returns:

the artis of the main line but not the one of the errorbars

Example:

Make a simple plot

>>> errorband(dat)

Make a plot in red with dashed lines

>>> errorband(dat, 'red', 'dashed')
pymule.plot.format_label_string_with_exponent(ax, axis='both')

Format the label string with the exponent from the ScalarFormatter

pymule.plot.kplot(sigma, labelx='$x_e$', labelsigma=None, labelknlo='$\\delta K^{(1)}$', labelknnlo='$\\delta K^{(2)}$', legend={'lo': '$\\rm LO$', 'nlo': '$\\rm NLO$', 'nnlo': '$\\rm NNLO$'}, legendopts={'loc': 'upper right', 'what': 'l'}, linestyle2=':', show=[0, -1], showk=[1, 2], nomule=False)

produces a K factor plot in line with McMule’s design, i.e. a two-panel plot showing in the upper panel the cross sections and in the lower panel the K factor defined as

\[K^{(i)} = d\sigma^{(i)} / d\sigma^{(i-1)}\]
Parameters:

  • sigma – dict; the data to plot, given as a dict with keys lo, nlo, and possibly nnlo. Only pass the corrections, not the full distribution

  • labelx – str, optional; label for the x axis (supports LaTeX maths)

  • labelsigma – str, optional; label for the upper y axis (supports LaTeX maths)

  • labelknlo – str, optional; the labels for the NLO K factor

  • labelknnlo – str, optional; the labels for the NNLO K factor

  • show – list, optional; a list which cross sections to show, 0 indicates the LO cross section, 1 the NLO etc. -1 indicates the last given cross section

  • showk – list, optional; a list which K factors to show, 0 indicates the LO cross section, 1 the NLO etc. -1 indicates the last given cross section

  • legend – dict, optional; a dict with the legend for lo, nlo, nnlo. The keys nlo2 and nnlo2 are optional and will be drawn dashed in the lower panel.

  • legendopts –

    dict, optional; a kwargs dict of options to be passed to legend(..) as well as the what key indicating whether the legend such be placed in the lower panel (l, default), upper panel (u), or as a figlegend (fig). Notable is the loc-key that places the legend inside the object specified by what. Possible values are (cf. legend)

    • upper right

    • upper left

    • lower left

    • lower right

    • right

    • center left

    • center right

    • lower center

    • upper center

    • center

  • nomule – bool, optional; if set to True, no mule will be printed

Returns:

the figure as well as all axis created

Example:

An NNLO K factor plot

>>> fig, (ax1, ax2, ax3) = kplot(
...   {
...     'lo':  lodata['thetae'],
...     'nlo': nlodata['thetae'],
...     'nnlo':nnlodata['thetae'],
...   },
...   labelx="$\theta_e\,/\,{\rm mrad}$",
...   labelsigma="$\D\sigma/\D\theta_e\ /\ {\rm\upmu b}$",
...   legend={
...     'lo': '$\sigma^{(0)}$',
...     'nlo': '$\sigma^{(1)}$',
...     'nnlo': '$\sigma^{(2)}$'
...   },
...   legendopts={'what': 'u', 'loc': 'lower right'}
... )
pymule.plot.setup_pgf()

setupf_pgf() ensures that Matplotlib exports PGF compatible plots.

pymule.plot.threepanel(labelx='', upleft=[], labupleft='', colupleft=['C0', 'C1', 'C2', 'C3', 'C4', 'C5', 'C6', 'C7', 'C8', 'C9'], middleleft=[], labmiddleleft='', colmiddleleft=['C0', 'C1', 'C2', 'C3', 'C4', 'C5', 'C6', 'C7', 'C8', 'C9'], downleft=[], labdownleft='', coldownleft=['C0', 'C1', 'C2', 'C3', 'C4', 'C5', 'C6', 'C7', 'C8', 'C9'])

creates three panel plot, accommodating at most three axes (upper, middle, lower). The x axis is naturally shared.

Parameters:

  • labelx – str, optional; label for the x axis

  • upleft – Nx3 numpy matrix or list thereof, optional; data plotted in the upper-left axes

  • colupleft – colour for upper-left data, defaults to colour scheme defined in colours

  • labupleft – str, optional; the label for the upper-left data

  • midleft – Nx3 numpy matrix or list thereof, optional; data plotted in the middle-left axes

  • colmidleft – colour for middle-left data, defaults to colour scheme defined in colours

  • labmidleft – str, optional; the label for the middle-left data

  • downleft – Nx3 numpy matrix or list thereof, optional; data plotted in the lower-left axes

  • coldownleft – colour for lower-left data, defaults to colour scheme defined in colours

  • labdownleft – str, optional; the label for the lower-left data

Returns:

the figure and a list of all axes created

pymule.plot.twopanel(labelx='', upleft=[], labupleft='', colupleft=['C0', 'C1', 'C2', 'C3', 'C4', 'C5', 'C6', 'C7', 'C8', 'C9'], downleft=[], labdownleft='', coldownleft=['C0', 'C1', 'C2', 'C3', 'C4', 'C5', 'C6', 'C7', 'C8', 'C9'], upright=[], labupright='', colupright=['C0', 'C1', 'C2', 'C3', 'C4', 'C5', 'C6', 'C7', 'C8', 'C9'], downright=[], labdownright='', coldownright=['C0', 'C1', 'C2', 'C3', 'C4', 'C5', 'C6', 'C7', 'C8', 'C9'], upalign=[], downalign=[])

creates two panel plot, accommodating at most four axes (upper left, upper right, lower left, and lower right). The x axis is naturally shared.

Parameters:

  • labelx – str, optional; label for the x axis

  • upleft – Nx3 numpy matrix or list thereof, optional; data plotted in the upper-left axes

  • colupleft – colour for upper-left data, defaults to colour scheme defined in colours

  • labupleft – str, optional; the label for the upper-left data

  • upright – Nx3 numpy matrix or list thereof, optional; data plotted in the upper-right axes

  • colupright – colour for upper-right data, defaults to colour scheme defined in colours

  • labupright – str, optional; the label for the upper-right data

  • downleft – Nx3 numpy matrix or list thereof, optional; data plotted in the lower-left axes

  • coldownleft – colour for lower-left data, defaults to colour scheme defined in colours

  • labdownleft – str, optional; the label for the lower-left data

  • downright – Nx3 numpy matrix or list thereof, optional; data plotted in the lower-right axes

  • coldownright – colour for lower-right data, defaults to colour scheme defined in colours

  • labdownright – str, optional; the label for the lower-right data

  • upalign – list of two values, optional; align the first and second values of the left and right y axes in the upper panel

  • downalign – list of two values, optional; align the first and second values of the left and right y axes in the lower panel

Returns:

the figure and a list of all axes created

Example:

make a comparison plot between dat and dat_ref as a \(\mathrm{d}\sigma/\mathrm{d}\theta_e\)

>>> fig,(ax1,ax2)=twopanel(
...     r'$\theta_e\,/\,{\rm mrad}$',
...     upleft=[dat, dat_ref],
...     downleft=divideplots(dat, dat_ref),
...     labupleft=r"$\D\sigma/\D\theta_e\,/\,\upmu{\rm b}$",
...     labdownleft=r'$\rm rel. difference$'
... )
pymule.plot.watermark(fig, txt='PRELIMINARY', fontsize=60, rotation=20)

watermarks a figure

Parameters:

  • fig – the figure to watermark

  • txt – str, optional; the watermark text to use

  • fontsize – int, optional; the fontsize of the watermark

  • rotation – int, optional; the angle of the watermark in deg

Example:

Watermark a figure as preliminary

>>> fig = figure()
>>> ...
>>> watermark(fig)

Watermark a figure as incomplete

>>> fig = figure()
>>> ...
>>> watermark(fig, "INCOMPLETE")
pymule.colours.alpha_composite(bg: str, fg: str, alpha: float) → str

calculates the result of alpha-compositing two colours

Parameters:

  • bf – colour specifier for the background

  • fg – colour specifier for the foreground

  • alpha – float; alpha value

Result:

the resulting colour

pymule.mule.mulify(fig, delx=0, dely=0, col='lightgray', realpha=True)

adds the McMule logo to a figure

Parameters:

  • fig – figure to add the logo

  • delx – float, optional; shift the logo in x direction

  • dely – float, optional; shift the logo in x direction

  • col – colour specifier, optional; colour to use for the logo

  • realpha – bool, optional; whether to re-run the alpha channel

pymule.mpl_axes_aligner.yaxes(ax1, ax2, y1=1, y2=None)

yaxes(ax1, ax2, y=1) changes the limits of ax1 and ax2 to align the values of y on both axis.

yaxes(ax1, y1, ax2, y2) changes the limits of the axis ax1 and ax2 such that the value for y1 on ax1 is aligned to the value of y2 on ax2.

Useful other functions

pymule.compress.uncompress(b)

uncompress(“string”) recovers an expression from a compressed string representation generated by Mathematica’s Compress. Only lists, numbers, and strings are supported. Lists can be nested.

pymule.maths.Li2(x)

Li2(x) returns PolyLog[2, x] for x as a number, a list, or an np.ndarray.

pymule.maths.Li3(x)

Li3(x) returns PolyLog[3, x] for x as a number, a list, or an np.ndarray.

Validator

The validator (cf. Section Observable validator) allows the checking of user files. It is implemented as part of pymule.

Phase space generation

exception pymule.validator.psgen.NotEnoughEnergyError(min: float, ms: list[float])

An error that indicates that not enough energy is available to generate this process. This may depend on the kinematic point.

pymule.validator.psgen.psgen(ra: <MagicMock name='mock.typing.NDArray.__getitem__().__ror__()' id='140397261699856'>, Ms: <MagicMock name='mock.typing.NDArray.__getitem__().__ror__()' id='140397261699856'>) → <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>

General purpose momentum generator

Parameters:

  • ra – array of random numbers to generate the momentum from. see note in psgen_X_FKS() for details

  • Ms – array of masses including 1 initial state masses and (>=2) final state masses

pymule.validator.psgen.psgen_FKS(ra: <MagicMock name='mock.typing.NDArray.__getitem__().__ror__()' id='140397261699856'>, Ms: <MagicMock name='mock.typing.NDArray.__getitem__().__ror__()' id='140397261699856'>, NumPhotons: int) → <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>

General purpose FKS-photon momentum generator

Parameters:

  • ra – array of random numbers to generate the momentum from. see note in psgen_X_FKS() for details

  • Ms – array of masses including 1 initial state mass and (>=2) final state masses, NOT including FKS photons

pymule.validator.psgen.psgen_X(ra: <MagicMock name='mock.typing.NDArray.__getitem__().__ror__()' id='140397261699856'>, Ms: <MagicMock name='mock.typing.NDArray.__getitem__().__ror__()' id='140397261699856'>, scms: float) → <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>

General purpose cross-section momentum generator

Parameters:

  • ra – array of random numbers to generate the momentum from. see note in psgen_X_FKS for details

  • Ms – array of masses including 2 initial state masses and (>=2) final state masses

  • scms – center of momentum energy for the system, typically s in the literature.

pymule.validator.psgen.psgen_X_FKS(ra: <MagicMock name='mock.typing.NDArray.__getitem__().__ror__()' id='140397261699856'>, Ms: <MagicMock name='mock.typing.NDArray.__getitem__().__ror__()' id='140397261699856'>, NumPhotons: int, scms: float) → <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>

General purpose cross section and FKS-photon momentum generator

Parameters:

  • ra – array of random numbers to generate the momentum from.

  • Ms – array of masses including 2 initial state masses and (\(\ge 2\)) final state masses, NOT including FKS-photons

  • scms – centre of mass energy for the system

Note

Letting \(E_{\gamma}\) be a photon energy, \(n_{\gamma}\) be number of photons, \(E_{\rm inv}\) be decay energy, and \(n_{\rm masses}\) be number of masses, the random numbers are structured as:

\[[E_{\gamma}] \otimes n_{\gamma}, [\theta_{\gamma} \phi_{\gamma}] \otimes n_{\gamma}, [E_{\rm inv} \theta_{\rm masses} \phi_{\rm masses}] * (n_{\rm masses} - 4), [\theta_{\rm masses} \phi_{\rm masses}].\]

The - 4 comes from there being 2 masses in and 2 masses out in the X case, for non-cross section cases it is - 3 instead.

SO loader

class pymule.validator.importer.DistanceFunction

A distance function base class. This can be called which passes the momenta of two events and returns the distance.

classmethod build(func: ~typing.Callable[[<MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>], float]) → DistanceFunction

Used as a decorator to change a normal function into a “distance function”

class pymule.validator.importer.FortranDistanceFunction(fort_load: FortranLoader)
Parameters:

fort_load – An instance of FortranLoader which is shared between Measurement and Distance

Attribute distance_names:

List of names with which the distance function may be found

class pymule.validator.importer.FortranLoader(userfile: str)

A class which acts as a wrapper to the measurement function which is defined in the userfile.

Parameters:

userfile – name of the file to be loaded, including extension default: user.so

Attribute inituser_names:

List of names with which the inituser() function may be found. This function may be unavailable

Attribute filenamesuffixLen_names:

List of names with which the filenamesuffixLen variable may be found

Attribute nrq_names:

List of names with which the nr_q variable may be found

Attribute nrbins_names:

List of names with which the nr_bins variable may be found

class pymule.validator.importer.FortranMeasurementFunction(fort_load: FortranLoader)

A measurement function implemented as a shared library for use in McMule. Calling this will evaluate quant: and return the histogram values, :f:var:`pass_cut(), and names

Parameters:

fort_load – An instance of FortranLoader which is shared between Measurement and Distance

Attribute quant_names:

List of names with which the measurement function may be found

process_names() → None

takes the binary names from the userfile and converts them into a string, assuming UTF-8 encoding.

Returns:

a string of human-readable names.

class pymule.validator.importer.MeasurementFunction

A measurement function base class. This can be called which passes the momenta to the function and returns the histogram values, pass_cut, and names

pymule.validator.importer.names_t

alias of Any

Validation logic

class pymule.validator.validator.BaseValidator(npoints: int = 0, order_max: int = 100, t: int = 2, smoothness_rtol: float = 0.001, smoothness_ratol: float = 1.0, limit_rtol: float = 0.001, limit_atol: float = 0.001, allowed_to_fail: bool = False, fail_msg: str = '')

Base class for the testing of a MeasurementFunction (or function with the same signature)

Attribute npoints:

the number of Monte Carlo points to use. If zero, no MC tests are carried out

Attribute order_max:

the maximum order to use in Richardson extrapolation

Attribute t:

amount to reduce \(\xi_0\) by each time (\(\xi_0 / t^n\)) in the Richardson extrapolation

Attribute smoothness_rtol:

used to determine smoothness, if two consecutive results differ by r[1] / r[0] > smoothness_rtol then the function is not smooth

Attribute smoothness_atol:

similar to smoothness_rtol but provides an absolute tolerance relative to the largest result found, i.e. (r[1] - r[0]) > ratol * r[max] is not smooth

Attribute limit_rtol:

allowed relative tolerance between the limit at 0 and the final result or allowed relative tolerance between the modified and unmodified results

Attribute limit_atol:

allowed absolute tolerance between limit at 0 and final result or allowed absolute tolerance between the modified and unmodified results

Attribute allow_to_fail:

failing these test is not fatal (eg. RotationValidator)

graph(xi: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, results: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, limit_res: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, names: ~typing.Any, rtol: float | None = 0.001, atol: float | None = 0.001) → None
Parameters:

  • xi – x-axis, usually energy proportion of the soft photon

  • results – y-axis, list of results corresponding to each point

  • limit_res – horizontal line, result when \(\xi \to 0\)

  • names – legend, names which are used to identify each line produced by results

modmomenta(momenta: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → None
Parameters:

momenta – a nx4 numpy matrix of 4-momenta which are transformed in place

modpoint(points: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, **kwargs: ~typing.Any) → None
Parameters:

points – numpy array of random numbers to be modified in place

class pymule.validator.validator.BoundsValidator(npoints: int, meas_func: MeasurementFunction, dist_func: DistanceFunction | None = None, **kwargs: Any)

Tests the histogram bounds by Monte Carlo sampling the phase space

class pymule.validator.validator.EventSwappingValidator(dist_func: DistanceFunction, meas_func: MeasurementFunction | None = None, **kwargs: Any)

Validates the symmetry relation \(d(e_1,e_2) = d(e_2,e_1)\) for a distance function

class pymule.validator.validator.FlipValidator(axis: int, meas_func: MeasurementFunction, dist_func: DistanceFunction | None = None, **kwargs: Any)

Test spatial-symmetry by flipping across an axis. This test may fail

Parameters:

axis – int, axis to be flipped across 0 for \(x\), 1 for \(y\), 2 for \(z\)

modmomenta(momenta: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → None
Parameters:

momenta – a nx4 numpy matrix of 4-momenta which are transformed in place

class pymule.validator.validator.LorentzBoostValidator(frame_momentum: <MagicMock name='mock.typing.NDArray.__getitem__().__ror__()' id='140397261699856'>, meas_func: ~pymule.validator.importer.MeasurementFunction, dist_func: ~pymule.validator.importer.DistanceFunction | None = None, **kwargs: ~typing.Any)

Test the boost-invariance, usually only used for cross-sections

This test may fail

Parameters:

  • int (frame_momentum,) – boost into the restframe of this particle using boost_RF()

  • vector (frame_momentum, numpy) – the 4-momentum to boost the frame into using boost_back()

modmomenta(momenta: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → None
Parameters:

momenta – a nx4 numpy matrix of 4-momenta which are transformed in place

class pymule.validator.validator.MeasurementValidator

A mixin class for Validators that can test both measurement function and distance functions. It implements the test framework for discrete and continuous tests by creating a base point and modifying it.

Parameters:

  • meas_func – the measurement function callable

  • dist_func – the distance function callable

  • limit_rtol – relative tolerance

  • limit_atol – absolute tolerance

  • softPhotonIndex – the index or list of indices that correspond to soft photons

  • order_max – the maximum order to use in Richardson extrapolation

  • t – amount to reduce \(\xi_0\) by each time (\(\xi_0 / t^n\)) in the Richardson extrapolation :param t:

compare(base_result: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, base_pass_cut: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, mod_result: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, mod_pass_cut: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, rtol: float, atol: float | list[float]) → tuple[list[bool], list[bool]]

Compares a base point with a modified point.

Parameters:
Returns:

a tuple of 2 boolean lists containing if each number was close or not

class pymule.validator.validator.Piece(name: str, p_in: list[~pymule.utils.pid.Particle], p_out: list[~pymule.utils.pid.Particle], masses: <MagicMock name='mock.typing.NDArray.__getitem__().__ror__()' id='140397261699856'>, scms: float | None, randlen: int, psgenfunc: ~typing.Callable[[...], <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>], numphotons: int | None, soft_photons: int)

A which_piece that can be used for validation by BaseValidator. Normally this is constructed with Piece.from_whichpiece() or Piece.from_process()

Attribute name:

the name of the which_piece

Attribute p_in:

a list of incoming Particle

Attribute p_in:

a list of outgoing Particle

Attribute masses:

a list of masses

Attribute scms:

the centre-of-mass energy \(s\)

Attribute randlen:

the number of random numbers required

Attribute psgenfunc:

the phase space routine used by this

Attribute numphotons:

the number of photons in the process, not counting the first two particles

Attribute soft_photons:

the number of photons that actually can become soft

Can be used to call the phase space generator

static convert_masses(masses: list[str]) → list[float]
Parameters:

masses – list of strings which correspond to physical masses, eg mm for the muon mass \(m_\mu = 105\,\textrm{MeV}\)

Note

gets masses from pymule/constants.py

classmethod from_process(PIDTable_: PIDTable, process: dict[str, dict[str, list[str]]], order: str | None = None, scms: float | None = None) → Generator[Piece, None, None]

Yields all Piece required for a specific generic process up to a given order or all orders

Parameters:

  • PIDTable – the PIDTable

  • process – a dictionary encoding the process, cf. load_processes_from_json()

  • order – str, optional;

Returns:

a generator of Piece

classmethod from_whichpiece(PIDTable_: PIDTable, whichpiece: str, scms: float | None = None) → Piece

Generates a Piece from a which_piece and the PIDTable

Parameters:

  • PIDTable – the PIDTable

  • whichpiece – the which_piece, eg. em2em0

Returns:

a Piece

Warning

The number of photons numphotons is only the number of photons that can be generated using FKS. For processes such as ee2gg, these will not be listed

static get_soft_from_sigma(fxn: str) → int
Parameters:

fxn – sigma_n for \(n = 0..9\), trailing characters after n are ignored

Returns:

number of soft photons

Warning

currently doesn’t work for sigma_YFS (i.e. arbitrary photons)

static pick_psgen_func(p_in: List[Particle], p_out: List[Particle]) → Tuple[Callable[[...], Any], int | None]

return the appropriate psgen function based on particles in and out

Parameters:

  • p_in – particles into the system

  • p_out – particles out of the system

Returns:

a tuple of the phase space generator and the number of FKS photons (None if no FKS exist)

Note

there can be photons present which cannot be generated by FKS (e.g. ee2gg or eg2eg) this is handled by returning number of FKS photons not number of photons total

class pymule.validator.validator.RotationValidator(alpha: float, beta: float, gamma: float, meas_func: MeasurementFunction, dist_func: DistanceFunction | None = None, **kwargs: Any)

Test the rotational-invariance of the measurement function and distance function, e.g. rotating all the momenta by \(\pi/3\) around the \(x\)-axis shouldn’t change either result

This test may fail

Parameters:

  • alpha – Euler angle for rotation

  • beta – Euler angle for rotation

  • gamma – Euler angle for rotation

modmomenta(momenta: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → None
Parameters:

momenta – a nx4 numpy matrix of 4-momenta which are transformed in place

class pymule.validator.validator.SameZeroValidator(dist_func: DistanceFunction, meas_func: MeasurementFunction | None = None, **kwargs: Any)

Validates the relation \(d(e_1,e_1) = 0\) for a distance function

class pymule.validator.validator.SoftLimitValidator(softPhotonIndex: int | list[int] | tuple[int], meas_func: MeasurementFunction, dist_func: DistanceFunction | None = None, **kwargs: Any)

Test whether \(\lim_{\gamma \to 0} S(q, g) = S(q)\)

Parameters:

softPhotonIndex – either the index of the random number associated with the photon energy or a list of indices which should go to 0 with the same ratio

modpoint(points: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, t: int = 2, zero: bool = False) → None
Parameters:

points – numpy array of random numbers to be modified in place

class pymule.validator.validator.SoftSwapValidator(softPhotonIndex: int | list[int] | tuple[int], swapAIndex: int, swapBIndex: int, meas_func: MeasurementFunction, dist_func: DistanceFunction | None = None, **kwargs: Any)

Combination of the SoftLimitValidator and the SoftSwapValidator symmetry test, i.e.

\[\lim_{\gamma_1 \to 0} S(q, \gamma_1, \gamma_2) = S(q, \gamma_2) = \lim_{\gamma_2 \to 0} S(q, \gamma_2, \gamma_1) = S(q, \gamma_1)\]
Parameters:

  • softPhotonIndex – index of the random number associated with the photon energy

  • swapAIndex – index of the first particle to be swapped

  • swapBIndex – index of the second particle to be swapped

modmomenta(momenta: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → None
Parameters:

momenta – a nx4 numpy matrix of 4-momenta which are transformed in place

modpoint(points: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>, t: int = 2, zero: bool = False) → None
Parameters:

points – numpy array of random numbers to be modified in place

class pymule.validator.validator.SwapSymmetryValidator(swapAIndex: int, swapBIndex: int, meas_func: MeasurementFunction, dist_func: DistanceFunction | None = None, **kwargs: Any)

Test whether \(S(p_1, p_2, q) = S(p_2, p_1, q)\) for 2 identical particles \(p_1\), \(p_2\)

Parameters:

  • swapAIndex – index of the first particle to be swapped

  • swapBIndex – index of the second particle to be swapped

modmomenta(momenta: <MagicMock name='mock.typing.NDArray.__getitem__()' id='140397293428560'>) → None
Parameters:

momenta – a nx4 numpy matrix of 4-momenta which are transformed in place

class pymule.validator.validator.TriangleValidator(dist_func: DistanceFunction, meas_func: MeasurementFunction | None = None, **kwargs: Any)

Validates the triangle inequality \(d(e_1,e_2) \ge d(e_1,e_3) + d(e_3,e_2)\) for a distance function

class pymule.validator.validator.ValidatorBuilder(meas_func: MeasurementFunction | None, dist_func: DistanceFunction | None = None, npoints: int = 0, order_max: int = 100, t: int = 2, smoothness_rtol: float = 0.001, smoothness_ratol: float = 1.0, limit_rtol: float = 0.001, limit_atol: float = 0.001, allowed_to_fail: bool = False, fail_msg: str = '')

Builder class to build tests from BaseValidator with a common set of initialisations

Attribute npoints:

the number of Monte Carlo points to use. If zero, no MC tests are carried out

Attribute order_max:

the maximum order to use in Richardson extrapolation

Attribute t:

amount to reduce \(\xi_0\) by each time (\(\xi_0 / t^n\)) in the Richardson extrapolation

Attribute smoothness_rtol:

used to determine smoothness, if two consecutive results differ by r[1] / r[0] > smoothness_rtol then the function is not smooth

Attribute smoothness_atol:

similar to smoothness_rtol but provides an absolute tolerance relative to the largest result found, i.e. (r[1] - r[0]) > ratol * r[max] is not smooth

Attribute limit_rtol:

allowed relative tolerance between the limit at 0 and the final result or allowed relative tolerance between the modified and unmodified results

Attribute limit_atol:

allowed absolute tolerance between limit at 0 and final result or allowed absolute tolerance between the modified and unmodified results

Attribute allow_to_fail:

failing these test is not fatal (eg. RotationValidator)

Attribute fail_msg:

the message being printed if the test fails

Example:
>>> # Looking at the m2ennR process where 1 photon is generated
>>> piece = Piece.from_whichpiece(table, "m2ennR")
>>> builder = ValidatorBuilder(
    meas_func=FortranMeasurementFunction(),
)
>>> # check if swapping the neutrinos causes a change
>>> swap_sym = builder.build_swap_symmetry(3, 4)
>>> swap_sym.run_test(piece)
build_all(piece: Piece) → Generator[BaseValidator, None, None]

Generates all tests that can be performed for this Piece with the given options

Returns:

an iterator of BaseValidator

build_boost(frame_momentum: <MagicMock name='mock.typing.NDArray.__getitem__().__ror__()' id='140397261699856'>) → LorentzBoostValidator

build an instance of LorentzBoostValidator with parameters from the builder

Parameters:

  • int (frame_momentum,) – boost into the restframe of this particle using boost_RF()

  • vector (frame_momentum, numpy) – the 4-momentum to boost the frame into using boost_back()

Returns:

instance of LorentzBoostValidator

build_bounds() → BoundsValidator

build an instance of BoundsValidator with parameters from the builder

Parameters:

npoints – int, number of Monte Carlo points to use

Returns:

instance of BoundsValidator

build_eventswap() → EventSwappingValidator

build an instance of EventSwappingValidator with parameters from the builder

Returns:

instance of EventSwappingValidator

build_flip(axis: int) → FlipValidator

build an instance of FlipValidator with parameters from the builder

Parameters:

axis – int, axis to be flipped across 0 for \(x\), 1 for \(y\), 2 for \(z\)

Returns:

instance of FlipValidator

build_rotation(alpha: float, beta: float, gamma: float) → RotationValidator

build an instance of RotationValidator with parameters from the builder

Parameters:

  • alpha – Euler angle for rotation

  • beta – Euler angle for rotation

  • gamma – Euler angle for rotation

Returns:

instance of RotationValidator

build_samezero() → SameZeroValidator

build an instance of SameZeroValidator with parameters from the builder

Returns:

instance of SameZeroValidator

build_soft_limit(softPhotonIndex: int | list[int] | tuple[int]) → SoftLimitValidator

build an instance of SoftLimitValidator with parameters from the builder

Parameters:

softPhotonIndex – index of the random number associated with the photon energy

Returns:

instance of SoftLimitValidator

build_soft_swap(softPhotonIndex: int | tuple[int] | list[int], swapAIndex: int, swapBIndex: int) → SoftSwapValidator

build an instance of SoftSwapValidator with parameters from the builder

Parameters:

  • softPhotonIndex – index of the random number associated with the photon energy

  • swapAIndex – index of the first particle to be swapped

  • swapBIndex – index of the second particle to be swapped

Returns:

instance of SoftSwapValidator

build_swap_symmetry(swapAIndex: int, swapBIndex: int) → SwapSymmetryValidator

build an instance of SwapSymmetryValidator with parameters from the builder

Parameters:

  • swapAIndex – index of the first particle to be swapped

  • swapBIndex – index of the second particle to be swapped

Returns:

instance of SwapSymmetryValidator

build_triangle() → TriangleValidator

build an instance of TriangleValidator with parameters from the builder

Returns:

instance of TriangleValidator

pymule.validator.validator.convert_to_PID(p_in: list[Particle], p_out: list[Particle]) → str

convert lists containing the incoming and outgoing :class:’Particle’s to the generic group they belong to

Example:
>>> p_in = [Muon(-1)]
>>> p_out = [Electron(-1), Neutrino(""), Neutrino("")]
>>> convert_to_PID(p_in, p_out)
"m2enn"
pymule.validator.validator.do_test_from_whichpiece(PIDTable_: PIDTable, whichpiece: str, measurement_func: MeasurementFunction, distance_function: DistanceFunction | None = None, scms: float | None = None, **kwargs: Any) → dict[BaseValidator, tuple[int, int, int]]

Performs all tests for one which_piece. See ValidatorBuilder for more options

Parameters:

  • PIDTable – instance of PIDTable in which whichpiece can be found

  • measurement_func – instance of MeasurementFunction or similar to be tested

  • whichpiece – the string which_piece to generate tests for

Example:
>>> PIDTable_ = PIDTable()
>>> M = MeasurementFunction()
>>> do_test_from_whichpiece(PIDTable_, M, "m2enng0")
...
pymule.validator.validator.do_tests_from_process(PIDTable_: PIDTable, generic_process: str, measurement_function: MeasurementFunction, distance_function: DistanceFunction | None = None, order: str | None = None, scms: float | None = None, **kwargs: Any) → dict[BaseValidator, tuple[int, int, int]]

Perform tests for all whichpieces under a generic process at a certain order (or all orders). See ValidatorBuilder for more options

Parameters:

  • PIDTable – instance of PIDTable in which the generic_process can be found

  • MeasurementFunction – instance of :class:’MeasurementFunction: to be tested

  • generic_process – the generic process to be tested with

  • order – str: LO, NLO, NNLO (sometimes split into el. NNLO and full NNLO) None: Test all orders available

Example:
>>> do_tests_from_process(..., generic_process="m2enn", order="NLO",...)
INFO:Testing whichpiece: 'm2ennV'...
INFO:Testing whichpiece: 'm2ennC'...
INFO:Testing whichpiece: 'm2ennR'...
Example:
>>> do_tests_from_process(..., generic_process="m2enn", order=None,...)
INFO:Testing at order LO:
INFO:Testing whichpiece: 'm2enn0'...
INFO:Testing at order NLO:
INFO:Testing whichpiece: 'm2ennF'...
INFO:Testing whichpiece: 'm2ennR'...
INFO:Testing at order NNLO:
...
pymule.validator.validator.load_processes_from_json(filename: str, path: Path) → dict[str, dict[str, dict[str, list[str]]]]
Parameters:

  • filename – .json file to be loaded

  • path – path to the file

Returns:

dict-object containing the .json data