  The observation of double and triple scalar boson production at hadron
colliders could provide key information on the Higgs self couplings and the
potential. As for single Higgs production the largest rates for multiple Higgs
production come from gluon-gluon fusion processes mediated by a top-quark loop.
However, at variance with single Higgs production, top-quark mass and width
effects from the loops cannot be neglected. Computations including the exact
top-quark mass dependence are only available at the leading order, and
currently predictions at higher orders are obtained by means of approximations
based on the Higgs-gluon effective field theory (HEFT). In this work we present
a reweighting technique that, starting from events obtained via the MC@NLO
method in the HEFT, allows to exactly include the top-quark mass and width
effects coming from one- and two-loop amplitudes. We describe our approach and
apply it to double Higgs production at NLO in QCD, computing the needed
one-loop amplitudes and using approximations for the unknown two-loop ones. The
results are compared to other approaches used in the literature, arguing that
they provide more accurate predictions for distributions and for total rates as
well. As a novel application of our procedure we present predictions at NLO in
QCD for triple Higgs production at hadron colliders.
