Update from Sandun
[u/mrichter/AliRoot.git] / JETAN / AliAnalysisTaskDiJets.cxx
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 /* $Id$ */
17 /*
18  Analysis Task 
19  for Dijet Analysis
20  based on AOD
21 */
22
23 #include "AliAnalysisTaskDiJets.h"
24 #include "AliAODEvent.h"
25 #include "AliAODJet.h"
26 #include "AliAODDiJet.h"
27 #include <TClonesArray.h>
28 #include <TLorentzVector.h>
29 #include <TStyle.h>
30
31
32 ClassImp(AliAnalysisTaskDiJets)
33
34 ////////////////////////////////////////////////////////////////////////
35
36 AliAnalysisTaskDiJets::AliAnalysisTaskDiJets():
37     AliAnalysisTaskSE(),
38     fDiJets(0),
39     fDiJetsIn(0),
40     fUseAODInput(kFALSE),
41     fFillAOD(kFALSE),
42     fJetBranch("jets"),
43     fAOD(0),
44     fHistList(0),
45     fH1DeltaPt(0),
46     fH1DeltaPhi(0),
47     fH1PhiImbal(0),
48     fH1Asym(0),
49     fH2Pt2vsPt1(0),
50     fH2DifvsSum(0)
51 {
52   // Default constructor
53 }
54
55 //----------------------------------------------------------------------
56 AliAnalysisTaskDiJets::AliAnalysisTaskDiJets(const char* name):
57     AliAnalysisTaskSE(name),
58     fDiJets(0),
59     fDiJetsIn(0),
60     fUseAODInput(kFALSE),
61     fFillAOD(kFALSE),
62     fJetBranch("jets"),
63     fAOD(0),
64     fHistList(0),
65     fH1DeltaPt(0),
66     fH1DeltaPhi(0),
67     fH1PhiImbal(0),
68     fH1Asym(0),
69     fH2Pt2vsPt1(0),
70     fH2DifvsSum(0)
71 {
72   // Default constructor
73     DefineOutput(1, TList::Class());  
74 }
75
76 //----------------------------------------------------------------------
77 void AliAnalysisTaskDiJets::UserCreateOutputObjects()
78 {
79 // Create the output container
80 //
81     if (fDebug) printf("AnalysisTaskDiJets::CreateOutPutData() \n");
82     fDiJets = new TClonesArray("AliAODDiJet", 0);
83     if (fFillAOD){
84       fDiJets->SetName(Form("dijets_%s",fJetBranch.Data()));
85       AddAODBranch("TClonesArray", &fDiJets);
86         }
87
88     if (!fHistList) fHistList = new TList();
89     fHistList->SetOwner();
90     Float_t pi=TMath::Pi();
91     gStyle->SetPalette(1);
92
93     fH1DeltaPt  = new TH1F("DeltaPt","Difference between the jets' Pt;#Deltap_{T} (GeV/c);Entries",150,0.,150.);
94     fH1DeltaPt->SetMarkerSize(0.6);
95     fH1DeltaPt->SetMarkerColor(4);
96     fH1DeltaPt->SetMarkerStyle(21);
97     fH1DeltaPt->SetOption("E");
98
99     fH1DeltaPhi = new TH1F("DeltaPhi","Difference in the azimuthal angle;#Delta#phi;Entries",100,0.,pi);
100     fH1DeltaPhi->SetMarkerSize(0.6);
101     fH1DeltaPhi->SetMarkerColor(4);
102     fH1DeltaPhi->SetMarkerStyle(21);
103     fH1DeltaPhi->SetOption("E");
104
105     fH1PhiImbal = new TH1F("PhiImb","Phi imbalance;#phi;Entries",100,-pi,pi);
106     fH1PhiImbal->SetMarkerSize(0.6);
107     fH1PhiImbal->SetMarkerColor(4);
108     fH1PhiImbal->SetMarkerStyle(21);
109     fH1PhiImbal->SetOption("E");
110
111     fH1Asym     = new TH1F("Asym","Pt asymmetry;#Deltap_{T}/(p_{T,1}+p_{T,2});Entries",50,0.,1.);
112     fH1Asym->SetMarkerSize(0.6);
113     fH1Asym->SetMarkerColor(4);
114     fH1Asym->SetMarkerStyle(21);
115     fH1Asym->SetOption("E");
116
117     fH2Pt2vsPt1 = new TH2F("Pt2vsPt1","Pt2 versus Pt1;p_{T,1} (GeV/c);p_{T,2} (GeV/c)",250,0.,250.,250,0.,250.);
118     fH2Pt2vsPt1->SetOption("cont0");
119
120     fH2DifvsSum = new TH2F("DifvsSum","Pt difference vs Pt sum;p_{T,1}+p_{T,2} (GeV/c);#Deltap_{T} (GeV/c)",400,0.,400.,150,0.,150.);
121     fH2DifvsSum->SetOption("cont0");
122
123     fHistList->Add(fH1DeltaPt);
124     fHistList->Add(fH1DeltaPhi);
125     fHistList->Add(fH1PhiImbal);
126     fHistList->Add(fH1Asym);
127     fHistList->Add(fH2Pt2vsPt1);
128     fHistList->Add(fH2DifvsSum);
129 }
130
131 //----------------------------------------------------------------------
132 void AliAnalysisTaskDiJets::Init()
133 {
134     // Initialization
135     if (fDebug) printf("AnalysisTaskDiJets::Init() \n");
136 }
137
138 //----------------------------------------------------------------------
139 void AliAnalysisTaskDiJets::UserExec(Option_t */*option*/)
140 {
141 // Execute analysis for current event
142 //
143     if (fDiJets) fDiJets->Delete();
144
145     if(fUseAODInput){
146       fAOD = dynamic_cast<AliAODEvent*> (InputEvent());
147       if(!fAOD){
148         // We do not have an input AOD, look in the output
149         if (fDebug) printf("%s:%d No AOD event in the input\n",(char*)__FILE__,__LINE__);
150         return;
151       }
152     } else {
153       fAOD = AODEvent();
154       if(!fAOD){
155         if (fDebug) printf("%s:%d AODEvent not found in the Output",(char*)__FILE__,__LINE__);
156         return;
157       }
158     }
159
160     TClonesArray* jets = (TClonesArray*) fAOD->FindListObject(fJetBranch.Data());
161     // N.B. if we take the aod from the output this is always
162     // empty and since it is the same as fDiJets 
163     fDiJetsIn = (TClonesArray*) (fAOD->GetList()->FindObject("dijets"));
164
165     if (fDiJetsIn) {
166       if (fDebug) printf("Found %d dijets in old list \n", fDiJetsIn->GetEntries());
167       AliAODJet* jj1, *jj2;
168       AliAODDiJet* testJ;
169
170       if (fDiJetsIn->GetEntries() > 0) {
171         testJ = (AliAODDiJet*) (fDiJetsIn->At(0));
172         jj1 = testJ->Jet(0);
173         jj1->Print("");
174         jj2 = testJ->Jet(1);
175         jj2->Print("");
176       }
177     }
178
179     Int_t nj = jets->GetEntriesFast();
180     if (fDebug) printf("There are %5d jets in the event \n", nj);
181
182     if (nj < 2){
183       PostData(1, fHistList);
184       return;
185     }
186     AliAODJet* jet1 = (AliAODJet*) (jets->At(0));
187     TLorentzVector v1 = *(jet1->MomentumVector());
188     AliAODJet* jet2 = (AliAODJet*) (jets->At(1));
189     TLorentzVector v2 = *(jet2->MomentumVector());
190     TLorentzVector v = v1 + v2;
191     if (fDiJets) {
192         Int_t ndi = fDiJets->GetEntriesFast();
193         TClonesArray &lref = *fDiJets;
194         new(lref[ndi]) AliAODDiJet(v);
195         AliAODDiJet* dijet = (AliAODDiJet*) (fDiJets->At(ndi));
196         dijet->SetJetRefs(jet1, jet2);
197         fH1DeltaPhi->Fill(dijet->DeltaPhi());
198         fH1PhiImbal->Fill(dijet->PhiImbalance());
199
200     }
201     
202     fH1DeltaPt->Fill(jet1->Pt()-jet2->Pt());
203     fH1Asym->Fill((jet1->Pt()-jet2->Pt())/(jet1->Pt()+jet2->Pt()));
204     fH2Pt2vsPt1->Fill(jet1->Pt(),jet2->Pt());
205     fH2DifvsSum->Fill(jet1->Pt()+jet2->Pt(),jet1->Pt()-jet2->Pt());
206
207     PostData(1, fHistList);
208     return;
209 }
210
211 //----------------------------------------------------------------------
212 void AliAnalysisTaskDiJets::Terminate(Option_t */*option*/)
213 {
214 // Terminate analysis
215 //
216     if (fDebug) printf("AnalysisDiJets: Terminate() \n");
217 }
218