[u/mrichter/AliRoot.git] / PWG3 / AliAnalysisTaskMuonAODfromGeneral.cxx

#define AliAnalysisTaskMuonAODfromGeneral_cxx

// 19 Nov 2007
// Class implementation for the specific dimuon AOD generation from a general AOD

#include "TTree.h"
#include "TROOT.h"
#include "TH1.h"
#include "TCanvas.h"
#include "TSystem.h"
#include "TRandom.h"

#include "AliAODEvent.h"
#include "AliAnalysisTaskMuonAODfromGeneral.h"
#include "AliAODHandler.h"
#include "AliAnalysisManager.h"
#include "AliAODEventInfo.h"
#include "AliAODDimuon.h"

ClassImp(AliAnalysisTaskMuonAODfromGeneral)

//________________________________________________________________________
AliAnalysisTaskMuonAODfromGeneral::AliAnalysisTaskMuonAODfromGeneral(const char *name, Double_t BeamEnergy):
	AliAnalysisTask(name, "AnalysisTaskMuonAODfromGeneral"),
	fInfos(0),
	fDimuons(0),
	fOrgAOD(0), 
	fNewAOD(0), 
	ft(0),
	fBeamEnergy(0)
{  // Constructor.
  // Input slot #0 works with an Ntuple
  SetBeamEnergy(BeamEnergy);
  DefineInput(0, TChain::Class());
  // Output slot #0 writes into a TTree container
  DefineOutput(0, TTree::Class());
}

//________________________________________________________________________
void AliAnalysisTaskMuonAODfromGeneral::ConnectInputData(Option_t *) {
  printf("   ConnectInputData %s\n", GetName());
  // New connection to the AOD
  fChain = (TChain*)GetInputData(0);
  fOrgAOD = new AliAODEvent();
  fOrgAOD->ReadFromTree(fChain); // This checks also if the branch address is already set
}

//________________________________________________________________________
void AliAnalysisTaskMuonAODfromGeneral::CreateOutputObjects() {
  printf("Creating output objects\n");
  if(!fNewAOD)fNewAOD=new AliAODEvent();
  fNewAOD->CreateStdContent();
 
  if (!ft) {
    ft = new TTree("AOD","Muon AOD tree",0);
    ft->BranchRef();
    ft->Branch(fNewAOD->GetList());
    //ft->Print();
 }
 
  //Add new Objects: MuonInfos & Dimuons
   fInfos = new AliAODEventInfo();
   fNewAOD->AddObject(fInfos);
   const char *nameInfos = "MuonInfos";
   fInfos->SetName(nameInfos);
   ft->Branch(nameInfos, &fInfos);
   fDimuons = new TClonesArray("AliAODDimuon",0);
   fNewAOD->AddObject(fDimuons);
   const char *nameDimuons = "Dimuons";
   fDimuons->SetName(nameDimuons);
   ft->Branch(nameDimuons, &fDimuons);
      
}

