d493cfac47830d43e079e58e20f49cc399c4432e
[u/mrichter/AliRoot.git] / PHOS / AliPHOSTrackSegmentMakerv1.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 /* $Id$ */
16 //_________________________________________________________________________
17 // Implementation version 1 of algorithm class to construct PHOS track segments
18 // Track segment for PHOS is list of 
19 //        EMC RecPoint + (possibly) CPV RecPoint + (possibly) PPSD RecPoint
20 // To find TrackSegments we do the following: for each EMC RecPoints we look at
21 // CPV/PPSD RecPoints in the radious fR0. If there is such a CPV RecPoint, 
22 // we make "Link" it is just indexes of EMC and CPV/PPSD RecPoint and distance
23 // between them in the PHOS plane. Then we sort "Links" and starting from the 
24 // least "Link" pointing to the unassined EMC and CPV RecPoints assing them to 
25 // new TrackSegment. If there is no CPV/PPSD RecPoint we make TrackSegment 
26 // consisting from EMC along. There is no TrackSegments without EMC RecPoint.
27 //
28 // In principle this class should be called from AliPHOSReconstructioner, but 
29 // one can use it as well in standalong mode.
30 // User case:
31 // root [0] AliPHOSTrackSegmentMakerv1 * t = new AliPHOSTrackSegmentMaker("galice.root")
32 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
33 // root [1] t->ExecuteTask()
34 // root [2] t->SetMaxEmcPpsdDistance(5)
35 // root [3] t->SetTrackSegmentsBranch("max distance 5 cm")
36 // root [4] t->ExecuteTask("deb all time") 
37 //                 
38 //*-- Author: Dmitri Peressounko (RRC Ki & SUBATECH)
39 //
40
41 // --- ROOT system ---
42 #include "TROOT.h"
43 #include "TFile.h"
44 #include "TTree.h"
45 #include "TSystem.h"
46 #include "TBenchmark.h"
47 // --- Standard library ---
48
49 #include <iostream.h>
50 #include <iomanip.h>
51
52 // --- AliRoot header files ---
53
54 #include "AliPHOSTrackSegmentMakerv1.h"
55 #include "AliPHOSClusterizerv1.h"
56 #include "AliPHOSTrackSegment.h"
57 #include "AliPHOSCpvRecPoint.h"
58 #include "AliPHOSPpsdRecPoint.h"
59 #include "AliPHOSLink.h"
60 #include "AliPHOSv0.h"
61 #include "AliRun.h"
62
63 ClassImp( AliPHOSTrackSegmentMakerv1) 
64
65
66 //____________________________________________________________________________
67  AliPHOSTrackSegmentMakerv1::AliPHOSTrackSegmentMakerv1() : AliPHOSTrackSegmentMaker()
68 {
69   // ctor
70   SetTitle("version 1") ;
71   SetName("AliPHOSTrackSegmentMaker") ;
72   fR0 = 10. ;   
73   fEmcFirst = 0 ;    
74   fEmcLast  = 0 ;   
75   fCpvFirst = 0 ;   
76   fCpvLast  = 0 ;   
77   fPpsdFirst= 0 ;   
78   fPpsdLast = 0 ;   
79   fLinkLowArray = 0 ;
80   fLinkUpArray  = 0 ;
81   fIsInitialized = kFALSE ;
82 }
83 //____________________________________________________________________________
84  AliPHOSTrackSegmentMakerv1::AliPHOSTrackSegmentMakerv1(const char* headerFile, const char* branchTitle): 
85 AliPHOSTrackSegmentMaker()
86 {
87   // ctor
88   SetTitle("version 1") ;
89   SetName("AliPHOSTrackSegmentMaker") ;
90   fR0 = 10. ;   
91   fEmcFirst = 0 ;    
92   fEmcLast  = 0 ;   
93   fCpvFirst = 0 ;   
94   fCpvLast  = 0 ;   
95   fPpsdFirst= 0 ;   
96   fPpsdLast = 0 ;   
97
98   fHeaderFileName = headerFile ;
99   fRecPointsBranchTitle = branchTitle ;
100     
101   TFile * file = (TFile*) gROOT->GetFile(fHeaderFileName.Data() ) ;
102   
103   if(file == 0){
104     file = new TFile(fHeaderFileName.Data(),"update") ;
105     gAlice = (AliRun *) file->Get("gAlice") ;
106   }
107   
108   AliPHOS * phos = (AliPHOS *) gAlice->GetDetector("PHOS") ;    
109   fGeom  = AliPHOSGeometry::GetInstance(phos->GetGeometry()->GetName(),phos->GetGeometry()->GetTitle() );
110   
111   fEmcRecPoints = new TObjArray(200) ;
112   fCpvRecPoints = new TObjArray(200) ;
113   fClusterizer  = new AliPHOSClusterizerv1() ;
114   
115   fTrackSegments = new TClonesArray("AliPHOSTrackSegment",200) ;
116   
117   fLinkLowArray = new TClonesArray("AliPHOSLink", 1000);
118   fLinkUpArray  = new TClonesArray("AliPHOSLink", 1000); 
119   
120   fIsInitialized = kTRUE ;
121
122 }
123 //____________________________________________________________________________
124 void  AliPHOSTrackSegmentMakerv1::Init(){
125   //Make all memory allokations not possible in default constructor
126
127   if(!fIsInitialized){
128     if(fHeaderFileName.IsNull())
129       fHeaderFileName = "galice.root" ;
130     
131     
132     TFile * file = (TFile*) gROOT->GetFile(fHeaderFileName.Data() ) ;
133     
134     if(file == 0){
135       file = new TFile(fHeaderFileName.Data(),"update") ;
136       gAlice = (AliRun *) file->Get("gAlice") ;
137     }
138     
139     AliPHOS * phos = (AliPHOS *) gAlice->GetDetector("PHOS") ;    
140     fGeom  = AliPHOSGeometry::GetInstance(phos->GetGeometry()->GetName(),phos->GetGeometry()->GetTitle() );
141
142
143     fEmcRecPoints = new TObjArray(200) ;
144     fCpvRecPoints = new TObjArray(200) ;
145     fClusterizer  = new AliPHOSClusterizerv1() ;
146
147     
148     fTrackSegments = new TClonesArray("AliPHOSTrackSegment",200) ;
149
150     fLinkLowArray = new TClonesArray("AliPHOSLink", 1000);
151     fLinkUpArray  = new TClonesArray("AliPHOSLink", 1000); 
152     
153     fIsInitialized = kTRUE ;
154    }
155 }
156
157 //____________________________________________________________________________
158  AliPHOSTrackSegmentMakerv1::~AliPHOSTrackSegmentMakerv1()
159
160   // dtor
161   if(fLinkLowArray) delete fLinkLowArray ;
162   if(fLinkUpArray)  delete fLinkUpArray  ;
163 }
164
165 //____________________________________________________________________________
166 void  AliPHOSTrackSegmentMakerv1::FillOneModule()
167 {
168   // Finds first and last indexes between which 
169   // clusters from one PHOS module are
170  
171
172   //First EMC clusters
173   Int_t totalEmc = fEmcRecPoints->GetEntriesFast() ;
174   for(fEmcFirst = fEmcLast; (fEmcLast < totalEmc) &&  
175         (((AliPHOSRecPoint *) fEmcRecPoints->At(fEmcLast))->GetPHOSMod() == fModule ); 
176       fEmcLast ++)  ;
177   
178   
179   //Now CPV clusters
180   Int_t totalCpv = fCpvRecPoints->GetEntriesFast() ;
181
182   if(fModule <= fGeom->GetNCPVModules()){ // in CPV geometry
183     
184     for(fCpvFirst = fCpvLast; (fCpvLast < totalCpv) && 
185           (((AliPHOSRecPoint *) fCpvRecPoints->At(fCpvLast))->GetPHOSMod() == fModule ); 
186         fCpvLast ++) ;
187     
188     fPpsdFirst = fCpvLast ; //To avoid scanning RecPoints between fPpsdFirst and fPpsdLast
189     fPpsdLast  = fCpvLast ; //and to be ready to switch to mixed geometry 
190   }
191   else{  //in PPSD geometry    
192     fCpvLast = fPpsdLast ;
193     //Upper layer first
194     for(fCpvFirst = fCpvLast; (fCpvLast < totalCpv) &&  
195           (((AliPHOSPpsdRecPoint *) fCpvRecPoints->At(fCpvLast))->GetPHOSMod() == fModule ) &&
196           (((AliPHOSPpsdRecPoint *) fCpvRecPoints->At(fCpvLast))->GetUp()) ; 
197         fCpvLast ++)  ;
198     
199     fPpsdLast= fCpvLast ;
200     for(fPpsdFirst = fPpsdLast; (fPpsdLast < totalCpv)  &&
201           (((AliPHOSPpsdRecPoint *) fCpvRecPoints->At(fPpsdLast))->GetPHOSMod() == fModule ) &&
202           (!((AliPHOSPpsdRecPoint *) fCpvRecPoints->At(fPpsdLast))->GetUp()) ; 
203         fPpsdLast ++) ;
204   }
205     
206 }
207 //____________________________________________________________________________
208 Float_t  AliPHOSTrackSegmentMakerv1::GetDistanceInPHOSPlane(AliPHOSEmcRecPoint * emcClu,AliPHOSRecPoint * cpvClu, Bool_t &toofar)
209 {
210   // Calculates the distance between the EMC RecPoint and the PPSD RecPoint
211   //clusters are sorted in "rows" and "columns" of width 1 cm
212
213   Float_t delta = 1 ;  // Width of the rows in sorting of RecPoints (in cm)
214                        // if you change this value, change it as well in xxxRecPoint::Compare()
215   Float_t r = fR0 ;
216  
217   TVector3 vecEmc ;
218   TVector3 vecCpv ;
219   
220   emcClu->GetLocalPosition(vecEmc) ;
221   cpvClu->GetLocalPosition(vecCpv)  ; 
222
223   if(emcClu->GetPHOSMod() == cpvClu->GetPHOSMod()){ 
224     if(vecCpv.X() <= vecEmc.X() + fR0 + 2*delta ){ 
225
226       vecCpv = vecCpv  - vecEmc ; 
227       r = vecCpv.Mag() ;
228       toofar = kFALSE ;
229
230     } // if  xPpsd >= xEmc + ...
231     else 
232       toofar = kTRUE ;
233   } 
234   else 
235     toofar = kTRUE ;
236
237   //toofar = kFALSE ;
238  
239   
240   return r ;
241 }
242
243 //____________________________________________________________________________
244 void  AliPHOSTrackSegmentMakerv1::MakeLinks()
245
246   // Finds distances (links) between all EMC and PPSD clusters, 
247   // which are not further apart from each other than fR0 
248   // and sort them in accordance with this distance
249   
250   fLinkUpArray->Clear() ;    
251   fLinkLowArray->Clear() ;
252
253   AliPHOSRecPoint * ppsd ; 
254   AliPHOSRecPoint * cpv ;
255   AliPHOSEmcRecPoint * emcclu ;
256
257   Int_t iLinkLow = 0 ;
258   Int_t iLinkUp  = 0 ;
259   
260   Int_t iEmcRP;
261   for(iEmcRP = fEmcFirst; iEmcRP < fEmcLast; iEmcRP++ ) {
262     emcclu = (AliPHOSEmcRecPoint *) fEmcRecPoints->At(iEmcRP) ;
263
264     Bool_t toofar ;    
265     Int_t iPpsd ;
266     for(iPpsd = fPpsdFirst; iPpsd < fPpsdLast;iPpsd++ ) {
267       
268       ppsd = (AliPHOSRecPoint *) fCpvRecPoints->At(iPpsd) ;
269       Float_t r = GetDistanceInPHOSPlane(emcclu, ppsd, toofar) ;
270
271       if(toofar) 
272         break ;  
273       if(r < fR0)
274         new ((*fLinkLowArray)[iLinkLow++])  AliPHOSLink(r, iEmcRP, iPpsd) ;
275     }
276     
277     Int_t iCpv = 0 ;    
278     for(iCpv = fCpvFirst; iCpv < fCpvLast;iCpv++ ) { 
279       
280       cpv = (AliPHOSRecPoint *) fCpvRecPoints->At(iCpv) ;
281       Float_t r = GetDistanceInPHOSPlane(emcclu, cpv, toofar) ;
282       
283       if(toofar)
284         break ;  
285       if(r < fR0) { 
286         new ((*fLinkUpArray)[iLinkUp++])  AliPHOSLink(r, iEmcRP, iCpv) ;
287       }      
288     }
289   } 
290   
291   fLinkLowArray->Sort() ; //first links with smallest distances
292   fLinkUpArray->Sort() ;
293 }
294
295 //____________________________________________________________________________
296 void  AliPHOSTrackSegmentMakerv1::MakePairs()
297
298   // Using the previously made list of "links", we found the smallest link - i.e. 
299   // link with the least distance betwing EMC and CPV and pointing to still 
300   // unassigned RecParticles. We assign these RecPoints to TrackSegment and 
301   // remove them from the list of "unassigned". 
302   
303   //Make arrays to mark clusters already chousen
304   Int_t * emcExist = 0;
305   if(fEmcLast > fEmcFirst)
306     emcExist = new Int_t[fEmcLast-fEmcFirst] ;
307   
308   Int_t index;
309   for(index = 0; index <fEmcLast-fEmcFirst; index ++)
310     emcExist[index] = 1 ;
311   
312   Bool_t * cpvExist = 0;
313   if(fCpvLast > fCpvFirst)
314     cpvExist = new Bool_t[fCpvLast-fCpvFirst] ;
315   for(index = 0; index <fCpvLast-fCpvFirst; index ++)
316     cpvExist[index] = kTRUE ;
317   
318   Bool_t * ppsdExist = 0;
319   if(fPpsdLast > fPpsdFirst)
320     ppsdExist = new Bool_t[fPpsdLast-fPpsdFirst] ;
321   for(index = 0; index <fPpsdLast-fPpsdFirst; index ++)
322     ppsdExist[index] = kTRUE ;
323   
324   // Finds the smallest links and makes pairs of CPV and EMC clusters with smallest distance 
325   TIter nextLow(fLinkLowArray) ;
326   TIter nextUp(fLinkUpArray) ;
327   
328   AliPHOSLink * linkLow ;
329   AliPHOSLink * linkUp ;
330
331
332   AliPHOSRecPoint * nullpointer = 0 ;
333
334   while ( (linkLow =  (AliPHOSLink *)nextLow() ) ){
335   
336     if( (emcExist[linkLow->GetEmc()-fEmcFirst]> 0) && 
337         ppsdExist[linkLow->GetPpsd()-fPpsdFirst]  ){ // RecPoints not removed yet 
338       new ((*fTrackSegments)[fNTrackSegments]) AliPHOSTrackSegment((AliPHOSEmcRecPoint *) fEmcRecPoints->At(linkLow->GetEmc()), 
339                                                  nullpointer, 
340                                                 (AliPHOSRecPoint *)fCpvRecPoints->At(linkLow->GetPpsd()) ) ;
341          
342       ((AliPHOSTrackSegment* )fTrackSegments->At(fNTrackSegments))->SetIndexInList(fNTrackSegments);    
343       //replace index of emc to negative and shifted index of TS      
344       emcExist[linkLow->GetEmc()-fEmcFirst] = -2 - fNTrackSegments ;  
345       //mark ppsd as used
346       ppsdExist[linkLow->GetPpsd()-fPpsdFirst] = kFALSE ; 
347       fNTrackSegments++ ;
348     } 
349   } 
350          
351
352   while ( (linkUp =  (AliPHOSLink *)nextUp() ) ){  
353     if(emcExist[linkUp->GetEmc()-fEmcFirst] != -1){ //without ppsd Up yet 
354
355       if(cpvExist[linkUp->GetPpsd()-fCpvFirst]){ //CPV still exist
356         
357         if(emcExist[linkUp->GetEmc()-fEmcFirst] > 0){ //without ppsd Low => create new TS
358
359           new ((* fTrackSegments)[fNTrackSegments]) AliPHOSTrackSegment((AliPHOSEmcRecPoint *) fEmcRecPoints->At(linkUp->GetEmc()) , 
360                                                                       (AliPHOSRecPoint *)fCpvRecPoints->At(linkUp->GetPpsd()), 
361                                                                       nullpointer) ;
362           ((AliPHOSTrackSegment *) fTrackSegments->At(fNTrackSegments))->SetIndexInList(fNTrackSegments);
363           fNTrackSegments++ ;
364         }
365         else{ // append ppsd Up to existing TS
366           ((AliPHOSTrackSegment *)fTrackSegments->At(-2-emcExist[linkUp->GetEmc()-fEmcFirst]))->SetCpvRecPoint((AliPHOSCpvRecPoint *)fCpvRecPoints->At(linkUp->GetPpsd()));
367         }
368
369         emcExist[linkUp->GetEmc()-fEmcFirst] = -1 ; //Mark emc  that Cpv was found 
370         //mark CPV recpoint as already used 
371         cpvExist[linkUp->GetPpsd()-fCpvFirst] = kFALSE ;
372       } //if ppsdUp still exist
373     } 
374   }      
375
376   //look through emc recPoints left without CPV/PPSD
377   if(emcExist){ //if there is emc rec point
378     Int_t iEmcRP ;
379     for(iEmcRP = 0; iEmcRP < fEmcLast-fEmcFirst  ; iEmcRP++ ){
380       if(emcExist[iEmcRP] > 0 ){
381         new ((*fTrackSegments)[fNTrackSegments])  AliPHOSTrackSegment((AliPHOSEmcRecPoint *)fEmcRecPoints->At(iEmcRP+fEmcFirst), 
382                                                                     nullpointer, 
383                                                                     nullpointer ) ;
384         ((AliPHOSTrackSegment *) fTrackSegments->At(fNTrackSegments))->SetIndexInList(fNTrackSegments);
385         fNTrackSegments++;    
386       } 
387     }
388   }
389   
390 }
391
392 //____________________________________________________________________________
393 void  AliPHOSTrackSegmentMakerv1::Exec(Option_t * option)
394 {
395   //STEERing method
396
397   if(! fIsInitialized) Init() ;
398
399   if(strstr(option,"tim"))
400     gBenchmark->Start("PHOSTSMaker");  
401
402   Int_t nEvents = (Int_t) gAlice->TreeE()->GetEntries() ;
403   
404   for(fEvent = 0;fEvent< nEvents; fEvent++){
405     if(!ReadRecPoints())  //reads RecPoints for event fEvent
406       return;
407     
408     for(fModule = 1; fModule <= fGeom->GetNModules() ; fModule++ ){
409       
410       FillOneModule() ; 
411       
412       MakeLinks() ;
413       
414       MakePairs() ;
415       
416     }
417
418     WriteTrackSegments() ;
419     if(strstr(option,"deb"))
420       PrintTrackSegments(option) ;
421   }
422
423   if(strstr(option,"tim")){
424     gBenchmark->Stop("PHOSTSMaker");
425     cout << "AliPHOSTSMaker:" << endl ;
426     cout << "  took " << gBenchmark->GetCpuTime("PHOSTSMaker") << " seconds for making TS " 
427          <<  gBenchmark->GetCpuTime("PHOSTSMaker")/nEvents << " seconds per event " << endl ;
428     cout << endl ;
429   }
430
431
432 }
433 //____________________________________________________________________________
434 void AliPHOSTrackSegmentMakerv1::Print(Option_t * option)const {
435   if(fIsInitialized){
436     cout <<  "======== AliPHOSTrackSegmentMakerv1 ========" << endl ;
437     cout <<  "Making Track segments "<< endl ;
438     cout <<  "    Headers file:                   " << fHeaderFileName.Data() << endl ;
439     cout <<  "    RecPoints branch file name:     " << fRecPointsBranchTitle.Data() << endl ;
440     cout <<  "    TrackSegments Branch file name: " << fTSBranchTitle.Data() << endl ;
441     cout <<  "with parameters: " << endl ;
442     cout <<  "    Maximal EMC - CPV (PPSD) distance (cm)" << fR0 << endl ;
443     cout <<  "============================================" << endl ;
444   }
445   else
446     cout << "AliPHOSTrackSegmentMakerv1 not initialized " << endl ;
447 }
448 //____________________________________________________________________________
449 Bool_t AliPHOSTrackSegmentMakerv1::ReadRecPoints(){
450   // Reads Emc and CPV recPoints with given title (fRecPointsBranchTitle) 
451   // made priveously with Clusterizer.
452
453
454   //Make some initializations 
455   fEmcRecPoints->Clear() ;
456   fCpvRecPoints->Clear() ;
457   fTrackSegments->Clear() ;
458   fNTrackSegments = 0 ;
459   fEmcFirst = 0 ;    
460   fEmcLast  = 0 ;   
461   fCpvFirst = 0 ;   
462   fCpvLast  = 0 ;   
463   fPpsdFirst= 0 ;   
464   fPpsdLast = 0 ;   
465
466
467   gAlice->GetEvent(fEvent) ;
468
469   // Get TreeR header from file
470   if(gAlice->TreeR()==0){
471     char treeName[20]; 
472     sprintf(treeName,"TreeR%d",fEvent);
473     cout << "Error in AliPHOSTrackSegmentMakerv1 : no "<<treeName << endl  ;
474     cout << "   Do nothing " << endl ;
475     return kFALSE ;
476   }
477
478
479   //Find RecPoints with title fRecPointsBranchTitle
480   TBranch * emcBranch = 0;
481   TBranch * cpvBranch = 0;
482   TBranch * clusterizerBranch = 0;
483
484   TObjArray * branches = gAlice->TreeR()->GetListOfBranches() ;
485   Int_t ibranch;
486   Bool_t emcNotFound = kTRUE ;
487   Bool_t cpvNotFound = kTRUE ;  
488   Bool_t clusterizerNotFound = kTRUE ;
489   
490   for(ibranch = 0;ibranch <branches->GetEntries();ibranch++){
491
492     if(emcNotFound){
493       emcBranch=(TBranch *) branches->At(ibranch) ;
494       if( fRecPointsBranchTitle.CompareTo(emcBranch->GetTitle())==0 )
495         if( strcmp(emcBranch->GetName(),"PHOSEmcRP") == 0) {
496           emcNotFound = kFALSE ;
497         }
498     }
499     
500     if(cpvNotFound){
501       cpvBranch=(TBranch *) branches->At(ibranch) ;
502       if( fRecPointsBranchTitle.CompareTo(cpvBranch->GetTitle())==0 )
503         if( strcmp(cpvBranch->GetName(),"PHOSCpvRP") == 0) 
504           cpvNotFound = kFALSE ;
505     }
506     
507     if(clusterizerNotFound){
508       clusterizerBranch = (TBranch *) branches->At(ibranch) ;
509       if( fRecPointsBranchTitle.CompareTo(clusterizerBranch->GetTitle()) == 0)
510         if( strcmp(clusterizerBranch->GetName(),"AliPHOSClusterizer") == 0) 
511           clusterizerNotFound = kFALSE ;
512     }
513     
514   }
515
516   if(clusterizerNotFound || emcNotFound || cpvNotFound){
517     cout << "AliPHOSTrackSegmentMakerv1: " << endl ;
518     cout << "    Can't find Branch with RecPoints or Clusterizer " ;
519     cout << "    Do nothing" <<endl  ;
520     return kFALSE ;
521   }
522   
523   emcBranch->SetAddress(&fEmcRecPoints) ;
524   cpvBranch->SetAddress(&fCpvRecPoints) ;
525   clusterizerBranch->SetAddress(&fClusterizer) ;
526   
527   gAlice->TreeR()->GetEvent(0) ;
528   
529   return kTRUE ;
530   
531 }
532 //____________________________________________________________________________
533 void AliPHOSTrackSegmentMakerv1::WriteTrackSegments(){
534   // Writes found TrackSegments to TreeR. Creates branches 
535   // "PHOSTS" and "AliPHOSTrackSegmentMaker" with the same title.
536   // In the former branch found TrackSegments are stored, while 
537   // in the latter all parameters, with which TS were made. 
538   // ROOT does not allow overwriting existing branches, therefore
539   // first we chesk, if branches with the same title alredy exist.
540   // If yes - exits without writing.
541   
542   //First, check, if branches already exist
543   TBranch * tsMakerBranch = 0;
544   TBranch * tsBranch = 0;
545   
546   TObjArray * branches = gAlice->TreeR()->GetListOfBranches() ;
547   Int_t ibranch;
548   Bool_t tsMakerNotFound = kTRUE ;
549   Bool_t tsNotFound = kTRUE ;
550   
551   for(ibranch = 0;(ibranch <branches->GetEntries())&&(tsMakerNotFound||tsNotFound);ibranch++){
552     if(tsMakerNotFound){
553       tsMakerBranch=(TBranch *) branches->At(ibranch) ;
554       if( (strcmp(tsMakerBranch->GetName(),"AliPHOSTrackSegmentMaker") == 0) &&
555           (fTSBranchTitle.CompareTo( tsMakerBranch->GetTitle())==0 ))
556         tsMakerNotFound = kFALSE ;
557     }
558     if(tsNotFound){
559       tsBranch=(TBranch *) branches->At(ibranch) ;
560       if( (strcmp(tsBranch->GetName(),"PHOSTS") == 0)  &&
561           (fTSBranchTitle.CompareTo( tsBranch->GetTitle())==0 ))
562         tsNotFound = kFALSE ;
563     }
564   }
565
566   if(!(tsMakerNotFound && tsNotFound )){ 
567     cout << "AliPHOSTrackSegmentMakerv1 error:"<< endl ;
568     cout << "       Branches PHOSTS and AliPHOSTrackSegementMaker " << endl ;
569     cout << "       with title '"<<fTSBranchTitle.Data() << "' already exist " << endl ;
570     cout << "       can not overwrite " << endl ;
571     return ;
572   }
573
574   //Make branch in TreeR for TrackSegments 
575   char * filename = 0;
576   if(gSystem->Getenv("CONFIG_SPLIT_FILE")!=0){   //generating file name
577     filename = new char[strlen(gAlice->GetBaseFile())+20] ;
578     sprintf(filename,"%s/PHOS.Reco.root",gAlice->GetBaseFile()) ; 
579   }
580
581   TDirectory *cwd = gDirectory;
582   
583   //First TS
584   Int_t bufferSize = 32000 ;    
585   tsBranch = gAlice->TreeR()->Branch("PHOSTS",&fTrackSegments,bufferSize);
586   tsBranch->SetTitle(fTSBranchTitle.Data());
587   if (filename) {
588     tsBranch->SetFile(filename);
589     TIter next( tsBranch->GetListOfBranches());
590     while ((tsBranch=(TBranch*)next())) {
591       tsBranch->SetFile(filename);
592     }   
593     cwd->cd();
594   } 
595   
596   //Second -TSMaker
597   Int_t splitlevel = 0 ;
598   AliPHOSTrackSegmentMakerv1 * ts = this ;
599   tsMakerBranch = gAlice->TreeR()->Branch("AliPHOSTrackSegmentMaker","AliPHOSTrackSegmentMakerv1",
600                                           &ts,bufferSize,splitlevel);
601   tsMakerBranch->SetTitle(fTSBranchTitle.Data());
602   if (filename) {
603     tsMakerBranch->SetFile(filename);
604     TIter next( tsMakerBranch->GetListOfBranches());
605     while ((tsMakerBranch=(TBranch*)next())) {
606       tsMakerBranch->SetFile(filename);
607     }   
608     cwd->cd();
609   } 
610   
611   gAlice->TreeR()->Fill() ;    
612   gAlice->TreeR()->Write(0,kOverwrite) ;  
613   
614 }
615
616
617 //____________________________________________________________________________
618 void AliPHOSTrackSegmentMakerv1::PrintTrackSegments(Option_t * option){
619   // option deb - prints # of found TrackSegments
620   // option deb all - prints as well indexed of found RecParticles assigned to the TS
621
622   
623   cout << "AliPHOSTrackSegmentMakerv1: " << endl ;
624   cout << "       Found " << fTrackSegments->GetEntriesFast() << "  trackSegments " << endl ;
625   
626   if(strstr(option,"all")) {  // printing found TS
627     cout << "TrackSegment # " << "    EMC RP#    " << "    CPV RP#    " << "     PPSD RP#" << endl ; 
628     
629     Int_t index;
630     for (index = 0 ; index <fTrackSegments->GetEntriesFast() ; index++) {
631       AliPHOSTrackSegment * ts = (AliPHOSTrackSegment * )fTrackSegments->At(index) ; 
632       cout<<"   "<< setw(4) << ts->GetIndexInList() << "            " 
633           <<setw(4) << ts->GetEmcIndex()<< "            " 
634           <<setw(4) << ts->GetCpvIndex()<< "            " 
635           <<setw(4) << ts->GetPpsdIndex()<< endl ;
636     }   
637     
638     cout << "-------------------------------------------------------"<< endl ;
639   }
640 }