Skip to content

madspace MLM: follow madevent's setclscales/rewgt and add a clustering-measure switch - #156

Open
oliviermattelaer wants to merge 32 commits into
feat-mlmfrom
claude/madspace-madevent-scales
Open

oliviermattelaer wants to merge 32 commits into
feat-mlmfrom
claude/madspace-madevent-scales

Conversation

@oliviermattelaer

@oliviermattelaer oliviermattelaer commented Sep 15, 2026 •

Copy link
Copy Markdown
Contributor

Brings madspace's MLM merging (dynamical_scale_choice = "mlm") in line with madevent's ickkw=1 scale setting and reweighting, validated against madevent on W+jets, tt~+jets and HEFT h(->tau tau)+jets.

Branch state. Both this branch and feat-mlm were merged with the latest main (aa9d133), so the diff is only this work: 17 files, +2647/-108, almost all in madspace/ (kernel src/kernels/mlm.hpp, src/phasespace/{mlm_clustering,integrand,scale}.cpp, src/driver/lhe_output.cpp, tests) plus the mg7 run card, banner.py and the launcher. The invariant-mass cut, auto_ptj_mjj and the combine fix listed below already reached main through their own PRs.

What changes

  • Scales: beam.scale_scheme gives madevent's mu_R/mu_F definition from setclscales. beam.parton_line_scheme adds its goodjet bookkeeping. Per-beam factorisation scales, mu_R and mu_F go into the LHE <scales> tag. beam.min_scale sets a floor, and scales are bounded to the PDF grid's Q range.
  • Reweighting: alpha_s is reweighted per clustering vertex (beam.alphas_reweighting), at the scale setclscales leaves behind. Beam densities are reweighted along the clustering ladder (beam.pdf_reweighting, madevent's hidden pdfwgt).
  • Merging jets and the merging cut: legs count as merging jets the way madevent's iqjets decides, and the xqcut cut applies only at a vertex that made a jet.
  • Fixes to follow setclscales / cluster.f: mt2last at an s-channel QCD root; each end of a t-channel chain gets its own flavour and mass; the root is measured in the frame cluster.f leaves it in.
  • Generation-level cuts: phasespace.auto_ptj_mjj (madevent's automatic ptj/mjj adjustment when xqcut > 0) and a two-particle invariant-mass cut that also bounds the integration.
  • New switch beam.clustering_measure: fxfx (default, unchanged) or madevent (LO DJ of kin_functions.f, whose massless-massive pair ties an initial-state clustering exactly, settled as cluster.f does). The FxFx score of q* > q W / g* > g h pairs is often below the jet pt, which left W + 1 jet and HEFT h + 1 jet low against madevent.

Validation (madspace vs madevent, same cards)

  • madspace test suite after the merge with main: 1648 passed (excluding test_pdf.py/test_systematics.py, which need a working lhapdf, and the torch-only test_flow.py/test_mlp.py); 103 of them in test_mlm.py.
  • tt~+jets, per multiplicity: 0-2j agree within 1%, 3j is still 4-8% high.
  • HEFT p p > h +0-3j HIG<=1 HIW=0 QED=1, h > ta+ ta-, madevent measure, xqcut 15: merged cross section 1.072 at qcut 15 and 1.007 at qcut 75; 2j matching efficiency identical in both codes.
  • Pure VBF p p > h j j QCD=0 + p p > h j j j QCD=1: parton level 1.011; showered 2j 0.99/0.98 and 3j 0.92/0.97 at qcut 20/60.

Known limitation (follow-up)

madevent chooses the clustering history partly from the integration channel: it reclusters along the channel graph when the jet count differs from the channel's first event. madspace takes the minimal-measure history over all diagrams. As a result:

  • HEFT with QED unrestricted clusters VBF events as gluon fusion (showered 2j at large qcut about 0.5 of madevent);
  • HEFT gg > h gg with a g > g h attachment is 5-15% high at 2j/3j;
  • in pure VBF, some 3j events get different merging-jet tags.

Selecting the history by per-diagram |A_i|^2 is being worked on in claude/mlm-diagram-history, which will target this branch.

🤖 Generated with Claude Code

oliviermattelaer and others added 30 commits September 7, 2026 22:53
SCALUP holds one number and the two codes disagree about which: unwgt.f
writes sqrt(max(q2fact)), i.e. mu_F, while madspace writes ren_scale.
Comparing that field between them is apples to oranges, which is how the
earlier claim that the mu_R definition explains the high-multiplicity
cross-section difference came to rest on nothing.

Emit both explicitly as the standard LHEF3 attributes - mur, muf, and
muf1/muf2 since madspace can carry a different factorisation scale per
beam - so a comparison does not have to guess. Nothing that was written
before changes.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
NOT VALIDATED - see below. Committed so the machinery is not lost, but
this should not be merged as it stands.

madevent reads its scales off the clustering history by following each
beam's parton line and taking

    mu_R   = (s[jlast1] s[jcentral1] s[jlast2] s[jcentral2])^(1/4)
    mu_F,b = sqrt(s[jlast_b] s[jcentral_b])

where jlast is the last initial-state clustering while that beam's line
is still a jet and jcentral the last one while it is still coloured.
madspace instead takes the geometric mean of every clustering scale and
a single factorisation scale. beam.scale_scheme selects between them,
defaulting to the existing behaviour because unlike jet_scale_scheme
this one moves the cross section.

The walk needs two things the clustering already had after the tracing
work: beam j's line always sits in slot j, because a clustering keeps
the lower of the two indices, so the slot index is madevent's ibeam(j)
for free; and the third state-machine word had spare bits for the
mother's flavour, which is what says where a line stops being a jet and
where it stops being coloured.

Why it is not validated: the check was going to be the mean scale per
jet multiplicity, measured from the LHE. That comparison is not
meaningful as it stands, because the two write different quantities into
SCALUP - unwgt.f writes sqrt(max(q2fact)) there, i.e. mu_F, while
madspace writes ren_scale. Comparing them was apples to oranges, so
neither the earlier claim that the mu_R definition explains the
high-multiplicity cross-section difference, nor this implementation,
rests on evidence yet.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Below the lowest Q of a PDF grid the densities are undefined, and a scale
of a few MeV comes back as a NaN that no later cut can remove: the veto
zeroes the weight, but the pdf has already been evaluated and NaN times
zero is still NaN. So the floor both vetoes the event, through a new
scale_weight output, and clamps the scales so that what is computed on
the way stays finite.

madevent applies the same floor to mu_F at 2 GeV, in reweight.f, as part
of setclscales. Putting it in EnergyScale instead of in the MLM kernel
keeps it out of the merging path and makes it apply to every dynamical
scale choice, which is where it belongs.

Verified directly: with min_scale = 2 on a Drell-Yan + 2 jet topology
whose madevent-scheme mu_F reaches 0.15 GeV, fact_scale1 comes back
floored at exactly 2 and scale_weight is 0 for the events below it.

Also extends the madevent scale scheme to take the final 2 -> 1
clustering step. madevent's loop runs one step further than this kernel,
and that step - which absorbs what is left of the final state into a
beam - is where jcentral usually lands, which is why its factorisation
scale comes out as the transverse mass of the lepton pair rather than a
jet scale. Rather than restructuring the state machine, the step is
taken in the kernel: after the walk exactly one final-state slot is
still alive and the clustering is fixed by which beam it goes to.

The madevent scheme is NOT working: it makes some channel of
p p > e+ e- + 0,1,2,3 jets return a NaN weight, where clustering_mean on
the same directory with the same floor integrates fine. Not diagnosed
yet, so the scheme should not be used and this should not be merged as
it stands.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
beam.min_scale kept mu_R and mu_F off the bottom of the grid, where a
density of a few MeV is undefined. The top end has the same hole, and the
madevent MLM scheme reaches it: it multiplies four clustering scales under
a fourth root and produced mu_F up to 12.3 TeV against NNPDF23's Qmax of
10 TeV. The density came back NaN, and since the PDF is evaluated before
the veto weight multiplies it, no later cut could remove it - one such
point aborted the whole channel.

So the floor becomes a range. The upper bound is not a new run-card knob:
it is read off the grid itself in launch.py, which is where the right
number actually lives. Events outside are both vetoed and clamped, the
clamp only so what is computed on the way stays finite.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
They were only recorded under partial_weights, so for an ordinary run the
event buffer had no fact_scale1/2, ChannelEventGenerator never set the
beam flags, and fill_lhe_event fell back to writing mu_R into all three of
mur, muf1 and muf2. The PDF itself was always evaluated at the right scale
- only the reporting was wrong - but it made every mu_F comparison against
madevent silently a mu_R comparison.

A dynamical scale choice can give the two beams different scales, and
madevent's MLM definition always does, so these are worth the two doubles
per event whether or not the run is being reweighted.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Five separate divergences from Template/LO/SubProcesses/reweight.f, each
found by comparing mu_F per jet multiplicity for Z+0,1,2,3j against a
madevent reference and each measurable on its own:

  - mt2ij was never recorded. Every initial-state clustering keeps the mT
    of its final-state daughter (cluster.f:841) and pt2ijcl(jcentral) is
    overwritten by it (reweight.f:1059).

  - The 2 -> 1 root lies on both beam lines, and reweight.f's do j=1,2
    walks it once per beam. The extra step here was attached to a single
    "winning" beam, so jlast never reached the root. Largest effect.

  - pt2ijcl(jlast) is raised in place, so when jlast and jcentral are the
    same vertex - the usual case, both lines running QCD to the root - the
    central read sees the raised value and mu_F is that one scale rather
    than a geometric mean of two.

  - The root's scale is a lab-frame transverse mass: cluster.f rotates and
    boosts back before taking it. Computing it in the frame the walk ended
    up in agrees at 0 jets, where no boost has happened yet, and drifts
    from 1 jet on - which is exactly how it was found.

  - The chain handling a beam line that is never coloured was missing.

mu_F per multiplicity, against madevent with pdfwgt off:

           0j     1j     2j     3j
  madevent 78.1  121.0  167.2  218.6
  madspace 78.1  138.1  184.8  220.4

What is left is goodjet: madevent propagates a per-line "every clustering
so far was a jet one" flag, while this kernel reads each daughter's static
flavour. That makes the parton line survive to the root too often, which
shows up as too many events sitting exactly at mT(ll) at 2 and 3 jets.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
reweight.f decides whether a beam line carries on past a vertex from
goodjet: a propagated flag saying every clustering that line has been
through so far was a jet vertex, so one non-jet vertex anywhere in its
history stops it for good. The kernel was reading each daughter's static
flavour instead, which cannot see that history.

Both are available through beam.parton_line_scheme, defaulting to the
faithful "goodjet". The flag only matters under scale_scheme = "madevent".

Two smaller divergences went with it. goodjet is cleared on a beam line
once it stops being a parton line, and at the 2 -> 1 root the mother seen
from one beam is the *other* beam's line as it stands there - ipdgcl(imo)
with imo = idacl(n,3-i) - not the external parton it started as, so its
flavour and colour are now tracked through the clustering rather than read
off beam_flags.

Measured, this changes almost nothing: for Z+0,1,2,3j the two schemes give
mu_F within the run-to-run scatter at every multiplicity. That is a real
result rather than a dead switch - forcing the emitted object to count as
a non-jet moves 3 jets from 219 to 180 GeV, so the code path is live, and
the schemes coincide because every object a beam emits in this process is
an external parton whose goodjet was never cleared. It is therefore not
the explanation for the residual gap against madevent at 2 and 3 jets:

           0j     1j     2j     3j
  madevent 78.1  121.0  167.2  218.6
  madspace 78.1  137.9  182.7  218.9

I could not build a unit test that separates the two schemes: it needs a
beam to emit an object formed at a non-jet vertex and then cluster again,
and neither test fixture produces that on a winning history, with or
without photons substituted into the final state.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
madevent adjusts the jet cuts whenever xqcut is on, in
Template/LO/SubProcesses/setcuts.f: it raises ptj to xqcut, which costs
no phase space that survives the merging anyway and buys a lot of
integration efficiency, and it drops the jet dR cuts, which xqcut
supersedes and which would otherwise carve a hole out of exactly the
region the parton shower is meant to fill. Without it the two codes were
not being compared on the same phase space: madspace kept ptj at 20 GeV
and a jet dR of 0.4 that madevent had already removed.

Only the ptj half is applied. mjj has no counterpart here - "mass" is a
single-particle observable in madspace, so a jet-mass cut asks each
massless parton to be heavier than xqcut and rejects every event. In
madevent mjj is an efficiency device rather than a physics cut, so
nothing is lost but some sampling speed.

Measured on p p > e+ ve + 0..3j, cross section per jet multiplicity in pb
against madevent at the same xqcut:

  xqcut = 15      1j     2j     3j
    madevent    2941   1682    907
    before      1527    544    202
    after       2116    844    333

The remaining deficit at high multiplicity is a separate question, as is
the 3-jet channel's poor convergence at tight xqcut.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
…e-madevent-scales

# Conflicts:
#	madgraph/iolibs/template_files/mg7/launch.py
#	madspace/src/phasespace/topology.cpp
madevent's auto_ptj_mjj sets both ptj and mmjj to xqcut. Only the ptj half
was ported, because "mass" was a single-particle observable here and there
was no way to ask for the invariant mass of a jet pair. The merge brings in
obs_m_inv, which is that cut and which also raises the propagator's e_min,
so the jet-pair invariants a composite s-channel node samples are bounded
the way madevent bounds them.

This matters for more than tidiness: ptj only reaches the phase space for
jets that are leaf children of the root. A jet inside a composite node has
its pt_min reset to zero in PhaseSpaceMapping, and mjj is what bounds that
node instead.

It moves the cross section by about 3% at a loose cut (9853 -> 9549 pb for
p p > e+ ve + 0..3j at xqcut = 15).

It does NOT fix the instability of the highest multiplicity at a tight cut,
which I had expected it to. Suppression between xqcut 15 and 60, per
multiplicity, against madevent:

              0j     1j     2j     3j
  with mjj   1.1x  10.8x  24.5x   0.8x
  ptj only   1.0x  10.3x  23.6x   1.0x
  madevent   1.0x   9.6x  21.2x  47.8x

That is a separate defect, and the measurement above says it is not the
merging cut missing from the phase space. It is overweight replication: at
xqcut = 60 the 3-jet sample was 1635 events drawn from 30 distinct
phase-space points, 54.5 copies of each, while every other multiplicity
sits at exactly one copy per point. The integral stays right because the
replication conserves it; only the event sample is wrong.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Each channel's share of the combined sample comes from its apportioned
target, not from how many unweighted events it actually wrote, so a channel
whose acceptance collapses can run out. EventFile::read then returns false
and leaves the buffer untouched, but read_and_combine ignored that and reset
buffer_index to 0 anyway, so the previous buffer was served again and its
events were emitted repeatedly.

For p p > e+ ve + 0..3j at xqcut = 60 GeV the 3-jet channel starves, and the
combined sample contained 1635 three-jet events drawn from 38 distinct
records, one of them repeated 120 times, byte for byte. The emission order
showed the replay directly: the first 32 were distinct, then it cycled back
through what was left in the buffer. Every other multiplicity was unaffected
at exactly one event per point.

The integral was never wrong - the requested number of events was still
written with the right weights - but the multiplicity composition was, so
jet rates, DJR and multiplicity fractions taken off the sample were
meaningless wherever a channel starved.

Now a channel that cannot supply its share is retired from the draw and the
event is taken from the others. If nothing is left the batch is short, which
ends the combine loop; a sample smaller than requested is the honest outcome
when the events do not exist.

Two things had to be handled to make that safe. Retiring a channel leaves it
with zero width in the cumulative distribution, but random_index == 0 still
resolves to the first entry, so a retired channel could still be drawn and
decrementing its count wrapped the unsigned counter - which is what made my
first two attempts at this spin. And retiring reduces the number of events
still available below the batch that was planned, so the draw has to stop
rather than sample from an empty distribution.

Measured, p p > e+ ve + 0..3j, sigma per multiplicity in pb, and the
suppression each multiplicity gets when xqcut goes from 15 to 60 GeV:

                   0j      1j     2j     3j
  before        5597.1   215.2   52.6  355.2      (3j suppression 1.9x)
  after         5942.1   235.1   54.6   11.2      (3j suppression 59.6x)
  madevent      5925.6   298.6   77.4   18.3      (3j suppression 47.8x)

Suppression has to grow with multiplicity, since every extra jet must pass
the same cut; before this it did not.

No new failures in the madspace suite against main.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
The helper says "as run(), plus the trailing xqcut_weight" but returned
whatever the clustering happened to produce, so every caller unpacked a
fixed count from a growing tuple. Adding an output to the kernel broke six
tests that had nothing to do with it.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
A merged event is generated at one renormalisation scale, but the ladder it
represents was not: each emission belongs at the scale of its own splitting.
madevent applies that correction in the rewgt loop of
Template/LO/SubProcesses/reweight.f, multiplying the weight by
alphas(pt_i)/alphas(mu_R) once per clustering vertex that produced a parton.
madspace kept alphas(mu_R)^n, so it was short by one factor per emission and
the deficit compounded with multiplicity.

The kernel now returns the scale of every clustering vertex, mu_R itself for
the vertices it does not reweight so that their ratio is one and the consumer
needs no mask, and the integrand turns those into the product. A vertex at or
below 2 GeV, where the coupling is not to be trusted, kills the event as it
does in reweight.f; that is its own weight rather than part of xqcut_weight,
which is the merging cut and nothing else.

p p > w+ 0-3j, xqcut 15, 150k events each side, as a ratio to madevent per
jet multiplicity:

              0j     1j     2j     3j   total
  before    1.000  0.764  0.600  0.461  0.839
  after     1.000  0.921  0.864  0.819  0.946

0 jets is untouched, as it must be: there is nothing to cluster.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Three ways to evaluate the coupling of a merged event, since the per-vertex
one is a choice and not the only defensible one:

  per_vertex      alphas(pt_i) at each clustering vertex, what madevent does
  geometric_mean  one coupling at the geometric mean of the pt_i, raised to
                  the number of vertices: one scale for the whole ladder
                  rather than one per rung
  none            a single coupling at mu_R, the behaviour before it existed

Default per_vertex, and only read under dynamical_scale_choice = "mlm".
"none" is there to measure what the reweighting is worth, which is how the
numbers in the previous commit were obtained.

The two reweighting schemes are not distinguishable in practice: for
p p > w+ 0-3j at xqcut 15, 150k events each, they give 10746.9 and 10760.4
pb, 0.1% apart, and their per-multiplicity ratios to madevent agree to
within half a percent at every multiplicity:

                     0j     1j     2j     3j   total
  geometric_mean   1.000  0.921  0.860  0.811  0.944
  per_vertex       1.000  0.921  0.864  0.819  0.946

That is worth knowing rather than assuming: the ladder is short enough here
that the spread of the individual pt_i does not matter, so the cheaper scheme
loses nothing. It may well not hold once the ladder is longer.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
The TOML card is rendered from a template of explicit %(section.key)s
placeholders, not from the parameter list, so a parameter added with
add_toml_param and not added here exists but never appears in a generated
card. parton_line_scheme and alphas_reweighting were both in that state:
readable if a user typed the line themselves, invisible otherwise.

Found while adding a third one and wondering why it was not in the card.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
The second thing madevent does to a merged event that a single scale does not
describe. A merged event's beam density belongs at the scale of the emission
that pulled the parton out of the beam, not at the scale of the hard process,
and madevent gets there in two moves (the pdf half of rewgt in
Template/LO/SubProcesses/reweight.f, under the hidden pdfwgt flag, which is on
by default there): it evaluates the density at min(pt(jfirst), mu_F), low on
the clustering ladder, then walks each beam line back up, multiplying by
f(x z, Q_i) / f(x z, Q_i-1) at every further clustering that line takes part
in, with the momentum fraction rescaled by that clustering's z as it goes.

The kernel decides the scales, the momentum fractions and the beam each step
sits on, and records the clustering z it did not need before. What it cannot
decide is the flavour: one clustering serves every flavour channel of a
subprocess, external_pdg_ids is a single representative with the signs
stripped, and the flavour is in any case sampled from the very densities this
is correcting. So the kernel names a *class* - the gluon, one of the two
beams' own flavours, or a flavour the diagram pinned down - and the integrand
resolves it once the sampling has happened. Gluon-versus-quark and an
unchanged quark line are exact; a line whose flavour changed along the way,
which for W+jets means past the W vertex, falls back to the representative.
That is an approximation, but only on lines the representative was already
guessing at, and a density ratio between two nearby scales is far more
sensitive to quark-versus-gluon than to which quark.

Under beam.pdf_reweighting, default true to match madevent, and only live
under dynamical_scale_choice = "mlm" with scale_scheme = "madevent", since
without that scheme there is no beam parton line to walk.

p p > w+ 0-3j, xqcut 15, 150k events each side, NNPDF23_lo_as_0130_qed, as a
ratio to madevent per jet multiplicity:

                      0j     1j     2j     3j   total
  madevent pdfwgt   on/off, for reference:
                  1.000  0.929  0.858  0.797  0.9455
  madspace off      1.000  0.991  1.005  1.028  1.0001   (vs madevent off)
  madspace on       1.000  0.999  1.005  0.981  0.9987   (vs madevent on)

10760.0 -> 10761.1 pb with both off, 11380.0 -> 11365.7 pb with both on. The
8-22% deficit that was left after the alpha_s reweighting was all of it this.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Same template gap as the previous commit, spotted while checking whether it
was on during a madevent comparison: declared in banner.py, read by launch.py,
absent from the rendered card. Its default is what a merged run wants, so
nothing was behaving wrongly - it just could not be seen or turned off.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
A merging cut exists to separate matrix-element radiation from what the parton
shower would have produced. It has no business rejecting an event whose jet is
already above it - and that is what this was doing, to 19% of a one-jet sample
at xqcut = 60.

The route in is a vertex that produced no radiation at all. In a quark-jet
topology there is a diagram in which the outgoing quark and the quark line
combine into the W:

  g u > e+ ve d   [['i1','i0','p0'], ['o0','o1','p1'], ['p0','o2','p1']]   p1 = W

so once the clustering has formed the W from the leptons, (W, jet) is a valid
clustering. Its measure is a kt between the W and the jet, which has nothing to
do with the jet's own transverse momentum and is frequently below xqcut, and
the jet was rejected there. A gluon has no W-gluon vertex, can only ever pair
with a beam, and its measure is then mT = pt, so gluon jets never saw this.

madevent gates the same cut on iqjets, which is only ever set for a leg emitted
at a jet vertex; a leg it does not set is exempt. The comment this replaces said
the cut "does not ask whether the vertex was a QCD one", which read reweight.f's
xqcut block without following iqjets back to where it is filled.

p p > e+ ve + 0-3j, 150k events each side, as a ratio to madevent:

              xqcut      0j     1j     2j     3j   total
  before         15   1.001  1.001  0.999  0.978  0.9988
  after          15   1.000  0.997  1.001  0.978  0.9976
  before         30   0.999  0.963  0.920  0.916  0.9882
  after          30   1.000  0.961  0.996  0.971  0.9938
  before         60   0.998  0.801  0.616  0.696  0.9833
  after          60   1.000  0.972  1.005  1.085  0.9993

and on the isolated case, g u > e+ ve d with ptj = 60, the cut goes from
removing 19% to removing nothing: 54.903 -> 68.328 against 67.751 with the cut
switched off.

This also supplies the demotion madspace was missing - madevent's iqjets is
worth 6.2% at xqcut = 15 - which is why the loose merging cut does not move
while the tight one does. Dropping the resonance preference instead fixes the
same veto but overshoots to 1.0168 at xqcut = 15, because it removes only one
half of that cancellation.

What is left is the 1-jet rate, 3-4% low at xqcut 30 and 60 and correct at 15.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
The previous commit gated the merging cut on the vertex being a QCD one, which
is a proxy for what madevent actually asks. This is the thing itself: iqjets,
built in the same walk that already produces jfirst/jlast/jcentral.

A leg is tagged 1 where it was emitted while the beam line was still a parton
line, or where it is a gluon - a gluon emission has a soft singularity the
shower has to be left free to fill. Otherwise it is tagged with jcode, a count
of how many times the walk has crossed between jet and non-jet vertices, and
resolved at the end: once either beam line has stopped being a parton line,
everything tagged only provisionally and before the last crossing was not a
merging emission after all. Final-state vertices tag their bare jet daughters
through isjetvx, and a vertex that made no jet demotes the legs whose line runs
into it. Only a leg this leaves above zero is checked against xqcut.

Three flavour facts the state machine did not carry are needed for that and are
added: whether the mother is an octet, whether it repeats either daughter's
flavour (isjetvx reads that as "emitted off a line" rather than "split into
something else"), and whether all three lines are colourless, which is what
keeps a W W Z or h h h vertex out of the demotion. Octet-ness per external leg
rides in the spare high half of jet_leg_mask.

Against the proxy it replaces, the two agree within statistics everywhere, so
this changes no number here - it is the rule that will still be right away from
W+jets, where "QCD vertex" and "vertex that made a jet" part company.

  madspace / madevent          0j     1j     2j     3j   total
    xqcut 15  proxy         1.000  0.997  1.001  0.978  0.9976
    xqcut 15  iqjets        0.999  0.999  1.003  0.981  0.9981
    xqcut 30  proxy         1.000  0.961  0.996  0.971  0.9938
    xqcut 30  iqjets        1.001  0.956  0.994  0.981  0.9936
    xqcut 60  proxy         1.000  0.972  1.005  1.085  0.9993
    xqcut 60  iqjets        0.999  0.965  1.012  1.110  0.9983

and on g u > e+ ve d with ptj = 60 the merging cut still removes nothing:
68.328 against 67.751 with it switched off.

What is left is the 1-jet rate, 3-4% low at xqcut 30 and 60 and correct at 15.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
rewgt reads pt2ijcl(n) after setclscales has already overwritten it - at the
central vertex with the emitted object's transverse mass, and at the last one to
keep it above the first - so the coupling never sees the raw clustering measure.
madspace was reweighting at the raw measure.

Measured on p p > e+ ve + 0-3j and p p > t t~ + 0-3j, this changes nothing at
all: byte-identical at xqcut 15 and 30 for W+jets, and 2904.6 -> 2903.2 pb and
1742.6 -> 1743.4 pb for t t~ at xqcut 20 and 30, both well inside the errors.
The reason is that the two overrides land on jcentral and jlast, which for these
processes are the 2 -> 1 root, and the alpha_s loop already skips the root.

Committed anyway because it is what reweight.f reads, and a process whose
central or last vertex is not the root would see it. It costs one array of
cluster_max doubles.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
When the last clustering joins two coloured final-state lines into a coloured
one - q q~ > g > t t~ is the case that matters - and both beam lines are still
parton lines at the 2 -> 1 root, setclscales replaces the root's transverse
mass by mt2last = sqrt(djb(d1) djb(d2)), the geometric mean of the pair's
transverse masses (cluster.f sets it, reweight.f uses it). madspace had no
counterpart, so the root kept the mT of the whole final state, which for
q q~ > t t~ is sqrt(shat) - 2.3 times mT(t) on average.

The kernel now records that mean at the walk's last clustering when it is a
final-state one, and uses it as the root's mT under madevent's conditions:
above 2 GeV (mt2last > 4 GeV^2), both jlast at the root, and the last clustering
a QCD vertex, which is exactly what bit 27 already encodes (both daughters and
the mother coloured).

p p > t t~ alone, 100k events per arm, MLM scale without reweighting:

                    madevent   madspace before   madspace after
  gg  > t t~          448.0        448.3              448.3
  qq~ > t t~           65.23        53.46              65.13
  total               513.20       501.72             513.41   (ratio 1.0004)

The qq~ deficit was closed beforehand by moving madspace's own events from
sqrt(shat) to mT(t) (53.46 -> 64.87 pb), and the scale is now mT(t) on every
event in both channels in both codes. W and Z + 0 jets never needed it: there
the root pair is colourless and sqrt(shat) equals the boson's mT exactly.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
find_t_vertices collects the t-channel propagators on its way back from beam 1,
so t_propagator_masses and t_propagator_pdg_ids start at the propagator next to
that beam. The state-machine compiler read them the other way round: the mask
built from beam 1 took at(t_count - 1 - k) and the one from beam 2 at(k). For a
chain with a single propagator the two coincide, which is why W + 1 jet and
every 0-jet process never showed it.

With two or more the flavour and mass of the two ends were swapped. In
g g > t t~ g, beam 0 emitting the jet leaves a gluon line, and the compiled
transition said massive and not a jet - the top's properties - while beam 1
taking back the t~ said massless jet. So the beam's parton line stopped at the
jet, jlast landed on the jet's own vertex, and mu_R became
(pt(jet) mu_F^3)^(1/4) instead of the central scale: in t t~ + 1j at xqcut 30
that was 76% of madspace's events, against mu_R == mu_F in 94.8% of madevent's.
The same swap also fed the wrong colour and pdf-reweighting flavour class into
every such initial-state clustering.

The new test decodes the compiled transitions of that chain; it fails on the
old ordering with (massive, not a jet) for beam 0 + jet.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
cluster.f boosts after an initial-state clustering only when nleft > 4 ("Don't
boost if boost vector too lightlike or last vertex"), i.e. never at the last
recorded clustering, and takes the 2 -> 1 root's scale as djb of what is left in
whatever frame that leaves: the lab after a final-state last clustering, which
it boosts back first, and the current frame after an initial-state one.
madspace boosted after every initial-state clustering and always took the root
in the lab.

In t t~ + 1j with the jet taken into a beam first, the frame madevent ends in is
the one where t and t~ balance, so the top's vertex and the root carry the same
scale and mu_R == mu_F. Taking the root in the lab left madspace's mu_R a few
per cent below mu_F even with the chain flavours right. W + 1 jet is unchanged
by construction: there is no earlier boost, and madevent's mu_F equals the lab
mT(e nu) in 84.8% of its events.

After a final-state last clustering cluster.f undoes only the most recent boost;
the lab copy undoes all of them, which is the same as long as at most one
boosting initial-state clustering came before.

The new test checks mu_R == mu_F on the jet-into-beam histories of
g g > t t~ g; it passed on no point with either this or the chain ordering
wrong.

Together with the previous commit, p p > t t~ + 0-3j against madevent, 50k
events each, madevent samples unchanged, ratio before -> after:

  xqcut     0j            1j            2j            3j          total
    20   0.990->0.990  1.029->0.984  1.033->0.976  1.103->1.041  1.039->0.995
    30   0.995->0.995  1.047->1.008  1.056->1.005  1.111->1.062  1.041->1.010
    45   0.995->0.997  1.028->0.996  1.073->1.025  1.124->1.078  1.024->1.005
    60   0.994->0.994  1.039->1.015  1.046->1.005  1.041->1.009  1.014->1.002

and p p > e+ ve + 0-3j moves only within statistics (totals within 0.1%, 0j/1j
unchanged). What is left in t t~ is the 3-jet rate, 4-8% high at xqcut 20-45.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
beam.clustering_measure picks how a candidate final-state clustering is
scored under dynamical_scale_choice = "mlm". "fxfx" (default, unchanged
behaviour) keeps the NLO cluster_scale; "madevent" uses the LO DJ of
kin_functions.f, where a massless-massive pair (q* > q W, g* > g h) scores
the massless leg's transverse mass times 1 + 1e-6, and settles the exact tie
with an initial-state clustering the way cluster.f does.

The FxFx score of those pairs is often below the jet's pt, which is what left
W + 1 jet and HEFT h + 1 jet low against madevent.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Conflicts: lhe_output (keep the <scales> tag with mu_R/mu_F and pt_clust, and main's <mgrwt>/<rwgt> blocks).

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
No content change: this branch already carried main and the same conflict resolutions (lhe_output keeps both the <scales> mu_R/mu_F attributes and main's <mgrwt>/<rwgt> blocks). The merge makes feat-mlm an ancestor, so the PR shows only this branch's work.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
The merging-cut adjustment still tightened a cut called jet-m_inv, the name the observable had when auto_ptj_mjj was written. main has since renamed it to pair_mass, so every mg7 run with phasespace.xqcut > 0 died building the cut: ms.Observable rejects observable='m_inv'.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
oliviermattelaer and others added 2 commits September 16, 2026 05:07
The merging-cut adjustment names a cut after the observable it wants, and parse_observable feeds that name straight to ms.Observable, but nothing checked the two agree: when madspace renamed m_inv to pair_mass, every merged run with xqcut > 0 died building its cuts. These call both methods unbound on a stand-in run card, so they need no process directory, and assert each cut written is one madspace accepts.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
They need a compiled madspace to check the launcher's cut names against ms.Observable, and the unittest jobs have none, so listed there they would only ever skip. This workflow already builds madspace once and every job restores it.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

This branch has not been deployed

No deployments
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant