[u/mrichter/AliRoot.git] / CORRFW / test / AliCFRsnTask.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

//-----------------------------------------------------------------------
// Author : R. Vernet, Consorzio Cometa - Catania (it)
//-----------------------------------------------------------------------


#include <TROOT.h>
#include <TInterpreter.h>

#include "AliCFRsnTask.h"
#include "TCanvas.h"
#include "AliStack.h"
#include "TParticle.h"
#include "TH1I.h"
#include "TChain.h"
#include "AliMCEventHandler.h"
#include "AliMCEvent.h"
#include "AliAnalysisManager.h"
#include "AliESDEvent.h"
#include "AliCFManager.h"
#include "AliCFCutBase.h"
#include "AliCFContainer.h"
#include "TChain.h"
#include "AliESDtrack.h"
#include "AliLog.h"
#include "AliRsnDaughter.h"
#include "AliCFPair.h"
#include <memory>

//__________________________________________________________________________
AliCFRsnTask::AliCFRsnTask() :
  fRsnPDG(313),
  fChain(0x0),
  fESD(0x0),
  fCFManager(0x0),
  fHistEventsProcessed(0x0)
{
  //Defual ctor
}
//___________________________________________________________________________
AliCFRsnTask::AliCFRsnTask(const Char_t* name) :
  AliAnalysisTask(name,"AliCFRsnTask"),
  fRsnPDG(313),
  fChain(0x0),
  fESD(0x0),
  fCFManager(0x0),
  fHistEventsProcessed(0x0)
{
  //
  // Constructor. Initialization of Inputs and Outputs
  //
  Info("AliCFRsnTask","Calling Constructor");
  DefineInput (0,TChain::Class());
  DefineOutput(0,TH1I::Class());
  DefineOutput(1,AliCFContainer::Class());
  //   DefineOutput(2,TList::Class());
}

//___________________________________________________________________________
AliCFRsnTask& AliCFRsnTask::operator=(const AliCFRsnTask& c) 
{
  //
  // Assignment operator
  //
  if (this!=&c) {
    AliAnalysisTask::operator=(c) ;
    fRsnPDG     = c.fRsnPDG;
    fChain      = c.fChain;
    fESD        = c.fESD;
    fCFManager  = c.fCFManager;
    fHistEventsProcessed = c.fHistEventsProcessed;
  }
  return *this;
}

//___________________________________________________________________________
AliCFRsnTask::AliCFRsnTask(const AliCFRsnTask& c) :
  AliAnalysisTask(c),
  fRsnPDG(c.fRsnPDG),
  fChain(c.fChain),
  fESD(c.fESD),
  fCFManager(c.fCFManager),
  fHistEventsProcessed(c.fHistEventsProcessed)
{
  //
  // Copy Constructor
  //
}

//___________________________________________________________________________
AliCFRsnTask::~AliCFRsnTask() {
  //
  //destructor
  //
  Info("~AliCFRsnTask","Calling Destructor");
  if (fChain)               delete fChain ;
  if (fESD)                 delete fESD ;
  if (fCFManager)           delete fCFManager ;
  if (fHistEventsProcessed) delete fHistEventsProcessed ;
}

//___________________________________________________________________________

