Resolve: add IRT output data structures

close #21

• add component eic::CherenkovPdgHypothesis
• add datatype eic::CherenkovParticleID

added priority label

assigned to @c.dilks

mentioned in merge request !61 (closed)

marked this merge request as ready

changed title from Draft: Irt init to Resolve: add IRT output data structures

changed the description

mentioned in merge request juggler!350 (closed)

mentioned in merge request !64 (closed)

Can we be more specific about the weight? What are the data model assumptions for this? Elsewhere we use weight as a quantity that can be used to weight events or tracks in histogram, with 1 as null effect value. Does this weight allow the same? Will the weights for all pdg assumptions add to 1?

The truth is that this is a placeholder. The exact meaning will refine as the algorithm gets more mature. As I see it now, the final evaluation can be done in a way where this weight has a meaning of chi^2 CCDF (so a value from 0 to 1, but then one more parameter like 'ndf' would not hurt), with the additive ingredients coming from (1) tracking cov.matrix, (2) either Chebyshev- of gaussian-based piece describing Cherenkov photon angle consistency with a mass hypothesis, (3) lambda chi^2 consistency with the expected number of detected photons per radiator, (4) consistency with the expected number of noise photons, (5) potentially a signature chi^2 along phi angle of the detected photons as well as a lambda chi^2 describing uniformity of the phi distribution. These are all relatively straightforward pieces to implement, but at present the weight just mimics the number of photons associated with a mass hypothesis. The numbers should not add to 1, this would be pointless. It would be a task of analysis guys to make sense of a particular combination of these weights (like "electron vs pion" or "kaon vs pion or proton"). One simply can not consider all the combinations of interest in a data model.

