]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG3/vertexingHF/AliMultiDimVector.cxx
Dplus task used AliAODPidHF via AliRDHFCutsDplustoKpipi (Renu, Francesco)
[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
16/* $Id: $ */
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
29
30#include "TH2.h"
31#include "AliMultiDimVector.h"
32#include "AliLog.h"
ac7636a0 33#include "TString.h"
34
35ClassImp(AliMultiDimVector)
36//___________________________________________________________________________
37AliMultiDimVector::AliMultiDimVector():TNamed("AliMultiDimVector","default"),
38fNVariables(0),
39fNPtBins(0),
40fVett(0),
41fNTotCells(0),
42fIsIntegrated(0)
43{
44 // default constructor
45}
46//___________________________________________________________________________
2350a1d9 47AliMultiDimVector::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 48fNVariables(npars),
49fNPtBins(nptbins),
50fVett(0),
51fNTotCells(0),
52fIsIntegrated(0){
53// standard constructor
54 ULong64_t ntot=1;
55 for(Int_t i=0;i<fNVariables;i++){
56 ntot*=nofcells[i];
57 fNCutSteps[i]=nofcells[i];
3aa6ce53 58 if(loosecuts[i] == tightcuts[i]){
59 if(loosecuts[i]!=0){
60 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]);
61 fMinLimits[i]=tightcuts[i]-0.1*tightcuts[i];
62 fMaxLimits[i]=tightcuts[i];
63 }else{
64 fMinLimits[i]=0;
65 fMaxLimits[i]=0.0001;
66 }
67 fGreaterThan[i]=kTRUE;
68 }
69 if(loosecuts[i] < tightcuts[i]){
ac7636a0 70 fMinLimits[i]=loosecuts[i];
71 fMaxLimits[i]=tightcuts[i];
72 fGreaterThan[i]=kTRUE;
73 }else{
74 fMinLimits[i]=tightcuts[i];
75 fMaxLimits[i]=loosecuts[i];
76 fGreaterThan[i]=kFALSE;
77 }
78 fAxisTitles[i]=axisTitles[i].Data();
79 }
bac542b0 80 fNTotCells=ntot*fNPtBins;
ac7636a0 81 fVett.Set(fNTotCells);
bac542b0 82 for(Int_t ipt=0;ipt<fNPtBins+1;ipt++) fPtLimits[ipt]=ptlimits[ipt];
83 for(Int_t ipt=fNPtBins+1;ipt<fgkMaxNPtBins+1;ipt++) fPtLimits[ipt]=999.;
ac7636a0 84 for (ULong64_t j=0;j<fNTotCells;j++) fVett.AddAt(0,j);
85}
86//___________________________________________________________________________
87AliMultiDimVector::AliMultiDimVector(const AliMultiDimVector &mv):TNamed(mv.GetName(),mv.GetTitle()),
88fNVariables(mv.fNVariables),
89fNPtBins(mv.fNPtBins),
90fVett(0),
91fNTotCells(mv.fNTotCells),
92fIsIntegrated(mv.fIsIntegrated)
93{
94 // copy constructor
95 for(Int_t i=0;i<fNVariables;i++){
96 fNCutSteps[i]=mv.GetNCutSteps(i);
97 fMinLimits[i]=mv.GetMinLimit(i);
98 fMaxLimits[i]=mv.GetMaxLimit(i);
99 fGreaterThan[i]=mv.GetGreaterThan(i);
100 fAxisTitles[i]=mv.GetAxisTitle(i);
101 }
102 fVett.Set(fNTotCells);
bac542b0 103
104 for(Int_t ipt=0;ipt<fNPtBins+1;ipt++) fPtLimits[ipt]=mv.GetPtLimit(ipt);
ac7636a0 105 for(ULong64_t i=0;i<fNTotCells;i++) fVett[i]=mv.GetElement(i);
106}
107//___________________________________________________________________________
108void AliMultiDimVector::CopyStructure(const AliMultiDimVector* mv){
109// Sets dimensions and limit from mv
110 fNVariables=mv->GetNVariables();
111 fNPtBins=mv->GetNPtBins();
112 fNTotCells=mv->GetNTotCells();
113 fIsIntegrated=mv->IsIntegrated();
114 for(Int_t i=0;i<fNVariables;i++){
115 fNCutSteps[i]=mv->GetNCutSteps(i);
116 fMinLimits[i]=mv->GetMinLimit(i);
117 fMaxLimits[i]=mv->GetMaxLimit(i);
118 fGreaterThan[i]=mv->GetGreaterThan(i);
119 fAxisTitles[i]=mv->GetAxisTitle(i);
120 }
bac542b0 121 for(Int_t ipt=0;ipt<fNPtBins+1;ipt++) fPtLimits[ipt]=mv->GetPtLimit(ipt);
122 fVett.Set(fNTotCells);
ac7636a0 123}
124//______________________________________________________________________
125Bool_t AliMultiDimVector::GetIndicesFromGlobalAddress(ULong64_t globadd, Int_t *ind, Int_t &ptbin) const {
126// returns matrix element indices and Pt bin from global index
127 if(globadd>=fNTotCells) return kFALSE;
128 ULong64_t r=globadd;
129 Int_t prod=1;
130 Int_t nOfCellsPlusLevel[fgkMaxNVariables+1];
131 for(Int_t k=0;k<fNVariables;k++) nOfCellsPlusLevel[k]=fNCutSteps[k];
132 nOfCellsPlusLevel[fNVariables]=fNPtBins;
133
134 for(Int_t i=0;i<fNVariables+1;i++) prod*=nOfCellsPlusLevel[i];
135 for(Int_t i=0;i<fNVariables+1;i++){
136 prod/=nOfCellsPlusLevel[i];
137 if(i<fNVariables) ind[i]=r/prod;
138 else ptbin=r/prod;
139 r=globadd%prod;
140 }
141 return kTRUE;
142}
143//______________________________________________________________________
144Bool_t AliMultiDimVector::GetCutValuesFromGlobalAddress(ULong64_t globadd, Float_t *cuts, Int_t &ptbin) const {
145 Int_t ind[fgkMaxNVariables];
146 Bool_t retcode=GetIndicesFromGlobalAddress(globadd,ind,ptbin);
147 if(!retcode) return kFALSE;
148 for(Int_t i=0;i<fNVariables;i++) cuts[i]=GetCutValue(i,ind[i]);
149 return kTRUE;
150}
151//______________________________________________________________________
2350a1d9 152ULong64_t AliMultiDimVector::GetGlobalAddressFromIndices(const Int_t *ind, Int_t ptbin) const {
ac7636a0 153 // Returns the global index of the cell in the matrix
154 Int_t prod=1;
155 ULong64_t elem=0;
156 Int_t indexPlusLevel[fgkMaxNVariables+1];
157 Int_t nOfCellsPlusLevel[fgkMaxNVariables+1];
158 for(Int_t i=0;i<fNVariables;i++){
159 indexPlusLevel[i]=ind[i];
160 nOfCellsPlusLevel[i]=fNCutSteps[i];
161 }
162 indexPlusLevel[fNVariables]=ptbin;
163 nOfCellsPlusLevel[fNVariables]=fNPtBins;
164
165 for(Int_t i=0;i<fNVariables+1;i++){
166 prod=indexPlusLevel[i];
167 if(i<fNVariables){
168 for(Int_t j=i+1;j<fNVariables+1;j++){
169 prod*=nOfCellsPlusLevel[j];
170 }
171 }
172 elem+=prod;
173 }
174 return elem;
175}
176//______________________________________________________________________
2350a1d9 177Bool_t AliMultiDimVector::GetIndicesFromValues(const Float_t *values, Int_t *ind) const {
ac7636a0 178 // Fills the array of matrix indices strating from variable values
179 for(Int_t i=0;i<fNVariables;i++){
180 if(fGreaterThan[i]){
181 if(values[i]<GetMinLimit(i)) return kFALSE;
182 ind[i]=(Int_t)((values[i]-fMinLimits[i])/GetCutStep(i));
183 if(ind[i]>=GetNCutSteps(i)) ind[i]=GetNCutSteps(i)-1;
184 }else{
185 if(values[i]>GetMaxLimit(i)) return kFALSE;
186 ind[i]=(Int_t)((fMaxLimits[i]-values[i])/GetCutStep(i));
187 if(ind[i]>=GetNCutSteps(i)) ind[i]=GetNCutSteps(i)-1;
188 }
189 }
190 return kTRUE;
191}
192//______________________________________________________________________
2350a1d9 193ULong64_t AliMultiDimVector::GetGlobalAddressFromValues(const Float_t *values, Int_t ptbin) const {
ac7636a0 194 // Returns the global index of the cell in the matrix
195 Int_t ind[fgkMaxNVariables];
196 Bool_t retcode=GetIndicesFromValues(values,ind);
197 if(retcode) return GetGlobalAddressFromIndices(ind,ptbin);
198 else{
199 AliError("Values out of range");
200 return fNTotCells+999;
201 }
202}
203//_____________________________________________________________________________
204void AliMultiDimVector::MultiplyBy(Float_t factor){
205 // multiply the AliMultiDimVector by a constant factor
206 for(ULong64_t i=0;i<fNTotCells;i++){
2350a1d9 207 if(fVett.At(i)>0.)
ac7636a0 208 fVett.AddAt(fVett.At(i)*factor,i);
209 else fVett.AddAt(-1,i);
210 }
211
212}
213//_____________________________________________________________________________
2350a1d9 214void AliMultiDimVector::Multiply(const AliMultiDimVector* mv,Float_t factor){
ac7636a0 215 // Sets AliMultiDimVector=mv*constant factor
216 for(ULong64_t i=0;i<fNTotCells;i++){
2350a1d9 217 if(mv->GetElement(i)>0.)
ac7636a0 218 fVett.AddAt(mv->GetElement(i)*factor,i);
219 else fVett.AddAt(-1,i);
220 }
221}
222//_____________________________________________________________________________
2350a1d9 223void AliMultiDimVector::Multiply(const AliMultiDimVector* mv1, const AliMultiDimVector* mv2){
ac7636a0 224 // Sets AliMultiDimVector=mv1*mv2
225 for(ULong64_t i=0;i<fNTotCells;i++){
2350a1d9 226 if(mv1->GetElement(i)>0. && mv2->GetElement(i)>0.)
ac7636a0 227 fVett.AddAt(mv1->GetElement(i)*mv2->GetElement(i),i);
228 else fVett.AddAt(-1,i);
229 }
230}
231//_____________________________________________________________________________
2350a1d9 232void AliMultiDimVector::Add(const AliMultiDimVector* mv){
ac7636a0 233 // Sums contents of mv to AliMultiDimVector
234 if (mv->GetNTotCells()!=fNTotCells){
235 AliError("Different dimension of the vectors!!");
236 }else{
237 for(ULong64_t i=0;i<fNTotCells;i++)
2350a1d9 238 if(mv->GetElement(i)>0. && fVett.At(i)>0.)
ac7636a0 239 fVett.AddAt(fVett.At(i)+mv->GetElement(i),i);
240 else fVett.AddAt(-1,i);
241 }
242}
243//_____________________________________________________________________________
2350a1d9 244void AliMultiDimVector::Sum(const AliMultiDimVector* mv1, const AliMultiDimVector* mv2){
ac7636a0 245 // Sets AliMultiDimVector=mv1+mv2
246 if (fNTotCells!=mv1->GetNTotCells()&&mv1->GetNTotCells()!=mv2->GetNTotCells()) {
247 AliError("Different dimension of the vectors!!");
248 }
249 else{
250 for(ULong64_t i=0;i<mv1->GetNTotCells();i++) {
2350a1d9 251 if(mv1->GetElement(i)>0. && mv2->GetElement(i)>0.)
ac7636a0 252 fVett.AddAt(mv1->GetElement(i)+mv2->GetElement(i),i);
253 else fVett.AddAt(-1,i);
254 }
255 }
256}
257//_____________________________________________________________________________
2350a1d9 258void AliMultiDimVector::LinearComb(const AliMultiDimVector* mv1, Float_t norm1, const AliMultiDimVector* mv2, Float_t norm2){
ac7636a0 259 // Sets AliMultiDimVector=n1*mv1+n2*mv2
260 if (fNTotCells!=mv1->GetNTotCells()&&mv1->GetNTotCells()!=mv2->GetNTotCells()) {
261 AliError("Different dimension of the vectors!!");
262 }
263 else{
264 for(ULong64_t i=0;i<mv1->GetNTotCells();i++) {
2350a1d9 265 if(mv1->GetElement(i)>0. && mv2->GetElement(i)>0.)
ac7636a0 266 fVett.AddAt(norm1*mv1->GetElement(i)+norm2*mv2->GetElement(i),i);
267 else fVett.AddAt(-1,i);
268 }
269 }
270}
271//_____________________________________________________________________________
2350a1d9 272void AliMultiDimVector::DivideBy(const AliMultiDimVector* mv){
ac7636a0 273 // Divide AliMulivector by mv
274 if (mv->GetNTotCells()!=fNTotCells) {
275 AliError("Different dimension of the vectors!!");
276 }
277 else{
278 for(ULong64_t i=0;i<fNTotCells;i++)
2350a1d9 279 if(mv->GetElement(i)!=0 &&mv->GetElement(i)>0. && fVett.At(i)>0.)
ac7636a0 280 fVett.AddAt(fVett.At(i)/mv->GetElement(i),i);
281 else fVett.AddAt(-1,i);
282 }
283
284}
285//_____________________________________________________________________________
2350a1d9 286void AliMultiDimVector::Divide(const AliMultiDimVector* mv1, const AliMultiDimVector* mv2){
ac7636a0 287 // Sets AliMultiDimVector=mv1/mv2
288 if (fNTotCells!=mv1->GetNTotCells()&&mv1->GetNTotCells()!=mv2->GetNTotCells()) {
289 AliError("Different dimension of the vectors!!");
290 }
291 else{
292 for(ULong64_t i=0;i<mv1->GetNTotCells();i++)
2350a1d9 293 if(mv2->GetElement(i)!=0&& mv2->GetElement(i)>0.&& mv1->GetElement(i)>0.)
ac7636a0 294 {
295 fVett.AddAt(mv1->GetElement(i)/mv2->GetElement(i),i);
296 }
297 else fVett.AddAt(-1,i);
298 }
299}
300//_____________________________________________________________________________
301void AliMultiDimVector::Sqrt(){
302 // Sqrt of elements of AliMultiDimVector
303 for(ULong64_t i=0;i<fNTotCells;i++) {
304 if(fVett.At(i)>=0) fVett.AddAt(TMath::Sqrt(fVett.At(i)),i);
305 else {
306 fVett.AddAt(-1,i);
307 }
308 }
309}
310//_____________________________________________________________________________
2350a1d9 311void AliMultiDimVector::Sqrt(const AliMultiDimVector* mv){
ac7636a0 312 // Sets AliMultiDimVector=sqrt(mv)
313 for(ULong64_t i=0;i<fNTotCells;i++)
314 if(mv->GetElement(i)>=0) fVett.AddAt(TMath::Sqrt(mv->GetElement(i)),i);
315 else fVett.AddAt(-1,i);
316}
317//_____________________________________________________________________________
318void AliMultiDimVector::FindMaximum(Float_t& maxValue, Int_t *ind , Int_t ptbin){
319 // finds the element with maximum contents
320 const ULong64_t nelem=fNTotCells/fNPtBins;
321 TArrayF vett;
322 vett.Set(nelem);
323 ULong64_t runningAddress;
324 for(ULong64_t i=0;i<nelem;i++){
325 runningAddress=ptbin+i*fNPtBins;
326 vett.AddAt(fVett[runningAddress],i);
327 }
328 maxValue=TMath::MaxElement(nelem,vett.GetArray());
329 ULong64_t maxAddress=TMath::LocMax(nelem,vett.GetArray());
330 ULong64_t maxGlobalAddress=ptbin+maxAddress*fNPtBins;
331 Int_t checkedptbin;
332 GetIndicesFromGlobalAddress(maxGlobalAddress,ind,checkedptbin);
333}
334
335//_____________________________________________________________________________
2350a1d9 336TH2F* AliMultiDimVector::Project(Int_t firstVar, Int_t secondVar, const Int_t* fixedVars, Int_t ptbin, Float_t norm){
ac7636a0 337 // Project the AliMultiDimVector on a 2D histogram
338
339 TString hisName=Form("hproj%s%dv%d",GetName(),secondVar,firstVar);
340 TString hisTit=Form("%s vs. %s",fAxisTitles[secondVar].Data(),fAxisTitles[firstVar].Data());
341 TH2F* h2=new TH2F(hisName.Data(),hisTit.Data(),fNCutSteps[firstVar],fMinLimits[firstVar],fMaxLimits[firstVar],fNCutSteps[secondVar],fMinLimits[secondVar],fMaxLimits[secondVar]);
342
343 Int_t index[fgkMaxNVariables];
344 for(Int_t i=0;i<fNVariables;i++){
345 index[i]=fixedVars[i];
346 }
347
348 for(Int_t i=0;i<fNCutSteps[firstVar];i++){
349 for(Int_t j=0;j<fNCutSteps[secondVar];j++){
350 index[firstVar]=i;
351 index[secondVar]=j;
352 Float_t cont=GetElement(index,ptbin)/norm;
353 Int_t bin1=i+1;
354 if(!fGreaterThan[firstVar]) bin1=fNCutSteps[firstVar]-i;
355 Int_t bin2=j+1;
356 if(!fGreaterThan[secondVar]) bin2=fNCutSteps[secondVar]-j;
357 h2->SetBinContent(bin1,bin2,cont);
358 }
359 }
360 return h2;
361}
362//_____________________________________________________________________________
363void AliMultiDimVector::GetIntegrationLimits(Int_t iVar, Int_t iCell, Int_t& minbin, Int_t& maxbin) const {
364 // computes bin limits for integrating the AliMultiDimVector
365 minbin=0;
366 maxbin=0;
367 if(iVar<fNVariables){
368 minbin=iCell;
369 maxbin=fNCutSteps[iVar]-1;
370 }
371}
372//_____________________________________________________________________________
373void AliMultiDimVector::GetFillRange(Int_t iVar, Int_t iCell, Int_t& minbin, Int_t& maxbin) const {
374 // computes range of cells passing the cuts for FillAndIntegrate
375 minbin=0;
376 maxbin=0;
377 if(iVar<fNVariables){
378 minbin=0; // bin 0 corresponds to loose cuts
379 maxbin=iCell;
380 }
381}
382//_____________________________________________________________________________
383void AliMultiDimVector::Integrate(){
384 // integrates the matrix
385 if(fIsIntegrated){
386 AliError("MultiDimVector already integrated");
387 return;
388 }
389 TArrayF integral(fNTotCells);
390 for(ULong64_t i=0;i<fNTotCells;i++) integral[i]=CountsAboveCell(i);
391 for(ULong64_t i=0;i<fNTotCells;i++) fVett[i]= integral[i];
392 fIsIntegrated=kTRUE;
bac542b0 393}//_____________________________________________________________________________
2350a1d9 394ULong64_t* AliMultiDimVector::GetGlobalAddressesAboveCuts(const Float_t *values, Int_t ptbin, Int_t& nVals) const{
bac542b0 395 // fills an array with global addresses of cells passing the cuts
396
397 Int_t ind[fgkMaxNVariables];
398 Bool_t retcode=GetIndicesFromValues(values,ind);
399 if(!retcode){
400 nVals=0;
401 return 0x0;
402 }
403 for(Int_t i=fNVariables; i<fgkMaxNVariables; i++) ind[i]=0;
404 Int_t mink[fgkMaxNVariables];
405 Int_t maxk[fgkMaxNVariables];
406 Int_t size=1;
407 for(Int_t i=0;i<fgkMaxNVariables;i++){
408 GetFillRange(i,ind[i],mink[i],maxk[i]);
409 size*=(maxk[i]-mink[i]+1);
410 }
411 ULong64_t* indexes=new ULong64_t[size];
412 nVals=0;
413 for(Int_t k0=mink[0]; k0<=maxk[0]; k0++){
414 for(Int_t k1=mink[1]; k1<=maxk[1]; k1++){
415 for(Int_t k2=mink[2]; k2<=maxk[2]; k2++){
416 for(Int_t k3=mink[3]; k3<=maxk[3]; k3++){
417 for(Int_t k4=mink[4]; k4<=maxk[4]; k4++){
418 for(Int_t k5=mink[5]; k5<=maxk[5]; k5++){
419 for(Int_t k6=mink[6]; k6<=maxk[6]; k6++){
420 for(Int_t k7=mink[7]; k7<=maxk[7]; k7++){
421 for(Int_t k8=mink[8]; k8<=maxk[8]; k8++){
422 for(Int_t k9=mink[9]; k9<=maxk[9]; k9++){
423 Int_t currentBin[fgkMaxNVariables]={k0,k1,k2,k3,k4,k5,k6,k7,k8,k9};
424 indexes[nVals++]=GetGlobalAddressFromIndices(currentBin,ptbin);
425 }
426 }
427 }
428 }
429 }
430 }
431 }
432 }
433 }
434 }
435 return indexes;
ac7636a0 436}
437//_____________________________________________________________________________
438Float_t AliMultiDimVector::CountsAboveCell(ULong64_t globadd) const{
439 // integrates the counts of cells above cell with address globadd
440 Int_t ind[fgkMaxNVariables];
441 Int_t ptbin;
442 GetIndicesFromGlobalAddress(globadd,ind,ptbin);
443 for(Int_t i=fNVariables; i<fgkMaxNVariables; i++) ind[i]=0;
444 Int_t mink[fgkMaxNVariables];
445 Int_t maxk[fgkMaxNVariables];
446 for(Int_t i=0;i<fgkMaxNVariables;i++){
447 GetIntegrationLimits(i,ind[i],mink[i],maxk[i]);
448 }
449 Float_t sumcont=0.;
450 for(Int_t k0=mink[0]; k0<=maxk[0]; k0++){
451 for(Int_t k1=mink[1]; k1<=maxk[1]; k1++){
452 for(Int_t k2=mink[2]; k2<=maxk[2]; k2++){
453 for(Int_t k3=mink[3]; k3<=maxk[3]; k3++){
454 for(Int_t k4=mink[4]; k4<=maxk[4]; k4++){
455 for(Int_t k5=mink[5]; k5<=maxk[5]; k5++){
456 for(Int_t k6=mink[6]; k6<=maxk[6]; k6++){
457 for(Int_t k7=mink[7]; k7<=maxk[7]; k7++){
458 for(Int_t k8=mink[8]; k8<=maxk[8]; k8++){
459 for(Int_t k9=mink[9]; k9<=maxk[9]; k9++){
460 Int_t currentBin[fgkMaxNVariables]={k0,k1,k2,k3,k4,k5,k6,k7,k8,k9};
461 sumcont+=GetElement(currentBin,ptbin);
462 }
463 }
464 }
465 }
466 }
467 }
468 }
469 }
470 }
471 }
472 return sumcont;
473}
474//_____________________________________________________________________________
475void AliMultiDimVector::Fill(Float_t* values, Int_t ptbin){
476 // fills the cells of AliMultiDimVector corresponding to values
477 if(fIsIntegrated){
478 AliError("MultiDimVector already integrated -- Use FillAndIntegrate");
479 return;
480 }
481 Int_t ind[fgkMaxNVariables];
482 Bool_t retcode=GetIndicesFromValues(values,ind);
483 for(Int_t i=fNVariables; i<fgkMaxNVariables; i++) ind[i]=0;
484 if(retcode) IncrementElement(ind,ptbin);
485}
486//_____________________________________________________________________________
487void AliMultiDimVector::FillAndIntegrate(Float_t* values, Int_t ptbin){
488 // fills the cells of AliMultiDimVector passing the cuts
489 // The number of nested loops must match fgkMaxNVariables!!!!
490 fIsIntegrated=kTRUE;
491 Int_t ind[fgkMaxNVariables];
492 Bool_t retcode=GetIndicesFromValues(values,ind);
493 if(!retcode) return;
494 for(Int_t i=fNVariables; i<fgkMaxNVariables; i++) ind[i]=0;
495 Int_t mink[fgkMaxNVariables];
496 Int_t maxk[fgkMaxNVariables];
497 for(Int_t i=0;i<fgkMaxNVariables;i++){
498 GetFillRange(i,ind[i],mink[i],maxk[i]);
499 }
500 for(Int_t k0=mink[0]; k0<=maxk[0]; k0++){
501 for(Int_t k1=mink[1]; k1<=maxk[1]; k1++){
502 for(Int_t k2=mink[2]; k2<=maxk[2]; k2++){
503 for(Int_t k3=mink[3]; k3<=maxk[3]; k3++){
504 for(Int_t k4=mink[4]; k4<=maxk[4]; k4++){
505 for(Int_t k5=mink[5]; k5<=maxk[5]; k5++){
506 for(Int_t k6=mink[6]; k6<=maxk[6]; k6++){
507 for(Int_t k7=mink[7]; k7<=maxk[7]; k7++){
508 for(Int_t k8=mink[8]; k8<=maxk[8]; k8++){
509 for(Int_t k9=mink[9]; k9<=maxk[9]; k9++){
510 Int_t currentBin[fgkMaxNVariables]={k0,k1,k2,k3,k4,k5,k6,k7,k8,k9};
511 IncrementElement(currentBin,ptbin);
512 }
513 }
514 }
515 }
516 }
517 }
518 }
519 }
520 }
521 }
522
523}
524//_____________________________________________________________________________
2350a1d9 525void AliMultiDimVector::SuppressZeroBKGEffect(const AliMultiDimVector* mvBKG){
ac7636a0 526 // Sets to zero elements for which mvBKG=0
527 for(ULong64_t i=0;i<fNTotCells;i++)
2350a1d9 528 if(mvBKG->GetElement(i)<0.00000001) fVett.AddAt(0,i);
ac7636a0 529}
530//_____________________________________________________________________________
531AliMultiDimVector* AliMultiDimVector:: ShrinkPtBins(Int_t firstBin, Int_t lastBin){
532 // sums the elements of pt bins between firstBin and lastBin
533 if(firstBin<0 || lastBin>=fNPtBins || firstBin>=lastBin){
534 AliError("Bad numbers of Pt bins to be shrinked");
535 return 0;
536 }
537 Int_t nofcells[fgkMaxNVariables];
538 Float_t loosecuts[fgkMaxNVariables];
539 Float_t tightcuts[fgkMaxNVariables];
540 TString axisTitles[fgkMaxNVariables];
541 for(Int_t i=0;i<fNVariables;i++){
542 nofcells[i]=fNCutSteps[i];
543 if(fGreaterThan[i]){
544 loosecuts[i]=fMinLimits[i];
545 tightcuts[i]=fMaxLimits[i];
546 }else{
547 loosecuts[i]=fMaxLimits[i];
548 tightcuts[i]=fMinLimits[i];
549 }
550 axisTitles[i]=fAxisTitles[i];
551 }
552 Int_t newNptbins=fNPtBins-(lastBin-firstBin);
bac542b0 553 Float_t ptlimits[fgkMaxNPtBins+1];
554 for(Int_t ipt=0; ipt<=firstBin;ipt++) ptlimits[ipt]=fPtLimits[ipt];
555 for(Int_t ipt=firstBin+1; ipt<newNptbins+1;ipt++) ptlimits[ipt]=fPtLimits[ipt+(lastBin-firstBin)];
556 AliMultiDimVector* shrinkedMV=new AliMultiDimVector(GetName(),GetTitle(),newNptbins,ptlimits,fNVariables,nofcells,loosecuts,tightcuts,axisTitles);
ac7636a0 557
558 ULong64_t nOfPointsPerPtbin=fNTotCells/fNPtBins;
559 ULong64_t addressOld,addressNew;
560 Int_t npb,opb;
561 for(npb=0;npb<firstBin;npb++){
562 opb=npb;
563 for(ULong64_t k=0;k<nOfPointsPerPtbin;k++){
564 addressOld=opb+k*fNPtBins;
565 addressNew=npb+k*newNptbins;
566 shrinkedMV->SetElement(addressNew,fVett[addressOld]);
567 }
568 }
569 npb=firstBin;
570 for(ULong64_t k=0;k<nOfPointsPerPtbin;k++){
571 Float_t summedValue=0.;
572 for(opb=firstBin;opb<=lastBin;opb++){
573 addressOld=opb+k*fNPtBins;
574 summedValue+=fVett[addressOld];
575 }
576 addressNew=npb+k*newNptbins;
577 shrinkedMV->SetElement(addressNew,summedValue);
578 }
579 for(npb=firstBin+1;npb<newNptbins;npb++){
580 opb=npb+(lastBin-firstBin);
581 for(ULong64_t k=0;k<nOfPointsPerPtbin;k++){
582 addressOld=opb+k*fNPtBins;
583 addressNew=npb+k*newNptbins;
584 shrinkedMV->SetElement(addressNew,fVett[addressOld]);
585 }
586 }
587 return shrinkedMV;
588}
bac542b0 589//_____________________________________________________________________________
3aa6ce53 590void AliMultiDimVector::SetNewLimits(Float_t* loose,Float_t* tight){
591 for(Int_t i=0;i<fNVariables;i++){
592 if(loose[i] < tight[i]){
593 fMinLimits[i]=loose[i];
594 fMaxLimits[i]=tight[i];
595 fGreaterThan[i]=kTRUE;
596 }else{
597 fMinLimits[i]=tight[i];
598 fMaxLimits[i]=loose[i];
599 fGreaterThan[i]=kFALSE;
600 }
601 }
602}
603//_____________________________________________________________________________
604void AliMultiDimVector::SwapLimits(Int_t ivar){
605 Float_t oldmin = fMinLimits[ivar];
606 fMinLimits[ivar] = fMaxLimits[ivar];
607 fMaxLimits[ivar] = oldmin;
608 if(fGreaterThan[ivar])fGreaterThan[ivar]=kFALSE;
609 else fGreaterThan[ivar]=kTRUE;
610}
611//_____________________________________________________________________________
bac542b0 612void AliMultiDimVector::PrintStatus(){
613 //
614 printf("Number of Pt bins = %d\n",fNPtBins);
615 printf("Limits of Pt bins = ");
616 for(Int_t ib=0;ib<fNPtBins+1;ib++) printf("%6.2f ",fPtLimits[ib]);
617 printf("\n");
618 printf("Number of cut variables = %d\n",fNVariables);
619 for(Int_t iv=0;iv<fNVariables;iv++){
620 printf("- Variable %d: %s\n",iv,fAxisTitles[iv].Data());
621 printf(" Nsteps= %d Rage = %6.2f %6.2f\n",
622 fNCutSteps[iv],fMinLimits[iv],fMaxLimits[iv]);
623 }
624}