1 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * Author: The ALICE Off-line Project. *
4 * Contributors are mentioned in the code where appropriate. *
6 * Permission to use, copy, modify and distribute this software and its *
7 * documentation strictly for non-commercial purposes is hereby granted *
8 * without fee, provided that the above copyright notice appears in all *
9 * copies and that both the copyright notice and this permission notice *
10 * appear in the supporting documentation. The authors make no claims *
11 * about the suitability of this software for any purpose. It is *
12 * provided "as is" without express or implied warranty. *
13 **************************************************************************/
18 // An event cut class for the flow framework
20 // origin: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
30 #include "AliVVertex.h"
31 #include "AliVEvent.h"
32 #include "AliESDEvent.h"
33 #include "AliAODEvent.h"
34 #include "AliAODHeader.h"
35 #include "AliCentrality.h"
36 #include "AliESDVZERO.h"
37 #include "AliMultiplicity.h"
38 #include "AliMCEvent.h"
39 #include "AliFlowEventCuts.h"
40 #include "AliFlowTrackCuts.h"
41 #include "AliTriggerAnalysis.h"
42 #include "AliCollisionGeometry.h"
43 #include "AliGenEventHeader.h"
47 ClassImp(AliFlowEventCuts)
49 //-----------------------------------------------------------------------
50 AliFlowEventCuts::AliFlowEventCuts():
51 AliFlowEventSimpleCuts(),
53 fCutNumberOfTracks(kFALSE),
54 fNumberOfTracksMax(INT_MAX),
55 fNumberOfTracksMin(INT_MIN),
57 fRefMultMethod(kTPConly),
58 fUseAliESDtrackCutsRefMult(kFALSE),
59 fRefMultMethodAliESDtrackCuts(AliESDtrackCuts::kTrackletsITSTPC),
64 fStandardTPCcuts(NULL),
65 fStandardGlobalCuts(NULL),
66 fCutPrimaryVertexX(kFALSE),
67 fPrimaryVertexXmax(INT_MAX),
68 fPrimaryVertexXmin(INT_MIN),
69 fCutPrimaryVertexY(kFALSE),
70 fPrimaryVertexYmax(INT_MAX),
71 fPrimaryVertexYmin(INT_MIN),
72 fCutPrimaryVertexZ(kFALSE),
73 fPrimaryVertexZmax(INT_MAX),
74 fPrimaryVertexZmin(INT_MIN),
75 fCutNContributors(kFALSE),
76 fNContributorsMax(INT_MAX),
77 fNContributorsMin(INT_MIN),
81 fCutSPDvertexerAnomaly(kFALSE),
82 fCutSPDTRKVtxZ(kFALSE),
83 fCutTPCmultiplicityOutliers(kFALSE),
84 fCutTPCmultiplicityOutliersAOD(kFALSE),
85 fUseCentralityUnchecked(kFALSE),
86 fCentralityPercentileMethod(kTPConly),
87 fCutZDCtiming(kFALSE),
89 fCutImpactParameter(kFALSE),
90 fImpactParameterMin(0.0),
91 fImpactParameterMax(100.0),
92 fhistTPCvsGlobalMult(0),
98 //-----------------------------------------------------------------------
99 AliFlowEventCuts::AliFlowEventCuts(const char* name, const char* title):
100 AliFlowEventSimpleCuts(name, title),
102 fCutNumberOfTracks(kFALSE),
103 fNumberOfTracksMax(INT_MAX),
104 fNumberOfTracksMin(INT_MIN),
106 fRefMultMethod(kTPConly),
107 fUseAliESDtrackCutsRefMult(kFALSE),
108 fRefMultMethodAliESDtrackCuts(AliESDtrackCuts::kTrackletsITSTPC),
109 fRefMultMax(INT_MAX),
110 fRefMultMin(INT_MIN),
113 fStandardTPCcuts(AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010()),
114 fStandardGlobalCuts(AliFlowTrackCuts::GetStandardGlobalTrackCuts2010()),
115 fCutPrimaryVertexX(kFALSE),
116 fPrimaryVertexXmax(INT_MAX),
117 fPrimaryVertexXmin(INT_MIN),
118 fCutPrimaryVertexY(kFALSE),
119 fPrimaryVertexYmax(INT_MAX),
120 fPrimaryVertexYmin(INT_MIN),
121 fCutPrimaryVertexZ(kFALSE),
122 fPrimaryVertexZmax(INT_MAX),
123 fPrimaryVertexZmin(INT_MIN),
124 fCutNContributors(kFALSE),
125 fNContributorsMax(INT_MAX),
126 fNContributorsMin(INT_MIN),
128 fMeanPtMax(-DBL_MAX),
130 fCutSPDvertexerAnomaly(kFALSE),
131 fCutSPDTRKVtxZ(kFALSE),
132 fCutTPCmultiplicityOutliers(kFALSE),
133 fCutTPCmultiplicityOutliersAOD(kFALSE),
134 fUseCentralityUnchecked(kFALSE),
135 fCentralityPercentileMethod(kTPConly),
136 fCutZDCtiming(kFALSE),
138 fCutImpactParameter(kFALSE),
139 fImpactParameterMin(0.0),
140 fImpactParameterMax(100.0),
141 fhistTPCvsGlobalMult(0),
147 ////-----------------------------------------------------------------------
148 AliFlowEventCuts::AliFlowEventCuts(const AliFlowEventCuts& that):
149 AliFlowEventSimpleCuts(that),
151 fCutNumberOfTracks(that.fCutNumberOfTracks),
152 fNumberOfTracksMax(that.fNumberOfTracksMax),
153 fNumberOfTracksMin(that.fNumberOfTracksMin),
154 fCutRefMult(that.fCutRefMult),
155 fRefMultMethod(that.fRefMultMethod),
156 fUseAliESDtrackCutsRefMult(that.fUseAliESDtrackCutsRefMult),
157 fRefMultMethodAliESDtrackCuts(that.fRefMultMethodAliESDtrackCuts),
158 fRefMultMax(that.fRefMultMax),
159 fRefMultMin(that.fRefMultMin),
162 fStandardTPCcuts(NULL),
163 fStandardGlobalCuts(NULL),
164 fCutPrimaryVertexX(that.fCutPrimaryVertexX),
165 fPrimaryVertexXmax(that.fPrimaryVertexXmax),
166 fPrimaryVertexXmin(that.fPrimaryVertexXmin),
167 fCutPrimaryVertexY(that.fCutPrimaryVertexX),
168 fPrimaryVertexYmax(that.fPrimaryVertexYmax),
169 fPrimaryVertexYmin(that.fPrimaryVertexYmin),
170 fCutPrimaryVertexZ(that.fCutPrimaryVertexX),
171 fPrimaryVertexZmax(that.fPrimaryVertexZmax),
172 fPrimaryVertexZmin(that.fPrimaryVertexZmin),
173 fCutNContributors(that.fCutNContributors),
174 fNContributorsMax(that.fNContributorsMax),
175 fNContributorsMin(that.fNContributorsMin),
176 fCutMeanPt(that.fCutMeanPt),
177 fMeanPtMax(that.fMeanPtMax),
178 fMeanPtMin(that.fMeanPtMin),
179 fCutSPDvertexerAnomaly(that.fCutSPDvertexerAnomaly),
180 fCutSPDTRKVtxZ(that.fCutSPDTRKVtxZ),
181 fCutTPCmultiplicityOutliers(that.fCutTPCmultiplicityOutliers),
182 fCutTPCmultiplicityOutliersAOD(that.fCutTPCmultiplicityOutliersAOD),
183 fUseCentralityUnchecked(that.fUseCentralityUnchecked),
184 fCentralityPercentileMethod(that.fCentralityPercentileMethod),
185 fCutZDCtiming(that.fCutZDCtiming),
187 fCutImpactParameter(that.fCutImpactParameter),
188 fImpactParameterMin(that.fImpactParameterMin),
189 fImpactParameterMax(that.fImpactParameterMax),
190 fhistTPCvsGlobalMult(that.fhistTPCvsGlobalMult),
191 fData2011(that.fData2011)
193 if (that.fQA) DefineHistograms();
195 if (that.fRefMultCuts)
196 fRefMultCuts = new AliFlowTrackCuts(*(that.fRefMultCuts));
197 if (that.fMeanPtCuts)
198 fMeanPtCuts = new AliFlowTrackCuts(*(that.fMeanPtCuts));
199 fStandardTPCcuts = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010();
200 fStandardGlobalCuts = AliFlowTrackCuts::GetStandardGlobalTrackCuts2010();
203 ////-----------------------------------------------------------------------
204 AliFlowEventCuts::~AliFlowEventCuts()
209 delete fStandardGlobalCuts;
210 delete fStandardTPCcuts;
211 if (fQA) { fQA->SetOwner(); fQA->Delete(); delete fQA; }
214 ////-----------------------------------------------------------------------
215 AliFlowEventCuts& AliFlowEventCuts::operator=(const AliFlowEventCuts& that)
218 if (this==&that) return *this;
227 fQA = static_cast<TList*>(that.fQA->Clone());
236 fCutNumberOfTracks=that.fCutNumberOfTracks;
237 fNumberOfTracksMax=that.fNumberOfTracksMax;
238 fNumberOfTracksMin=that.fNumberOfTracksMin;
239 fCutRefMult=that.fCutRefMult;
240 fRefMultMethod=that.fRefMultMethod;
241 fUseAliESDtrackCutsRefMult=that.fUseAliESDtrackCutsRefMult;
242 fRefMultMethodAliESDtrackCuts=that.fRefMultMethodAliESDtrackCuts;
243 fRefMultMax=that.fRefMultMax;
244 fRefMultMin=that.fRefMultMin;
245 if (that.fRefMultCuts) *fRefMultCuts=*(that.fRefMultCuts);
246 if (that.fMeanPtCuts) *fMeanPtCuts=*(that.fMeanPtCuts);
247 fStandardTPCcuts = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010();
248 fStandardGlobalCuts = AliFlowTrackCuts::GetStandardGlobalTrackCuts2010();
249 fCutPrimaryVertexX=that.fCutPrimaryVertexX;
250 fPrimaryVertexXmax=that.fPrimaryVertexXmax;
251 fPrimaryVertexXmin=that.fPrimaryVertexXmin;
252 fCutPrimaryVertexY=that.fCutPrimaryVertexY;
253 fPrimaryVertexYmax=that.fPrimaryVertexYmax;
254 fPrimaryVertexYmin=that.fPrimaryVertexYmin;
255 fCutPrimaryVertexZ=that.fCutPrimaryVertexZ;
256 fPrimaryVertexZmax=that.fPrimaryVertexZmax;
257 fPrimaryVertexZmin=that.fPrimaryVertexZmin;
258 fCutNContributors=that.fCutNContributors;
259 fNContributorsMax=that.fNContributorsMax;
260 fNContributorsMin=that.fNContributorsMin;
261 fCutMeanPt=that.fCutMeanPt;
262 fMeanPtMax=that.fMeanPtMax;
263 fMeanPtMin=that.fMeanPtMin;
264 fCutSPDvertexerAnomaly=that.fCutSPDvertexerAnomaly;
265 fCutSPDTRKVtxZ=that.fCutSPDTRKVtxZ;
266 fCutTPCmultiplicityOutliers=that.fCutTPCmultiplicityOutliers;
267 fCutTPCmultiplicityOutliersAOD=that.fCutTPCmultiplicityOutliersAOD;
268 fUseCentralityUnchecked=that.fUseCentralityUnchecked;
269 fCentralityPercentileMethod=that.fCentralityPercentileMethod;
270 fCutZDCtiming=that.fCutZDCtiming;
271 fhistTPCvsGlobalMult=that.fhistTPCvsGlobalMult;
272 fData2011=that.fData2011;
276 //-----------------------------------------------------------------------
277 Bool_t AliFlowEventCuts::IsSelected(TObject* obj, TObject* objmc)
280 AliVEvent* vevent = dynamic_cast<AliVEvent*>(obj);
281 AliMCEvent* mcevent = dynamic_cast<AliMCEvent*>(objmc);
282 if (vevent) return PassesCuts(vevent,mcevent);;
283 return kFALSE; //when passed wrong type of object
285 //-----------------------------------------------------------------------
286 Bool_t AliFlowEventCuts::PassesCuts(AliVEvent *event, AliMCEvent *mcevent)
288 ///check if event passes cuts
289 const AliVVertex* pvtx=event->GetPrimaryVertex();
290 Double_t pvtxx = pvtx->GetX();
291 Double_t pvtxy = pvtx->GetY();
292 Double_t pvtxz = pvtx->GetZ();
293 Int_t ncontrib = pvtx->GetNContributors();
295 AliESDEvent* esdevent = dynamic_cast<AliESDEvent*>(event);
296 AliAODEvent* aodevent = dynamic_cast<AliAODEvent*>(event);
298 Int_t multGlobal = 0;
299 // these estimates only work for esd's
300 multTPC = fStandardTPCcuts->Count(event);
301 multGlobal = fStandardGlobalCuts->Count(event);
303 if ( fCutTPCmultiplicityOutliers && esdevent )
305 //this is pretty slow as we check the event track by track twice
306 //this cut will work for 2010 PbPb data and is dependent on
307 //TPC and ITS reco efficiency (e.g. geometry, calibration etc)
308 if (multTPC > ( 23+1.216*multGlobal)) {pass=kFALSE;}
309 if (multTPC < (-20+1.087*multGlobal)) {pass=kFALSE;}
312 if(fCutTPCmultiplicityOutliersAOD && aodevent) {
313 //similar (slow) cut for aod's. will work for both 2010 and 2010 pbpb data.
314 //this should be moved to AliFlowTrackCuts::Count()
315 //but at this moment the flow track cuts does not know the data that is passed
316 Int_t nTracks(aodevent->GetNumberOfTracks());
317 for(Int_t iTracks = 0; iTracks < nTracks; iTracks++) {
318 AliAODTrack* track = aodevent->GetTrack(iTracks);
320 if (!track || track->Pt() < .2 || track->Pt() > 5.0 || TMath::Abs(track->Eta()) > .8 || track->GetTPCNcls() < 70 || !track->GetDetPid() || track->GetDetPid()->GetTPCsignal() < 10.0) continue; // general quality cut
321 if (track->TestFilterBit(1) && track->Chi2perNDF() > 0.2) multTPC++;
322 if (!track->TestFilterBit(16) || track->Chi2perNDF() < 0.1) continue;
323 Double_t b[2] = {-99., -99.};
324 Double_t bCov[3] = {-99., -99., -99.};
325 AliAODTrack copy(*track);
326 if (copy.PropagateToDCA(event->GetPrimaryVertex(), event->GetMagneticField(), 100., b, bCov) && TMath::Abs(b[0]) < 0.3 && TMath::Abs(b[1]) < 0.3) multGlobal++;
328 if(!fData2011 && (multTPC < (-40.3+1.22*multGlobal) || multTPC > (32.1+1.59*multGlobal))) pass = kFALSE;
329 if(fData2011 && (multTPC < (-36.73 + 1.48*multGlobal) || multTPC > (62.87 + 1.78*multGlobal))) pass = kFALSE;
334 QAbefore(0)->Fill(pvtxz);
335 QAbefore(1)->Fill(multGlobal,multTPC);
338 if (fCutNContributors)
340 if (ncontrib < fNContributorsMin || ncontrib >= fNContributorsMax) pass=kFALSE;
342 if (fCutPrimaryVertexX)
344 if (pvtxx < fPrimaryVertexXmin || pvtxx >= fPrimaryVertexXmax) pass=kFALSE;
346 if (fCutPrimaryVertexY)
348 if (pvtxy < fPrimaryVertexYmin || pvtxy >= fPrimaryVertexYmax) pass=kFALSE;
350 if (fCutPrimaryVertexZ)
352 if (pvtxz < fPrimaryVertexZmin || pvtxz >= fPrimaryVertexZmax)
355 if (fCutCentralityPercentile&&esdevent)
357 AliCentrality* centr = esdevent->GetCentrality();
358 if (fUseCentralityUnchecked)
360 if (!centr->IsEventInCentralityClassUnchecked( fCentralityPercentileMin,
361 fCentralityPercentileMax,
362 CentrMethName(fCentralityPercentileMethod) ))
369 if (!centr->IsEventInCentralityClass( fCentralityPercentileMin,
370 fCentralityPercentileMax,
371 CentrMethName(fCentralityPercentileMethod) ))
377 if (fCutSPDvertexerAnomaly&&esdevent)
379 const AliESDVertex* sdpvertex = esdevent->GetPrimaryVertexSPD();
380 if (sdpvertex->GetNContributors()<1) pass=kFALSE;
381 if (sdpvertex->GetDispersion()>0.04) pass=kFALSE;
382 if (sdpvertex->GetZRes()>0.25) pass=kFALSE;
383 const AliESDVertex* tpcvertex = esdevent->GetPrimaryVertexTPC();
384 if (tpcvertex->GetNContributors()<1) pass=kFALSE;
385 const AliMultiplicity* tracklets = esdevent->GetMultiplicity();
386 if (tpcvertex->GetNContributors()<(-10.0+0.25*tracklets->GetNumberOfITSClusters(0)))
391 if (fCutZDCtiming&&esdevent)
393 if (!fTrigAna.ZDCTimeTrigger(esdevent))
398 if(fCutNumberOfTracks) {if ( event->GetNumberOfTracks() < fNumberOfTracksMin ||
399 event->GetNumberOfTracks() >= fNumberOfTracksMax ) pass=kFALSE;}
400 if((fCutRefMult&&mcevent)||(fCutRefMult&&esdevent))
402 //reference multiplicity still to be defined
403 Double_t refMult = RefMult(event,mcevent);
404 if (refMult < fRefMultMin || refMult >= fRefMultMax )
413 Double_t tVtxZ = aodevent->GetPrimaryVertex()->GetZ();
414 Double_t tSPDVtxZ = aodevent->GetPrimaryVertexSPD()->GetZ();
415 if( TMath::Abs(tVtxZ-tSPDVtxZ) > 0.5 ) pass = kFALSE;
417 AliCentrality* centr = aodevent->GetHeader()->GetCentralityP();
418 if(fCutTPCmultiplicityOutliers){
419 Double_t v0Centr = centr->GetCentralityPercentile("V0M");
420 Double_t trkCentr = centr->GetCentralityPercentile("TRK");
421 if(TMath::Abs(v0Centr-trkCentr) > 5) pass = kFALSE;
423 if (fCutCentralityPercentile) {
424 if (fUseCentralityUnchecked) {
425 if (!centr->IsEventInCentralityClassUnchecked( fCentralityPercentileMin,
426 fCentralityPercentileMax,
427 CentrMethName(fCentralityPercentileMethod) )) {
431 if (!centr->IsEventInCentralityClass( fCentralityPercentileMin,
432 fCentralityPercentileMax,
433 CentrMethName(fCentralityPercentileMethod) )) {
443 Int_t ntracks=event->GetNumberOfTracks();
445 for (Int_t i=0; i<ntracks; i++)
447 AliVParticle* track = event->GetTrack(i);
448 if (!track) continue;
449 Bool_t localpass=kTRUE;
450 if (fMeanPtCuts) localpass=fMeanPtCuts->IsSelected(track);
453 meanpt += track->Pt();
457 if (nselected) meanpt=meanpt/nselected;
458 if (meanpt<fMeanPtMin || meanpt >= fMeanPtMax) pass=kFALSE;
461 //impact parameter cut
462 if(fCutImpactParameter) {
463 Double_t gImpactParameter = 0.;
465 AliCollisionGeometry* headerH = dynamic_cast<AliCollisionGeometry*>(dynamic_cast<AliMCEvent*>(mcevent)->GenEventHeader());
467 gImpactParameter = headerH->ImpactParameter();
469 if ((gImpactParameter < fImpactParameterMin) || (gImpactParameter >= fImpactParameterMax ))
475 QAafter(1)->Fill(multGlobal,multTPC);
476 QAafter(0)->Fill(pvtxz);
481 //-----------------------------------------------------------------------
482 Float_t AliFlowEventCuts::GetCentrality(AliVEvent* event, AliMCEvent* /*mcEvent*/)
484 //get the centrality percentile of the event
485 AliESDEvent* esdEvent = dynamic_cast<AliESDEvent*>(event);
486 AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(event);
488 Float_t centrality=-1.;
490 AliCentrality* centr = NULL;
492 centr = esdEvent->GetCentrality();
494 centr = aodEvent->GetHeader()->GetCentralityP();
496 if (!centr) return -1.;
498 if (fUseCentralityUnchecked)
499 centrality=centr->GetCentralityPercentileUnchecked(CentrMethName(fCentralityPercentileMethod));
501 centrality=centr->GetCentralityPercentile(CentrMethName(fCentralityPercentileMethod));
506 //-----------------------------------------------------------------------
507 const char* AliFlowEventCuts::CentrMethName(refMultMethod method) const
509 //get the string for refmultmethod, for use with AliCentrality in
510 //the cut on centrality percentile
525 //-----------------------------------------------------------------------
526 AliFlowEventCuts* AliFlowEventCuts::StandardCuts()
528 //make a set of standard event cuts, caller becomes owner
529 AliFlowEventCuts* cuts = new AliFlowEventCuts();
533 //-----------------------------------------------------------------------
534 Int_t AliFlowEventCuts::RefMult(AliVEvent* event, AliMCEvent *mcEvent)
536 //calculate the reference multiplicity, if all fails return 0
537 AliESDVZERO* vzero = NULL;
538 AliESDEvent* esdevent = dynamic_cast<AliESDEvent*>(event);
540 if (fUseAliESDtrackCutsRefMult && esdevent)
542 //use the standard ALICE reference multiplicity with the default eta range
543 return AliESDtrackCuts::GetReferenceMultiplicity(esdevent, fRefMultMethodAliESDtrackCuts);
546 if (fRefMultMethod==kTPConly && !fRefMultCuts)
548 fRefMultCuts = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts();
549 fRefMultCuts->SetEtaRange(-0.8,0.8);
550 fRefMultCuts->SetPtMin(0.15);
552 else if (fRefMultMethod==kSPDtracklets && !fRefMultCuts)
554 fRefMultCuts = new AliFlowTrackCuts("tracklet refmult cuts");
555 fRefMultCuts->SetParamType(AliFlowTrackCuts::kSPDtracklet);
556 fRefMultCuts->SetEtaRange(-0.8,0.8);
558 else if (fRefMultMethod==kVZERO)
560 if (!esdevent) return 0;
561 vzero=esdevent->GetVZEROData();
562 if (!vzero) return 0;
563 return TMath::Nint(vzero->GetMTotV0A()+vzero->GetMTotV0C());
565 else if (fRefMultMethod==kSPD1clusters)
567 if (!esdevent) return 0;
568 const AliMultiplicity* mult = esdevent->GetMultiplicity();
570 return mult->GetNumberOfITSClusters(1);
574 fRefMultCuts->SetEvent(event,mcEvent);
575 Int_t numberOfInputObjects = fRefMultCuts->GetNumberOfInputObjects();
576 for (Int_t i=0; i<numberOfInputObjects; i++) {
577 if (fRefMultCuts->IsSelected(fRefMultCuts->GetInputObject(i),i))
582 //_____________________________________________________________________________
583 void AliFlowEventCuts::DefineHistograms()
588 Bool_t adddirstatus = TH1::AddDirectoryStatus();
589 TH1::AddDirectory(kFALSE);
590 fQA = new TList(); fQA->SetOwner();
591 fQA->SetName(Form("%s QA",GetName()));
592 TList* before = new TList(); before->SetOwner();
593 before->SetName("before");
594 TList* after = new TList(); after->SetOwner();
595 after->SetName("after");
598 before->Add(new TH1F("zvertex",";z;event cout",500,-15.,15.)); //0
599 after->Add(new TH1F("zvertex",";z;event cout",500,-15.,15.)); //0
600 before->Add(new TH2F("fTPCvsGlobalMult","TPC only vs Global track multiplicity;global;TPC only",500,0,2500,500,0,3500));//1
601 after->Add(new TH2F("fTPCvsGlobalMult","TPC only vs Global track multiplicity;global;TPC only",500,0,2500,500,0,3500));//1
602 TH1::AddDirectory(adddirstatus);
605 //---------------------------------------------------------------//
606 void AliFlowEventCuts::Browse(TBrowser* b)
608 //some browsing capabilities
609 if (fQA) b->Add(fQA);
612 //---------------------------------------------------------------//
613 Long64_t AliFlowEventCuts::Merge(TCollection* list)
617 AliFlowEventCuts* obj;
619 if (list->GetEntries()<1) return 0;
621 while ( (obj = dynamic_cast<AliFlowEventCuts*>(next())) )
623 if (obj==this) continue;
625 listwrapper.Add(obj->GetQA());
626 fQA->Merge(&listwrapper);