//________________________________________________________________________
void AliAnalysisTaskMuonAODfromGeneral::Exec(Option_t *) {
  static int ncall=0;
  // Task making a Muon AOD
  // Get input data
  TChain *chain = (TChain*)GetInputData(0);
  Long64_t ientry = chain->GetReadEntry();
  ientry=ientry;

  if (!fOrgAOD) return;
  
  Int_t nTracks=fOrgAOD->GetNumberOfTracks();
  Int_t nMuTracks=0;
  Int_t nPosTracks = 0;
  Int_t mutrNumb[10]; for (Int_t i=0; i<10; mutrNumb[i++]=0);
  for (Int_t iTrack=0; iTrack<nTracks; iTrack++){
       const Double_t *trackpid=fOrgAOD->GetTrack(iTrack)->PID();
       if (trackpid[AliAODTrack::kMuon]==1.) {
       		mutrNumb[nMuTracks]=iTrack;
		nMuTracks++;
		}
       if (fOrgAOD->GetTrack(iTrack)->Charge()> 0) nPosTracks++;
       }
  
  Bool_t ExistMuon=kFALSE;
  if (nMuTracks>0) ExistMuon=kTRUE;
  
  //--------------------------------------------------------------------
  // set arrays and pointers
  Double_t pos[3];
  Double_t covVtx[6];
  
    // Access to the header
    AliAODHeader *header = fNewAOD->GetHeader();

  // fill the header
  header->SetRunNumber       (fOrgAOD->GetRunNumber()       );
  header->SetBunchCrossNumber(fOrgAOD->GetBunchCrossNumber());
  header->SetOrbitNumber     (fOrgAOD->GetOrbitNumber()     );
  header->SetPeriodNumber    (fOrgAOD->GetPeriodNumber()    );
  header->SetTriggerMask     (fOrgAOD->GetTriggerMask()     ); 
  header->SetTriggerCluster  (fOrgAOD->GetTriggerCluster()  );
  header->SetEventType       (fOrgAOD->GetEventType()       );
  header->SetMagneticField   (fOrgAOD->GetMagneticField()   );
  header->SetZDCN1Energy     (fOrgAOD->GetZDCN1Energy()     );
  header->SetZDCP1Energy     (fOrgAOD->GetZDCP1Energy()     );
  header->SetZDCN2Energy     (fOrgAOD->GetZDCN2Energy()     );
  header->SetZDCP2Energy     (fOrgAOD->GetZDCP2Energy()     );
  header->SetZDCEMEnergy     (0,fOrgAOD->GetZDCEMEnergy(0)  );
  header->SetZDCEMEnergy     (1,fOrgAOD->GetZDCEMEnergy(1)  );
  header->SetRefMultiplicity   (nTracks);
  header->SetRefMultiplicityPos(nPosTracks);
  header->SetRefMultiplicityNeg(nTracks - nPosTracks);
  header->SetMuonMagFieldScale(-999.); // FIXME
  header->SetCentrality(-999.);        // FIXME

  Int_t nVertices = 0;
  const AliAODVertex *vtx = fOrgAOD->GetPrimaryVertex();
  if(vtx)nVertices = 1;

  Int_t nDimuons=nMuTracks*(nMuTracks-1)/2;
  
  // Access to the AOD container of vertices
  TClonesArray &vertices = *(fNewAOD->GetVertices());
  Int_t jVertices=0;

  // Access to the AOD container of tracks
   TClonesArray &tracks = *(fNewAOD->GetTracks());
   tracks.Delete();
   Int_t jTracks=0; 


  vtx->GetCovMatrix(covVtx); //covariance matrix

  AliAODVertex * primary = new(vertices[jVertices++])
      AliAODVertex(pos, covVtx, vtx->GetChi2perNDF(), NULL, AliAODVertex::kPrimary);

  static int ncal=0;
  static int numtracks=0;
  for (Int_t iTrack=0; iTrack<nMuTracks; iTrack++) {
	AliAODTrack *aodTrack;
	primary->AddDaughter( aodTrack=new(tracks[jTracks++])
	AliAODTrack((*(fOrgAOD->GetTrack(mutrNumb[iTrack])))));
        ncal++;
  }
  numtracks+=jTracks;

  fDimuons->Delete();
  fDimuons->Expand(nDimuons);
  TClonesArray &rDimuons = *fDimuons;
    
  Int_t jDimuons=0;
  for(Int_t i=0; i<nMuTracks; i++){
     for(Int_t j=i+1; j<nMuTracks; j++){
       new(rDimuons[jDimuons++]) AliAODDimuon(tracks[i],tracks[j],fInfos);
     } 
  }
  
  fInfos->SetBeamEnergy(fBeamEnergy);
  fInfos->ev=fNewAOD;
  fInfos->ei=fInfos;
  fInfos->he=header;
  fInfos->tr=fNewAOD->GetTracks();
  fInfos->di=fDimuons;
  fInfos->fMUON_Single_LPt_L0=0;
  fInfos->fMUON_Single_HPt_L0=1;
  fInfos->fMUON_Like_LPt_L0=2; 
  fInfos->fMUON_Like_HPt_L0=3;  
  fInfos->fMUON_Unlike_LPt_L0=4;
  fInfos->fMUON_Unlike_HPt_L0=5;
  if(ExistMuon) ft->Fill();
  ncall++;

  // Post final data. It will be written to a file with option "RECREATE"
  PostData(0, ft);
}      