void AliCFRsnTask::Init()
{

}
//_________________________________________________
void AliCFRsnTask::Exec(Option_t *)
{
  //
  // Main loop function
  //
  Info("Exec","") ;
  // Get the mc truth
  AliMCEventHandler* mcTruth = (AliMCEventHandler*)((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());

  if (!mcTruth) Error("Exec","NO MC INFO FOUND... EXITING");

  // transform possible old AliESD into AliESDEvent

  if (fESD->GetAliESDOld()) fESD->CopyFromOldESD(); //transition to new ESD format

  //pass the MC evt handler to the cuts that need it 

  fCFManager->SetEventInfo(mcTruth);

  // Get the MC event 
  AliMCEvent* mcEvent = mcTruth->MCEvent();
  AliStack*   stack   = mcEvent->Stack();

  // MC-event selection
  Double_t containerInput[2] ;
        
  //loop on the MC event
  Info("Exec","Looping on MC event");
  for (Int_t ipart=0; ipart<stack->GetNprimary(); ipart++) { 
    AliMCParticle *mcPart  = mcEvent->GetTrack(ipart);

    //check the MC-level cuts
    if (!fCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,mcPart)) continue;
    containerInput[0] = mcPart->Pt();
    containerInput[1] = mcPart->Eta() ;
    //fill the container for Gen-level selection
    fCFManager->GetParticleContainer()->Fill(containerInput,kStepGenerated);
    
    //check the Acceptance-level cuts
    if (!fCFManager->CheckParticleCuts(AliCFManager::kPartAccCuts,mcPart)) continue;
    //fill the container for Acceptance-level selection
    fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructible);
  }    


  //Now go to rec level
  Info("Exec","Looping on ESD event");

  //SET THE ESD AS EVENT INFO IN RECONSTRUCTION CUTS
  TObjArray* fCutsReco = fCFManager->GetParticleCutsList(AliCFManager::kPartRecCuts);
  TObjArray* fCutsPID = fCFManager->GetParticleCutsList(AliCFManager::kPartSelCuts);
  TObjArrayIter iter1(fCutsReco);
  TObjArrayIter iter2(fCutsPID);
  AliCFCutBase *cut = 0;
  while ( (cut = (AliCFCutBase*)iter1.Next()) ) {
    cut->SetEvtInfo(fESD);
  }
  while ( (cut = (AliCFCutBase*)iter2.Next()) ) {
    cut->SetEvtInfo(fESD);
  }
  
  for (Int_t iTrack1 = 0; iTrack1<fESD->GetNumberOfTracks(); iTrack1++) {
    AliESDtrack* esdTrack1 = fESD->GetTrack(iTrack1);
    //track1 is negative
    if (esdTrack1->Charge()>=0) continue;
    Int_t esdLabel1 = esdTrack1->GetLabel();
    if (esdLabel1<0) continue;

    for (Int_t iTrack2 = 0; iTrack2<fESD->GetNumberOfTracks(); iTrack2++) {
      AliESDtrack* esdTrack2 = fESD->GetTrack(iTrack2);
      //track2 is positive
      if (esdTrack2->Charge()<=0) continue;
      Int_t esdLabel2 = esdTrack2->GetLabel();
      if (esdLabel2<0) continue;
	
      // copy ESDtrack data into RsnDaughter (and make Bayesian PID)
      //if problem with copy constructor, use the following special command to destroy the pointer when exiting scope
      //std::auto_ptr<AliRsnDaughter> track1(AliRsnDaughter::Adopt(esdTrack1,iTrack1));

      AliRsnDaughter* tmp1=AliRsnDaughter::Adopt(esdTrack1,iTrack1);
      AliRsnDaughter* tmp2=AliRsnDaughter::Adopt(esdTrack2,iTrack2);
      AliRsnDaughter track1(*tmp1);
      AliRsnDaughter track2(*tmp2);
      delete tmp1;
      delete tmp2;

      TParticle *part1 = stack->Particle(esdLabel1);
      track1.SetTruePDG(part1->GetPdgCode());
      Int_t mother1 = part1->GetFirstMother();
      track1.SetMother(mother1);
      if (mother1 >= 0) {
	TParticle *mum = stack->Particle(mother1);
	track1.SetMotherPDG(mum->GetPdgCode());
      }
      
      TParticle *part2 = stack->Particle(esdLabel2);
      track2.SetTruePDG(part2->GetPdgCode());
      Int_t mother2 = part2->GetFirstMother();
      track2.SetMother(mother2);
      if (mother2 >= 0) {
	TParticle *mum = stack->Particle(mother2);
	track2.SetMotherPDG(mum->GetPdgCode());
      }
	
      //make a mother resonance from the 2 candidate daughters
      AliRsnDaughter rsn = AliRsnDaughter::Sum(track1,track2) ;
      AliCFPair pair(esdTrack1,esdTrack2);

      //check if true resonance
      if (rsn.GetMotherPDG() != fRsnPDG) continue;
      if (!fCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,&pair)) continue;

      //check if associated MC resonance passes the cuts
      Int_t motherLabel=rsn.GetLabel();
      if (motherLabel<0) continue ;
      AliMCParticle* mcRsn = mcEvent->GetTrack(motherLabel);
      if (!mcRsn) continue;
      if (!fCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,mcRsn)) continue; 
    
      //fill the container
      containerInput[0] = rsn.GetPt() ;
      containerInput[1] = GetRapidity(rsn.GetEnergy(),rsn.GetPz());
      fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructed) ;   

      if (!fCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,&pair)) continue ;
      fCFManager->GetParticleContainer()->Fill(containerInput,kStepSelected);
    }
  }
    
  fHistEventsProcessed->Fill(0);
  PostData(0,fHistEventsProcessed) ;
  PostData(1,fCFManager->GetParticleContainer()) ;
  
  //   TList * list = new TList();
  //   fCFManager->AddQAHistosToList(list);
  //   PostData(2,list) ;
}


