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