//________________________________________________________________________
void AliAnalysisTaskMuonAODfromGeneral::Terminate(Option_t *) {
  // Draw some histogram at the end.
  ft->GetUserInfo()->Add(fNewAOD);
  ft->Write("",TObject::kOverwrite);
  ft->ResetBranchAddresses();
  if (!gROOT->IsBatch()) {
    TCanvas *c1 = new TCanvas("c1","Pt",10,10,310,310);
    c1->SetFillColor(10);
    c1->SetHighLightColor(10);
    
    c1->cd(1)->SetLeftMargin(0.15);
    c1->cd(1)->SetBottomMargin(0.15);  
    c1->cd(1)->SetLogy();
    ft->Draw("tracks.Pt()");
    //ft->Draw("dimuons.M()");
  }
}
Commit	Line	Data
71b7d225	1	#define AliAnalysisTaskMuonAODfromGeneral_cxx
	2
	3	// 19 Nov 2007
	4	// Class implementation for the specific dimuon AOD generation from a general AOD
	5
	6	#include "TTree.h"
	7	#include "TROOT.h"
	8	#include "TH1.h"
	9	#include "TCanvas.h"
	10	#include "TSystem.h"
	11	#include "TRandom.h"
	12
	13	#include "AliAODEvent.h"
	14	#include "AliAnalysisTaskMuonAODfromGeneral.h"
	15	#include "AliAODHandler.h"
	16	#include "AliAnalysisManager.h"
	17	#include "AliAODEventInfo.h"
	18	#include "AliAODDimuon.h"
	19
	20	ClassImp(AliAnalysisTaskMuonAODfromGeneral)
	21
	22	//________________________________________________________________________
	23	AliAnalysisTaskMuonAODfromGeneral::AliAnalysisTaskMuonAODfromGeneral(const char *name, Double_t BeamEnergy):
	24	AliAnalysisTask(name, "AnalysisTaskMuonAODfromGeneral"),
	25	fInfos(0),
	26	fDimuons(0),
	27	fOrgAOD(0),
	28	fNewAOD(0),
	29	ft(0),
	30	fBeamEnergy(0)
	31	{ // Constructor.
	32	// Input slot #0 works with an Ntuple
	33	SetBeamEnergy(BeamEnergy);
	34	DefineInput(0, TChain::Class());
	35	// Output slot #0 writes into a TTree container
	36	DefineOutput(0, TTree::Class());
	37	}
	38
	39	//________________________________________________________________________
	40	void AliAnalysisTaskMuonAODfromGeneral::ConnectInputData(Option_t *) {
	41	printf(" ConnectInputData %s\n", GetName());
	42	// New connection to the AOD
	43	fChain = (TChain*)GetInputData(0);
	44	fOrgAOD = new AliAODEvent();
	45	fOrgAOD->ReadFromTree(fChain); // This checks also if the branch address is already set
	46	}
	47
	48	//________________________________________________________________________
	49	void AliAnalysisTaskMuonAODfromGeneral::CreateOutputObjects() {
	50	printf("Creating output objects\n");
	51	if(!fNewAOD)fNewAOD=new AliAODEvent();
	52	fNewAOD->CreateStdContent();
	53
	54	if (!ft) {
	55	ft = new TTree("AOD","Muon AOD tree",0);
	56	ft->BranchRef();
	57	ft->Branch(fNewAOD->GetList());
	58	//ft->Print();
	59	}
	60
	61	//Add new Objects: MuonInfos & Dimuons
	62	fInfos = new AliAODEventInfo();
	63	fNewAOD->AddObject(fInfos);
	64	const char *nameInfos = "MuonInfos";
65	fInfos->SetName(nameInfos);
66	ft->Branch(nameInfos, &fInfos);
67	fDimuons = new TClonesArray("AliAODDimuon",0);
68	fNewAOD->AddObject(fDimuons);
69	const char *nameDimuons = "Dimuons";
70	fDimuons->SetName(nameDimuons);
71	ft->Branch(nameDimuons, &fDimuons);
72
73	}
74
75	//________________________________________________________________________
76	void AliAnalysisTaskMuonAODfromGeneral::Exec(Option_t *) {
77	static int ncall=0;
78	// Task making a Muon AOD
79	// Get input data
80	TChain chain = (TChain)GetInputData(0);
81	Long64_t ientry = chain->GetReadEntry();
82	ientry=ientry;
83
84	if (!fOrgAOD) return;
85
86	Int_t nTracks=fOrgAOD->GetNumberOfTracks();
87	Int_t nMuTracks=0;
88	Int_t nPosTracks = 0;
89	Int_t mutrNumb[10]; for (Int_t i=0; i<10; mutrNumb[i++]=0);
90	for (Int_t iTrack=0; iTrack<nTracks; iTrack++){
91	const Double_t *trackpid=fOrgAOD->GetTrack(iTrack)->PID();
92	if (trackpid[AliAODTrack::kMuon]==1.) {
93	mutrNumb[nMuTracks]=iTrack;
94	nMuTracks++;
95	}
96	if (fOrgAOD->GetTrack(iTrack)->Charge()> 0) nPosTracks++;
97	}
98
99	Bool_t ExistMuon=kFALSE;
100	if (nMuTracks>0) ExistMuon=kTRUE;
101
102	//--------------------------------------------------------------------
103	// set arrays and pointers
104	Double_t pos[3];
105	Double_t covVtx[6];
106
107	// Access to the header
108	AliAODHeader *header = fNewAOD->GetHeader();
109
110	// fill the header
111	header->SetRunNumber (fOrgAOD->GetRunNumber() );
112	header->SetBunchCrossNumber(fOrgAOD->GetBunchCrossNumber());
113	header->SetOrbitNumber (fOrgAOD->GetOrbitNumber() );
114	header->SetPeriodNumber (fOrgAOD->GetPeriodNumber() );
115	header->SetTriggerMask (fOrgAOD->GetTriggerMask() );
116	header->SetTriggerCluster (fOrgAOD->GetTriggerCluster() );
117	header->SetEventType (fOrgAOD->GetEventType() );
118	header->SetMagneticField (fOrgAOD->GetMagneticField() );
119	header->SetZDCN1Energy (fOrgAOD->GetZDCN1Energy() );
120	header->SetZDCP1Energy (fOrgAOD->GetZDCP1Energy() );
121	header->SetZDCN2Energy (fOrgAOD->GetZDCN2Energy() );
122	header->SetZDCP2Energy (fOrgAOD->GetZDCP2Energy() );
123	header->SetZDCEMEnergy (0,fOrgAOD->GetZDCEMEnergy(0) );
124	header->SetZDCEMEnergy (1,fOrgAOD->GetZDCEMEnergy(1) );
125	header->SetRefMultiplicity (nTracks);
126	header->SetRefMultiplicityPos(nPosTracks);
127	header->SetRefMultiplicityNeg(nTracks - nPosTracks);
128	header->SetMuonMagFieldScale(-999.); // FIXME
129	header->SetCentrality(-999.); // FIXME
130
131	Int_t nVertices = 0;
132	const AliAODVertex *vtx = fOrgAOD->GetPrimaryVertex();
133	if(vtx)nVertices = 1;
134
135	Int_t nDimuons=nMuTracks*(nMuTracks-1)/2;
136
137	// Access to the AOD container of vertices
138	TClonesArray &vertices = *(fNewAOD->GetVertices());
139	Int_t jVertices=0;
140
141	// Access to the AOD container of tracks
142	TClonesArray &tracks = *(fNewAOD->GetTracks());
143	tracks.Delete();
144	Int_t jTracks=0;
145
146
147	vtx->GetCovMatrix(covVtx); //covariance matrix
148
149	AliAODVertex * primary = new(vertices[jVertices++])
150	AliAODVertex(pos, covVtx, vtx->GetChi2perNDF(), NULL, AliAODVertex::kPrimary);
151
152	static int ncal=0;
153	static int numtracks=0;
154	for (Int_t iTrack=0; iTrack<nMuTracks; iTrack++) {
155	AliAODTrack *aodTrack;
156	primary->AddDaughter( aodTrack=new(tracks[jTracks++])
157	AliAODTrack((*(fOrgAOD->GetTrack(mutrNumb[iTrack])))));
158	ncal++;
159	}
160	numtracks+=jTracks;
161
162	fDimuons->Delete();
163	fDimuons->Expand(nDimuons);
164	TClonesArray &rDimuons = *fDimuons;
165
166	Int_t jDimuons=0;
167	for(Int_t i=0; i<nMuTracks; i++){
168	for(Int_t j=i+1; j<nMuTracks; j++){
169	new(rDimuons[jDimuons++]) AliAODDimuon(tracks[i],tracks[j],fInfos);
170	}
171	}
172
173	fInfos->SetBeamEnergy(fBeamEnergy);
174	fInfos->ev=fNewAOD;
175	fInfos->ei=fInfos;
176	fInfos->he=header;
177	fInfos->tr=fNewAOD->GetTracks();
178	fInfos->di=fDimuons;
179	fInfos->fMUON_Single_LPt_L0=0;
180	fInfos->fMUON_Single_HPt_L0=1;
181	fInfos->fMUON_Like_LPt_L0=2;
182	fInfos->fMUON_Like_HPt_L0=3;
183	fInfos->fMUON_Unlike_LPt_L0=4;
184	fInfos->fMUON_Unlike_HPt_L0=5;
185	if(ExistMuon) ft->Fill();
186	ncall++;
187
188	// Post final data. It will be written to a file with option "RECREATE"
189	PostData(0, ft);
190	}
191
192	//________________________________________________________________________
193	void AliAnalysisTaskMuonAODfromGeneral::Terminate(Option_t *) {
194	// Draw some histogram at the end.
195	ft->GetUserInfo()->Add(fNewAOD);
196	ft->Write("",TObject::kOverwrite);
197	ft->ResetBranchAddresses();
198	if (!gROOT->IsBatch()) {
199	TCanvas *c1 = new TCanvas("c1","Pt",10,10,310,310);
200	c1->SetFillColor(10);
201	c1->SetHighLightColor(10);
202
203	c1->cd(1)->SetLeftMargin(0.15);
204	c1->cd(1)->SetBottomMargin(0.15);
205	c1->cd(1)->SetLogy();
206	ft->Draw("tracks.Pt()");
207	//ft->Draw("dimuons.M()");
208	}
209	}