fix mem leak and compiler warning
[u/mrichter/AliRoot.git] / ANALYSIS / AliReaderAOD.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 //______________________________________________________________________________
19 ////////////////////////////////////////////////////////////////////////////////
20 //                                                                            //
21 // class AliReaderAOD                                                         //
22 //                                                                            //
23 // Reader and Writer for AOD format.                                          //
24 // AODs are stored in a tree named by the variable fgkTreeName.               //
25 // There is stored 1 or 2 branches. Each of them stores AOD objects           //
26 // First branch is named by the variable fgkReconstructedDataBranchName       //
27 // ("reconstructed.") and keeps reconstructed data.                           //
28 // Second branch is called by the variable fgkSimulatedDataBranchName         //
29 // ("simulated.") and stores Monte carlo truth. If both branches are present  //
30 // AODs are parallel, i.e. nth particle in one branch corresponds to the nth  //
31 // particle in the other one.                                                 //
32 //                                                                            //
33 // Since we accept different formats of particles that are stored in AODs     //
34 // reader must take care of that fact: clean buffer if the next file contains //
35 // different particle type.                                                   //
36 //                                                                            //
37 // Piotr.Skowronski@cern.ch                                                   //
38 //                                                                            //
39 ////////////////////////////////////////////////////////////////////////////////
40
41
42 #include <TError.h>
43 #include <TFile.h>
44 #include <TTree.h>
45 #include <TH1.h>
46
47 #include "AliAOD.h"
48 #include "AliLog.h"
49 #include "AliReaderAOD.h"
50
51 const TString AliReaderAOD::fgkTreeName("TAOD");
52 const TString AliReaderAOD::fgkReconstructedDataBranchName("reconstructed.");
53 const TString AliReaderAOD::fgkSimulatedDataBranchName("simulated.");
54
55 ClassImp(AliReaderAOD)
56
57 AliReaderAOD::AliReaderAOD(const Char_t* aodfilename):
58  fFileName(aodfilename),
59  fReadSim(kFALSE),
60  fReadRec(kTRUE),
61  fTree(0x0),
62  fFile(0x0),
63  fSimBuffer(0x0),
64  fRecBuffer(0x0)
65 {
66   //ctor
67 }
68 /********************************************************************/
69
70 AliReaderAOD::~AliReaderAOD()
71 {
72 //dtor
73   if (fEventSim == fSimBuffer )
74    {
75     fEventSim = 0x0;
76     fEventRec = 0x0;
77    }
78   delete fSimBuffer;
79   delete fRecBuffer;
80   
81   delete fTree;
82   delete fFile;
83 }
84 /********************************************************************/
85
86 void AliReaderAOD::Rewind()
87 {
88 //Rewinds reading
89   delete fTree;
90   fTree = 0x0;
91   delete fFile;
92   fFile = 0x0;
93   fCurrentDir = 0;
94   fNEventsRead= 0;
95 }
96 /********************************************************************/
97 Int_t AliReaderAOD::ReadNext()
98 {
99 //Reads next event
100   
101   Info("ReadNext","Entered");
102   do  //do{}while; is OK even if 0 dirs specified. In that case we try to read from "./"
103     {
104       if (fFile == 0x0)
105        {
106          Int_t openfailed = OpenFile(fCurrentDir);//rl is opened here
107          if (openfailed)
108            {
109              //Error("ReadNext","Error Occured while opening directory number %d",fCurrentDir);
110              fCurrentDir++;
111              continue;
112            }
113          fCurrentEvent = 0;
114        }
115       //Tree must exist because OpenFile would reuturn error in the other case
116       if ( fCurrentEvent >= fTree->GetEntries() )
117        {
118          delete fTree;
119          fTree = 0x0;
120          delete fFile;
121          fFile = 0x0;
122          
123          delete fSimBuffer;
124          delete fRecBuffer;
125          
126          fSimBuffer = 0x0;
127          fRecBuffer = 0x0;
128          fCurrentDir++;
129          continue;
130        }
131       
132       Info("ReadNext","Getting event %d",fCurrentEvent);
133       fTree->GetEvent(fCurrentEvent);
134       Info("ReadNext","Getting event %d Done",fCurrentEvent);
135       
136       Int_t retval = 0;
137       if (fReadRec && fReadSim)
138        {
139          retval = ReadRecAndSim();
140        }
141       else
142        {
143          if (fReadRec) retval = ReadRec();
144          if (fReadSim) retval = ReadSim();
145        } 
146
147       fCurrentEvent++;
148       if (retval != 0) 
149         {
150           //something wrong has happend during reading this event, take next
151           continue;
152         }
153
154       fNEventsRead++;
155       return retval;//success -> read one event
156       
157     }while(fCurrentDir < GetNumberOfDirs());//end of loop over directories specified in fDirs Obj Array  
158    
159   return 1; //no more directories to read
160   
161   
162 }
163 /********************************************************************/
164
165 Int_t AliReaderAOD::ReadRecAndSim()
166 {
167 //Reads raconstructed and simulated data 
168
169  Info("ReadRecAndSim","Found %d reconstructed tracks and %d simulated particles",
170        fRecBuffer->GetNumberOfParticles(),fSimBuffer->GetNumberOfParticles());
171
172  if (fCuts->GetEntriesFast() == 0x0)
173   {//if there is no cuts we return pointer to the buffer
174     if (fEventRec != fRecBuffer)
175      {
176        delete fEventRec;
177        delete fEventSim;
178      }
179     fEventRec = fRecBuffer;//fEventRec is the pointer that the user gets when he asks about an event
180     fEventSim = fSimBuffer;
181   }
182  else
183   {//if there are cuts specified
184     if ( (fEventRec == 0x0) || (fEventRec == fRecBuffer) )
185      {//we need to create a new event, if it is not existing  or it is the same as branch buffer
186        fEventRec = new AliAOD();
187        fEventSim = new AliAOD();
188
189        fEventRec->SetParticleClass( fRecBuffer->GetParticleClass() );
190        fEventSim->SetParticleClass( fSimBuffer->GetParticleClass() );
191      } 
192     else
193      {//or simply reset it in case it already exists
194        fEventRec->Reset();
195        fEventSim->Reset();
196      }
197
198     Int_t npart = fRecBuffer->GetNumberOfParticles();
199     
200     if (npart != fSimBuffer->GetNumberOfParticles())
201      {
202        Error("ReadRecAndSim","There is different number of simulated and reconstructed particles!",
203                               fSimBuffer->GetNumberOfParticles(),npart);
204        return 1;
205      } 
206     for (Int_t i = 0; i < npart; i++)
207      {
208        AliVAODParticle* prec = fRecBuffer->GetParticle(i);
209        AliVAODParticle* psim = fSimBuffer->GetParticle(i);
210        
211        if (prec == 0x0)
212         {
213           Error("ReadRecAndSim","Reconstructed Particle is NULL !!!");
214           continue;
215         }
216        if (psim == 0x0)
217         {
218           Error("ReadRecAndSim","Simulated Particle is NULL !!!");
219           continue;
220         }
221        
222        if (Rejected(prec)) continue;//we make cuts only on reconstructed data
223        
224        fEventRec->AddParticle(prec);
225        fEventSim->AddParticle( fSimBuffer->GetParticle(i));
226      }
227   }
228
229  Info("ReadRecAndSim","Read %d reconstructed tracks and %d simulated particles",
230        fEventRec->GetNumberOfParticles(),fEventSim->GetNumberOfParticles());
231  
232  fTrackCounter->Fill(fEventRec->GetNumberOfParticles());
233  
234  return 0;
235 }
236 /********************************************************************/
237
238 Int_t AliReaderAOD::ReadRec()
239 {
240 //Reads reconstructed data only
241
242  Info("ReadRec","Found %d reconstructed tracks",fRecBuffer->GetNumberOfParticles());
243
244  if (fCuts->GetEntriesFast() == 0x0)
245   {//if there is no cuts we return pointer to the buffer
246     if (fEventRec != fRecBuffer)
247      {
248        delete fEventRec;
249      }
250     fEventRec = fRecBuffer;//fEventRec is the pointer that the user gets when he asks about an event
251   }
252  else
253   {//if there are cuts specified
254     if ( (fEventRec == 0x0) || (fEventRec == fRecBuffer) )
255      {//we need to create a new event, if it is not existing  or it is the same as branch buffer
256        fEventRec = new AliAOD();
257
258        fEventRec->SetParticleClass( fRecBuffer->GetParticleClass() );
259      } 
260     else
261      {//or simply reset it in case it already exists
262        fEventRec->Reset();
263      }
264
265     Int_t npart = fRecBuffer->GetNumberOfParticles();
266     for (Int_t i = 0; i < npart; i++)
267      {
268        AliVAODParticle* prec = fRecBuffer->GetParticle(i);
269        if (Rejected(prec)) continue;//we make cuts only on simulated data
270
271        fEventRec->AddParticle(prec);
272      }
273   }
274
275  Info("ReadRec","Read %d reconstructed tracks",fEventRec->GetNumberOfParticles());
276  fTrackCounter->Fill(fEventRec->GetNumberOfParticles());
277
278  return 0;
279 }
280 /********************************************************************/
281
282 Int_t AliReaderAOD::ReadSim()
283 {
284 //Reads simulated data only
285
286  Info("ReadSim","Found %d simulated particles",fSimBuffer->GetNumberOfParticles());
287
288  if (fCuts->GetEntriesFast() == 0x0)
289   {//if there is no cuts we return pointer to the buffer
290     if (fEventSim != fSimBuffer)
291      {
292        delete fEventSim;
293      }
294     fEventSim = fSimBuffer;
295   }
296  else
297   {//if there are cuts specified
298     if ( (fEventSim == 0x0) || (fEventSim == fSimBuffer) )
299      {//we need to create a new event, if it is not existing  or it is the same as branch buffer
300        fEventSim = new AliAOD();
301
302        fEventSim->SetParticleClass( fSimBuffer->GetParticleClass() );
303      } 
304     else
305      {//or simply reset it in case it already exists
306        fEventSim->Reset();
307      }
308
309     Int_t npart = fSimBuffer->GetNumberOfParticles();
310     for (Int_t i = 0; i < npart; i++)
311      {
312        AliVAODParticle* prec = fSimBuffer->GetParticle(i);
313        if (Rejected(prec)) continue;//we make cuts only on simulated data
314        fEventSim->AddParticle(prec);
315      }
316   }
317
318  Info("ReadSim","Read %d simulated particles",fEventSim->GetNumberOfParticles());
319  fTrackCounter->Fill(fEventSim->GetNumberOfParticles());
320
321
322  return 0;
323 }
324 /********************************************************************/
325
326 Int_t AliReaderAOD::OpenFile(Int_t n)
327 {
328 //opens fFile with  tree
329
330 // Info("ReadNext","Opening File %d",n);
331  const TString dirname = GetDirName(n);
332  if (dirname == "")
333   {
334     AliDebug(3,"Got empty string as a directory name."); 
335     return 1;
336   }
337  
338  TString filename = dirname +"/"+ fFileName;
339  fFile = TFile::Open(filename.Data()); 
340  if ( fFile == 0x0)
341   {
342     Error("OpenFile","Can't open fFile %s",filename.Data());
343     return 2;
344   }
345  if (!fFile->IsOpen())
346   {
347     Error("OpenFile","Can't open fFile  %s",filename.Data());
348     delete fFile;
349     fFile = 0x0;
350     return 3;
351   }
352
353  Info("ReadNext","File %s Is Opened, Getting the TREE",filename.Data());
354  
355  fTree = dynamic_cast<TTree*>(fFile->Get(fgkTreeName));
356  if (fTree == 0x0)
357   {
358     AliDebug(3,Form("Can not find TTree object named %s",fgkTreeName.Data()));
359     delete fFile;
360     fFile = 0x0;
361     return 4;
362   }
363
364 //  Info("ReadNext","Got TREE, Setting branch addresses");
365
366   if (fReadRec)
367    {
368      TBranch* branch = fTree->GetBranch(fgkReconstructedDataBranchName);
369      if (branch == 0x0)
370       {
371         Error("OpenFile","Can not find branch %s in file %s",
372               fgkReconstructedDataBranchName.Data(),filename.Data());
373         
374         delete fTree;
375         fTree = 0x0;
376         delete fFile;
377         fFile = 0x0;
378         return 5;
379       }
380      fTree->SetBranchAddress(fgkReconstructedDataBranchName,&fRecBuffer);
381    }
382   
383
384   if (fReadSim)
385    {
386      TBranch* branch = fTree->GetBranch(fgkSimulatedDataBranchName);
387      if (branch == 0x0)
388       {
389         Error("OpenFile","Can not find branch %s in file %s",
390               fgkSimulatedDataBranchName.Data(),filename.Data());
391         
392         delete fTree;
393         fTree = 0x0;
394         delete fFile;
395         fFile = 0x0;
396         return 6;
397       }
398      fTree->SetBranchAddress(fgkSimulatedDataBranchName,&fSimBuffer);
399    }
400 //  Info("ReadNext","Got TREE, Addresses are set.");
401 //  Info("ReadNext","Quitting the method.");
402   
403   return 0;
404  
405 }
406 /********************************************************************/
407
408 Int_t AliReaderAOD::WriteAOD(AliReader* reader, const char* outfilename, const char* pclassname,  Bool_t /*multcheck*/)
409 {
410 //reads tracks from runs and writes them to file
411   ::Info("AliReaderAOD::Write","________________________________________________________");
412   ::Info("AliReaderAOD::Write","________________________________________________________");
413   ::Info("AliReaderAOD::Write","________________________________________________________");
414   
415   if (reader == 0x0)
416    {
417      ::Error("AliReaderAOD::Write","Input Reader is NULL");
418      return -1;
419    }
420   TFile *outfile = TFile::Open(outfilename,"recreate");
421   if (outfile == 0x0)
422    {
423      ::Error("AliReaderAOD::Write","Can not open output file %s",outfilename);
424      return -1;
425    }
426
427   TTree *tree = new TTree(fgkTreeName,"Tree with tracks");
428   
429   TBranch *recbranch = 0x0, *simbranch = 0x0;
430   
431   AliAOD* eventrec = new AliAOD();//must be created before Branch is called. Otherwise clones array is not splitted
432   AliAOD* eventsim = new AliAOD();//AOD together with fParticles clones array knowing exact type of particles
433   
434   eventrec->SetParticleClassName(pclassname);
435   eventsim->SetParticleClassName(pclassname);
436  
437   AliAOD* recbuffer = eventrec;
438   AliAOD* simbuffer = eventsim;
439   
440   if (reader->ReadsRec()) recbranch = tree->Branch(fgkReconstructedDataBranchName,"AliAOD",&recbuffer,32000,99);
441   if (reader->ReadsSim()) simbranch = tree->Branch(fgkSimulatedDataBranchName,"AliAOD",&simbuffer,32000,99);
442
443   reader->Rewind();
444   while (reader->Next() == kFALSE)
445    {
446      
447      if (reader->ReadsRec())
448       {//here we can get AOD that has different particle type
449         AliAOD* event = reader->GetEventRec();
450         if ( eventrec->GetParticleClass() != event->GetParticleClass() )
451          {//if class type is not what what we whant we copy particles
452            eventrec->CopyData(event);
453            recbuffer = eventrec;
454          }
455         else
456          {//else just pointer to event from input reader is passed
457            recbuffer = event;
458          } 
459       }
460
461      if (reader->ReadsSim())
462       {
463         AliAOD* event = reader->GetEventSim();
464         if ( eventsim->GetParticleClass() != event->GetParticleClass() )
465          {//if class type is not what what we whant we copy particles
466            eventsim->CopyData(event);
467            simbuffer = eventrec;
468          }
469         else
470          {//else just pointer to event from input reader is passed
471            simbuffer = event;
472          } 
473       }
474      tree->Fill();
475    }
476   
477   ::Info("AliReaderAOD::Write","Written %d events",tree->GetEntries());
478   outfile->cd();
479   tree->Write();
480
481   delete eventsim;
482   delete eventrec;
483   
484   delete tree;
485   delete outfile;
486   return 0; 
487 }
488