DAS::Unfolding::ExtractHistogram Namespace Reference

Classes
struct	Dist
struct	Matrix

Functions
template<typename T >
void	nodeLoop (TDirectory *dOut, TUnfoldBinning *binning, T t, ostream &cout, TString name, int depth)
template<typename T >
void	axisLoop (TDirectory *dOut, TUnfoldBinning *node, T t, ostream &cout, vector< int > ibins, int shift, TString name="", TString title="")

Detailed Description

Two functors and two functions to extract the distributions and matrices with physical axes from bin IDs.

The two functions allow to loop over

• the nodes of a global binning
• and the axes, regardless of the dimensionality.

The two functors are similar:

• the first one just handles 1D distribution, and is a bit trivial, but necessary to allow the second one;
• the second one can keep in memory the result of the loops on the first axis and triggers a secound round of loops for the second axis.

Function Documentation

axisLoop()

void DAS::Unfolding::ExtractHistogram::axisLoop	(	TDirectory *	dOut,
		TUnfoldBinning *	node,
		T	t,
		ostream &	cout,
		vector< int >	ibins,
		int	shift,
		TString	name = "",
		TString	title = ""
	)

Loop on the axes of the node of a given binning.

< either Dist or Matrix

Parameters

dOut	output directory
node	child node
t	expecting either a `Dist` or a `Matrix` object
cout	standard output
ibins	global bin IDs of present iteration
shift	shift to continue the bin index despite the change of node/binning

{
    int iaxis = count_if(next(begin(ibins)), end(ibins),
                         [](int ibin) { return ibin == 0; });
    int& ibin = ibins[iaxis];
    assert(ibin == 0);
    if (iaxis > 0) {
        TString axisname = node->GetDistributionAxisLabel(iaxis);
        cout << "\t\taxis `" << axisname << "`" << endl;
        name += "_" + axisname + "bin";
        const TVectorD * edges = node->GetDistributionBinning(iaxis);
        for (ibin = 1; ibin < edges->GetNrows(); ++ibin) {
            TString nameExt = name + (ibin + shift),
                    titleExt = title + Form("[%f, %f] ", (*edges)[ibin-1], (*edges)[ibin]);
            axisLoop(dOut, node, t, cout, ibins, shift, nameExt, titleExt);
        }
    }
    else {
        cout << "\t\t\tpreparing `" << name << "`" << endl;
        size_t naxes = ibins.size();
        vector<double> x(naxes);
        for (size_t jaxis = 1; jaxis < naxes; ++jaxis) {
            const TVectorD * edges = node->GetDistributionBinning(jaxis);
            int ibin = ibins[jaxis];
            assert(ibin > 0);
            x[jaxis] = (*edges)[ibin-1]; // vector component (hist bin) index starts from 0 (1)
        }
 
        if constexpr (is_same_v<T,Dist>)
            t(dOut, name, node, x, title);
        else if constexpr (is_same_v<T,Matrix>) {
            if (t.node1 == nullptr) {
                t.name1 = name;
                t.node1 = node;
                t.ibins1 = ibins;
                t.x1 = x;
                t.title1 = title;
                nodeLoop<T>(dOut, t.binning2, t, cout);
            }
            else
                t(dOut, name, node, x, title);
        }
        else {
            cerr << red << "Unexpected type. Skipping.\n" << def;
            return;
        }
    }
}

nodeLoop()

void nodeLoop	(	TDirectory *	dOut,
		TUnfoldBinning *	binning,
		T	t,
		ostream &	cout,
		TString	name,
		int	depth
	)

Loop on the nodes of a given binning.

< either Dist or Matrix

Todo:

d > 2?

Todo:

d > 2?

Parameters

dOut	output directory
binning	global binning or parent node
t	expecting either a `Dist` or a `Matrix` object
cout	standard output
name	name of output directory
depth	trick to know the number of recursive calls

{
    int shift = 0;
    for (auto node = binning->GetChildNode(); node != nullptr;
              node = node->GetNextNode()) {
        TString dirname = name + node->GetName();
        if (node->GetDistributionNumberOfBins() == 0) {
            nodeLoop(dOut, node, t, cout, dirname, depth+1);
            continue;
        }
        cout << "\tnode `" << dirname << "`" << endl;
        auto ddOut = dOut->mkdir(dirname, node->GetTitle(), true /* = `-p` */);
        ddOut->cd();
 
        int dim = node->GetDistributionDimension();
        TString axisname = node->GetDistributionAxisLabel(0);
        axisLoop<T>(ddOut, node, t, cout, vector<int>(dim,0), shift, axisname);
        if (depth > 1) { // NOTE: this logic is tuned for dijet mass, with different resolution in central and forward regions
            const TVectorD * edges = node->GetDistributionBinning(1 );
            shift += edges->GetNrows()-1;
        }
    } // end of loop on nodes
}