]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG3/vertexingHF/AliMultiDimVector.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Update (Francesco, Gian Michele)
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliMultiDimVector.cxx
CommitLineData
ac7636a0 1/**************************************************************************
2 * Copyright(c) 2007-2009, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
7b45817b 16/* $Id$ */
ac7636a0 17
18///////////////////////////////////////////////////////////////////
19// //
20// Implementation of the class to store the number of signal //
21// and background events in bins of the cut values //
22// Origin: Elena Bruna (bruna@to.infn.it) //
23// Updated: Sergey Senyukov (senyukov@to.infn.it) //
3aa6ce53 24// Francesco Prino (prino@to.infn.it) //
25// Last Updated: Giacomo Ortona (ortona@to.infn.it) //
ac7636a0 26// //
27///////////////////////////////////////////////////////////////////
28
0adeb69c 29#include <fstream>
30#include <Riostream.h>
ac7636a0 31#include "TH2.h"
32#include "AliMultiDimVector.h"
33#include "AliLog.h"
ac7636a0 34#include "TString.h"
35
36ClassImp(AliMultiDimVector)
37//___________________________________________________________________________
38AliMultiDimVector::AliMultiDimVector():TNamed("AliMultiDimVector","default"),
39fNVariables(0),
40fNPtBins(0),
41fVett(0),
42fNTotCells(0),
43fIsIntegrated(0)
44{
45 // default constructor
46}
47//___________________________________________________________________________
2350a1d9 48AliMultiDimVector::AliMultiDimVector(const char *name,const char *title, const Int_t nptbins, const Float_t* ptlimits, const Int_t npars, const Int_t *nofcells,const Float_t *loosecuts, const Float_t *tightcuts, const TString *axisTitles):TNamed(name,title),
ac7636a0 49fNVariables(npars),
50fNPtBins(nptbins),
51fVett(0),
52fNTotCells(0),
53fIsIntegrated(0){
54// standard constructor
55 ULong64_t ntot=1;
56 for(Int_t i=0;i<fNVariables;i++){
57 ntot*=nofcells[i];
58 fNCutSteps[i]=nofcells[i];
3aa6ce53 59 if(loosecuts[i] == tightcuts[i]){
60 if(loosecuts[i]!=0){
61 printf("AliMultiDimVector::AliMultiDimVector: WARNING! same tight/loose variable for variable number %d. AliMultiDimVector with run with the following values: loose: %f; tight: %f\n",i,tightcuts[i]-0.1*tightcuts[i],tightcuts[i]);
62 fMinLimits[i]=tightcuts[i]-0.1*tightcuts[i];
63 fMaxLimits[i]=tightcuts[i];
64 }else{
65 fMinLimits[i]=0;
66 fMaxLimits[i]=0.0001;
67 }
68 fGreaterThan[i]=kTRUE;
69 }
70 if(loosecuts[i] < tightcuts[i]){
ac7636a0 71 fMinLimits[i]=loosecuts[i];
72 fMaxLimits[i]=tightcuts[i];
73 fGreaterThan[i]=kTRUE;
74 }else{
75 fMinLimits[i]=tightcuts[i];
76 fMaxLimits[i]=loosecuts[i];
77 fGreaterThan[i]=kFALSE;
78 }
79 fAxisTitles[i]=axisTitles[i].Data();
80 }
bac542b0 81 fNTotCells=ntot*fNPtBins;
ac7636a0 82 fVett.Set(fNTotCells);
bac542b0 83 for(Int_t ipt=0;ipt<fNPtBins+1;ipt++) fPtLimits[ipt]=ptlimits[ipt];
84 for(Int_t ipt=fNPtBins+1;ipt<fgkMaxNPtBins+1;ipt++) fPtLimits[ipt]=999.;
ac7636a0 85 for (ULong64_t j=0;j<fNTotCells;j++) fVett.AddAt(0,j);
86}
87//___________________________________________________________________________
88AliMultiDimVector::AliMultiDimVector(const AliMultiDimVector &mv):TNamed(mv.GetName(),mv.GetTitle()),
89fNVariables(mv.fNVariables),
90fNPtBins(mv.fNPtBins),
91fVett(0),
92fNTotCells(mv.fNTotCells),
93fIsIntegrated(mv.fIsIntegrated)
94{
95 // copy constructor
96 for(Int_t i=0;i<fNVariables;i++){
97 fNCutSteps[i]=mv.GetNCutSteps(i);
98 fMinLimits[i]=mv.GetMinLimit(i);
99 fMaxLimits[i]=mv.GetMaxLimit(i);
100 fGreaterThan[i]=mv.GetGreaterThan(i);
101 fAxisTitles[i]=mv.GetAxisTitle(i);
102 }
103 fVett.Set(fNTotCells);
bac542b0 104
105 for(Int_t ipt=0;ipt<fNPtBins+1;ipt++) fPtLimits[ipt]=mv.GetPtLimit(ipt);
ac7636a0 106 for(ULong64_t i=0;i<fNTotCells;i++) fVett[i]=mv.GetElement(i);
107}
108//___________________________________________________________________________
109void AliMultiDimVector::CopyStructure(const AliMultiDimVector* mv){
110// Sets dimensions and limit from mv
111 fNVariables=mv->GetNVariables();
112 fNPtBins=mv->GetNPtBins();
113 fNTotCells=mv->GetNTotCells();
114 fIsIntegrated=mv->IsIntegrated();
115 for(Int_t i=0;i<fNVariables;i++){
116 fNCutSteps[i]=mv->GetNCutSteps(i);
117 fMinLimits[i]=mv->GetMinLimit(i);
118 fMaxLimits[i]=mv->GetMaxLimit(i);
119 fGreaterThan[i]=mv->GetGreaterThan(i);
120 fAxisTitles[i]=mv->GetAxisTitle(i);
121 }
bac542b0 122 for(Int_t ipt=0;ipt<fNPtBins+1;ipt++) fPtLimits[ipt]=mv->GetPtLimit(ipt);
123 fVett.Set(fNTotCells);
ac7636a0 124}
125//______________________________________________________________________
126Bool_t AliMultiDimVector::GetIndicesFromGlobalAddress(ULong64_t globadd, Int_t *ind, Int_t &ptbin) const {
127// returns matrix element indices and Pt bin from global index
128 if(globadd>=fNTotCells) return kFALSE;
129 ULong64_t r=globadd;
130 Int_t prod=1;
131 Int_t nOfCellsPlusLevel[fgkMaxNVariables+1];
132 for(Int_t k=0;k<fNVariables;k++) nOfCellsPlusLevel[k]=fNCutSteps[k];
133 nOfCellsPlusLevel[fNVariables]=fNPtBins;
134
135 for(Int_t i=0;i<fNVariables+1;i++) prod*=nOfCellsPlusLevel[i];
136 for(Int_t i=0;i<fNVariables+1;i++){
137 prod/=nOfCellsPlusLevel[i];
138 if(i<fNVariables) ind[i]=r/prod;
139 else ptbin=r/prod;
140 r=globadd%prod;
141 }
142 return kTRUE;
143}
144//______________________________________________________________________
145Bool_t AliMultiDimVector::GetCutValuesFromGlobalAddress(ULong64_t globadd, Float_t *cuts, Int_t &ptbin) const {
146 Int_t ind[fgkMaxNVariables];
147 Bool_t retcode=GetIndicesFromGlobalAddress(globadd,ind,ptbin);
148 if(!retcode) return kFALSE;
149 for(Int_t i=0;i<fNVariables;i++) cuts[i]=GetCutValue(i,ind[i]);
150 return kTRUE;
151}
152//______________________________________________________________________
2350a1d9 153ULong64_t AliMultiDimVector::GetGlobalAddressFromIndices(const Int_t *ind, Int_t ptbin) const {
ac7636a0 154 // Returns the global index of the cell in the matrix
155 Int_t prod=1;
156 ULong64_t elem=0;
157 Int_t indexPlusLevel[fgkMaxNVariables+1];
158 Int_t nOfCellsPlusLevel[fgkMaxNVariables+1];
159 for(Int_t i=0;i<fNVariables;i++){
160 indexPlusLevel[i]=ind[i];
161 nOfCellsPlusLevel[i]=fNCutSteps[i];
162 }
163 indexPlusLevel[fNVariables]=ptbin;
164 nOfCellsPlusLevel[fNVariables]=fNPtBins;
165
166 for(Int_t i=0;i<fNVariables+1;i++){
167 prod=indexPlusLevel[i];
168 if(i<fNVariables){
169 for(Int_t j=i+1;j<fNVariables+1;j++){
170 prod*=nOfCellsPlusLevel[j];
171 }
172 }
173 elem+=prod;
174 }
175 return elem;
176}
177//______________________________________________________________________
2350a1d9 178Bool_t AliMultiDimVector::GetIndicesFromValues(const Float_t *values, Int_t *ind) const {
ac7636a0 179 // Fills the array of matrix indices strating from variable values
180 for(Int_t i=0;i<fNVariables;i++){
181 if(fGreaterThan[i]){
182 if(values[i]<GetMinLimit(i)) return kFALSE;
183 ind[i]=(Int_t)((values[i]-fMinLimits[i])/GetCutStep(i));
184 if(ind[i]>=GetNCutSteps(i)) ind[i]=GetNCutSteps(i)-1;
185 }else{
186 if(values[i]>GetMaxLimit(i)) return kFALSE;
187 ind[i]=(Int_t)((fMaxLimits[i]-values[i])/GetCutStep(i));
188 if(ind[i]>=GetNCutSteps(i)) ind[i]=GetNCutSteps(i)-1;
189 }
190 }
191 return kTRUE;
192}
193//______________________________________________________________________
2350a1d9 194ULong64_t AliMultiDimVector::GetGlobalAddressFromValues(const Float_t *values, Int_t ptbin) const {
ac7636a0 195 // Returns the global index of the cell in the matrix
196 Int_t ind[fgkMaxNVariables];
197 Bool_t retcode=GetIndicesFromValues(values,ind);
198 if(retcode) return GetGlobalAddressFromIndices(ind,ptbin);
199 else{
200 AliError("Values out of range");
201 return fNTotCells+999;
202 }
203}
204//_____________________________________________________________________________
205void AliMultiDimVector::MultiplyBy(Float_t factor){
206 // multiply the AliMultiDimVector by a constant factor
207 for(ULong64_t i=0;i<fNTotCells;i++){
2350a1d9 208 if(fVett.At(i)>0.)
ac7636a0 209 fVett.AddAt(fVett.At(i)*factor,i);
210 else fVett.AddAt(-1,i);
211 }
212
213}
214//_____________________________________________________________________________
2350a1d9 215void AliMultiDimVector::Multiply(const AliMultiDimVector* mv,Float_t factor){
ac7636a0 216 // Sets AliMultiDimVector=mv*constant factor
217 for(ULong64_t i=0;i<fNTotCells;i++){
2350a1d9 218 if(mv->GetElement(i)>0.)
ac7636a0 219 fVett.AddAt(mv->GetElement(i)*factor,i);
220 else fVett.AddAt(-1,i);
221 }
222}
223//_____________________________________________________________________________
2350a1d9 224void AliMultiDimVector::Multiply(const AliMultiDimVector* mv1, const AliMultiDimVector* mv2){
ac7636a0 225 // Sets AliMultiDimVector=mv1*mv2
226 for(ULong64_t i=0;i<fNTotCells;i++){
2350a1d9 227 if(mv1->GetElement(i)>0. && mv2->GetElement(i)>0.)
ac7636a0 228 fVett.AddAt(mv1->GetElement(i)*mv2->GetElement(i),i);
229 else fVett.AddAt(-1,i);
230 }
231}
232//_____________________________________________________________________________
2350a1d9 233void AliMultiDimVector::Add(const AliMultiDimVector* mv){
ac7636a0 234 // Sums contents of mv to AliMultiDimVector
235 if (mv->GetNTotCells()!=fNTotCells){
236 AliError("Different dimension of the vectors!!");
237 }else{
238 for(ULong64_t i=0;i<fNTotCells;i++)
2350a1d9 239 if(mv->GetElement(i)>0. && fVett.At(i)>0.)
ac7636a0 240 fVett.AddAt(fVett.At(i)+mv->GetElement(i),i);
241 else fVett.AddAt(-1,i);
242 }
243}
244//_____________________________________________________________________________
2350a1d9 245void AliMultiDimVector::Sum(const AliMultiDimVector* mv1, const AliMultiDimVector* mv2){
ac7636a0 246 // Sets AliMultiDimVector=mv1+mv2
247 if (fNTotCells!=mv1->GetNTotCells()&&mv1->GetNTotCells()!=mv2->GetNTotCells()) {
248 AliError("Different dimension of the vectors!!");
249 }
250 else{
251 for(ULong64_t i=0;i<mv1->GetNTotCells();i++) {
2350a1d9 252 if(mv1->GetElement(i)>0. && mv2->GetElement(i)>0.)
ac7636a0 253 fVett.AddAt(mv1->GetElement(i)+mv2->GetElement(i),i);
254 else fVett.AddAt(-1,i);
255 }
256 }
257}
258//_____________________________________________________________________________
2350a1d9 259void AliMultiDimVector::LinearComb(const AliMultiDimVector* mv1, Float_t norm1, const AliMultiDimVector* mv2, Float_t norm2){
ac7636a0 260 // Sets AliMultiDimVector=n1*mv1+n2*mv2
261 if (fNTotCells!=mv1->GetNTotCells()&&mv1->GetNTotCells()!=mv2->GetNTotCells()) {
262 AliError("Different dimension of the vectors!!");
263 }
264 else{
265 for(ULong64_t i=0;i<mv1->GetNTotCells();i++) {
2350a1d9 266 if(mv1->GetElement(i)>0. && mv2->GetElement(i)>0.)
ac7636a0 267 fVett.AddAt(norm1*mv1->GetElement(i)+norm2*mv2->GetElement(i),i);
268 else fVett.AddAt(-1,i);
269 }
270 }
271}
272//_____________________________________________________________________________
2350a1d9 273void AliMultiDimVector::DivideBy(const AliMultiDimVector* mv){
ac7636a0 274 // Divide AliMulivector by mv
275 if (mv->GetNTotCells()!=fNTotCells) {
276 AliError("Different dimension of the vectors!!");
277 }
278 else{
279 for(ULong64_t i=0;i<fNTotCells;i++)
2350a1d9 280 if(mv->GetElement(i)!=0 &&mv->GetElement(i)>0. && fVett.At(i)>0.)
ac7636a0 281 fVett.AddAt(fVett.At(i)/mv->GetElement(i),i);
282 else fVett.AddAt(-1,i);
283 }
284
285}
286//_____________________________________________________________________________
2350a1d9 287void AliMultiDimVector::Divide(const AliMultiDimVector* mv1, const AliMultiDimVector* mv2){
ac7636a0 288 // Sets AliMultiDimVector=mv1/mv2
289 if (fNTotCells!=mv1->GetNTotCells()&&mv1->GetNTotCells()!=mv2->GetNTotCells()) {
290 AliError("Different dimension of the vectors!!");
291 }
292 else{
293 for(ULong64_t i=0;i<mv1->GetNTotCells();i++)
2350a1d9 294 if(mv2->GetElement(i)!=0&& mv2->GetElement(i)>0.&& mv1->GetElement(i)>0.)
ac7636a0 295 {
296 fVett.AddAt(mv1->GetElement(i)/mv2->GetElement(i),i);
297 }
298 else fVett.AddAt(-1,i);
299 }
300}
301//_____________________________________________________________________________
302void AliMultiDimVector::Sqrt(){
303 // Sqrt of elements of AliMultiDimVector
304 for(ULong64_t i=0;i<fNTotCells;i++) {
305 if(fVett.At(i)>=0) fVett.AddAt(TMath::Sqrt(fVett.At(i)),i);
306 else {
307 fVett.AddAt(-1,i);
308 }
309 }
310}
311//_____________________________________________________________________________
2350a1d9 312void AliMultiDimVector::Sqrt(const AliMultiDimVector* mv){
ac7636a0 313 // Sets AliMultiDimVector=sqrt(mv)
314 for(ULong64_t i=0;i<fNTotCells;i++)
315 if(mv->GetElement(i)>=0) fVett.AddAt(TMath::Sqrt(mv->GetElement(i)),i);
316 else fVett.AddAt(-1,i);
317}
318//_____________________________________________________________________________
319void AliMultiDimVector::FindMaximum(Float_t& maxValue, Int_t *ind , Int_t ptbin){
320 // finds the element with maximum contents
321 const ULong64_t nelem=fNTotCells/fNPtBins;
322 TArrayF vett;
323 vett.Set(nelem);
324 ULong64_t runningAddress;
325 for(ULong64_t i=0;i<nelem;i++){
326 runningAddress=ptbin+i*fNPtBins;
327 vett.AddAt(fVett[runningAddress],i);
328 }
329 maxValue=TMath::MaxElement(nelem,vett.GetArray());
330 ULong64_t maxAddress=TMath::LocMax(nelem,vett.GetArray());
331 ULong64_t maxGlobalAddress=ptbin+maxAddress*fNPtBins;
332 Int_t checkedptbin;
333 GetIndicesFromGlobalAddress(maxGlobalAddress,ind,checkedptbin);
334}
335
0adeb69c 336//_____________________________________________________________________________
337//Int_t* AliMultiDimVector::FindLocalMaximum(Float_t& maxValue, Bool_t *isFree,Int_t* indFixed, Int_t ptbin){
338Int_t* AliMultiDimVector::FindLocalMaximum(Float_t& maxValue, Int_t *numFixed,Int_t* indFixed, Int_t nfixed,Int_t ptbin){
339 //return the elements with maximum content (maxValue) given fixed step for not free variables
340 //numFixed[nfixed] is the indices of the fixed variables in the cuts array [fNVariables]={kTRUE,kTRUE,...,kFALSE,...,kFALSE,...,kTRUE}
341 //indFixed[nfixed]={1,2} //nfixed is the number of false in isFree; indFixed contains the step for the i-th variable
342 //!!take care of deleting the array of index returned!!
343
344 // Int_t nfixed=0,nfree=0;
345 //Int_t indtmp[fNVariables];
346 if(nfixed>fNVariables)cout<<"AliMultiDimVector::FindLocalMaximum:ERROR! too many variables"<<endl;
347 ULong64_t nelem=1;
348 Int_t* indMax=new Int_t[fNVariables];
349 //Get the number of fixed vars
350 /*
351 for (Int_t iv=0;iv<fNVariables;iv++){
352 if(isFree[iv]){
353 nfree++;
354 nelem*=fNCutSteps[iv];
355 indMax[iv]=0;
356 }
357 else {
358 indMax[iv]=indFixed[nfixed];
359 if(indFixed[nfixed]>=GetNCutSteps(iv)){
360 indMax[iv]=0;
361 cout<<"AliMultiDimVector::FindLocalMaximum:ERROR! called fixed ind "<< indFixed[nfixed]<<" but "<<iv<<" var has only "<<GetNCutSteps(iv)<<" steps"<<endl;
362 }
363 nfixed++;
364 }
365 }
366 */
367 for (Int_t iv=0;iv<fNVariables;iv++)indMax[iv]=0;
368 for(Int_t i=0;i<nfixed;i++){
369 indMax[numFixed[i]]=indFixed[i];
370 }
371 //Get position of fixed vars
372 /*
373 Int_t fixedIndexes[nfixed];
374 Int_t iforfixed=0;
375 for (Int_t iv=0;iv<fNVariables;iv++){
376 if(!isFree[iv]){
377 fixedIndexes[iforfixed]=iv;
378 iforfixed++;
379 }
380 }
381 */
382 TArrayF vett;
383 vett.Set(nelem);
384
385 ULong64_t first=fNTotCells/fNPtBins*ptbin;
386 ULong64_t last=first+fNTotCells/fNPtBins;
387 Int_t dummyptbin;
388
389 maxValue=fVett[GetGlobalAddressFromIndices(indMax,ptbin)];
390 Int_t tmpInd[fNVariables];
391
392 //loop on multidimvector global addresses
393 for(ULong64_t iga=first;iga<last;iga++){
394 GetIndicesFromGlobalAddress(iga,tmpInd,dummyptbin);
395 Bool_t goodCell=kTRUE;
396 for(Int_t ifix=0;ifix<nfixed&&goodCell;ifix++){
397 // if(indFixed[ifix]!=indMax[fixedIndexes[ifix]])goodCell=kFALSE;
398 if(indFixed[ifix]!=tmpInd[numFixed[ifix]])goodCell=kFALSE;
399 }
400 if(goodCell){
401 if(fVett[iga]>maxValue){
402 maxValue=fVett[iga];
403 // GetIndicesFromGlobalAddress(iga,indMax,dummyptbin);
404 for(Int_t inv=0;inv<fNVariables;inv++)indMax[inv]=tmpInd[inv];
405 }
406 }
407 }
408
409 return indMax;
410}
411
ac7636a0 412//_____________________________________________________________________________
2350a1d9 413TH2F* AliMultiDimVector::Project(Int_t firstVar, Int_t secondVar, const Int_t* fixedVars, Int_t ptbin, Float_t norm){
ac7636a0 414 // Project the AliMultiDimVector on a 2D histogram
415
416 TString hisName=Form("hproj%s%dv%d",GetName(),secondVar,firstVar);
417 TString hisTit=Form("%s vs. %s",fAxisTitles[secondVar].Data(),fAxisTitles[firstVar].Data());
418 TH2F* h2=new TH2F(hisName.Data(),hisTit.Data(),fNCutSteps[firstVar],fMinLimits[firstVar],fMaxLimits[firstVar],fNCutSteps[secondVar],fMinLimits[secondVar],fMaxLimits[secondVar]);
419
420 Int_t index[fgkMaxNVariables];
421 for(Int_t i=0;i<fNVariables;i++){
422 index[i]=fixedVars[i];
423 }
424
425 for(Int_t i=0;i<fNCutSteps[firstVar];i++){
426 for(Int_t j=0;j<fNCutSteps[secondVar];j++){
427 index[firstVar]=i;
428 index[secondVar]=j;
429 Float_t cont=GetElement(index,ptbin)/norm;
430 Int_t bin1=i+1;
431 if(!fGreaterThan[firstVar]) bin1=fNCutSteps[firstVar]-i;
432 Int_t bin2=j+1;
433 if(!fGreaterThan[secondVar]) bin2=fNCutSteps[secondVar]-j;
434 h2->SetBinContent(bin1,bin2,cont);
435 }
436 }
437 return h2;
438}
439//_____________________________________________________________________________
440void AliMultiDimVector::GetIntegrationLimits(Int_t iVar, Int_t iCell, Int_t& minbin, Int_t& maxbin) const {
441 // computes bin limits for integrating the AliMultiDimVector
442 minbin=0;
443 maxbin=0;
444 if(iVar<fNVariables){
445 minbin=iCell;
446 maxbin=fNCutSteps[iVar]-1;
447 }
448}
449//_____________________________________________________________________________
450void AliMultiDimVector::GetFillRange(Int_t iVar, Int_t iCell, Int_t& minbin, Int_t& maxbin) const {
451 // computes range of cells passing the cuts for FillAndIntegrate
452 minbin=0;
453 maxbin=0;
454 if(iVar<fNVariables){
455 minbin=0; // bin 0 corresponds to loose cuts
456 maxbin=iCell;
457 }
458}
459//_____________________________________________________________________________
460void AliMultiDimVector::Integrate(){
461 // integrates the matrix
462 if(fIsIntegrated){
463 AliError("MultiDimVector already integrated");
464 return;
465 }
466 TArrayF integral(fNTotCells);
467 for(ULong64_t i=0;i<fNTotCells;i++) integral[i]=CountsAboveCell(i);
468 for(ULong64_t i=0;i<fNTotCells;i++) fVett[i]= integral[i];
469 fIsIntegrated=kTRUE;
bac542b0 470}//_____________________________________________________________________________
2350a1d9 471ULong64_t* AliMultiDimVector::GetGlobalAddressesAboveCuts(const Float_t *values, Int_t ptbin, Int_t& nVals) const{
bac542b0 472 // fills an array with global addresses of cells passing the cuts
473
474 Int_t ind[fgkMaxNVariables];
475 Bool_t retcode=GetIndicesFromValues(values,ind);
476 if(!retcode){
477 nVals=0;
478 return 0x0;
479 }
480 for(Int_t i=fNVariables; i<fgkMaxNVariables; i++) ind[i]=0;
481 Int_t mink[fgkMaxNVariables];
482 Int_t maxk[fgkMaxNVariables];
483 Int_t size=1;
484 for(Int_t i=0;i<fgkMaxNVariables;i++){
485 GetFillRange(i,ind[i],mink[i],maxk[i]);
486 size*=(maxk[i]-mink[i]+1);
487 }
488 ULong64_t* indexes=new ULong64_t[size];
489 nVals=0;
490 for(Int_t k0=mink[0]; k0<=maxk[0]; k0++){
491 for(Int_t k1=mink[1]; k1<=maxk[1]; k1++){
492 for(Int_t k2=mink[2]; k2<=maxk[2]; k2++){
493 for(Int_t k3=mink[3]; k3<=maxk[3]; k3++){
494 for(Int_t k4=mink[4]; k4<=maxk[4]; k4++){
495 for(Int_t k5=mink[5]; k5<=maxk[5]; k5++){
496 for(Int_t k6=mink[6]; k6<=maxk[6]; k6++){
497 for(Int_t k7=mink[7]; k7<=maxk[7]; k7++){
498 for(Int_t k8=mink[8]; k8<=maxk[8]; k8++){
499 for(Int_t k9=mink[9]; k9<=maxk[9]; k9++){
500 Int_t currentBin[fgkMaxNVariables]={k0,k1,k2,k3,k4,k5,k6,k7,k8,k9};
501 indexes[nVals++]=GetGlobalAddressFromIndices(currentBin,ptbin);
502 }
503 }
504 }
505 }
506 }
507 }
508 }
509 }
510 }
511 }
512 return indexes;
ac7636a0 513}
514//_____________________________________________________________________________
515Float_t AliMultiDimVector::CountsAboveCell(ULong64_t globadd) const{
516 // integrates the counts of cells above cell with address globadd
517 Int_t ind[fgkMaxNVariables];
518 Int_t ptbin;
519 GetIndicesFromGlobalAddress(globadd,ind,ptbin);
520 for(Int_t i=fNVariables; i<fgkMaxNVariables; i++) ind[i]=0;
521 Int_t mink[fgkMaxNVariables];
522 Int_t maxk[fgkMaxNVariables];
523 for(Int_t i=0;i<fgkMaxNVariables;i++){
524 GetIntegrationLimits(i,ind[i],mink[i],maxk[i]);
525 }
526 Float_t sumcont=0.;
527 for(Int_t k0=mink[0]; k0<=maxk[0]; k0++){
528 for(Int_t k1=mink[1]; k1<=maxk[1]; k1++){
529 for(Int_t k2=mink[2]; k2<=maxk[2]; k2++){
530 for(Int_t k3=mink[3]; k3<=maxk[3]; k3++){
531 for(Int_t k4=mink[4]; k4<=maxk[4]; k4++){
532 for(Int_t k5=mink[5]; k5<=maxk[5]; k5++){
533 for(Int_t k6=mink[6]; k6<=maxk[6]; k6++){
534 for(Int_t k7=mink[7]; k7<=maxk[7]; k7++){
535 for(Int_t k8=mink[8]; k8<=maxk[8]; k8++){
536 for(Int_t k9=mink[9]; k9<=maxk[9]; k9++){
537 Int_t currentBin[fgkMaxNVariables]={k0,k1,k2,k3,k4,k5,k6,k7,k8,k9};
538 sumcont+=GetElement(currentBin,ptbin);
539 }
540 }
541 }
542 }
543 }
544 }
545 }
546 }
547 }
548 }
549 return sumcont;
550}
551//_____________________________________________________________________________
552void AliMultiDimVector::Fill(Float_t* values, Int_t ptbin){
553 // fills the cells of AliMultiDimVector corresponding to values
554 if(fIsIntegrated){
555 AliError("MultiDimVector already integrated -- Use FillAndIntegrate");
556 return;
557 }
558 Int_t ind[fgkMaxNVariables];
559 Bool_t retcode=GetIndicesFromValues(values,ind);
560 for(Int_t i=fNVariables; i<fgkMaxNVariables; i++) ind[i]=0;
561 if(retcode) IncrementElement(ind,ptbin);
562}
563//_____________________________________________________________________________
564void AliMultiDimVector::FillAndIntegrate(Float_t* values, Int_t ptbin){
565 // fills the cells of AliMultiDimVector passing the cuts
566 // The number of nested loops must match fgkMaxNVariables!!!!
567 fIsIntegrated=kTRUE;
568 Int_t ind[fgkMaxNVariables];
569 Bool_t retcode=GetIndicesFromValues(values,ind);
570 if(!retcode) return;
571 for(Int_t i=fNVariables; i<fgkMaxNVariables; i++) ind[i]=0;
572 Int_t mink[fgkMaxNVariables];
573 Int_t maxk[fgkMaxNVariables];
574 for(Int_t i=0;i<fgkMaxNVariables;i++){
575 GetFillRange(i,ind[i],mink[i],maxk[i]);
576 }
577 for(Int_t k0=mink[0]; k0<=maxk[0]; k0++){
578 for(Int_t k1=mink[1]; k1<=maxk[1]; k1++){
579 for(Int_t k2=mink[2]; k2<=maxk[2]; k2++){
580 for(Int_t k3=mink[3]; k3<=maxk[3]; k3++){
581 for(Int_t k4=mink[4]; k4<=maxk[4]; k4++){
582 for(Int_t k5=mink[5]; k5<=maxk[5]; k5++){
583 for(Int_t k6=mink[6]; k6<=maxk[6]; k6++){
584 for(Int_t k7=mink[7]; k7<=maxk[7]; k7++){
585 for(Int_t k8=mink[8]; k8<=maxk[8]; k8++){
586 for(Int_t k9=mink[9]; k9<=maxk[9]; k9++){
587 Int_t currentBin[fgkMaxNVariables]={k0,k1,k2,k3,k4,k5,k6,k7,k8,k9};
588 IncrementElement(currentBin,ptbin);
589 }
590 }
591 }
592 }
593 }
594 }
595 }
596 }
597 }
598 }
599
600}
601//_____________________________________________________________________________
2350a1d9 602void AliMultiDimVector::SuppressZeroBKGEffect(const AliMultiDimVector* mvBKG){
ac7636a0 603 // Sets to zero elements for which mvBKG=0
604 for(ULong64_t i=0;i<fNTotCells;i++)
2350a1d9 605 if(mvBKG->GetElement(i)<0.00000001) fVett.AddAt(0,i);
ac7636a0 606}
607//_____________________________________________________________________________
608AliMultiDimVector* AliMultiDimVector:: ShrinkPtBins(Int_t firstBin, Int_t lastBin){
609 // sums the elements of pt bins between firstBin and lastBin
610 if(firstBin<0 || lastBin>=fNPtBins || firstBin>=lastBin){
611 AliError("Bad numbers of Pt bins to be shrinked");
612 return 0;
613 }
614 Int_t nofcells[fgkMaxNVariables];
615 Float_t loosecuts[fgkMaxNVariables];
616 Float_t tightcuts[fgkMaxNVariables];
617 TString axisTitles[fgkMaxNVariables];
e11ae259 618 for(Int_t j=0;j<fgkMaxNVariables;j++) {
619 nofcells[j]=0;
620 loosecuts[j]=0.;
621 tightcuts[j]=0.;
622 axisTitles[j]="";
623 }
ac7636a0 624 for(Int_t i=0;i<fNVariables;i++){
625 nofcells[i]=fNCutSteps[i];
626 if(fGreaterThan[i]){
627 loosecuts[i]=fMinLimits[i];
628 tightcuts[i]=fMaxLimits[i];
629 }else{
630 loosecuts[i]=fMaxLimits[i];
631 tightcuts[i]=fMinLimits[i];
632 }
633 axisTitles[i]=fAxisTitles[i];
634 }
635 Int_t newNptbins=fNPtBins-(lastBin-firstBin);
bac542b0 636 Float_t ptlimits[fgkMaxNPtBins+1];
637 for(Int_t ipt=0; ipt<=firstBin;ipt++) ptlimits[ipt]=fPtLimits[ipt];
638 for(Int_t ipt=firstBin+1; ipt<newNptbins+1;ipt++) ptlimits[ipt]=fPtLimits[ipt+(lastBin-firstBin)];
639 AliMultiDimVector* shrinkedMV=new AliMultiDimVector(GetName(),GetTitle(),newNptbins,ptlimits,fNVariables,nofcells,loosecuts,tightcuts,axisTitles);
ac7636a0 640
641 ULong64_t nOfPointsPerPtbin=fNTotCells/fNPtBins;
642 ULong64_t addressOld,addressNew;
643 Int_t npb,opb;
644 for(npb=0;npb<firstBin;npb++){
645 opb=npb;
646 for(ULong64_t k=0;k<nOfPointsPerPtbin;k++){
647 addressOld=opb+k*fNPtBins;
648 addressNew=npb+k*newNptbins;
649 shrinkedMV->SetElement(addressNew,fVett[addressOld]);
650 }
651 }
652 npb=firstBin;
653 for(ULong64_t k=0;k<nOfPointsPerPtbin;k++){
654 Float_t summedValue=0.;
655 for(opb=firstBin;opb<=lastBin;opb++){
656 addressOld=opb+k*fNPtBins;
657 summedValue+=fVett[addressOld];
658 }
659 addressNew=npb+k*newNptbins;
660 shrinkedMV->SetElement(addressNew,summedValue);
661 }
662 for(npb=firstBin+1;npb<newNptbins;npb++){
663 opb=npb+(lastBin-firstBin);
664 for(ULong64_t k=0;k<nOfPointsPerPtbin;k++){
665 addressOld=opb+k*fNPtBins;
666 addressNew=npb+k*newNptbins;
667 shrinkedMV->SetElement(addressNew,fVett[addressOld]);
668 }
669 }
670 return shrinkedMV;
671}
bac542b0 672//_____________________________________________________________________________
3aa6ce53 673void AliMultiDimVector::SetNewLimits(Float_t* loose,Float_t* tight){
674 for(Int_t i=0;i<fNVariables;i++){
675 if(loose[i] < tight[i]){
676 fMinLimits[i]=loose[i];
677 fMaxLimits[i]=tight[i];
678 fGreaterThan[i]=kTRUE;
679 }else{
680 fMinLimits[i]=tight[i];
681 fMaxLimits[i]=loose[i];
682 fGreaterThan[i]=kFALSE;
683 }
684 }
685}
686//_____________________________________________________________________________
687void AliMultiDimVector::SwapLimits(Int_t ivar){
688 Float_t oldmin = fMinLimits[ivar];
689 fMinLimits[ivar] = fMaxLimits[ivar];
690 fMaxLimits[ivar] = oldmin;
691 if(fGreaterThan[ivar])fGreaterThan[ivar]=kFALSE;
692 else fGreaterThan[ivar]=kTRUE;
693}
694//_____________________________________________________________________________
bac542b0 695void AliMultiDimVector::PrintStatus(){
696 //
697 printf("Number of Pt bins = %d\n",fNPtBins);
698 printf("Limits of Pt bins = ");
699 for(Int_t ib=0;ib<fNPtBins+1;ib++) printf("%6.2f ",fPtLimits[ib]);
700 printf("\n");
701 printf("Number of cut variables = %d\n",fNVariables);
702 for(Int_t iv=0;iv<fNVariables;iv++){
703 printf("- Variable %d: %s\n",iv,fAxisTitles[iv].Data());
704 printf(" Nsteps= %d Rage = %6.2f %6.2f\n",
705 fNCutSteps[iv],fMinLimits[iv],fMaxLimits[iv]);
706 }
707}