85316aa13e4572a6ee30215b60d7e87f16a2e984
[u/mrichter/AliRoot.git] / JETAN / AliJetFinder.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 //---------------------------------------------------------------------
17 // Jet finder base class
18 // manages the search for jets 
19 // Authors: jgcn@mda.cinvestav.mx
20 //          andreas.morsch@cern.ch
21 //---------------------------------------------------------------------
22
23 #include <Riostream.h>
24 #include <TFile.h>
25 #include "AliJetFinder.h"
26 #include "AliJet.h"
27 #include "AliJetReader.h"
28 #include "AliJetReaderHeader.h"
29 #include "AliJetControlPlots.h"
30 #include "AliLeading.h"
31
32 ClassImp(AliJetFinder)
33
34 AliJetFinder::AliJetFinder():
35     fPlotMode(kFALSE),
36     fJets(0),
37     fGenJets(0),
38     fLeading(0),
39     fReader(0x0),
40     fPlots(0x0),
41     fOut(0x0)
42     
43 {
44   // Constructor
45   fJets    = new AliJet();
46   fGenJets = new AliJet();
47   fLeading = new AliLeading();
48 }
49
50 ////////////////////////////////////////////////////////////////////////
51
52 AliJetFinder::~AliJetFinder()
53 {
54   // destructor
55   // here reset and delete jets
56   fJets->ClearJets();
57   delete fJets;
58   fGenJets->ClearJets();
59   delete fGenJets;
60   // close file
61   if (fOut) {
62     fOut->Close();
63     fOut->Delete();
64   }
65   delete fOut;
66   // reset and delete control plots
67   if (fPlots) delete fPlots;
68 }
69
70 ////////////////////////////////////////////////////////////////////////
71
72 void AliJetFinder::SetOutputFile(const char *name)
73 {
74   // opens output file 
75   fOut = new TFile(name,"recreate");
76 }
77
78 ////////////////////////////////////////////////////////////////////////
79
80 void AliJetFinder::PrintJets()
81 {
82   // Print jet information
83   cout << " Jets found with jet algorithm:" << endl;
84   fJets->PrintJets();
85   cout << " Jets found by pythia:" << endl;
86   fGenJets->PrintJets();
87 }
88
89 ////////////////////////////////////////////////////////////////////////
90
91 void AliJetFinder::SetPlotMode(Bool_t b)
92 {
93   // Sets the plotting mode
94   fPlotMode=b;
95   if (b && !fPlots) fPlots = new AliJetControlPlots(); 
96 }
97
98 ////////////////////////////////////////////////////////////////////////
99
100 void AliJetFinder::WriteJetsToFile(Int_t i)
101 {
102   // Writes the jets to file
103   fOut->cd();
104   char hname[30];
105   sprintf(hname,"TreeJ%d",i);
106   TTree* jetT = new TTree(hname,"AliJet");
107   jetT->Branch("FoundJet",&fJets,1000);
108   jetT->Branch("GenJet",&fGenJets,1000);
109   jetT->Branch("LeadingPart",&fLeading,1000);
110   jetT->Fill();
111   jetT->Write(hname);
112   delete jetT;
113 }
114
115 ////////////////////////////////////////////////////////////////////////
116
117 void AliJetFinder::WriteRHeaderToFile()
118 {
119   // write reader header
120   fOut->cd();
121   AliJetReaderHeader *rh = fReader->GetReaderHeader();
122   rh->Write();
123 }
124
125 ////////////////////////////////////////////////////////////////////////
126
127
128 void AliJetFinder::Run()
129 {
130   // Do some initialization
131   Init();
132   // connect files
133   fReader->OpenInputFiles();
134
135   // write headers
136   WriteHeaders();
137   // loop over events
138   Int_t nFirst,nLast;
139   nFirst = fReader->GetReaderHeader()->GetFirstEvent();
140   nLast = fReader->GetReaderHeader()->GetLastEvent();
141   // loop over events
142   for (Int_t i=nFirst;i<nLast;i++) {
143       fReader->FillMomentumArray(i);
144       fLeading->FindLeading(fReader);
145       fReader->GetGenJets(fGenJets);
146       FindJets();
147       if (fOut) WriteJetsToFile(i);
148       if (fPlots) fPlots->FillHistos(fJets);
149       fLeading->Reset();
150       fGenJets->ClearJets();
151       Reset();
152   } 
153   // write out
154   if (fPlots) {
155       fPlots->Normalize();
156       fPlots->PlotHistos();
157   }
158   if (fOut) {
159       fOut->cd();
160       fPlots->Write();
161       fOut->Close();
162   }
163 }
164
165
166 //
167 // The following methods have been added to allow for event steering from the outside
168 //
169
170 void AliJetFinder::ConnectTree(TTree* tree)
171 {
172     // Connect the input file
173     fReader->ConnectTree(tree);
174 }
175
176 void AliJetFinder::WriteHeaders()
177 {
178     // Write the Headers
179     if (fOut) {
180         fOut->cd();
181         WriteRHeaderToFile();
182         WriteJHeaderToFile();
183     }
184 }
185
186
187 Bool_t AliJetFinder::ProcessEvent(Long64_t entry)
188 {
189 //
190 // Process one event
191 //
192     printf("<<<<< Processing Event %5d >>>>> \n", (Int_t) entry);
193     Bool_t ok = fReader->FillMomentumArray(entry);
194     if (!ok) return kFALSE;
195     fLeading->FindLeading(fReader);
196     FindJets();
197     if (fOut)   WriteJetsToFile(entry);
198     if (fPlots) fPlots->FillHistos(fJets);
199     fLeading->Reset();
200     fGenJets->ClearJets();
201     Reset();  
202     return kTRUE;
203 }
204
205 void AliJetFinder::FinishRun()
206 {
207     // Finish a run
208   if (fPlots) {
209       fPlots->Normalize();
210       fPlots->PlotHistos();
211   }
212   if (fOut) {
213       fOut->cd();
214       fPlots->Write();
215       fOut->Close();
216   }   
217 }
218