- //
- // Build binning
- //
- // 1. calculate number of nodes
- // 2. calculate first terminal row
- // 3. initialize index array
- // 4. non recursive building of the binary tree
- //
- //1.
- fNnodes = fNpoints/fBucketSize-1;
- if (fNpoints%fBucketSize) fNnodes++;
-
- //2.
- fRowT0=0;
- for (;(fNnodes+1)>(1<<fRowT0);fRowT0++);
- fRowT0-=1;
- // 2 = 2**0 + 1
- // 3 = 2**1 + 1
- // 4 = 2**1 + 2
- // 5 = 2**2 + 1
- // 6 = 2**2 + 2
- // 7 = 2**2 + 3
- // 8 = 2**2 + 4
-
- //3.
- // allocate space for boundaries
- fRange = new Value[2*fNDim];
- fIndPoints= new Index[fNpoints];
- for (Index i=0; i<fNpoints; i++) fIndPoints[i] = i;
- fAxis = new UChar_t[fNnodes];
- fValue = new Value[fNnodes];
- //
- fCrossNode = (1<<(fRowT0+1))-1;
- if (fCrossNode<fNnodes) fCrossNode = 2*fCrossNode+1;
- //
- // fOffset = (((fNnodes+1)-(1<<fRowT0)))*2;
- Int_t over = (fNnodes+1)-(1<<fRowT0);
- Int_t filled = ((1<<fRowT0)-over)*fBucketSize;
- fOffset = fNpoints-filled;
- //
-// printf("Row0 %d\n", fRowT0);
-// printf("CrossNode %d\n", fCrossNode);
-// printf("Offset %d\n", fOffset);
- //
- //
- //4.
- // stack for non recursive build - size 128 bytes enough
- Int_t rowStack[128];
- Int_t nodeStack[128];
- Int_t npointStack[128];
- Int_t posStack[128];
- Int_t currentIndex = 0;
- Int_t iter =0;
- rowStack[0] = 0;
- nodeStack[0] = 0;
- npointStack[0] = fNpoints;
- posStack[0] = 0;
- //
- Int_t nbucketsall =0;
- while (currentIndex>=0){
- iter++;
- //
- Int_t npoints = npointStack[currentIndex];
- if (npoints<=fBucketSize) {
- //printf("terminal node : index %d iter %d\n", currentIndex, iter);
- currentIndex--;
- nbucketsall++;
- continue; // terminal node
- }
- Int_t crow = rowStack[currentIndex];
- Int_t cpos = posStack[currentIndex];
- Int_t cnode = nodeStack[currentIndex];
- //printf("currentIndex %d npoints %d node %d\n", currentIndex, npoints, cnode);
- //
- // divide points
- Int_t nbuckets0 = npoints/fBucketSize; //current number of buckets
- if (npoints%fBucketSize) nbuckets0++; //
- Int_t restRows = fRowT0-rowStack[currentIndex]; // rest of fully occupied node row
- if (restRows<0) restRows =0;
- for (;nbuckets0>(2<<restRows); restRows++);
- Int_t nfull = 1<<restRows;
- Int_t nrest = nbuckets0-nfull;
- Int_t nleft =0, nright =0;
- //
- if (nrest>(nfull/2)){
- nleft = nfull*fBucketSize;
- nright = npoints-nleft;
- }else{
- nright = nfull*fBucketSize/2;
- nleft = npoints-nright;
- }
-
- //
- //find the axis with biggest spread
- Value maxspread=0;
- Value tempspread, min, max;
- Index axspread=0;
- Value *array;
- for (Int_t idim=0; idim<fNDim; idim++){
- array = fData[idim];
- Spread(npoints, array, fIndPoints+cpos, min, max);
- tempspread = max - min;
- if (maxspread < tempspread) {
- maxspread=tempspread;
- axspread = idim;
- }
- if(cnode) continue;
- //printf("set %d %6.3f %6.3f\n", idim, min, max);
- fRange[2*idim] = min; fRange[2*idim+1] = max;
- }
- array = fData[axspread];
- KOrdStat(npoints, array, nleft, fIndPoints+cpos);
- fAxis[cnode] = axspread;
- fValue[cnode] = array[fIndPoints[cpos+nleft]];
- //printf("Set node %d : ax %d val %f\n", cnode, node->fAxis, node->fValue);
- //
- //
- npointStack[currentIndex] = nleft;
- rowStack[currentIndex] = crow+1;
- posStack[currentIndex] = cpos;
- nodeStack[currentIndex] = cnode*2+1;
- currentIndex++;
- npointStack[currentIndex] = nright;
- rowStack[currentIndex] = crow+1;
- posStack[currentIndex] = cpos+nleft;
- nodeStack[currentIndex] = (cnode*2)+2;
- //
- if (0){
- // consistency check
- Info("Build()", Form("points %d left %d right %d", npoints, nleft, nright));
- if (nleft<nright) Warning("Build", "Problem Left-Right");
- if (nleft<0 || nright<0) Warning("Build()", "Problem Negative number");
- }
- }
-
- //printf("NBuckets\t%d\n", nbucketsall);
- //fData = 0;
+ //
+ // Build binning
+ //
+ // 1. calculate number of nodes
+ // 2. calculate first terminal row
+ // 3. initialize index array
+ // 4. non recursive building of the binary tree
+
+ //
+ // Teh tree is recursivaly divided. The division point is choosen in such a way, that the left side
+ // alway has 2^k points, while at the same time trying
+ //
+ //1.
+ fNnodes = fNpoints/fBucketSize-1;
+ if (fNpoints%fBucketSize) fNnodes++;
+
+ //2.
+ fRowT0=0;
+ for (;(fNnodes+1)>(1<<fRowT0);fRowT0++);
+ fRowT0-=1;
+ // 2 = 2**0 + 1
+ // 3 = 2**1 + 1
+ // 4 = 2**1 + 2
+ // 5 = 2**2 + 1
+ // 6 = 2**2 + 2
+ // 7 = 2**2 + 3
+ // 8 = 2**2 + 4
+
+ //3.
+ // allocate space for boundaries
+ fRange = new Value[2*fNDim];
+ fIndPoints= new Index[fNpoints];
+ for (Index i=0; i<fNpoints; i++) fIndPoints[i] = i;
+ fAxis = new UChar_t[fNnodes];
+ fValue = new Value[fNnodes];
+ //
+ fCrossNode = (1<<(fRowT0+1))-1;
+ if (fCrossNode<fNnodes) fCrossNode = 2*fCrossNode+1;
+ //
+ // fOffset = (((fNnodes+1)-(1<<fRowT0)))*2;
+ Int_t over = (fNnodes+1)-(1<<fRowT0);
+ Int_t filled = ((1<<fRowT0)-over)*fBucketSize;
+ fOffset = fNpoints-filled;
+ //
+ // printf("Row0 %d\n", fRowT0);
+ // printf("CrossNode %d\n", fCrossNode);
+ // printf("Offset %d\n", fOffset);
+ //
+ //
+ //4.
+ // stack for non recursive build - size 128 bytes enough
+ Int_t rowStack[128];
+ Int_t nodeStack[128];
+ Int_t npointStack[128];
+ Int_t posStack[128];
+ Int_t currentIndex = 0;
+ Int_t iter =0;
+ rowStack[0] = 0;
+ nodeStack[0] = 0;
+ npointStack[0] = fNpoints;
+ posStack[0] = 0;
+ //
+ Int_t nbucketsall =0;
+ while (currentIndex>=0){
+ iter++;
+ //
+ Int_t npoints = npointStack[currentIndex];
+ if (npoints<=fBucketSize) {
+ //printf("terminal node : index %d iter %d\n", currentIndex, iter);
+ currentIndex--;
+ nbucketsall++;
+ continue; // terminal node
+ }
+ Int_t crow = rowStack[currentIndex];
+ Int_t cpos = posStack[currentIndex];
+ Int_t cnode = nodeStack[currentIndex];
+ //printf("currentIndex %d npoints %d node %d\n", currentIndex, npoints, cnode);
+ //
+ // divide points
+ Int_t nbuckets0 = npoints/fBucketSize; //current number of buckets
+ if (npoints%fBucketSize) nbuckets0++; //
+ Int_t restRows = fRowT0-rowStack[currentIndex]; // rest of fully occupied node row
+ if (restRows<0) restRows =0;
+ for (;nbuckets0>(2<<restRows); restRows++);
+ Int_t nfull = 1<<restRows;
+ Int_t nrest = nbuckets0-nfull;
+ Int_t nleft =0, nright =0;
+ //
+ if (nrest>(nfull/2)){
+ nleft = nfull*fBucketSize;
+ nright = npoints-nleft;
+ }else{
+ nright = nfull*fBucketSize/2;
+ nleft = npoints-nright;
+ }
+
+ //
+ //find the axis with biggest spread
+ Value maxspread=0;
+ Value tempspread, min, max;
+ Index axspread=0;
+ Value *array;
+ for (Int_t idim=0; idim<fNDim; idim++){
+ array = fData[idim];
+ Spread(npoints, array, fIndPoints+cpos, min, max);
+ tempspread = max - min;
+ if (maxspread < tempspread) {
+ maxspread=tempspread;
+ axspread = idim;
+ }
+ if(cnode) continue;
+ //printf("set %d %6.3f %6.3f\n", idim, min, max);
+ fRange[2*idim] = min; fRange[2*idim+1] = max;
+ }
+ array = fData[axspread];
+ KOrdStat(npoints, array, nleft, fIndPoints+cpos);
+ fAxis[cnode] = axspread;
+ fValue[cnode] = array[fIndPoints[cpos+nleft]];
+ //printf("Set node %d : ax %d val %f\n", cnode, node->fAxis, node->fValue);
+ //
+ //
+ npointStack[currentIndex] = nleft;
+ rowStack[currentIndex] = crow+1;
+ posStack[currentIndex] = cpos;
+ nodeStack[currentIndex] = cnode*2+1;
+ currentIndex++;
+ npointStack[currentIndex] = nright;
+ rowStack[currentIndex] = crow+1;
+ posStack[currentIndex] = cpos+nleft;
+ nodeStack[currentIndex] = (cnode*2)+2;
+ //
+ if (0){
+ // consistency check
+ Info("Build()", Form("points %d left %d right %d", npoints, nleft, nright));
+ if (nleft<nright) Warning("Build", "Problem Left-Right");
+ if (nleft<0 || nright<0) Warning("Build()", "Problem Negative number");
+ }
+ }