1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 //_________________________________________________________________________
20 //*-- Yves Schutz (SUBATECH)
21 // Reconstruction class. Redesigned from the old AliReconstructionner class and
22 // derived from STEER/AliReconstructor.
24 // --- ROOT system ---
26 // --- Standard library ---
28 // --- AliRoot header files ---
29 #include "AliESDEvent.h"
30 #include "AliESDCaloCluster.h"
31 #include "AliPHOSReconstructor.h"
32 #include "AliPHOSClusterizerv1.h"
33 #include "AliPHOSTrackSegmentMakerv1.h"
34 #include "AliPHOSPIDv1.h"
35 #include "AliPHOSGetter.h"
36 #include "AliPHOSTracker.h"
37 #include "AliRawReader.h"
38 #include "AliPHOSTrigger.h"
39 #include "AliPHOSGeometry.h"
40 #include "AliPHOSRecoParamEmc.h"
41 #include "AliPHOSRecoParamCpv.h"
43 ClassImp(AliPHOSReconstructor)
45 Bool_t AliPHOSReconstructor::fgDebug = kFALSE ;
46 AliPHOSRecoParam* AliPHOSReconstructor::fgkRecoParamEmc =0; // EMC rec. parameters
47 AliPHOSRecoParam* AliPHOSReconstructor::fgkRecoParamCpv =0; // CPV rec. parameters
49 //____________________________________________________________________________
50 AliPHOSReconstructor::AliPHOSReconstructor()
54 if (!fgkRecoParamEmc) {
55 AliWarning("The Reconstruction parameters for EMC nonitialized - Used default one");
56 fgkRecoParamEmc = AliPHOSRecoParamEmc::GetEmcDefaultParameters();
59 if (!fgkRecoParamCpv) {
60 AliWarning("The Reconstruction parameters for CPV nonitialized - Used default one");
61 fgkRecoParamCpv = AliPHOSRecoParamCpv::GetCpvDefaultParameters();
66 //____________________________________________________________________________
67 AliPHOSReconstructor::~AliPHOSReconstructor()
73 //____________________________________________________________________________
74 void AliPHOSReconstructor::Reconstruct(AliRunLoader* runLoader) const
76 // method called by AliReconstruction;
77 // Only the clusterization is performed,; the rest of the reconstruction is done in FillESD because the track
78 // segment maker needs access to the AliESDEvent object to retrieve the tracks reconstructed by
79 // the global tracking.
81 TString headerFile(runLoader->GetFileName()) ;
82 TString branchName(runLoader->GetEventFolder()->GetName()) ;
84 AliPHOSClusterizerv1 clu(headerFile, branchName);
85 clu.SetEventRange(0, -1) ; // do all the events
87 clu.ExecuteTask("deb all") ;
93 //____________________________________________________________________________
94 void AliPHOSReconstructor::Reconstruct(AliRunLoader* runLoader, AliRawReader* rawreader) const
96 // method called by AliReconstruction;
97 // Only the clusterization is performed,; the rest of the reconstruction is done in FillESD because the track
98 // segment maker needs access to the AliESDEvent object to retrieve the tracks reconstructed by
99 // the global tracking.
100 // Here we reconstruct from Raw Data
103 TString headerFile(runLoader->GetFileName()) ;
104 TString branchName(runLoader->GetEventFolder()->GetName()) ;
106 AliPHOSClusterizerv1 clu(headerFile, branchName);
107 clu.SetEventRange(0, -1) ; // do all the events
108 clu.SetRawReader(rawreader);
110 TString option = GetOption();
111 if (option.Contains("OldRCUFormat"))
112 clu.SetOldRCUFormat(kTRUE);
115 clu.ExecuteTask("deb all") ;
117 clu.ExecuteTask("") ;
121 //____________________________________________________________________________
122 void AliPHOSReconstructor::FillESD(AliRunLoader* runLoader, AliESDEvent* esd) const
124 // This function creates AliESDtracks from AliPHOSRecParticles
126 // writes them to the ESD
128 Int_t eventNumber = runLoader->GetEventNumber() ;
131 AliPHOSGetter *gime = AliPHOSGetter::Instance();
132 gime->Event(eventNumber, "DRTP") ;
133 TClonesArray *recParticles = gime->RecParticles();
134 Int_t nOfRecParticles = recParticles->GetEntries();
137 esd->SetNumberOfPHOSClusters(nOfRecParticles) ;
138 esd->SetFirstPHOSCluster(esd->GetNumberOfCaloClusters()) ;
140 Int_t maxClu = esd->GetNumberOfPHOSClusters() ;
142 AliDebug(2,Form("%d digits and %d rec. particles in event %d, option %s",gime->Digits()->GetEntries(),nOfRecParticles,eventNumber,GetOption()));
145 //#########Calculate trigger and set trigger info###########
148 // tr.SetPatchSize(1);//create 4x4 patches
151 Float_t maxAmp2x2 = tr.Get2x2MaxAmplitude();
152 Float_t maxAmpnxn = tr.GetnxnMaxAmplitude();
153 Float_t ampOutOfPatch2x2 = tr.Get2x2AmpOutOfPatch() ;
154 Float_t ampOutOfPatchnxn = tr.GetnxnAmpOutOfPatch() ;
156 AliPHOSGeometry * geom = gime->PHOSGeometry();
158 Int_t iSM2x2 = tr.Get2x2SuperModule();
159 Int_t iSMnxn = tr.GetnxnSuperModule();
160 Int_t iCrystalPhi2x2 = tr.Get2x2CrystalPhi();
161 Int_t iCrystalPhinxn = tr.GetnxnCrystalPhi();
162 Int_t iCrystalEta2x2 = tr.Get2x2CrystalEta();
163 Int_t iCrystalEtanxn = tr.GetnxnCrystalEta();
165 AliDebug(2, Form("Trigger 2x2 max amp %f, out amp %f, SM %d, iphi %d ieta %d", maxAmp2x2, ampOutOfPatch2x2, iSM2x2,iCrystalPhi2x2, iCrystalEta2x2));
166 AliDebug(2, Form("Trigger 4x4 max amp %f , out amp %f, SM %d, iphi %d, ieta %d", maxAmpnxn, ampOutOfPatchnxn, iSMnxn,iCrystalPhinxn, iCrystalEtanxn));
168 Int_t iRelId2x2 []= {iSM2x2+1,0,iCrystalPhi2x2,iCrystalEta2x2};// PHOS modules in order to calculate AbsId need to be 1-5 not 0-4 as returns trigger.
170 Int_t iRelIdnxn []= {iSMnxn+1,0,iCrystalPhinxn,iCrystalEtanxn};// PHOS modules in order to calculate AbsId need to be 1-5 not 0-4 as returns trigger.
172 TVector3 pos2x2(-1,-1,-1);
173 TVector3 posnxn(-1,-1,-1);
174 geom->RelToAbsNumbering(iRelId2x2, iAbsId2x2);
175 geom->RelToAbsNumbering(iRelIdnxn, iAbsIdnxn);
176 geom->RelPosInAlice(iAbsId2x2, pos2x2);
177 geom->RelPosInAlice(iAbsIdnxn, posnxn);
179 TArrayF triggerPosition(6);
180 triggerPosition[0] = pos2x2(0) ;
181 triggerPosition[1] = pos2x2(1) ;
182 triggerPosition[2] = pos2x2(2) ;
183 triggerPosition[3] = posnxn(0) ;
184 triggerPosition[4] = posnxn(1) ;
185 triggerPosition[5] = posnxn(2) ;
187 TArrayF triggerAmplitudes(4);
188 triggerAmplitudes[0] = maxAmp2x2 ;
189 triggerAmplitudes[1] = ampOutOfPatch2x2 ;
190 triggerAmplitudes[2] = maxAmpnxn ;
191 triggerAmplitudes[3] = ampOutOfPatchnxn ;
193 //esd->SetPHOSTriggerCells(triggerPosition);
194 esd->AddPHOSTriggerPosition(triggerPosition);
195 esd->AddPHOSTriggerAmplitudes(triggerAmplitudes);
197 //######################################
200 const Float_t kBigShort = std::numeric_limits<short int>::max() - 1;
201 const Float_t nsec100 = 1e9*100.; // units of 0.01 ns
202 const Float_t gev500 = 500.; // units of GeV/500
204 for (Int_t recpart = 0 ; recpart < nOfRecParticles ; recpart++) {
205 AliPHOSRecParticle * rp = dynamic_cast<AliPHOSRecParticle*>(recParticles->At(recpart));
208 // Get track segment and EMC rec.point associated with this rec.particle
209 AliPHOSTrackSegment *ts = gime->TrackSegment(rp->GetPHOSTSIndex());
210 AliPHOSEmcRecPoint *emcRP = gime->EmcRecPoint(ts->GetEmcIndex());
211 AliESDCaloCluster *ec = new AliESDCaloCluster() ;
214 for (Int_t ixyz=0; ixyz<3; ixyz++)
215 xyz[ixyz] = rp->GetPos()[ixyz];
217 AliDebug(2,Form("Global position xyz=(%f,%f,%f)",xyz[0],xyz[1],xyz[2]));
219 //Create digits lists
220 Int_t digitMult = emcRP->GetDigitsMultiplicity();
221 Int_t *digitsList = emcRP->GetDigitsList();
222 Short_t *amplList = new Short_t[digitMult];
223 Short_t *timeList = new Short_t[digitMult];
224 Short_t *digiList = new Short_t[digitMult];
226 // Convert Float_t* and Int_t* to Short_t* to save memory
227 for (Int_t iDigit=0; iDigit<digitMult; iDigit++) {
228 AliPHOSDigit *digit = gime->Digit(digitsList[iDigit]);
230 (Short_t)(TMath::Min(digit->GetEnergy()*gev500,kBigShort)); // Energy in units of GeV/500
232 (Short_t)(TMath::Min(digit->GetTime()*nsec100,kBigShort)); // time in units of 0.01 ns
233 digiList[iDigit] = (Short_t)(digit->GetId());
238 Int_t *primInts = emcRP->GetPrimaries(primMult);
239 Short_t *primList = new Short_t[primMult];
240 for (Int_t ipr=0; ipr<primMult; ipr++)
241 primList[ipr] = (Short_t)(primInts[ipr]);
243 // fills the ESDCaloCluster
246 ec->SetPosition(xyz); //rec.point position in MARS
247 ec->SetE(rp->Energy()); //total particle energy
248 ec->SetClusterDisp(emcRP->GetDispersion()); //cluster dispersion
249 ec->SetPid (rp->GetPID()) ; //array of particle identification
250 ec->SetM02(emcRP->GetM2x()) ; //second moment M2x
251 ec->SetM20(emcRP->GetM2z()) ; //second moment M2z
252 ec->SetNExMax(emcRP->GetNExMax()); //number of local maxima
253 ec->SetEmcCpvDistance(-1); //not yet implemented
254 ec->SetClusterChi2(-1); //not yet implemented
255 ec->SetM11(-1) ; //not yet implemented
258 TArrayS arrayAmpList(digitMult,amplList);
259 TArrayS arrayTimeList(digitMult,timeList);
260 TArrayS arrayIndexList(digitMult,digiList);
261 ec->AddDigitAmplitude(arrayAmpList);
262 ec->AddDigitTime(arrayTimeList);
263 ec->AddDigitIndex(arrayIndexList);
265 //Distance to the nearest bad crystal
266 ec->SetDistanceToBadChannel(emcRP->GetDistanceToBadCrystal());
268 //Array of MC indeces
269 TArrayS arrayPrim(primMult,primList);
270 ec->AddLabels(arrayPrim);
272 //Array of tracks uncomment when available in future
273 //TArrayS arrayTrackMatched(1);// Only one track, temporal solution.
274 //arrayTrackMatched[0]= (Short_t)(matchedTrack[iClust]);
275 //ec->AddTracksMatched(arrayTrackMatched);
277 // add the track to the esd object
278 esd->AddCaloCluster(ec);
283 //____________________________________________________________________________
284 void AliPHOSReconstructor::FillESD(AliRunLoader* runLoader,
285 AliRawReader* rawReader, AliESDEvent* esd) const
287 //This function creates AliESDtracks from AliPHOSRecParticles
288 //and writes them to the ESD in the case of raw data reconstruction.
290 Int_t eventNumber = runLoader->GetEventNumber() ;
292 AliPHOSGetter *gime = AliPHOSGetter::Instance();
293 gime->Event(eventNumber, "DRTP") ;
295 TClonesArray *recParticles = gime->RecParticles();
296 Int_t nOfRecParticles = recParticles->GetEntries();
298 esd->SetNumberOfPHOSClusters(nOfRecParticles) ;
299 esd->SetFirstPHOSCluster(esd->GetNumberOfCaloClusters()) ;
301 AliDebug(2,Form("%d digits and %d rec. particles in event %d, option %s",gime->Digits()->GetEntries(),nOfRecParticles,eventNumber,GetOption()));
303 for (Int_t recpart = 0 ; recpart < nOfRecParticles ; recpart++) {
304 AliPHOSRecParticle * rp = dynamic_cast<AliPHOSRecParticle*>(recParticles->At(recpart));
308 for (Int_t ixyz=0; ixyz<3; ixyz++)
309 xyz[ixyz] = rp->GetPos()[ixyz];
311 AliDebug(2,Form("Global position xyz=(%f,%f,%f)",xyz[0],xyz[1],xyz[2]));
313 AliPHOSTrackSegment *ts = gime->TrackSegment(rp->GetPHOSTSIndex());
314 AliPHOSEmcRecPoint *emcRP = gime->EmcRecPoint(ts->GetEmcIndex());
315 AliESDCaloCluster *ec = new AliESDCaloCluster() ;
317 Int_t digitMult = emcRP->GetDigitsMultiplicity();
318 Int_t *digitsList = emcRP->GetDigitsList();
319 Short_t *amplList = new Short_t[digitMult];
320 Short_t *digiList = new Short_t[digitMult];
322 // Convert Float_t* and Int_t* to UShort_t* to save memory
323 for (Int_t iDigit=0; iDigit<digitMult; iDigit++) {
324 AliPHOSDigit *digit = gime->Digit(digitsList[iDigit]);
326 AliFatal(Form("Digit not found at the expected position %d!",iDigit));
329 amplList[iDigit] = (Short_t)digit->GetEnergy();
330 digiList[iDigit] = (Short_t)(digit->GetId());
331 //timeList[iDigit] = (Short_t)(digit->GetTime());
336 ec->SetPosition(xyz); //rec.point position in MARS
337 ec->SetE(rp->Energy()); //total particle energy
338 ec->SetClusterDisp(emcRP->GetDispersion()); //cluster dispersion
339 ec->SetPid (rp->GetPID()) ; //array of particle identification
340 ec->SetM02(emcRP->GetM2x()) ; //second moment M2x
341 ec->SetM20(emcRP->GetM2z()) ; //second moment M2z
342 ec->SetNExMax(emcRP->GetNExMax()); //number of local maxima
344 ec->SetEmcCpvDistance(-1); //not yet implemented
345 ec->SetClusterChi2(-1); //not yet implemented
346 ec->SetM11(-1) ; //not yet implemented
347 // TArrayS arrayAmpList(digitMult,amplList);
348 // TArrayS arrayTimeList(digitMult,timeList);
349 // TArrayS arrayIndexList(digitMult,digiList);
350 // ec->AddDigitAmplitude(arrayAmpList);
351 // ec->AddDigitTime(arrayTimeList);
352 // ec->AddDigitIndex(arrayIndexList);
353 //Distance to the nearest bad crystal
354 ec->SetDistanceToBadChannel(emcRP->GetDistanceToBadCrystal());
356 // add the track to the esd object
357 esd->AddCaloCluster(ec);
366 AliTracker* AliPHOSReconstructor::CreateTracker(AliRunLoader* runLoader) const
368 // creates the PHOS tracker
369 if (!runLoader) return NULL;
370 return new AliPHOSTracker(runLoader);