<html>
<head>
<meta content="text/html; charset=ISO-8859-1"
http-equiv="Content-Type">
</head>
<body bgcolor="#FFFFFF" text="#000000">
Hi Andy, <br>
<br>
I looked at a handful events with the E-scheme and you are right,
the mismatch between jets from Rivet and FastJet is gone. I see,
it's the rescaling of particles to make them massless in FastJet
which causes the difference. I think that Rivet should follow
FastJet in this respect. <br>
<br>
I can understand that there is a good reason for having the
Rivet::Jet interface. However, I'm not sure why you have to
recalculate the jet four-vector after FastJet has already done this
work for you. When I use jets, I expect that the four-vector of the
jet is calculated according to the prescription of the recombination
scheme which I pass to FastJet in the constructor:<br>
addProjection(FastJets(dfs, fastjet::kt_algorithm,
fastjet::pt_scheme,
1.), "FastJet");
<br>
This is exactly what happens in the data analysis, so it's best to
use FastJet's calculation. Also I think it's safe for you to assume
that any analysis that uses FastJet, will use its implementation of
the recombination schemes. By the way, the recombination scheme has
to be defined *before* the jet finding is performed, because it
influences the jet finding itself. <br>
<br>
Concerning the pt-weighted jet properties in Rivet, I can't help
you. I'm not sure where these are used and why they have been
implemented. It's of course always useful to have additional methods
available, but to be honest, I find the methods from Rivet::Jet a
bit confusing. For example, what's the difference between<br>
Jet::ptWeightedEta() and Jet::eta()? When using the pt-recombination
scheme there should not be a difference. However, we found out that
there is a difference because of a different treatment of the
particle masses, this can be confusing. Anyways, maybe you would
still like to keep both methods for backward compatibility. I saw
that there is only one analysis which uses the Jet::ptSum() method:
<a class="el"
href="http://projects.hepforge.org/rivet/code/dev/a00315_source.html#l00086">CDF_2001_S4751469</a>.
Do you know which recombination scheme this analysis exactly used?<br>
<br>
In any event, I think the important point is that the jet
four-vector should be filled by FastJet.<br>
<br>
Cheers, <br>
Roman<br>
<br>
<br>
<br>
<br>
<br>
<br>
On 11/18/11 0:40, Andy Buckley wrote:
<blockquote cite="mid:4EC59B73.6080108@ed.ac.uk" type="cite">Hi
Roman,
<br>
<br>
I just had a bit of a look at this -- to check that I've
understood this correctly, is this mismatch issue only seen when
using a jet recombination scheme other than E (i.e. pure 4-vector
addition)?
<br>
<br>
Cross-checking FastJet's manual re. the pT scheme and Rivet's
pT-weighted jet properties, there seems to be a difference in that
the FastJet one rescales momenta to be massless before combining.
Rivet doesn't do that: we just do this
<br>
<br>
foreach (const FourMomentum& p, momenta()) {
<br>
double pt = p.pT();
<br>
ptsum += pt;
<br>
ptwetasum += pt * p.pseudorapidity();
<br>
ptwdphisum += pt * mapAngleMPiToPi(phi0 -
p.azimuthalAngle());
<br>
}
<br>
_ptWeightedEta = ptwetasum/ptsum;
<br>
[...]
<br>
<br>
Maybe that's definitely wrong -- if so we can fix it. Or just
remove it / base it on the accumulated 4-vector. There is probably
work that can be done to improve the Rivet::Jet interface, which
was not very carefully designed and doesn't provide options for
recombination other than the E scheme (but does have other reasons
to exist).
<br>
<br>
Anyway, if I've understood correctly this is an area where I'd
like to hear suggestions about how we can do better, but I'm
breathing a sign of relief because I think all the "official"
analyses that I know are using E-recombination and should be
unaffected.
<br>
<br>
Andy
<br>
<br>
<br>
On 14/11/11 14:12, Roman Kogler wrote:
<br>
<blockquote type="cite">Dear Rivet developers,
<br>
<br>
I believe I found a bug in the jet definition of Rivet. It's a
quite
<br>
annoying feature and it took me quite some time to figure it
out.
<br>
<br>
First, I couldn't reproduce the results that come out of my jet
analysis
<br>
in H1. To find the difference I put some Sherpa events into the
format
<br>
we use in H1 for our analyses. Then I started to compare the
output
<br>
even-by-event, using the same events in my H1 analysis and the
Rivet
<br>
analysis I wrote. The boost of the HFS to the Breit frame is
identical,
<br>
as well as the number of jets found. Even the jet constituents
of each
<br>
jet are the same in both analyses. So finally I tracked the
difference
<br>
to the jet four-vector. It seems that Rivet overwrites the jet
<br>
four-vector as calculated by FastJet. To clarify that I looked
into the
<br>
FastJets projection, Projections/FastJets.C. I put some output
<br>
statements into the method _pseudojetsToJets, the output is
labelled
<br>
"from FastJet" below (all code can be found at the end of this
email).
<br>
Also, I put some output into my analysis using the standard
<br>
jet::momentum() method ("from Rivet"). I then recalculate the
jet
<br>
four-vector with the pt-recombination scheme ("my calculation").
Here
<br>
are the results for one event:
<br>
<br>
jet 0: from FastJet, Pt = 14.6163 phi = 6.27012 eta = 0.406647
<br>
jet 1: from FastJet, Pt = 14.3426 phi = 3.08307 eta = 0.240164
<br>
<br>
jet 0: from Rivet, Pt = 14.0764 phi = 6.27134 eta = 0.443595, pt
<br>
weighted phi = -0.0106254 pt weighted eta = 0.406647
<br>
jet 1: from Rivet, Pt = 13.5777 phi = 3.08787 eta = 0.259143, pt
<br>
weighted phi = 3.09267 pt weighted eta = 0.240164
<br>
<br>
jet 0: my calculation: pt = 14.6163 phi = 6.27012 eta = 0.406647
<br>
jet 1: my calculation: pt = 14.3426 phi = 3.08307 eta = 0.240164
<br>
<br>
The jet four-vectors from FastJet agree with the ones using my
<br>
implementation of the pt-recombination scheme. Also, these
numbers are
<br>
identical with the numbers I get using the H1 analysis framework
(with
<br>
its own interface to FastJet). But the jets stored by Rivet have
<br>
different pt! How could this happen? It is really quite a
dangerous bug,
<br>
I believe. Please have a look at it soon.
<br>
<br>
Cheers,
<br>
Roman
<br>
<br>
<br>
P.S. Here is the code I used for the output above:
<br>
<br>
<br>
-----------------------------------------------
<br>
Projections/FastJets.C:
<br>
<br>
Jets FastJets::_pseudojetsToJets(const PseudoJets& pjets)
const {
<br>
Jets rtn;
<br>
int i=0;
<br>
foreach (const fastjet::PseudoJet& pj, pjets) {
<br>
cout << "jet " << i++ << ": ";
<br>
cout << "from FastJet, Pt = " << pj.perp() <<
" phi = " << pj.phi() << "
<br>
eta = " << pj.rap() << endl;
<br>
...
<br>
}
<br>
<br>
<br>
-----------------------------------------------
<br>
H1 jet analysis:
<br>
<br>
void init() {
<br>
...
<br>
addProjection(FastJets(dfs, fastjet::kt_algorithm,
fastjet::pt_scheme,
<br>
1.), "FastJet");
<br>
...
<br>
}
<br>
<br>
void analyze(const Event& event) {
<br>
<br>
const FastJets & fj = applyProjection<FastJets>(event,
"FastJet");
<br>
Jets jets = fj.jetsByPt(3., 80.);
<br>
<br>
vector<FourMomentum> jet_vecs;
<br>
int njet = 0;
<br>
<br>
foreach (const Jet& j, jets){
<br>
<br>
cout << "jet " << njet++ << ": ";
<br>
cout << "from Rivet, Pt = " << j.momentum().pT()
<< " phi = " <<
<br>
j.momentum().phi() << " eta = " <<
j.momentum().eta()
<br>
<< " pt weighted phi = " << j.ptWeightedPhi()
<< " pt weighted eta = "
<br>
<< j.ptWeightedEta() << endl;
<br>
<br>
// recalculate the jet four-vector with my own implementation of
the
<br>
pt-recombination scheme
<br>
const vector<Particle> parts = j.particles();
<br>
FourMomentum jet_vec(0, 0, 0, 0);
<br>
<br>
foreach(const Particle& jp, parts){
<br>
<br>
// scale particle momtenta to make them massless
<br>
double pp = sqrt(jp.momentum().px()*jp.momentum().px() +
<br>
jp.momentum().py()*jp.momentum().py() +
<br>
jp.momentum().pz()*jp.momentum().pz() );
<br>
FourMomentum pv(pp, jp.momentum().px(), jp.momentum().py(),
<br>
jp.momentum().pz());
<br>
if (ij==0){
<br>
jet_vec = pv;
<br>
} else {
<br>
jet_vec = MergeParts(pv, jet_vec);
<br>
}
<br>
++ij;
<br>
}
<br>
<br>
cout << "my calculation: pt = " << jet_vec.pT()
<< " phi = " <<
<br>
jet_vec.phi() << " eta = " << jet_vec.eta() <<
endl;
<br>
}
<br>
<br>
<br>
FourMomentum MergeParts(FourMomentum a, FourMomentum b)
<br>
{
<br>
// merge momenta a and b with pt-recombination scheme
<br>
<br>
double pt = a.pT() + b.pT();
<br>
double phi_a = a.phi();
<br>
double phi_b = b.phi();
<br>
if (phi_a - phi_b > fastjet::pi) phi_b += fastjet::twopi;
<br>
if (phi_a - phi_b < -fastjet::pi) phi_b -= fastjet::twopi;
<br>
double phi_sum = a.pT()*phi_a + b.pT()*phi_b;
<br>
double phi = phi_sum / pt;
<br>
double eta_sum = a.pT()*a.pseudorapidity() +
b.pT()*b.pseudorapidity();
<br>
double eta = eta_sum / pt;
<br>
<br>
FourMomentum v(pt*cosh(eta), pt*cos(phi), pt*sin(phi),
pt*sinh(eta));
<br>
return v;
<br>
<br>
}
<br>
<br>
<br>
</blockquote>
<br>
<br>
</blockquote>
</body>
</html>