Corral/src/ZvtxAverage.h

0001 #ifndef ZVTXAVERAGE_H
0002 #define ZVTXAVERAGE_H
0003 //
0004 // Shared by CalcRm::Calculate and corral::Finalize (README_Finalize.md), so that the per-run and the
0005 // full-statistics Zvtx averages use exactly the same rules. Moved verbatim from CalcRm.cxx.
0006 //
0007 #include <cmath>
0008 #include "TH1.h"
0009
0010 //------------------------------------------------------------
0011 // Zvtx average of a CF-type histogram, bin by bin, over the slices where that bin is VALID
0012 // (README_Crossing sec 8): slice ok (computed), weight>0, and mask>0 (mask = rho2(M) of the slice,
0013 // or the crossing-corrected validity map); mask==0 -> valid iff the slice's bin error>0 (used for
0014 // the 1D projections, which are only filled from valid bins). Weights are renormalized per bin:
0015 // v = sum(w v)/sum(w), e = sqrt(sum(w^2 e^2))/sum(w). A bin valid in no slice gets `sentinel`, error 0.
0016 // Per-run errors only; Finalize replaces them with subgroup errors.
0017 inline void ZvtxAverageValid(TH1* avg, TH1** sl, TH1** mask, const double* w, const bool* ok, int n, double sentinel){
0018     avg->Reset();
0019     double nent = 0;
0020     for (int z=0;z<n;z++) if (ok[z] && w[z]>0.) nent += w[z]*sl[z]->GetEntries();
0021     for (int ib=0;ib<avg->GetNcells();ib++){
0022         double sw=0, swv=0, sw2e2=0;
0023         for (int z=0;z<n;z++){
0024             if (!ok[z] || w[z]<=0.) continue;
0025             bool valid  = mask ? (mask[z]->GetBinContent(ib)>0.) : (sl[z]->GetBinError(ib)>0.);
0026             if (!valid) continue;
0027             double e    = sl[z]->GetBinError(ib);
0028             sw      += w[z];
0029             swv     += w[z]*sl[z]->GetBinContent(ib);
0030             sw2e2   += w[z]*w[z]*e*e;
0031         }
0032         if (sw>0.){ avg->SetBinContent(ib,swv/sw);   avg->SetBinError(ib,sqrt(sw2e2)/sw); }
0033         else      { avg->SetBinContent(ib,sentinel); avg->SetBinError(ib,0.); }
0034     }
0035     avg->SetEntries(nent);
0036 }
0037
0038 #endif