[Rivet] Bug in Rivet?

[Andy Buckley]

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>> 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