remove dependency to aliroot libraries, access of ESDEvent object through abstract...
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALReconstructor.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 //-- Yves Schutz (SUBATECH) 
21 // Reconstruction class. Redesigned from the old AliReconstructionner class and 
22 // derived from STEER/AliReconstructor. 
23 // 
24 //-- Aleksei Pavlinov : added staf for EMCAL jet trigger 9Apr 25, 2008)
25 //                    : fgDigitsArr should read just once at event
26
27 // --- ROOT system ---
28 #include <TList.h>
29 #include <TClonesArray.h>
30 #include <TH2.h>
31 #include "TGeoManager.h"
32 #include "TGeoMatrix.h"
33
34 // --- Standard library ---
35
36 // --- AliRoot header files ---
37 #include "AliEMCALReconstructor.h"
38
39 #include "AliCodeTimer.h"
40 #include "AliCaloCalibPedestal.h"
41 #include "AliEMCALCalibData.h"
42 #include "AliESDEvent.h"
43 #include "AliESDCaloCluster.h"
44 #include "AliESDCaloCells.h"
45 #include "AliESDtrack.h"
46 #include "AliEMCALLoader.h"
47 #include "AliEMCALRawUtils.h"
48 #include "AliEMCALDigit.h"
49 #include "AliEMCALClusterizerv1.h"
50 #include "AliEMCALClusterizerNxN.h"
51 #include "AliEMCALRecPoint.h"
52 #include "AliEMCALPID.h"
53 #include "AliEMCALTrigger.h"
54 #include "AliRawReader.h"
55 #include "AliCDBEntry.h"
56 #include "AliCDBManager.h"
57 #include "AliEMCALGeometry.h"
58 #include "AliEMCAL.h"
59 #include "AliESDVZERO.h"
60 #include "AliCDBManager.h"
61 #include "AliRunLoader.h"
62 #include "AliRun.h"
63 #include "AliEMCALTriggerData.h"
64 #include "AliEMCALTriggerElectronics.h"
65 #include "AliEMCALTriggerDCSConfigDB.h"
66 #include "AliEMCALTriggerDCSConfig.h"
67 #include "AliEMCALTriggerData.h"
68 #include "AliEMCALTriggerRawDigit.h"
69 #include "AliEMCALTriggerPatch.h"
70 #include "AliEMCALTriggerTypes.h"
71
72 ClassImp(AliEMCALReconstructor) 
73
74 const AliEMCALRecParam*     AliEMCALReconstructor::fgkRecParam        = 0;   // EMCAL rec. parameters
75 AliEMCALRawUtils*           AliEMCALReconstructor::fgRawUtils         = 0;   // EMCAL raw utilities class
76 AliEMCALClusterizer*        AliEMCALReconstructor::fgClusterizer      = 0;   // EMCAL clusterizer class
77 TClonesArray*               AliEMCALReconstructor::fgDigitsArr        = 0;   // list of digits, to be used multiple times
78 TObjArray*                  AliEMCALReconstructor::fgClustersArr      = 0;   // list of clusters, to be used multiple times
79 TClonesArray*               AliEMCALReconstructor::fgTriggerDigits    = 0;   // list of trigger digits, to be used multiple times
80 AliEMCALTriggerElectronics* AliEMCALReconstructor::fgTriggerProcessor = 0x0;
81 //____________________________________________________________________________
82 AliEMCALReconstructor::AliEMCALReconstructor() 
83   : fDebug(kFALSE), fList(0), fGeom(0),fCalibData(0),fPedestalData(0),fTriggerData(0x0) 
84 {
85   // ctor
86
87   fgRawUtils = new AliEMCALRawUtils;
88
89   //To make sure we match with the geometry in a simulation file,
90   //let's try to get it first.  If not, take the default geometry
91   AliRunLoader *rl = AliRunLoader::Instance();
92   if (rl->GetAliRun()){
93     AliEMCAL * emcal = dynamic_cast<AliEMCAL*>(rl->GetAliRun()->GetDetector("EMCAL"));
94     if(emcal) fGeom = emcal->GetGeometry();
95   }
96   
97   if(!fGeom) {
98     AliInfo(Form("Using default geometry in reconstruction"));
99     fGeom =  AliEMCALGeometry::GetInstance(AliEMCALGeometry::GetDefaultGeometryName());
100   }
101
102   //Get calibration parameters  
103   if(!fCalibData)
104     {
105                 AliCDBEntry *entry = (AliCDBEntry*) 
106                 AliCDBManager::Instance()->Get("EMCAL/Calib/Data");
107                 if (entry) fCalibData =  (AliEMCALCalibData*) entry->GetObject();
108     }
109         
110   if(!fCalibData)
111                 AliFatal("Calibration parameters not found in CDB!");
112         
113   //Get calibration parameters  
114   if(!fPedestalData)
115     {
116       AliCDBEntry *entry = (AliCDBEntry*) 
117       AliCDBManager::Instance()->Get("EMCAL/Calib/Pedestals");
118       if (entry) fPedestalData =  (AliCaloCalibPedestal*) entry->GetObject();
119     }
120   
121   if(!fPedestalData)
122     AliFatal("Dead map not found in CDB!");
123         
124   if(!fGeom) AliFatal(Form("Could not get geometry!"));
125
126   AliEMCALTriggerDCSConfigDB* dcsConfigDB = AliEMCALTriggerDCSConfigDB::Instance();
127
128   const AliEMCALTriggerDCSConfig* dcsConfig = dcsConfigDB->GetTriggerDCSConfig();
129
130   if (!dcsConfig) AliFatal("No Trigger DCS Configuration from OCDB!");
131   fgTriggerProcessor = new AliEMCALTriggerElectronics( dcsConfig );
132         
133   fTriggerData = new AliEMCALTriggerData();
134
135   //Init temporary list of digits
136   fgDigitsArr     = new TClonesArray("AliEMCALDigit",1000);
137   fgClustersArr   = new TObjArray(1000);
138   fgTriggerDigits = new TClonesArray("AliEMCALTriggerRawDigit",1000);   
139
140
141 //____________________________________________________________________________
142 AliEMCALReconstructor::~AliEMCALReconstructor()
143 {
144   // dtor
145
146   if(fGeom)              delete fGeom;
147   
148   //No need to delete, recovered from OCDB
149   //if(fCalibData)         delete fCalibData;
150   //if(fPedestalData)      delete fPedestalData;
151   
152   if(fgDigitsArr){
153     fgDigitsArr->Clear("C");
154     delete fgDigitsArr; 
155   }
156   
157   if(fgClustersArr){
158     fgClustersArr->Clear();
159     delete fgClustersArr; 
160   }
161   
162   if(fgTriggerDigits){
163     fgTriggerDigits->Clear();
164     delete fgTriggerDigits; 
165   }
166   
167   if(fgRawUtils)         delete fgRawUtils;
168   if(fgClusterizer)      delete fgClusterizer;
169   if(fgTriggerProcessor) delete fgTriggerProcessor;
170   
171   AliCodeTimer::Instance()->Print();
172
173
174 // //____________________________________________________________________________
175 // void AliEMCALReconstructor::Init()
176 // {
177 //   // Trigger hists - Oct 24, 2007
178 //    fList = AliEMCALHistoUtilities::GetTriggersListOfHists(kTRUE);
179 // }
180
181 //____________________________________________________________________________                                  
182 void AliEMCALReconstructor::InitClusterizer() const
183 {
184   //Init the clusterizer with geometry and calibration pointers, avoid doing it twice.                          
185   Int_t clusterizerType = -1;
186   Int_t eventType = -1;
187   if(GetRecParam()) {
188     clusterizerType = GetRecParam()->GetClusterizerFlag();
189     eventType       = GetRecParam()->GetEventSpecie();
190   }
191   else{
192     AliCDBEntry *entry = (AliCDBEntry*)
193     AliCDBManager::Instance()->Get("EMCAL/Calib/RecoParam");
194     //Get The reco param for the default event specie                                                           
195     if (entry) {
196       AliEMCALRecParam *recParam  = (AliEMCALRecParam*)((TObjArray *) entry->GetObject())->At(0);
197       if(recParam) clusterizerType = recParam->GetClusterizerFlag(); 
198     }
199   }
200   
201   //Check if clusterizer previously set corresponds to what is needed for this event type                       
202   if(fgClusterizer){
203     if(eventType!=AliRecoParam::kCalib){
204       //printf("ReCreate clusterizer? Clusterizer set <%d>, Clusterizer in use <%s>\n",
205         //     clusterizerType, fgClusterizer->Version());
206       
207       if     (clusterizerType == AliEMCALRecParam::kClusterizerv1 && !strcmp(fgClusterizer->Version(),"clu-v1")) return;
208       
209       else if(clusterizerType == AliEMCALRecParam::kClusterizerNxN && !strcmp(fgClusterizer->Version(),"clu-NxN")) return;
210       
211       //Need to create new clusterizer, the one set previously is not the correct one     
212       delete fgClusterizer;
213     }
214     else return;
215   }
216
217   if (clusterizerType  == AliEMCALRecParam::kClusterizerv1)
218   {
219     fgClusterizer = new AliEMCALClusterizerv1(fGeom, fCalibData,fPedestalData);
220   }
221   else
222   {
223     fgClusterizer = new AliEMCALClusterizerNxN(fGeom, fCalibData,fPedestalData);
224   }
225 }
226
227 //____________________________________________________________________________
228 void AliEMCALReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const
229 {
230   // method called by AliReconstruction; 
231   // Only the clusterization is performed,; the rest of the reconstruction is done in FillESD because the track
232   // segment maker needs access to the AliESD object to retrieve the tracks reconstructed by 
233   // the global tracking.
234   // Works on the current event.
235
236   AliCodeTimerAuto("",0)
237
238   ReadDigitsArrayFromTree(digitsTree);
239
240   InitClusterizer();
241   
242   fgClusterizer->InitParameters();
243   fgClusterizer->SetOutput(clustersTree);
244   
245   //Skip clusterization of LED events
246   if (GetRecParam()->GetEventSpecie()!=AliRecoParam::kCalib){
247         
248                 if(fgDigitsArr && fgDigitsArr->GetEntries()) {
249
250                   fgClusterizer->SetInput(digitsTree);
251     
252                   if(Debug())
253                           fgClusterizer->Digits2Clusters("deb all") ;
254                   else
255                           fgClusterizer->Digits2Clusters("");
256     
257                   fgClusterizer->Clear();
258
259           }//digits array exists and has somethind
260   }//not a LED event
261         
262   clustersTree->Fill(); 
263 }
264
265 //____________________________________________________________________________
266 void AliEMCALReconstructor::ConvertDigits(AliRawReader* rawReader, TTree* digitsTree) const
267
268 {
269   // Conversion from raw data to
270   // EMCAL digits.
271   // Works on a single-event basis
272
273   rawReader->Reset() ; 
274
275   fTriggerData->SetMode(1);     
276
277   if(fgDigitsArr) fgDigitsArr->Clear("C");
278
279   TClonesArray *digitsTrg = new TClonesArray("AliEMCALTriggerRawDigit", 32 * 96);
280
281   Int_t bufsize = 32000;
282   digitsTree->Branch("EMCAL", &fgDigitsArr, bufsize);
283   digitsTree->Branch("EMTRG", &digitsTrg, bufsize);
284         
285   //Skip calibration events do the rest
286   Bool_t doFit = kTRUE;
287   if ( !(GetRecParam()->FitLEDEvents()) && GetRecParam()->GetEventSpecie()==AliRecoParam::kCalib) doFit = kFALSE;
288   if (doFit){
289           //must be done here because, in constructor, option is not yet known
290           fgRawUtils->SetOption(GetOption());
291
292           fgRawUtils->SetRawFormatHighLowGainFactor(GetRecParam()->GetHighLowGainFactor());
293           fgRawUtils->SetRawFormatOrder(GetRecParam()->GetOrderParameter());
294           fgRawUtils->SetRawFormatTau(GetRecParam()->GetTau());
295           fgRawUtils->SetNoiseThreshold(GetRecParam()->GetNoiseThreshold());
296           fgRawUtils->SetNPedSamples(GetRecParam()->GetNPedSamples());
297           fgRawUtils->SetRemoveBadChannels(GetRecParam()->GetRemoveBadChannels());
298           fgRawUtils->SetFittingAlgorithm(GetRecParam()->GetFittingAlgorithm());
299           fgRawUtils->SetFALTROUsage(GetRecParam()->UseFALTRO());
300           fgRawUtils->SetTimeMin(GetRecParam()->GetTimeMin());
301           fgRawUtils->SetTimeMax(GetRecParam()->GetTimeMax());
302         
303           fgRawUtils->Raw2Digits(rawReader,fgDigitsArr,fPedestalData,digitsTrg,fTriggerData);
304   }//skip calibration event
305   else{
306         AliDebug(1," Calibration Event, skip!");
307   }
308         
309   digitsTree->Fill();
310   digitsTrg->Delete();
311   delete digitsTrg;
312
313 }
314
315
316 //____________________________________________________________________________
317 void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, 
318                                     AliESDEvent* esd) const
319 {
320   // Called by AliReconstruct after Reconstruct() and global tracking and vertexing 
321   // and V0 
322   // Works on the current event
323   // printf(" ## AliEMCALReconstructor::FillESD() is started ### \n ");
324   //return;
325
326   //########################################
327   // Trigger
328   //########################################
329
330    Int_t v0M[2] = {0, 0};
331         
332    AliESDVZERO* esdV0 = esd->GetVZEROData();
333
334    if (esdV0) 
335    {
336           for (Int_t i = 0; i < 32; i++)
337           {
338                   v0M[0] += (Int_t)esdV0->GetAdcV0C(i);
339                   v0M[1] += (Int_t)esdV0->GetAdcV0A(i);
340           }
341    }
342    else
343    {
344           AliWarning("Cannot retrieve V0 ESD! Run w/ null V0 charges");
345    }
346
347    if (fgTriggerDigits) fgTriggerDigits->Clear();
348
349    TBranch *branchtrg = digitsTree->GetBranch("EMTRG");
350         
351    if (!branchtrg) 
352    { 
353           AliError("Can't get the branch with the EMCAL trigger digits!");
354           return;
355    }
356         
357    branchtrg->SetAddress(&fgTriggerDigits);
358    branchtrg->GetEntry(0);
359   
360    // Note: fgTriggerProcessor reset done at the end of this method
361    fgTriggerProcessor->Digits2Trigger(fgTriggerDigits, v0M, fTriggerData);
362
363    // Fill ESD
364    AliESDCaloTrigger* trgESD = esd->GetCaloTrigger("EMCAL");
365   
366    if (trgESD)
367    {
368           trgESD->Allocate(fgTriggerDigits->GetEntriesFast());
369           
370           for (Int_t i = 0; i < fgTriggerDigits->GetEntriesFast(); i++)
371           {       
372                   AliEMCALTriggerRawDigit* rdig = (AliEMCALTriggerRawDigit*)fgTriggerDigits->At(i);
373                   
374                   Int_t px, py;
375                   if (fGeom->GetPositionInEMCALFromAbsFastORIndex(rdig->GetId(), px, py))
376                   {
377                           Int_t a = -1, t = -1, times[10]; 
378                           
379                           rdig->GetMaximum(a, t);
380                           rdig->GetL0Times(times);
381                                   
382                           trgESD->Add(px, py, a, t, times, rdig->GetNL0Times(), rdig->GetL1TimeSum(), rdig->GetTriggerBits());
383                   }
384           }
385           
386           trgESD->SetL1Threshold(0, fTriggerData->GetL1GammaThreshold());
387           
388           trgESD->SetL1Threshold(1, fTriggerData->GetL1JetThreshold()  );
389    }
390
391    // Resetting
392    fTriggerData->Reset();
393   
394   //########################################
395   //##############Fill CaloCells###############
396   //########################################
397   ReadDigitsArrayFromTree(digitsTree);
398
399 //   TClonesArray *digits = new TClonesArray("AliEMCALDigit",1000);
400 //   TBranch *branchdig = digitsTree->GetBranch("EMCAL");
401 //   if (!branchdig) { 
402 //     AliError("can't get the branch with the EMCAL digits !");
403 //     return;
404 //   }
405 //   branchdig->SetAddress(&digits);
406 //   digitsTree->GetEvent(0);
407   Int_t nDigits = fgDigitsArr->GetEntries(), idignew = 0 ;
408   AliDebug(1,Form("%d digits",nDigits));
409
410   AliESDCaloCells &emcCells = *(esd->GetEMCALCells());
411   emcCells.CreateContainer(nDigits);
412   emcCells.SetType(AliVCaloCells::kEMCALCell);
413   Float_t energy = 0;
414   for (Int_t idig = 0 ; idig < nDigits ; idig++) {
415     const AliEMCALDigit * dig = (const AliEMCALDigit*)fgDigitsArr->At(idig);
416     if(dig->GetAmplitude() > 0 ){
417           energy = fgClusterizer->Calibrate(dig->GetAmplitude(),dig->GetTime(),dig->GetId()); //TimeR or Time?
418           if(energy > 0){ //Digits tagged as bad (dead, hot, not alive) are set to 0 in calibrate, remove them  
419                   emcCells.SetCell(idignew,dig->GetId(),energy, dig->GetTime());   
420                   idignew++;
421           }
422     }
423   }
424   emcCells.SetNumberOfCells(idignew);
425   emcCells.Sort();
426
427   //------------------------------------------------------------
428   //-----------------CLUSTERS-----------------------------
429   //------------------------------------------------------------
430   clustersTree->SetBranchStatus("*",0); //disable all branches
431   clustersTree->SetBranchStatus("EMCALECARP",1); //Enable only the branch we need
432   if(fgClustersArr) fgClustersArr->Clear();
433   TBranch *branch = clustersTree->GetBranch("EMCALECARP");
434   branch->SetAddress(&fgClustersArr);
435   branch->GetEntry(0);
436   //clustersTree->GetEvent(0);
437
438   Int_t nClusters = fgClustersArr->GetEntries(),  nClustersNew=0;
439   AliDebug(1,Form("%d clusters",nClusters));
440
441   //######################################################
442   //#######################TRACK MATCHING###############
443   //######################################################
444   //Fill list of integers, each one is index of track to which the cluster belongs.
445
446   // step 1 - initialize array of matched track indexes
447   Int_t *matchedTrack = new Int_t[nClusters];
448   for (Int_t iclus = 0; iclus < nClusters; iclus++)
449     matchedTrack[iclus] = -1;  // neg. index --> no matched track
450   
451   // step 2, change the flag for all matched clusters found in tracks
452   Int_t iemcalMatch = -1;
453   Int_t endtpc = esd->GetNumberOfTracks();
454   for (Int_t itrack = 0; itrack < endtpc; itrack++) {
455     AliESDtrack * track = esd->GetTrack(itrack) ; // retrieve track
456     iemcalMatch = track->GetEMCALcluster();
457     if(iemcalMatch >= 0) matchedTrack[iemcalMatch] = itrack;
458   } 
459   
460   //########################################
461   //##############Fill CaloClusters#############
462   //########################################
463   for (Int_t iClust = 0 ; iClust < nClusters ; iClust++) {
464     const AliEMCALRecPoint * clust = (const AliEMCALRecPoint*)fgClustersArr->At(iClust);
465     //if(clust->GetClusterType()== AliVCluster::kEMCALClusterv1) nRP++; else nPC++;
466     if (Debug()) clust->Print();
467     // Get information from EMCAL reconstruction points
468     Float_t xyz[3];
469     TVector3 gpos;
470     clust->GetGlobalPosition(gpos);
471     for (Int_t ixyz=0; ixyz<3; ixyz++)
472       xyz[ixyz] = gpos[ixyz];
473     Float_t elipAxis[2];
474     clust->GetElipsAxis(elipAxis);
475        //Create digits lists
476     Int_t cellMult = clust->GetMultiplicity();
477     //TArrayS digiList(digitMult);
478     Float_t *amplFloat = clust->GetEnergiesList();
479     Int_t   *digitInts = clust->GetAbsId();
480     TArrayS absIdList(cellMult);
481     TArrayD fracList(cellMult);
482
483     Int_t newCellMult = 0;
484     for (Int_t iCell=0; iCell<cellMult; iCell++) {
485       if (amplFloat[iCell] > 0) {
486       absIdList[newCellMult] = (UShort_t)(digitInts[iCell]);
487       //Calculate Fraction
488       if(emcCells.GetCellAmplitude(digitInts[iCell])>0 && GetRecParam()->GetUnfold())
489         fracList[newCellMult] = amplFloat[iCell]/(emcCells.GetCellAmplitude(digitInts[iCell]));//get cell calibration value 
490       else
491         fracList[newCellMult] = 0; 
492       newCellMult++;
493       }
494     }
495
496     absIdList.Set(newCellMult);
497     fracList.Set(newCellMult);
498     
499     if(newCellMult > 0) { // accept cluster if it has some digit
500       nClustersNew++;
501       //Primaries
502       Int_t  parentMult  = 0;
503       Int_t *parentList =  clust->GetParents(parentMult);
504       // fills the ESDCaloCluster
505       AliESDCaloCluster * ec = new AliESDCaloCluster() ;
506       ec->SetType(AliVCluster::kEMCALClusterv1);
507       ec->SetPosition(xyz);
508       ec->SetE(clust->GetEnergy());
509                 
510       //Distance to the nearest bad crystal
511       ec->SetDistanceToBadChannel(clust->GetDistanceToBadTower()); 
512
513       ec->SetNCells(newCellMult);
514       //Change type of list from short to ushort
515       UShort_t *newAbsIdList  = new UShort_t[newCellMult];
516       Double_t *newFracList   = new Double_t[newCellMult];
517       for(Int_t i = 0; i < newCellMult ; i++) {
518         newAbsIdList[i]=absIdList[i];
519         newFracList[i] =fracList[i];
520       }
521       ec->SetCellsAbsId(newAbsIdList);
522       ec->SetCellsAmplitudeFraction(newFracList);
523       ec->SetDispersion(clust->GetDispersion());
524       ec->SetChi2(-1); //not yet implemented
525       ec->SetM02(elipAxis[0]*elipAxis[0]) ;
526       ec->SetM20(elipAxis[1]*elipAxis[1]) ;
527       ec->SetTOF(clust->GetTime()) ; //time-of-fligh
528       ec->SetNExMax(clust->GetNExMax());          //number of local maxima
529       TArrayI arrayTrackMatched(1);// Only one track, temporal solution.
530       arrayTrackMatched[0]= matchedTrack[iClust];
531       ec->AddTracksMatched(arrayTrackMatched);
532
533       TArrayI arrayParents(parentMult,parentList);
534       ec->AddLabels(arrayParents);
535
536       // add the cluster to the esd object
537       esd->AddCaloCluster(ec);
538       delete ec;
539       delete [] newAbsIdList ;
540       delete [] newFracList ;
541    }
542  } // cycle on clusters
543
544  delete [] matchedTrack;
545
546  //Fill ESDCaloCluster with PID weights
547  AliEMCALPID *pid = new AliEMCALPID;
548  //pid->SetPrintInfo(kTRUE);
549  pid->SetReconstructor(kTRUE);
550  pid->RunPID(esd);
551  delete pid;
552     
553  //Store EMCAL misalignment matrixes
554  FillMisalMatrixes(esd) ;
555
556 }
557
558 //==================================================================================
559 void AliEMCALReconstructor::FillMisalMatrixes(AliESDEvent* esd)const{
560         //Store EMCAL matrixes in ESD Header
561         
562         //Check, if matrixes was already stored
563         for(Int_t sm = 0 ; sm < fGeom->GetNumberOfSuperModules(); sm++){
564                 if(esd->GetEMCALMatrix(sm)!=0)
565                         return ;
566         }
567         
568         //Create and store matrixes
569         if(!gGeoManager){
570                 AliError("Can not store misal. matrixes: no gGeoManager! \n") ;
571                 return ;
572         }
573         //Note, that owner of copied marixes will be header
574   const Int_t bufsize = 255;
575         char path[bufsize] ;
576         TGeoHMatrix * m = 0x0;
577         for(Int_t sm = 0; sm < fGeom->GetNumberOfSuperModules(); sm++){
578                 snprintf(path,bufsize,"/ALIC_1/XEN1_1/SMOD_%d",sm+1) ; //In Geometry modules numbered 1,2,.,5
579                 if(sm >= 10) snprintf(path,bufsize,"/ALIC_1/XEN1_1/SM10_%d",sm-10+1) ;
580                 
581                 if (gGeoManager->CheckPath(path)){
582                         gGeoManager->cd(path);
583                         m = gGeoManager->GetCurrentMatrix() ;
584 //                      printf("================================================= \n");
585 //                      printf("AliEMCALReconstructor::FixMisalMatrixes(), sm %d, \n",sm);
586 //                      m->Print("");
587                         esd->SetEMCALMatrix(new TGeoHMatrix(*m),sm) ;
588 //                      printf("================================================= \n");
589                 }
590                 else{
591                         esd->SetEMCALMatrix(NULL,sm) ;
592                 }
593         }
594 }
595
596
597
598 //__________________________________________________________________________
599 void AliEMCALReconstructor::ReadDigitsArrayFromTree(TTree *digitsTree) const
600 {
601   // See AliEMCALClusterizer::SetInput(TTree *digitsTree);
602   if(fgDigitsArr) {
603     // Clear previous digits 
604     fgDigitsArr->Clear("C");
605     //delete fgDigitsArr;
606   }
607   // Read the digits from the input tree
608   TBranch *branch = digitsTree->GetBranch("EMCAL");
609   if (!branch) { 
610     AliError("can't get the branch with the EMCAL digits !");
611     return;
612   }
613   fgDigitsArr = new TClonesArray("AliEMCALDigit",100);
614   branch->SetAddress(&fgDigitsArr);
615   branch->GetEntry(0);
616 }
617
618