Efficiency Object: linear representation of N-dimensional matrix

The full code is here: Efficiency.h.
It compiles and works in CMSSW, but the only dependencies are from boost. The class can be easily exported and used in a standalone C++ executable as needed.

What follows is an explaination of the crucial parts.

Compute the index of a linear representation of a matrix with a recursive function.

Given a N dimensional matrix where each dimension has its own number of elements it is useful to have a linear representation of it so that the number of dimensions can be dynamically changed.
One way of doing this is by using recursion.

The total number of elements in the matrix is simply the product of the sizes of all dimensions: S = Prod(s_i) with i=1...N

If we denote with “i” the index running from 1 to N and given “x_i” as the i-th index in the matrix, the index S of the element given by the vector of x_i is:

I = x_i + s_i*(x_(i-1) + size(x_(i-1))*(...))

This index can be computed with a recursive function:

unsigned int linearIndex(unsigned int i, unsigned int previous) {
  if( i< N ) return( linearIndex(i+1, x_i + size(x_i)*previous) );
  return (x_i + size(x_i)*previous);
}

by calling it with the first element:

linearIndex(1, 0)


it will return the I of the linearized representation.

This defines the order of the variables. The i=1 is the most external variable, this means that it must be the first also in any loop.

Initialization

Fill the array (or vector) from 0 to S.

N-k dimensional projections

The N dimensions are good to store the data, but sometimes there can be a need to visualize them or to reduce some of the dimensions.
Selecting the k dimensions that must be kept, a new array of size Prod(s_i) is created, where the product is now only on the k dimensions, or if you want, the other dimensions have sizes set to 0.

The values of the (N-k) variables must be collapsed into the remaining k.

The important thing is to compute indexes from the index “L” in the linear representation.

Starting from the first variable its index can be computed as:

int(L/size(N-1))


where size(N-1) is the product of the sizes for all the variables minus the first one.
Then the remainder of this division can be passed for the computation of the next variable:

return(L%size(N-1))


and the indexes can be computed recursively stopping at the biggest x among the k selected.
Note that when the value returned is 0 the computation can be stopped immediately and all the rest of the indexes set to 0.

Filling

When a set of values are passed the corresponding bins need to be found and the counters increased.
Finding the bins is just a matter of computing:

index = (value-min)/bins


and

if( value < min ) index = 0;

Note that this is not automatic when using unsigned int, see here: http://stackoverflow.com/questions/50605/signed-to-unsigned-conversion-in-c-is-it-always-safe.

if( index > bins ) index = bins-1;

Once the array of indexes is built it can be fed to the linearIndex function to get the index in the linear representation and use it to fill the value.