Support for HFE analysis added (Matthias)
[u/mrichter/AliRoot.git] / PWGHF / correlationHF / AliHFAssociatedTrackCuts.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-2010, 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 // Base class for cuts on Associated tracks for HF Correlation analysis
21 //
22 // Author: S.Bjelogrlic (Utrecht) sandro.bjelogrlic@cern.ch
23 ////////////////////////////////////////////////////////////////////////
24 #include <Riostream.h>
25 #include "AliHFAssociatedTrackCuts.h"
26 #include "AliAODPidHF.h"
27 #include "AliESDtrackCuts.h"
28 #include "AliESDtrack.h"
29 #include "AliAODv0.h"
30 #include "AliAODVertex.h"
31 #include "AliAODMCParticle.h"
32 #include "AliAnalysisManager.h"
33 #include "AliInputEventHandler.h"
34 #include "TString.h"
35
36 using std::cout;
37 using std::endl;
38
39 ClassImp(AliHFAssociatedTrackCuts)
40
41 //--------------------------------------------------------------------------
42 AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts():
43 AliAnalysisCuts(),
44 fESDTrackCuts(0),
45 fPidObj(0),
46
47 fPoolMaxNEvents(0), 
48 fPoolMinNTracks(0), 
49 fMinEventsToMix(0),
50 fNzVtxBins(0), 
51 fNzVtxBinsDim(0), 
52 fZvtxBins(0), 
53 fNCentBins(0), 
54 fNCentBinsDim(0), 
55 fCentBins(0),
56
57 fNofMCEventType(0),
58 fMCEventType(0),
59
60 fNTrackCuts(0),
61 fAODTrackCuts(0),
62 fTrackCutsNames(0),
63 fNvZeroCuts(0),
64 fAODvZeroCuts(0),
65 fvZeroCutsNames(0),
66 fBit(0),
67 fCharge(0),
68 fDescription("")
69
70 {
71         //
72         //default constructor
73         //
74         //
75         //default constructor
76         //
77         
78 }
79
80 //--------------------------------------------------------------------------
81 AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(const char* name, const char* title):
82 AliAnalysisCuts(name,title),
83 fESDTrackCuts(0),
84 fPidObj(0),
85
86 fPoolMaxNEvents(0), 
87 fPoolMinNTracks(0), 
88 fMinEventsToMix(0),
89 fNzVtxBins(0), 
90 fNzVtxBinsDim(0), 
91 fZvtxBins(0), 
92 fNCentBins(0), 
93 fNCentBinsDim(0), 
94 fCentBins(0),
95
96 fNofMCEventType(0),
97 fMCEventType(0),
98
99 fNTrackCuts(0),
100 fAODTrackCuts(0),
101 fTrackCutsNames(0),
102 fNvZeroCuts(0),
103 fAODvZeroCuts(0),
104 fvZeroCutsNames(0),
105 fBit(0),
106 fCharge(0),
107 fDescription("")
108
109 {
110         //
111         //default constructor
112         //
113         
114 }
115 //--------------------------------------------------------------------------
116 AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(const AliHFAssociatedTrackCuts &source) :
117 AliAnalysisCuts(source),
118 fESDTrackCuts(source.fESDTrackCuts),
119 fPidObj(source.fPidObj),
120
121 fPoolMaxNEvents(source.fPoolMaxNEvents), 
122 fPoolMinNTracks(source.fPoolMinNTracks), 
123 fMinEventsToMix(source.fMinEventsToMix),
124 fNzVtxBins(source.fNzVtxBins), 
125 fNzVtxBinsDim(source.fNzVtxBinsDim), 
126 fZvtxBins(source.fZvtxBins), 
127 fNCentBins(source.fNCentBins), 
128 fNCentBinsDim(source.fNCentBinsDim), 
129 fCentBins(source.fCentBins),
130
131 fNofMCEventType(source.fNofMCEventType),
132 fMCEventType(source.fMCEventType),
133
134 fNTrackCuts(source.fNTrackCuts),
135 fAODTrackCuts(source.fAODTrackCuts),
136 fTrackCutsNames(source.fTrackCutsNames),
137 fNvZeroCuts(source.fNvZeroCuts),
138 fAODvZeroCuts(source.fAODvZeroCuts),
139 fvZeroCutsNames(source.fvZeroCutsNames),
140 fBit(source.fBit),
141 fCharge(source.fCharge),
142 fDescription(source.fDescription)
143 {
144         //
145         // copy constructor
146         //
147         
148
149         AliInfo("AliHFAssociatedTrackCuts::Copy constructor ");
150         if(source.fESDTrackCuts) AddTrackCuts(source.fESDTrackCuts);
151         if(source.fAODTrackCuts) SetAODTrackCuts(source.fAODTrackCuts);
152         if(source.fAODvZeroCuts) SetAODvZeroCuts(source.fAODvZeroCuts);
153         if(source.fPidObj) SetPidHF(source.fPidObj);
154 }
155 //--------------------------------------------------------------------------
156 AliHFAssociatedTrackCuts &AliHFAssociatedTrackCuts::operator=(const AliHFAssociatedTrackCuts &source)
157 {
158         //
159         // assignment operator
160         //      
161         if(&source == this) return *this;
162         
163         AliAnalysisCuts::operator=(source);
164         fESDTrackCuts=source.fESDTrackCuts;
165         fPidObj=source.fPidObj;
166         fNTrackCuts=source.fNTrackCuts;
167         fAODTrackCuts=source.fAODTrackCuts;
168         fTrackCutsNames=source.fTrackCutsNames;
169         fNvZeroCuts=source.fNvZeroCuts;
170         fAODvZeroCuts=source.fAODvZeroCuts;
171         fvZeroCutsNames=source.fvZeroCutsNames;
172         fBit=source.fBit;
173         fCharge=source.fCharge;
174         
175         return *this;
176
177 }
178
179
180 //--------------------------------------------------------------------------
181 AliHFAssociatedTrackCuts::~AliHFAssociatedTrackCuts()
182 {
183         if(fESDTrackCuts) {delete fESDTrackCuts; fESDTrackCuts = 0;}
184         if(fPidObj) {delete fPidObj; fPidObj = 0;}
185         if(fZvtxBins) {delete[] fZvtxBins; fZvtxBins=0;} 
186         if(fCentBins) {delete[] fCentBins; fCentBins=0;}
187         if(fAODTrackCuts) {delete[] fAODTrackCuts; fAODTrackCuts=0;}
188         if(fTrackCutsNames) {delete[] fTrackCutsNames; fTrackCutsNames=0;}
189         if(fAODvZeroCuts){delete[] fAODvZeroCuts; fAODvZeroCuts=0;}
190         if(fvZeroCutsNames) {delete[] fvZeroCutsNames; fvZeroCutsNames=0;}
191
192         
193 }
194 //--------------------------------------------------------------------------
195 Bool_t AliHFAssociatedTrackCuts::IsInAcceptance()
196 {
197         printf("Careful: method AliHFAssociatedTrackCuts::IsInAcceptance is not implemented yet \n");
198         return kFALSE;
199 }
200 //--------------------------------------------------------------------------
201 Bool_t AliHFAssociatedTrackCuts::IsHadronSelected(AliAODTrack * track)
202 {
203         AliESDtrack esdtrack(track);
204         if(!fESDTrackCuts->IsSelected(&esdtrack)) return kFALSE;
205         
206         if(fBit && !track->TestFilterBit(fBit)) return kFALSE; // check the filter bit
207  
208         return kTRUE;
209         
210 }
211
212 //--------------------------------------------------------------------------
213 Bool_t AliHFAssociatedTrackCuts::CheckHadronKinematic(Double_t pt, Double_t d0) 
214 {
215         
216         
217         
218         if(pt < fAODTrackCuts[0]) return kFALSE;
219         if(pt > fAODTrackCuts[1]) return kFALSE;
220         if(d0 < fAODTrackCuts[2]) return kFALSE;
221         if(d0 > fAODTrackCuts[3]) return kFALSE;
222         
223         return kTRUE;
224
225         
226 }
227 //--------------------------------------------------------------------------
228
229 Bool_t AliHFAssociatedTrackCuts::Charge(Short_t charge, AliAODTrack* track) 
230 {// charge is the charge to compare to (for example, a daughter of a D meson)
231         
232         if(!fCharge) return kTRUE; // if fCharge is set to 0 (no selection on the charge), returns always true
233         if(track->Charge()!= fCharge*charge) return kFALSE;
234         return kTRUE;
235 }
236
237 //--------------------------------------------------------------------------
238 Bool_t AliHFAssociatedTrackCuts::CheckKaonCompatibility(AliAODTrack * track, Bool_t useMc, TClonesArray* mcArray, Int_t method)
239 {
240         Bool_t isKaon = kFALSE;
241         
242         if(useMc) { // on MC
243                 Int_t hadLabel = track->GetLabel();
244                 if(hadLabel < 0) return kFALSE;
245                 AliAODMCParticle* hadron = dynamic_cast<AliAODMCParticle*>(mcArray->At(hadLabel)); 
246                 if(hadron){
247                   Int_t pdg = TMath::Abs(hadron->GetPdgCode()); 
248                   if (pdg == 321) isKaon = kTRUE;
249                 }
250         }
251         
252         if(!useMc) { // on DATA
253           switch(method) {
254           case(1): {
255                 Bool_t isKTPC=kFALSE;
256                 Bool_t isPiTPC=kFALSE;
257                 Bool_t isPTPC=kFALSE;
258                 Bool_t isKTOF=kFALSE;
259                 Bool_t isPiTOF=kFALSE;
260                 Bool_t isPTOF=kFALSE;
261                 
262                 Bool_t KaonHyp = kFALSE;
263                 Bool_t PionHyp = kFALSE;
264                 Bool_t ProtonHyp = kFALSE;
265                 
266                 if(fPidObj->CheckStatus(track,"TOF")) {
267                         isKTOF=fPidObj->IsKaonRaw(track,"TOF");
268                         isPiTOF=fPidObj->IsPionRaw(track,"TOF");
269                         isPTOF=fPidObj->IsProtonRaw(track,"TOF");
270                 }
271                 if(fPidObj->CheckStatus(track,"TPC")){
272                         isKTPC=fPidObj->IsKaonRaw(track,"TPC");
273                         isPiTPC=fPidObj->IsPionRaw(track,"TPC");
274                         isPTPC=fPidObj->IsProtonRaw(track,"TPC");               
275                 }
276                 
277                 if (isKTOF && isKTPC) KaonHyp = kTRUE;
278                 if (isPiTOF && isPiTPC) PionHyp = kTRUE;
279                 if (isPTOF && isPTPC) ProtonHyp = kTRUE;
280                 
281                 if(KaonHyp && !PionHyp && !ProtonHyp) isKaon = kTRUE; 
282                 break;
283               }
284       case(2): {
285                 if(fPidObj->MakeRawPid(track,3)>=1) isKaon = kTRUE;
286                 break;
287               }
288       }
289         }
290         
291         return isKaon;
292         
293 }
294 //--------------------------------------------------------------------------
295 Bool_t AliHFAssociatedTrackCuts::IsKZeroSelected(AliAODv0 *vzero, AliAODVertex *vtx1)
296 {
297         
298         if(vzero->DcaV0Daughters()>fAODvZeroCuts[0]) return kFALSE;
299         if(vzero->Chi2V0()>fAODvZeroCuts[1]) return kFALSE;
300         if(vzero->DecayLength(vtx1) < fAODvZeroCuts[2]) return kFALSE;
301         if(vzero->DecayLength(vtx1) > fAODvZeroCuts[3]) return kFALSE;
302         if(vzero->OpenAngleV0() > fAODvZeroCuts[4]) return kFALSE;
303         if(vzero->Pt() < fAODvZeroCuts[5]) return kFALSE;
304         if(TMath::Abs(vzero->Eta()) > fAODvZeroCuts[6]) return kFALSE;
305
306         
307         return kTRUE;
308 }
309 //--------------------------------------------------------------------------
310 Bool_t *AliHFAssociatedTrackCuts::IsMCpartFromHF(Int_t label, TClonesArray*mcArray){
311   // Check origin in MC
312
313   AliAODMCParticle* mcParticle;
314   Int_t pdgCode = -1;
315         
316   Bool_t isCharmy = kFALSE;
317   Bool_t isBeauty = kFALSE;
318   Bool_t isD = kFALSE;
319   Bool_t isB = kFALSE;
320     
321      Bool_t *originvect = new Bool_t[4];
322     
323     originvect[0] = kFALSE;
324         originvect[1] = kFALSE;
325         originvect[2] = kFALSE;
326         originvect[3] = kFALSE;
327
328         if (label<0) return originvect;
329   
330         while(pdgCode!=2212){ // loops back to the collision to check the particle source
331
332     mcParticle = dynamic_cast<AliAODMCParticle*>(mcArray->At(label));
333     if(!mcParticle) {AliError("NO MC PARTICLE"); break;}
334     pdgCode =  TMath::Abs(mcParticle->GetPdgCode());
335
336     label = mcParticle->GetMother();
337
338
339     if((pdgCode>=400 && pdgCode <500) || (pdgCode>=4000 && pdgCode<5000 )) isD = kTRUE;
340     if((pdgCode>=500 && pdgCode <600) || (pdgCode>=5000 && pdgCode<6000 )) {isD = kFALSE; isB = kTRUE;}
341
342
343     if(pdgCode == 4) isCharmy = kTRUE;
344     if(pdgCode == 5) {isBeauty = kTRUE; isCharmy = kFALSE;}
345     if(label<0) break;
346
347   }
348
349         
350         originvect[0] = isCharmy;
351         originvect[1] = isBeauty;
352         originvect[2] = isD;
353         originvect[3] = isB;
354  
355     
356   return originvect;
357 }
358
359 //--------------------------------------------------------------------------
360 Bool_t AliHFAssociatedTrackCuts::InvMassDstarRejection(AliAODRecoDecayHF2Prong* d, AliAODTrack *track, Int_t hypD0) const {
361         //
362         // Calculates invmass of track+D0 and rejects if compatible with D*
363         // (to remove pions from D*)
364         // 
365         Double_t nsigma = 3.;
366         
367         Double_t mD0, mD0bar;
368         d->InvMassD0(mD0,mD0bar);
369         
370         Double_t invmassDstar1 = 0, invmassDstar2 = 0; 
371         Double_t e1Pi = d->EProng(0,211), e2K = d->EProng(1,321); //hyp 1 (pi,K) - D0
372         Double_t e1K = d->EProng(0,321), e2Pi = d->EProng(1,211); //hyp 2 (K,pi) - D0bar
373         Double_t psum2 = (d->Px()+track->Px())*(d->Px()+track->Px())
374         +(d->Py()+track->Py())*(d->Py()+track->Py())
375         +(d->Pz()+track->Pz())*(d->Pz()+track->Pz());
376         
377         switch(hypD0) {
378                 case 1:
379                         invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2);
380                         if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;
381                         break;
382                 case 2:
383                         invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2);
384                         if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;
385                         break;
386                 case 3:
387                         invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2);
388                         invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2);
389                         if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;
390                         if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;
391                         break;
392         }
393         
394         return kTRUE;
395 }
396 //________________________________________________________
397 void AliHFAssociatedTrackCuts::SetMCEventTypes(Int_t *MCEventTypeArray)
398 // set the array of event types you want to process in MonteCarlo (gluon splitting, pair production etc.)
399 {
400         if(!fMCEventType) fMCEventType = new Int_t[fNofMCEventType];
401         
402         for(Int_t k=0; k<fNofMCEventType; k++){
403                 fMCEventType[k] = MCEventTypeArray[k];
404         }
405         return; 
406 }
407
408 //________________________________________________________
409 void AliHFAssociatedTrackCuts::SetAODTrackCuts(Float_t *cutsarray)
410 {
411         if(!fAODTrackCuts) fAODTrackCuts = new Float_t[fNTrackCuts];
412         for(Int_t i =0; i<fNTrackCuts; i++){
413                 fAODTrackCuts[i] = cutsarray[i];
414         }
415         SetTrackCutsNames();
416         return;
417 }
418 //________________________________________________________
419 void AliHFAssociatedTrackCuts::SetTrackCutsNames(/*TString *namearray*/){
420         
421         fTrackCutsNames = new TString[4];
422         fTrackCutsNames[0]= "associated track:: pt min [GeV/c]................: ";
423         fTrackCutsNames[1]= "associated track:: pt max [GeV/c]................: ";
424         fTrackCutsNames[2]= "associated track:: d0 min [cm]...................: ";
425         fTrackCutsNames[3]= "associated track:: d0 max [cm]...................: ";
426         
427
428         
429         return;
430 }
431 //--------------------------------------------------------------------------
432 void AliHFAssociatedTrackCuts::SetAODvZeroCuts(Float_t *cutsarray)
433 {
434         
435
436         if(!fAODvZeroCuts) fAODvZeroCuts = new Float_t[fNvZeroCuts];
437         for(Int_t i =0; i<fNvZeroCuts; i++){
438                 fAODvZeroCuts[i] = cutsarray[i];
439         }
440         SetvZeroCutsNames();
441         return;
442 }
443 //--------------------------------------------------------------------------
444 void AliHFAssociatedTrackCuts::SetvZeroCutsNames(/*TString *namearray*/){
445         
446         fvZeroCutsNames = new TString[7];
447         fvZeroCutsNames[0] = "vZero:: max DCA between two daughters [cm].......: ";
448         fvZeroCutsNames[1] = "vZero:: max fit Chi Square between two daughters.: ";
449         fvZeroCutsNames[2] = "vZero:: min decay length [cm]....................: ";
450         fvZeroCutsNames[3] = "vZero:: max decay length [cm]....................: ";
451         fvZeroCutsNames[4] = "vZero:: max opening angle between daughters [rad]: ";
452         fvZeroCutsNames[5] = "vZero:: pt min [Gev/c]...........................: ";
453         fvZeroCutsNames[6] = "vZero:: |Eta| range <............................: ";
454         
455         
456         return;
457 }
458
459 //--------------------------------------------------------------------------
460 void AliHFAssociatedTrackCuts::SetPidAssociated()
461 {
462   //setting PidResponse
463   if(fPidObj->GetOldPid()==kFALSE && fPidObj->GetPidResponse()==0x0){
464     AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
465     AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
466     AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
467     fPidObj->SetPidResponse(pidResp);
468   }
469 }
470
471 //--------------------------------------------------------------------------
472 void AliHFAssociatedTrackCuts::PrintAll()
473 {
474         
475         printf("\n=================================================");
476         if(fESDTrackCuts){
477           printf("\nCuts for the associated track: \n \n");
478           
479           printf("ITS Refit........................................: %s\n",fESDTrackCuts->GetRequireITSRefit() ? "Yes" : "No");
480           printf("TPC Refit........................................: %s\n",fESDTrackCuts->GetRequireTPCRefit() ? "Yes" : "No");
481           printf("ITS SA...........................................: %s\n",fESDTrackCuts->GetRequireITSStandAlone() ? "Yes" : "No");
482           printf("TPC SA...........................................: %s\n",fESDTrackCuts->GetRequireTPCStandAlone() ? "Yes" : "No");
483           printf("Min number of ITS clusters.......................: %d\n",fESDTrackCuts->GetMinNClustersITS());
484           printf("Min number of TPC clusters.......................: %d\n",fESDTrackCuts->GetMinNClusterTPC());
485           Int_t spd = fESDTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD);
486           if(spd==0) std::cout <<  "SPD..............................................: kOff"  << std::endl;
487           if(spd==1) std::cout <<  "SPD..............................................: kNone"  << std::endl;
488           if(spd==2) std::cout <<  "SPD..............................................: kAny"  << std::endl;
489           if(spd==3) std::cout <<  "SPD..............................................: kFirst"  << std::endl;
490           if(spd==4) std::cout <<  "SPD..............................................: kOnlyFirst"  << std::endl;
491           if(spd==5) std::cout <<  "SPD..............................................: kSecond"  << std::endl;
492           if(spd==6) std::cout <<  "SPD..............................................: kOnlySecond"  << std::endl;
493           if(spd==7) std::cout <<  "SPD..............................................: kBoth"  << std::endl;
494         }
495         else printf("\nNo Cuts for Associated Tracks\n");
496         std::cout <<  "Filter Bit.......................................: " << fBit  << std::endl;
497         std::cout <<  "Charge...........................................: " << fCharge  << std::endl;
498         
499         if(fAODTrackCuts){
500           for(Int_t j=0;j<fNTrackCuts;j++){
501             std::cout << fTrackCutsNames[j] << fAODTrackCuts[j] << std::endl;
502           }
503         }
504
505         if(fAODvZeroCuts){
506           printf("\n");
507           printf("=================================================");
508           printf("\nCuts for the K0 candidates: \n \n");
509           for(Int_t k=0;k<fNvZeroCuts;k++){
510             std::cout << fvZeroCutsNames[k] <<  fAODvZeroCuts[k] << std::endl;
511           }
512         }
513         else printf("\nNo Cuts for the K0 candidates\n");
514         std::cout << " " << std::endl;
515         PrintPoolParameters();
516         PrintSelectedMCevents();
517         if(fDescription){
518           printf("=================================================");
519           printf("\nAdditional description\n");
520           std::cout << fDescription << std::endl;
521           printf("\n");
522         }
523
524 }
525
526 //--------------------------------------------------------------------------
527 void AliHFAssociatedTrackCuts::PrintPoolParameters()
528 {   
529         printf("=================================================");
530         printf("\nEvent Pool settings: \n \n");
531         
532         printf("Number of zVtx Bins: %d\n", fNzVtxBins);
533         printf("\nzVtx Bins:\n");
534         //Double_t zVtxbinLims[fNzVtxBins+1] = fNzVtxBins;
535         for(Int_t k=0; k<fNzVtxBins; k++){
536                 printf("Bin %d..............................................: %.1f - %.1f cm\n", k, fZvtxBins[k], fZvtxBins[k+1]);      
537         }
538         printf("\n");
539         printf("\nNumber of Centrality(multiplicity) Bins: %d\n", fNCentBins);
540         printf("\nCentrality(multiplicity) Bins:\n");
541         for(Int_t k=0; k<fNCentBins; k++){
542                 printf("Bin %d..............................................: %.1f - %.1f\n", k, fCentBins[k], fCentBins[k+1]);
543         }
544
545         
546         
547 }
548
549 //--------------------------------------------------------------------------
550 void AliHFAssociatedTrackCuts::PrintSelectedMCevents()
551 {
552         printf("\n=================================================");
553         
554         printf("\nSelected MC events: \n \n");
555         printf("Number of selected events: %d\n",fNofMCEventType);
556         
557         for(Int_t k=0; k<fNofMCEventType; k++){
558         if(fMCEventType[k]==28) printf("=> Flavour excitation \n");     
559         if(fMCEventType[k]==53) printf("=> Pair creation \n");  
560         if(fMCEventType[k]==68) printf("=> Gluon splitting \n");        
561         }
562         
563         printf("\n");
564         for(Int_t k=0; k<fNofMCEventType; k++){
565                 printf("MC process code %d \n",fMCEventType[k]);                
566         }
567         
568         printf("\n");
569         
570         
571                 
572         
573 }
574
575