//___________________________________________________________________________
void AliCFRsnTask::Terminate(Option_t*)
{
  // The Terminate() function is the last function to be called during
  // a query. It always runs on the client, it can be used to present
  // the results graphically or save the results to file.

  Info("Terminate","");
  AliAnalysisTask::Terminate();


  Double_t max1 = fCFManager->GetParticleContainer()->ShowProjection(0,0)->GetMaximum();
  Double_t max2 = fCFManager->GetParticleContainer()->ShowProjection(1,0)->GetMaximum();

  fCFManager->GetParticleContainer()->ShowProjection(0,0)->GetYaxis()->SetRangeUser(0,max1*1.2);
  fCFManager->GetParticleContainer()->ShowProjection(0,1)->GetYaxis()->SetRangeUser(0,max1*1.2);
  fCFManager->GetParticleContainer()->ShowProjection(0,2)->GetYaxis()->SetRangeUser(0,max1*1.2);
  fCFManager->GetParticleContainer()->ShowProjection(0,3)->GetYaxis()->SetRangeUser(0,max1*1.2);

  fCFManager->GetParticleContainer()->ShowProjection(1,0)->GetYaxis()->SetRangeUser(0,max2*1.2);
  fCFManager->GetParticleContainer()->ShowProjection(1,1)->GetYaxis()->SetRangeUser(0,max2*1.2);
  fCFManager->GetParticleContainer()->ShowProjection(1,2)->GetYaxis()->SetRangeUser(0,max2*1.2);
  fCFManager->GetParticleContainer()->ShowProjection(1,3)->GetYaxis()->SetRangeUser(0,max2*1.2);

  fCFManager->GetParticleContainer()->ShowProjection(0,0)->SetMarkerStyle(23) ;
  fCFManager->GetParticleContainer()->ShowProjection(0,1)->SetMarkerStyle(24) ;
  fCFManager->GetParticleContainer()->ShowProjection(0,2)->SetMarkerStyle(25) ;
  fCFManager->GetParticleContainer()->ShowProjection(0,3)->SetMarkerStyle(26) ;

  fCFManager->GetParticleContainer()->ShowProjection(1,0)->SetMarkerStyle(23) ;
  fCFManager->GetParticleContainer()->ShowProjection(1,1)->SetMarkerStyle(24) ;
  fCFManager->GetParticleContainer()->ShowProjection(1,2)->SetMarkerStyle(25) ;
  fCFManager->GetParticleContainer()->ShowProjection(1,3)->SetMarkerStyle(26) ;

  TCanvas * c =new TCanvas("c","",1400,800);
  c->Divide(4,2);

  c->cd(1);
  fCFManager->GetParticleContainer()->ShowProjection(0,0)->Draw("p");
  c->cd(2);
  fCFManager->GetParticleContainer()->ShowProjection(0,1)->Draw("p");
  c->cd(3);
  fCFManager->GetParticleContainer()->ShowProjection(0,2)->Draw("p");
  c->cd(4);
  fCFManager->GetParticleContainer()->ShowProjection(0,3)->Draw("p");
  c->cd(5);
  fCFManager->GetParticleContainer()->ShowProjection(1,0)->Draw("p");
  c->cd(6);
  fCFManager->GetParticleContainer()->ShowProjection(1,1)->Draw("p");
  c->cd(7);
  fCFManager->GetParticleContainer()->ShowProjection(1,2)->Draw("p");
  c->cd(8);
  fCFManager->GetParticleContainer()->ShowProjection(1,3)->Draw("p");

  c->SaveAs("plots.eps");

  delete fHistEventsProcessed ;
}