I guess your comments bring me back to the formulation I was going to use at first: what algebra does the weight follow? Is it multiplicative or additive? The distinction between (kaon) vs. (pion or proton) must imply something about how to treat the 'weights' of kaons, pions and protons. I fear that if we call it 'weight' it will invoke certain ideas about how it should be treated (in particular that the sum of the weights of pions and protons would correspond to the weight when you don't care which one it is). Other words for this quantity could be 'confidence', 'likelihood', 'odds', and each of these have very different but well-defined connotations. 'weight' in a classification problem would bring to mind logistic regression outcomes, not number of photons. In that sense, I don't like 'weight' for something that will often be > 1, as it sounds like from your description.

Anyway, maybe this is overthinking it. I'll let @sly2j make final comments.

Well, compared to the coding styles, this discussion I actually appreciate. As coded now (call it algorithm version v1), the weight is an artificial additive quantity. Given a known {x,p} parameterization from ideal tracking, each of the photon candidates after IRT pass has a PDF pho(theta) describing a probability of having been emitted under a certain polar angle 'theta'. For the pfRICH (no focusing optics) this would be a more or less flat distribution in a range [theta_min,theta_max] determined by the aerogel thickness. A given mass hypothesis either belongs to this range or not for each photon individually. This gives a total number of photons (npe), which a mass hypothesis can accommodate. The weight variable is only different in a sense it accounts for the actual PDF value at the nominal theta(mass) for a given {x,p} parameterization. A user can request an additional smearing when calculating these overlaps, either following a uniform +/-range or a gaussian distribution with some predefined sigma. The latter ones are supposed to mimic the effects of track uncertainties and refractive index variation. Since a PDF in basically a collection of uniform distributions, under the hood the overlaps are reduced to erf() calculation, in a worst case. This works reasonably well at the saturation angles, must be resilient to noise, does not care in which of the radiators photons were produced, and must be rather robust for multi-track configurations. It should also work out of the box in case of the magnetic field, since PDFs are sampled along the trajectory, no matter was it straight or bent. And it can be extended in a very straightforward way to account for an effective attenuation lengths (next thing on my agenda, since it will be needed for the v2 algorithm as well). Basically there is a pool of photons, and for every track each mass hypothesis tries to grab as many photon and photon tail distribution overlap integrals from this pool. The one which grabs more, wins. If none grabbed too many (Poisson), the estimate is most likely less reliable for the track to make it into a clean sample. The problem is that (1) it is not the most powerful estimate, (2) it can hardly work at threshold where photon count also matters, (3) it ignores track-per-track covariance matrices, and more. Comments are welcome. I agree this variable is not really a weight. One more half a page piece follows.

The v2, which I have in mind, will look like the following (comments are even more welcome, since it is not coded yet). Tracks will still be treated independently. Same as in v1 (where this generalization is presently optional), the common pool of photons is shared between the tracks, and any MC information about hit origin is ignored. For each track with a set of {x,p,cov} parameterizations along the trajectory, the group of {e,pi,K,p} mass hypotheses defines a subset of photons, which are loosely consistent with at least one of them. This defines N photons to be distributed between three categories: aerogel, gas, noise. Noise (if used) has to be injected into the data stream, and it's average rate is provided to the algorithm. Part of it will come from the Rayleigh scattering. Susceptibility to noise will be proportional to the emission angle for a gven mass hypothesis I believe. Refractive index distribution, which is of course known from the MC output, is provided to the plugin in a way similar to QE, and is applied to the individual photon PDFs, acting as effective smearing. I should mention that in principle one can probably also add this effect as one extra (diagonal) dimension to the tracking covariance matrix. Attenuation is applied to the individual photon PDFs at this stage as well, taking into account azimuthal emission angle and consequently varying path length in each radiator. From the technical point of view, all these effects will be added as sampling with appropriate binning, so math at the end does not change. One can optimize things later, but as long as CPU needs are bearable compared to calorimetry simulations, we are probably good. Each mass hypothesis tries to twist the track within its cov.matrix range, in a way to minimize a chi^2-like value, which is defined by a sum of the following ingredients (here for simplicity assume gaussian shape of each individual photon PDF, which is roughly true for dRICH, then the task has a correct mathematical formulation; in the limit of large photon count it may still be clean if one takes each PDF average and RMS as gaussian distribution equivalents): (1) track-based chi^2 term <dv|M|dv>, basically a track displacement dv in a 2D {theta_x,theta_y} track parameter space, which is optimal for this mass hypothesis, with a metric M defined by the covariance matrix reduced to 2D, (2) usual 1D chi^2 displacement terms dtheta^2/sigma^2 in theta Cherenkov space per every photon which this mass hypothesis accepts, (3) a Poissonic in nature "multinomial lambda chi^2" term in a sense of Baker and Cousins for a mismatch between a set of {N_aerogel, N_gas, N_noise} numbers, which minimize the respective overall chi^2 for this mass hypothesis and the expected average numbers, (4) perhaps also a "signature chi^2" along the emitted photon azimuthal angles for the accepted photons, to suppress situations which regular chi^2 does not capture, when say all ten accepted photons have smaller Cherenkov theta than the nominal one, would be typical for a wrong hypothesis, (5) presumably less important extra terms like one other "lambda chi^2" describing expected uniformity of Cherenkov azimuthal angle distribution. All in all one gets a canonic chi^2 distribution for a given number of degrees of freedom. Well, would be true in a large photon count limit. Not sure that ten is large anough, but we will see. The value which I called "weight" would then be a ROOT::Math::chisquared_cdf_(c?) (double x, double r), I always forget which one is which. Do we then prefer a "ccdf" name for it? This is a number, which should be uniformly distributed between 0 and 1 for the true mass hypothesis. For the false hypotheses, since they can not get rid of "wrong photons" without incurring a large penalty driven by the expected (small) noise level, the values would be sharply peaked around 0. It is a matter of basic histogramming to build the actual plots, see by how much they deviate from a flat distribution for the true hypothesis, and decide at which value like 0.01 to apply a cut to actually suppress the wrong ones. I believe this would suffice for our present purposes. Given what is written already, this is a relatively straightforward coding exercise. But I certainly need Sylvester's track parameterizations to play with.

close in favor of rebased !70 (closed)

closed

mentioned in merge request !70 (closed)