JEC to JSON - experience of first shot

[kirschen]

Hi @nsmith- , all,
I gave the JEC-conversion a shot and have a working "V0" version with some shortcomings. In particular, I wasn't sure how to pass parameters to the "Formula", so I transformed the formula max(0.0001,[0]+((x-[1])*([2]+((x-[1])*([3]+((x-[1])*[4])))))) to something messy like max(0.0001,1.7675362825393677+((x-21.433332443237305)*(-0.007086962461471558+((x-21.433332443237305)*(-0.006380382925271988+((x-21.433332443237305)*0.0013125771656632423)))))) with string replacement. FYI: The "V0" directly transforms what is parsed by JetCorrectorParameters (in CMSSW) to the JSON.

The string replacement and "V0" in general has some issues:

• It's stupid: Better to keep the formula and the parameters separate. Is that functionality already existing?
• The evaluator can not parse these formulae, if they contain numbers like "7.981487578945234e-05", i.e. it doesn't seem to comprehend the e-05 notation for floats
• Bin edges like "-5.191" become "-5.190999984741211" in the JSON. Is there an easy way to avoid that?
• Is there a CMSSW-version with Python3 as default? Now I source a ROOT 6.22-version after cmsenv to be able to access CMSSW-classes from Python3, but it gives numerous errors complaining about the version mismatch of libraries.

Thanks for your thoughts/suggestions :-)

[kirschen]

Also, in the first shot, the "clamping" for JEC is not implemented, yet, cf. also https://github.com/nsmith-/correctionlib/issues/5

[nsmith-]

Hi, thanks for the first pass at this. Regarding your points:

• I agree formulas and parameters are better kept separate. I wrote the evaluator with this in mind but it is not in v1 schema. It should be no problem to include it in v2 schema, and I wrote up my thoughts on it in Formulas should take parameters as well as input variables #22. I welcome your feedback on that, in particular with regards to formularef.
• The fact that the formula grammar doesn't support scientific notation is an oversight. Tracked in TFormula_v1 grammar doesn't support scientific notation #23
• This appears as though the value was declared as a single-precision float and then converted to double-precision before converting to string representation, e.g.

root [0] float f { -5.191 };
root [1] double d { -5.191 };
root [2] double fd { -5.190999984741211 };
root [3] f - d
(double) 1.5258789e-08
root [4] f - fd
(double) 0.0000000

In this library (and in json) every float is double-precision so I believe if you initially convert -5.191 from string to double, it should write out in the simpler form.

• CMSSW 11_2 and up can use ROOT from python3, but I don't think any release has python3 as a default yet.

[kirschen]

For the first problem, it looks like the clamping only has an effect at the edges of the pt range, correct? I assume you are using nested Binning nodes for the eta and pt dimensions (JEC are not rectangular binning right?) So one solution is, for each eta node, at the lowest and highest pt bin, split that bin into two, one using the regular formula, and the other using the value the formula would have taken, e.g. for the first line of the JEC, you have a pt bin [0.001, 8] with scalar value precomputed (whatever max(0.0001,[0]+((8-[1])*([2]+((8-[1])*([3]+((8-[1])*[4])))))) is) and another bin [8, 8.9532766] with the normal formula.

Yes, that would work reasonably well in this particular case, but I am not sure how well it generalizes to multi-dimensional formulae with headers like e.g. (L1-corrections)

{1 JetEta 3 Rho JetPt JetA max(0.0001,1-(z/y)*([1]*(x-[0])*(1.+[2]*log(y/15.)))) Correction L1FastJet}

One would have to introduce another three nested binnings for rho/JetPt/JetA and manipulate the formula each time at the boundaries (e.g. rename y->x, z->y for the x/Rho corner cases etc), right?

I was wondering whether instead one could add an optional ranges-parameter to the formula node that would change the behavior around https://github.com/nsmith-/correctionlib/blob/master/src/correction.cc#L279-L284 if needed.

[kirschen]

For the second issue, I guess it is up to you but personally I would sacrifice aesthetics for consistency. You should recover (or at least dramatically improve) agreement with the CMSSW implementation by casting the txtfile inputs with something like float(numpy.single("1.234")). Though if a factor of 1e-6 is OK for everyone (certainly is for me) and well-advertised (nobody likes surprises) then I can see why it would be nice to keep the string representations more compact.

I just tried reverting back to having the CMSSW C-float directly used as python float for the parameters, but it didn't really improve the consistency (maybe truncation when writing the JSON?)

['JetEta', 'JetPt']
[0.0, 9.0]
JSON: 1.1630146229835252; CMSSW: 1.1630146503448486; relative difference [%] (100*(J-C)/C): -2.352620701816858e-06

Personally, I am fine with 1e-6 agreement, but one would have to see whether it is max. 1e-6 or worse in some corners of phase space, I guess.

[nsmith-]

Another approach might be to embed such bounds in the formula, e.g. max(1, y) for JetPt or similar. I agree it would not be hard to add such bounds to the formula node, but for a technical reason (it increases the size of the Content std::variant which is already quite large) I'm trying to avoid that route. Another approach that might be useful in other contexts (e.g. #3) is to introduce a transform node that could rewrite the inputs according to some simple enumeration of rules like abs, clamp(min, max), etc.

[kirschen]

The node transformation approach sounds quite fitting to me, yepp.

For now, I have taken the "brute-force approach" of tampering with the formula expression and doing var--> max(min({var},[{upbound}]),[{lowbound}]) replacements, adding the up/low bounds to the parameter array.

That way, e.g.

max(0.0001,1-(z/y)*[4]*([1]*(x/[3]-[0])*(1.+[2]*log(y/15.))))
turns to
max(0.0001,1-(max(min(z,[5]),[4])/max(min(y,[3]),[2]))*[10]*([7]*(max(min(x,[1]),[0])/[9]-[6])*(1.+[8]*log(max(min(y,[3]),[2])/15.))))

Not very nice, but it gets the job done for now.

[nsmith-]

I'll have to check at some point, but its quite possible turning var into max(min(var, up), lo) leads to a more performant evaluator than rewriting the input because it saves a heap allocation to make the new vector. I did some profiling and found that the line at https://github.com/nsmith-/correctionlib/blob/c3710d48bb0dca6334bf03c41b79130f09c925fb/src/correction.cc#L281 eats up a lot of time and I made a patch to avoid that in #41. Nevertheless, I think this use case in addition to the abs(eta) makes #3 something worth having for schema v2. At least you only need to use the transform for overflow bins.