//___________________________________________________________________________
void AliCFRsnTask::ConnectInputData(Option_t *) {
  //
  // Initialize branches.
  //
  Info("ConnectInputData","ConnectInputData of task %s\n",GetName());

  fChain = (TChain*)GetInputData(0);
  fChain->SetBranchStatus("*FMD*",0);
  fChain->SetBranchStatus("*CaloClusters*",0);
  fESD = new AliESDEvent();
  fESD->ReadFromTree(fChain);
}

//___________________________________________________________________________
void AliCFRsnTask::CreateOutputObjects() {
  //HERE ONE CAN CREATE OUTPUT OBJECTS, IN PARTICULAR IF THE OBJECT PARAMETERS DON'T NEED
  //TO BE SET BEFORE THE EXECUTION OF THE TASK
  //
  Info("CreateOutputObjects","CreateOutputObjects of task %s\n", GetName());

  //slot #0
  fHistEventsProcessed = new TH1I("fHistEventsProcessed","",1,0,1) ;
}

//___________________________________________________________________________
Double_t AliCFRsnTask::GetRapidity(Double_t energy, Double_t pz) {
  if (energy == pz || energy == -pz) {
    printf("GetRapidity : ERROR : rapidity for 4-vector with E = Pz -- infinite result");
    return 999;
  }
  if (energy < pz) {
    printf("GetRapidity : ERROR : rapidity for 4-vector with E = Pz -- infinite result");
    return 999;
  }
  Double_t y = 0.5 * log((energy + pz) / (energy - pz));
  return y;
}
Commit	Line	Data
2fbc0b17	1	/**************************************************************************
	2	* Copyright(c) 1998-1999, 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	//-----------------------------------------------------------------------
	17	// Author : R. Vernet, Consorzio Cometa - Catania (it)
	18	//-----------------------------------------------------------------------
	19
	20
	21	#include <TROOT.h>
	22	#include <TInterpreter.h>
	23
	24	#include "AliCFRsnTask.h"
	25	#include "TCanvas.h"
	26	#include "AliStack.h"
	27	#include "TParticle.h"
	28	#include "TH1I.h"
	29	#include "TChain.h"
	30	#include "AliMCEventHandler.h"
	31	#include "AliMCEvent.h"
	32	#include "AliAnalysisManager.h"
	33	#include "AliESDEvent.h"
	34	#include "AliCFManager.h"
	35	#include "AliCFCutBase.h"
	36	#include "AliCFContainer.h"
	37	#include "TChain.h"
	38	#include "AliESDtrack.h"
	39	#include "AliLog.h"
	40	#include "AliRsnDaughter.h"
	41	#include "AliCFPair.h"
	42	#include <memory>
	43
	44	//__________________________________________________________________________
	45	AliCFRsnTask::AliCFRsnTask() :
	46	fRsnPDG(313),
	47	fChain(0x0),
	48	fESD(0x0),
	49	fCFManager(0x0),
	50	fHistEventsProcessed(0x0)
	51	{
	52	//Defual ctor
	53	}
	54	//___________________________________________________________________________
	55	AliCFRsnTask::AliCFRsnTask(const Char_t* name) :
	56	AliAnalysisTask(name,"AliCFRsnTask"),
	57	fRsnPDG(313),
	58	fChain(0x0),
	59	fESD(0x0),
	60	fCFManager(0x0),
	61	fHistEventsProcessed(0x0)
	62	{
	63	//
	64	// Constructor. Initialization of Inputs and Outputs
65	//
66	Info("AliCFRsnTask","Calling Constructor");
67	DefineInput (0,TChain::Class());
68	DefineOutput(0,TH1I::Class());
69	DefineOutput(1,AliCFContainer::Class());
70	// DefineOutput(2,TList::Class());
71	}
72
73	//___________________________________________________________________________
74	AliCFRsnTask& AliCFRsnTask::operator=(const AliCFRsnTask& c)
75	{
76	//
77	// Assignment operator
78	//
79	if (this!=&c) {
80	AliAnalysisTask::operator=(c) ;
81	fRsnPDG = c.fRsnPDG;
82	fChain = c.fChain;
83	fESD = c.fESD;
84	fCFManager = c.fCFManager;
85	fHistEventsProcessed = c.fHistEventsProcessed;
86	}
87	return *this;
88	}
89
90	//___________________________________________________________________________
91	AliCFRsnTask::AliCFRsnTask(const AliCFRsnTask& c) :
92	AliAnalysisTask(c),
93	fRsnPDG(c.fRsnPDG),
94	fChain(c.fChain),
95	fESD(c.fESD),
96	fCFManager(c.fCFManager),
97	fHistEventsProcessed(c.fHistEventsProcessed)
98	{
99	//
100	// Copy Constructor
101	//
102	}
103
104	//___________________________________________________________________________
105	AliCFRsnTask::~AliCFRsnTask() {
106	//
107	//destructor
108	//
109	Info("~AliCFRsnTask","Calling Destructor");
110	if (fChain) delete fChain ;
111	if (fESD) delete fESD ;
112	if (fCFManager) delete fCFManager ;
113	if (fHistEventsProcessed) delete fHistEventsProcessed ;
114	}
115
116	//___________________________________________________________________________
117
118	void AliCFRsnTask::Init()
119	{
120
121	}
122	//_________________________________________________
123	void AliCFRsnTask::Exec(Option_t *)
124	{
125	//
126	// Main loop function
127	//
128	Info("Exec","") ;
129	// Get the mc truth
130	AliMCEventHandler* mcTruth = (AliMCEventHandler*)((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());
131
132	if (!mcTruth) Error("Exec","NO MC INFO FOUND... EXITING");
133
134	// transform possible old AliESD into AliESDEvent
135
136	if (fESD->GetAliESDOld()) fESD->CopyFromOldESD(); //transition to new ESD format
137
138	//pass the MC evt handler to the cuts that need it
139
140	fCFManager->SetEventInfo(mcTruth);
141
142	// Get the MC event
143	AliMCEvent* mcEvent = mcTruth->MCEvent();
144	AliStack* stack = mcEvent->Stack();
145
146	// MC-event selection
147	Double_t containerInput[2] ;
148
149	//loop on the MC event
150	Info("Exec","Looping on MC event");
151	for (Int_t ipart=0; ipart<stack->GetNprimary(); ipart++) {
152	AliMCParticle *mcPart = mcEvent->GetTrack(ipart);
153
154	//check the MC-level cuts
155	if (!fCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,mcPart)) continue;
156	containerInput[0] = mcPart->Pt();
157	containerInput[1] = mcPart->Eta() ;
158	//fill the container for Gen-level selection
159	fCFManager->GetParticleContainer()->Fill(containerInput,kStepGenerated);
160
161	//check the Acceptance-level cuts
162	if (!fCFManager->CheckParticleCuts(AliCFManager::kPartAccCuts,mcPart)) continue;
163	//fill the container for Acceptance-level selection
164	fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructible);
165	}
166
167
168	//Now go to rec level
169	Info("Exec","Looping on ESD event");
170
171	//SET THE ESD AS EVENT INFO IN RECONSTRUCTION CUTS
172	TObjArray* fCutsReco = fCFManager->GetParticleCutsList(AliCFManager::kPartRecCuts);
173	TObjArray* fCutsPID = fCFManager->GetParticleCutsList(AliCFManager::kPartSelCuts);
174	TObjArrayIter iter1(fCutsReco);
175	TObjArrayIter iter2(fCutsPID);
176	AliCFCutBase *cut = 0;
177	while ( (cut = (AliCFCutBase*)iter1.Next()) ) {
178	cut->SetEvtInfo(fESD);
179	}
180	while ( (cut = (AliCFCutBase*)iter2.Next()) ) {
181	cut->SetEvtInfo(fESD);
182	}
183
184	for (Int_t iTrack1 = 0; iTrack1<fESD->GetNumberOfTracks(); iTrack1++) {
185	AliESDtrack* esdTrack1 = fESD->GetTrack(iTrack1);
186	//track1 is negative
187	if (esdTrack1->Charge()>=0) continue;
188	Int_t esdLabel1 = esdTrack1->GetLabel();
189	if (esdLabel1<0) continue;
190
191	for (Int_t iTrack2 = 0; iTrack2<fESD->GetNumberOfTracks(); iTrack2++) {
192	AliESDtrack* esdTrack2 = fESD->GetTrack(iTrack2);
193	//track2 is positive
194	if (esdTrack2->Charge()<=0) continue;
195	Int_t esdLabel2 = esdTrack2->GetLabel();
196	if (esdLabel2<0) continue;
197
198	// copy ESDtrack data into RsnDaughter (and make Bayesian PID)
199	//if problem with copy constructor, use the following special command to destroy the pointer when exiting scope
200	//std::auto_ptr<AliRsnDaughter> track1(AliRsnDaughter::Adopt(esdTrack1,iTrack1));
201
202	AliRsnDaughter* tmp1=AliRsnDaughter::Adopt(esdTrack1,iTrack1);
203	AliRsnDaughter* tmp2=AliRsnDaughter::Adopt(esdTrack2,iTrack2);
204	AliRsnDaughter track1(*tmp1);
205	AliRsnDaughter track2(*tmp2);
206	delete tmp1;
207	delete tmp2;
208
209	TParticle *part1 = stack->Particle(esdLabel1);
210	track1.SetTruePDG(part1->GetPdgCode());
211	Int_t mother1 = part1->GetFirstMother();
212	track1.SetMother(mother1);
213	if (mother1 >= 0) {
214	TParticle *mum = stack->Particle(mother1);
215	track1.SetMotherPDG(mum->GetPdgCode());
216	}
217
218	TParticle *part2 = stack->Particle(esdLabel2);
219	track2.SetTruePDG(part2->GetPdgCode());
220	Int_t mother2 = part2->GetFirstMother();
221	track2.SetMother(mother2);
222	if (mother2 >= 0) {
223	TParticle *mum = stack->Particle(mother2);
224	track2.SetMotherPDG(mum->GetPdgCode());
225	}
226
227	//make a mother resonance from the 2 candidate daughters
228	AliRsnDaughter rsn = AliRsnDaughter::Sum(track1,track2) ;
229	AliCFPair pair(esdTrack1,esdTrack2);
230
231	//check if true resonance
232	if (rsn.GetMotherPDG() != fRsnPDG) continue;
233	if (!fCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,&pair)) continue;
234
235	//check if associated MC resonance passes the cuts
236	Int_t motherLabel=rsn.GetLabel();
237	if (motherLabel<0) continue ;
238	AliMCParticle* mcRsn = mcEvent->GetTrack(motherLabel);
239	if (!mcRsn) continue;
240	if (!fCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,mcRsn)) continue;
241
242	//fill the container
243	containerInput[0] = rsn.GetPt() ;
244	containerInput[1] = GetRapidity(rsn.GetEnergy(),rsn.GetPz());
245	fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructed) ;
246
247	if (!fCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,&pair)) continue ;
248	fCFManager->GetParticleContainer()->Fill(containerInput,kStepSelected);
249	}
250	}
251
252	fHistEventsProcessed->Fill(0);
253	PostData(0,fHistEventsProcessed) ;
254	PostData(1,fCFManager->GetParticleContainer()) ;
255
256	// TList * list = new TList();
257	// fCFManager->AddQAHistosToList(list);
258	// PostData(2,list) ;
259	}
260
261
262	//___________________________________________________________________________
263	void AliCFRsnTask::Terminate(Option_t*)
264	{
265	// The Terminate() function is the last function to be called during
266	// a query. It always runs on the client, it can be used to present
267	// the results graphically or save the results to file.
268
269	Info("Terminate","");
270	AliAnalysisTask::Terminate();
271
272
273	Double_t max1 = fCFManager->GetParticleContainer()->ShowProjection(0,0)->GetMaximum();
274	Double_t max2 = fCFManager->GetParticleContainer()->ShowProjection(1,0)->GetMaximum();
275
276	fCFManager->GetParticleContainer()->ShowProjection(0,0)->GetYaxis()->SetRangeUser(0,max1*1.2);
277	fCFManager->GetParticleContainer()->ShowProjection(0,1)->GetYaxis()->SetRangeUser(0,max1*1.2);
278	fCFManager->GetParticleContainer()->ShowProjection(0,2)->GetYaxis()->SetRangeUser(0,max1*1.2);
279	fCFManager->GetParticleContainer()->ShowProjection(0,3)->GetYaxis()->SetRangeUser(0,max1*1.2);
280
281	fCFManager->GetParticleContainer()->ShowProjection(1,0)->GetYaxis()->SetRangeUser(0,max2*1.2);
282	fCFManager->GetParticleContainer()->ShowProjection(1,1)->GetYaxis()->SetRangeUser(0,max2*1.2);
283	fCFManager->GetParticleContainer()->ShowProjection(1,2)->GetYaxis()->SetRangeUser(0,max2*1.2);
284	fCFManager->GetParticleContainer()->ShowProjection(1,3)->GetYaxis()->SetRangeUser(0,max2*1.2);
285
286	fCFManager->GetParticleContainer()->ShowProjection(0,0)->SetMarkerStyle(23) ;
287	fCFManager->GetParticleContainer()->ShowProjection(0,1)->SetMarkerStyle(24) ;
288	fCFManager->GetParticleContainer()->ShowProjection(0,2)->SetMarkerStyle(25) ;
289	fCFManager->GetParticleContainer()->ShowProjection(0,3)->SetMarkerStyle(26) ;
290
291	fCFManager->GetParticleContainer()->ShowProjection(1,0)->SetMarkerStyle(23) ;
292	fCFManager->GetParticleContainer()->ShowProjection(1,1)->SetMarkerStyle(24) ;
293	fCFManager->GetParticleContainer()->ShowProjection(1,2)->SetMarkerStyle(25) ;
294	fCFManager->GetParticleContainer()->ShowProjection(1,3)->SetMarkerStyle(26) ;
295
296	TCanvas * c =new TCanvas("c","",1400,800);
297	c->Divide(4,2);
298
299	c->cd(1);
300	fCFManager->GetParticleContainer()->ShowProjection(0,0)->Draw("p");
301	c->cd(2);
302	fCFManager->GetParticleContainer()->ShowProjection(0,1)->Draw("p");
303	c->cd(3);
304	fCFManager->GetParticleContainer()->ShowProjection(0,2)->Draw("p");
305	c->cd(4);
306	fCFManager->GetParticleContainer()->ShowProjection(0,3)->Draw("p");
307	c->cd(5);
308	fCFManager->GetParticleContainer()->ShowProjection(1,0)->Draw("p");
309	c->cd(6);
310	fCFManager->GetParticleContainer()->ShowProjection(1,1)->Draw("p");
311	c->cd(7);
312	fCFManager->GetParticleContainer()->ShowProjection(1,2)->Draw("p");
313	c->cd(8);
314	fCFManager->GetParticleContainer()->ShowProjection(1,3)->Draw("p");
315
316	c->SaveAs("plots.eps");
317
318	delete fHistEventsProcessed ;
319	}
320
321	//___________________________________________________________________________
322	void AliCFRsnTask::ConnectInputData(Option_t *) {
323	//
324	// Initialize branches.
325	//
326	Info("ConnectInputData","ConnectInputData of task %s\n",GetName());
327
328	fChain = (TChain*)GetInputData(0);
329	fChain->SetBranchStatus("FMD",0);
330	fChain->SetBranchStatus("CaloClusters",0);
331	fESD = new AliESDEvent();
332	fESD->ReadFromTree(fChain);
333	}
334
335	//___________________________________________________________________________
336	void AliCFRsnTask::CreateOutputObjects() {
337	//HERE ONE CAN CREATE OUTPUT OBJECTS, IN PARTICULAR IF THE OBJECT PARAMETERS DON'T NEED
338	//TO BE SET BEFORE THE EXECUTION OF THE TASK
339	//
340	Info("CreateOutputObjects","CreateOutputObjects of task %s\n", GetName());
341
342	//slot #0
343	fHistEventsProcessed = new TH1I("fHistEventsProcessed","",1,0,1) ;
344	}
345
346	//___________________________________________________________________________
347	Double_t AliCFRsnTask::GetRapidity(Double_t energy, Double_t pz) {
348	if (energy == pz \|\| energy == -pz) {
349	printf("GetRapidity : ERROR : rapidity for 4-vector with E = Pz -- infinite result");
350	return 999;
351	}
352	if (energy < pz) {
353	printf("GetRapidity : ERROR : rapidity for 4-vector with E = Pz -- infinite result");
354	return 999;
355	}
356	Double_t y = 0.5 * log((energy + pz) / (energy - pz));
357	return y;
358	}
359
360