[Rivet] Bug in Rivet?

[Andy Buckley]

On 12/04/16 17:26, Dag Gillberg wrote:
> Thanks Andy!
>
> I will ask the HXSWG ppl about what they had in mind regarding the
> "dodgy physics" stuff.

Ok! Make sure not to quote me as as expert on this ;-) I'm interested to 
know if there is some unambiguous (or at least very good approximation) 
way to do this at arbitrary perturbative order.

> Regarding the B-mesons - I meant B-baryons ... and D-baryons apparently.
> I'll put a few print-outs below.
> This is for 13 TeV gg→H MC, disabling warnings for Lambda particles
> which are much more frequent than B baryons in this sample (there are no
> b-quarks at tree level).
>
> So 5232, 4232, 5122 ...

I tested these with the fixed definitions and they work fine. So I 
expect the patch and upcoming 2.4.2 release will also work fine for you. 
Let me know if you have trouble with the patch before I push out a 
release that doesn't actually fix the problem!

Andy


> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (2.2,-1.3,4.8,0), ID = -11
>    - partent -5232
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (1.7,-1.1,4.3,0), ID = 12
>    - partent -5232
> SimplifiedXSec::execute() INFO    Processed 0.3k events
>
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (1.2,-1.4,6.1,-0), ID = -16
>    - partent -15
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (0.3,-2.2,0.0,0), ID = 16
>    - partent -5232
> SimplifiedXSec::execute() INFO    Processed 0.4k events
>
> SimplifiedXSec::execute() INFO    Processed 0.5k events
>
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (0.3,-3.8,4.0,0), ID = -11
>    - partent 4232
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (0.1,-7.6,5.3,0), ID = 12
>    - partent 4232
> SimplifiedXSec::execute() INFO    Processed 0.6k events
>
> SimplifiedXSec::execute() INFO    Processed 0.7k events
>
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (0.2,-0.3,5.1,0), ID = 11
>    - partent -4232
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (1.7,0.1,3.0,0), ID = -12
>    - partent -4232
> SimplifiedXSec::execute() INFO    Processed 0.8k events
>
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (1.0,-0.9,5.8,0), ID = -11
>    - partent -5122
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (1.9,0.6,5.9,-0), ID = 12
>    - partent -5122
> SimplifiedXSec::execute() INFO    Processed 0.9k events
>
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (0.3,-1.9,5.9,0), ID = -13
>    - partent -5132
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (0.4,1.4,4.6,-0), ID = 14
>    - partent -5132
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (2.8,3.8,3.5,0), ID = 11
>    - partent 5122
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (0.8,3.3,2.3,0), ID = -12
>    - partent 5122
> SimplifiedXSec::execute() INFO    Processed 1.0k events
>
> SimplifiedXSec::execute() INFO    Processed 1.1k events
>
> SimplifiedXSec::execute() INFO    Processed 1.2k events
>
> SimplifiedXSec::execute() INFO    Processed 1.3k events
>
> SimplifiedXSec::execute() INFO    Processed 1.4k events
>
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (0.6,0.9,4.6,0), ID = 13
>    - partent 5122
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (2.1,1.5,5.8,-0), ID = -14
>    - partent 5122
> SimplifiedXSec::execute() INFO    Processed 1.5k events
>
> SimplifiedXSec::execute() INFO    Processed 1.6k events
>
> SimplifiedXSec::execute() INFO    Processed 1.7k events
>
> SimplifiedXSec::execute() INFO    Processed 1.8k events
>
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (1.2,3.4,4.8,0), ID = 11
>    - partent 11
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (1.6,3.5,5.0,0), ID = -11
>    - partent -11
> SimplifiedXSec::execute() INFO    Processed 1.9k events
>
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (0.7,3.3,6.2,0), ID = -13
>    - partent 4122
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (0.4,3.2,2.5,0), ID = 14
>    - partent 4122
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (0.6,-1.1,0.8,0), ID = 11
>    - partent -4232
> Warning: found unexpected particle!
>    (pT,eta,phi,m) = (0.3,-1.6,0.2,0), ID = -12
>    - partent -4232
>
>
>
>
>
> On Tue, Apr 12, 2016 at 11:52 AM, Andy Buckley <andy.buckley at cern.ch
> <mailto:andy.buckley at cern.ch>> wrote:
>
>     On 12/04/16 14:51, Dag Gillberg wrote:
>
>         Hi Andy,
>
>         Thanks!
>         I saw in some more tests that not only Lambda particles are
>         affected,
>         but also the B-mesons.
>         This can actually slightly have affected our old Rivet routine
>         (for the
>         Higgs differential combination) since we might claissify
>         electrons from
>         a semileptonic heavy flavor hadron decay as a non-hadronic..
>
>
>     Which b-mesons are causing trouble? I think my error was only for
>     baryons, and the digit ordering should be absolute for mesons.
>
>     I checked with an equivalent set of functions, and at least 511,
>     521, and their excited states do seem to work fine now, but I'm not
>     aware of having made any directly relevant change in the patch.
>
>     (BTW, in the replacement file I've completely removed the check for
>     baryons rather than enforce the exceptional conditions, since the
>     particle listing that I had handy in ASCII form includes some baryon
>     IDs which don't seem to fit into the published scheme at all.)
>
>         I'm currently using Jim's setup, to run Rivet within
>         AnalysisBase (i.e.
>         RootCore).
>         Hmm.. not sure how this works, but I see:
>
>         url = $ROOTCOREBIN/../HEPSoftware/Rivet-2.4.1.tar.gz
>
>         So I guess this is being downloaded and compiled somehow.
>
>
>     Sounds like you can make do for now, but I'd better get a new
>     version ready for validation and patch release.
>
>         Now, regarding the code.
>         We have two open questions, where we perhaps can pick your brain :)
>
>         First, note that this is not a pure fiducial/differential
>         measurement.
>         We are using kinematic information, but are allowed to "cheat"
>         and use
>         quarks etc.
>
>         1. We are supposed to classify
>             gg → ZH
>         vs
>             qq → ZH
>         Which, at NLO, I don't think is possible...
>
>         gu → ZHu
>
>         could be either of the above. Any ideas here?
>         Is this split simply impossible with NLO MC?
>
>
>     I think the split is physically quite dodgy in general! There's no
>     physical separation between gluon emissions and splittings that take
>     place in the initial state shower and which take place in the matrix
>     element... so at NLO where the ME can contain (sums of) such
>     splittings on the incoming legs without any "external" indication,
>     the distinction is false. The physical way would be to do
>     "approximate" subprocess classification by using kinematic info to
>     choose gg or qq (or whatever) enhanced regions... if there is a
>     kinematic distinction.
>
>     This is vague, armchair-physicist stuff from me, though, and
>     sometimes there are "loopholes" (like in single top where there
>     really is no interference between s- and t-channels because the
>     final states are sufficiently distinct).
>
>         2. Regarding the jets:
>         We want to take all stable particles, remove anything from the Higgs
>         (very easy with Rivet!).
>
>
>     As long as there's a Higgs in the event!
>
>         Also ignore particles from leptonic Z and W decays (leptons, e,
>         µ, and
>         FSR γ)
>
>         Is there a standard way to do this?
>
>
>     Well, this is pretty well-aligned with isPrompt(p). That will return
>     false if there are any hadrons in the ancestor chain, otherwise
>     true. This includes W and Z production, without accidentally
>     tripping up on the off-shell W's that Herwig++ sometimes likes to
>     insert into semileptonic hadron decays. But it will also (IIRC)
>     classify FSR photons from quarks as prompt -- and you probably don't
>     want to exclude those. The more explicit way would be to configure
>     WFinder and ZFinder instances and exclude particles which appear in
>     their particles() lists. We don't currently have any built-in
>     handling for Z and W decays via taus, though...
>
>     Let us know what works for you, and if there's any machinery we
>     could usefully add (or receive as a patch :-P )
>
>         I suppose we can just do exaclty what we did for the Higgs routine:
>
>         https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/HiggsPhys/HSG1/Analysis/trunk/GeneralTools/HggUnfoldingPaper2013/Rivet/Combination/HiggsDiffXSecCombination.cc#L110
>
>         which btw use the projections :)
>
>
>     :-)
>
>     Actually, I think I cleaned this one a bit recently, for our
>     forthcoming C++11-enabled 2.5.0 release:
>     https://rivet.hepforge.org/hg/rivet/file/66e20000b87f/src/Analyses/ATLAS_2015_I1364361.cc
>
>     Andy
>
>
>         On Tue, Apr 12, 2016 at 4:04 AM, Andy Buckley
>         <andy.buckley at cern.ch <mailto:andy.buckley at cern.ch>
>         <mailto:andy.buckley at cern.ch <mailto:andy.buckley at cern.ch>>> wrote:
>
>              On 12/04/16 03:20, Dag Gillberg wrote:
>
>                  Hi Andy,
>
>                  I'm writing Rivet code (yay! fun!).
>                  And I want to build jets out of stable particles that
>         are hadrons or
>                  from hadron decay.
>                  So I do:
>                        if ( PID::isHadron(p.pdgId()) || p.fromHadron() )
>
>                  Now, as a cross check, I printed out particles that
>         failed this,
>                  and I
>                  see a ton of strange baryons fail:
>                  Warning: found unexpected particle!
>                      (pT,eta,phi,m) = (0.9,-2.0,2.0,1), ID = 3122
>                  Warning: found unexpected particle!
>                      (pT,eta,phi,m) = (0.5,-4.3,1.9,1), ID = -3122
>                  Warning: found unexpected particle!
>                      (pT,eta,phi,m) = (0.9,-2.2,3.9,1), ID = 3122
>
>                  These are Lambda particles.
>                  They all fail this line in isBaryon:
>
>                  if (_digit(nq2,pid) < _digit(nq3,pid)) return false;
>
>                  since the second digit (1) is smaller than the third (2)
>                  ( not sure what this condition is based on though..)
>                  But a Lambda is certainly a Baryon..
>                  So this must be a bug, right?
>
>
>
>              Hi Dag,
>
>              Yes, you're right: I had tried to tighten up my
>         implementation of
>              PID checks in the MCUtils package, following the PDG scheme
>              definition, but missed that for baryons there are
>         exceptions to the
>              ordering.
>
>              I fixed this a while ago in MCUtils but had forgotten to
>         sync the
>              fixes into Rivet. I'm not sure how it managed to not show
>         up in our
>              release validation -- we should look into that.
>
>              Anyway I've done the sync now, and the fixed header file (from
>              include/Rivet/Tools) is attached -- are you able to patch,
>         or do you
>              rely on a pre-built copy of Rivet? Since it's a bug with real
>              physics implications, I'm tempted to get a patch release
>         out asap.
>
>              Thanks for the report!
>
>              Now a couple of unrelated comments on your code snippet ;-)
>
>                          Particles hadrons_notFromHiggs;
>                          FourMomentum sum(0,0,0,0);
>                          for ( const GenParticle *ptcl :
>                  Rivet::particles(event.genEvent()) ) {
>                            Particle p(ptcl);
>                            if (!p.isStable()) continue;
>
>
>              Since you only want stable particles, rather than loop over
>         HepMC
>              GenParticles, convert them to Rivet Particles, and ignore the
>              non-stable ones the "Rivet way" is to set up a FinalState
>         projection:
>
>              // In init()
>              addProjection(FinalState(), "FS")
>
>              // In apply()
>              const Particles& particles = applyProjection<FinalState>(event,
>              "FS").particles();
>              for (const Particle& p : particles) {
>                 ...
>
>              I am trying to move some core intelligence out of the
>         projections
>              and into standalone functions (or maybe Event methods) to
>         make it
>              easier to do many things without projections if wanted,
>         though --
>              the main benefit of projections is automatic caching, but
>         that's not
>              always the most important thing.
>
>              Cheers,
>              Andy
>
>              --
>              Dr Andy Buckley, Lecturer / Royal Society University
>         Research Fellow
>              Particle Physics Expt Group, University of Glasgow
>
>
>
>
>     --
>     Dr Andy Buckley, Lecturer / Royal Society University Research Fellow
>     Particle Physics Expt Group, University of Glasgow
>
>


-- 
Dr Andy Buckley, Lecturer / Royal Society University Research Fellow
Particle Physics Expt Group, University of Glasgow