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() ;
130 AliPHOSGetter *gime = AliPHOSGetter::Instance();
131 gime->Event(eventNumber, "DRTP") ;
132 TClonesArray *recParticles = gime->RecParticles();
133 Int_t nOfRecParticles = recParticles->GetEntries();
135 esd->SetNumberOfPHOSClusters(nOfRecParticles) ;
136 esd->SetFirstPHOSCluster(esd->GetNumberOfCaloClusters()) ;
138 AliDebug(2,Form("%d digits and %d rec. particles in event %d, option %s",gime->Digits()->GetEntries(),nOfRecParticles,eventNumber,GetOption()));
141 //#########Calculate trigger and set trigger info###########
144 // tr.SetPatchSize(1);//create 4x4 patches
147 Float_t maxAmp2x2 = tr.Get2x2MaxAmplitude();
148 Float_t maxAmpnxn = tr.GetnxnMaxAmplitude();
149 Float_t ampOutOfPatch2x2 = tr.Get2x2AmpOutOfPatch() ;
150 Float_t ampOutOfPatchnxn = tr.GetnxnAmpOutOfPatch() ;
152 AliPHOSGeometry * geom = gime->PHOSGeometry();
154 Int_t iSM2x2 = tr.Get2x2SuperModule();
155 Int_t iSMnxn = tr.GetnxnSuperModule();
156 Int_t iCrystalPhi2x2 = tr.Get2x2CrystalPhi();
157 Int_t iCrystalPhinxn = tr.GetnxnCrystalPhi();
158 Int_t iCrystalEta2x2 = tr.Get2x2CrystalEta();
159 Int_t iCrystalEtanxn = tr.GetnxnCrystalEta();
161 AliDebug(2, Form("Trigger 2x2 max amp %f, out amp %f, SM %d, iphi %d ieta %d", maxAmp2x2, ampOutOfPatch2x2, iSM2x2,iCrystalPhi2x2, iCrystalEta2x2));
162 AliDebug(2, Form("Trigger 4x4 max amp %f , out amp %f, SM %d, iphi %d, ieta %d", maxAmpnxn, ampOutOfPatchnxn, iSMnxn,iCrystalPhinxn, iCrystalEtanxn));
164 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.
166 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.
168 TVector3 pos2x2(-1,-1,-1);
169 TVector3 posnxn(-1,-1,-1);
170 geom->RelToAbsNumbering(iRelId2x2, iAbsId2x2);
171 geom->RelToAbsNumbering(iRelIdnxn, iAbsIdnxn);
172 geom->RelPosInAlice(iAbsId2x2, pos2x2);
173 geom->RelPosInAlice(iAbsIdnxn, posnxn);
175 TArrayF triggerPosition(6);
176 triggerPosition[0] = pos2x2(0) ;
177 triggerPosition[1] = pos2x2(1) ;
178 triggerPosition[2] = pos2x2(2) ;
179 triggerPosition[3] = posnxn(0) ;
180 triggerPosition[4] = posnxn(1) ;
181 triggerPosition[5] = posnxn(2) ;
183 TArrayF triggerAmplitudes(4);
184 triggerAmplitudes[0] = maxAmp2x2 ;
185 triggerAmplitudes[1] = ampOutOfPatch2x2 ;
186 triggerAmplitudes[2] = maxAmpnxn ;
187 triggerAmplitudes[3] = ampOutOfPatchnxn ;
189 //esd->SetPHOSTriggerCells(triggerPosition);
190 esd->AddPHOSTriggerPosition(triggerPosition);
191 esd->AddPHOSTriggerAmplitudes(triggerAmplitudes);
193 //######################################
196 const Float_t kBigShort = std::numeric_limits<short int>::max() - 1;
197 const Float_t nsec100 = 1e9*100.; // units of 0.01 ns
198 const Float_t gev500 = 500.; // units of GeV/500
200 for (Int_t recpart = 0 ; recpart < nOfRecParticles ; recpart++) {
201 AliPHOSRecParticle * rp = dynamic_cast<AliPHOSRecParticle*>(recParticles->At(recpart));
204 // Get track segment and EMC rec.point associated with this rec.particle
205 AliPHOSTrackSegment *ts = gime->TrackSegment(rp->GetPHOSTSIndex());
206 AliPHOSEmcRecPoint *emcRP = gime->EmcRecPoint(ts->GetEmcIndex());
207 AliESDCaloCluster *ec = new AliESDCaloCluster() ;
210 for (Int_t ixyz=0; ixyz<3; ixyz++)
211 xyz[ixyz] = rp->GetPos()[ixyz];
213 AliDebug(2,Form("Global position xyz=(%f,%f,%f)",xyz[0],xyz[1],xyz[2]));
215 //Create digits lists
216 Int_t digitMult = emcRP->GetDigitsMultiplicity();
217 Int_t *digitsList = emcRP->GetDigitsList();
218 Short_t *amplList = new Short_t[digitMult];
219 Short_t *timeList = new Short_t[digitMult];
220 Short_t *digiList = new Short_t[digitMult];
222 // Convert Float_t* and Int_t* to Short_t* to save memory
223 for (Int_t iDigit=0; iDigit<digitMult; iDigit++) {
224 AliPHOSDigit *digit = gime->Digit(digitsList[iDigit]);
226 (Short_t)(TMath::Min(digit->GetEnergy()*gev500,kBigShort)); // Energy in units of GeV/500
228 (Short_t)(TMath::Min(digit->GetTime()*nsec100,kBigShort)); // time in units of 0.01 ns
229 digiList[iDigit] = (Short_t)(digit->GetId());
234 Int_t *primInts = emcRP->GetPrimaries(primMult);
235 Short_t *primList = new Short_t[primMult];
236 for (Int_t ipr=0; ipr<primMult; ipr++)
237 primList[ipr] = (Short_t)(primInts[ipr]);
239 // fills the ESDCaloCluster
242 ec->SetPosition(xyz); //rec.point position in MARS
243 ec->SetE(rp->Energy()); //total particle energy
244 ec->SetClusterDisp(emcRP->GetDispersion()); //cluster dispersion
245 ec->SetPid (rp->GetPID()) ; //array of particle identification
246 ec->SetM02(emcRP->GetM2x()) ; //second moment M2x
247 ec->SetM20(emcRP->GetM2z()) ; //second moment M2z
248 ec->SetNExMax(emcRP->GetNExMax()); //number of local maxima
249 ec->SetEmcCpvDistance(-1); //not yet implemented
250 ec->SetClusterChi2(-1); //not yet implemented
251 ec->SetM11(-1) ; //not yet implemented
254 TArrayS arrayAmpList(digitMult,amplList);
255 TArrayS arrayTimeList(digitMult,timeList);
256 TArrayS arrayIndexList(digitMult,digiList);
257 ec->AddDigitAmplitude(arrayAmpList);
258 ec->AddDigitTime(arrayTimeList);
259 ec->AddDigitIndex(arrayIndexList);
261 //Distance to the nearest bad crystal
262 ec->SetDistanceToBadChannel(emcRP->GetDistanceToBadCrystal());
264 //Array of MC indeces
265 TArrayS arrayPrim(primMult,primList);
266 ec->AddLabels(arrayPrim);
268 //Array of tracks uncomment when available in future
269 //TArrayS arrayTrackMatched(1);// Only one track, temporal solution.
270 //arrayTrackMatched[0]= (Short_t)(matchedTrack[iClust]);
271 //ec->AddTracksMatched(arrayTrackMatched);
273 // add the track to the esd object
274 esd->AddCaloCluster(ec);
279 //____________________________________________________________________________
280 void AliPHOSReconstructor::FillESD(AliRunLoader* runLoader,
281 AliRawReader* rawReader, AliESDEvent* esd) const
283 //This function creates AliESDtracks from AliPHOSRecParticles
284 //and writes them to the ESD in the case of raw data reconstruction.
286 Int_t eventNumber = runLoader->GetEventNumber() ;
288 AliPHOSGetter *gime = AliPHOSGetter::Instance();
289 gime->Event(eventNumber, "DRTP") ;
291 TClonesArray *recParticles = gime->RecParticles();
292 Int_t nOfRecParticles = recParticles->GetEntries();
294 esd->SetNumberOfPHOSClusters(nOfRecParticles) ;
295 esd->SetFirstPHOSCluster(esd->GetNumberOfCaloClusters()) ;
297 AliDebug(2,Form("%d digits and %d rec. particles in event %d, option %s",gime->Digits()->GetEntries(),nOfRecParticles,eventNumber,GetOption()));
299 for (Int_t recpart = 0 ; recpart < nOfRecParticles ; recpart++) {
300 AliPHOSRecParticle * rp = dynamic_cast<AliPHOSRecParticle*>(recParticles->At(recpart));
304 for (Int_t ixyz=0; ixyz<3; ixyz++)
305 xyz[ixyz] = rp->GetPos()[ixyz];
307 AliDebug(2,Form("Global position xyz=(%f,%f,%f)",xyz[0],xyz[1],xyz[2]));
309 AliPHOSTrackSegment *ts = gime->TrackSegment(rp->GetPHOSTSIndex());
310 AliPHOSEmcRecPoint *emcRP = gime->EmcRecPoint(ts->GetEmcIndex());
311 AliESDCaloCluster *ec = new AliESDCaloCluster() ;
313 Int_t digitMult = emcRP->GetDigitsMultiplicity();
314 Int_t *digitsList = emcRP->GetDigitsList();
315 Short_t *amplList = new Short_t[digitMult];
316 Short_t *digiList = new Short_t[digitMult];
318 // Convert Float_t* and Int_t* to UShort_t* to save memory
319 for (Int_t iDigit=0; iDigit<digitMult; iDigit++) {
320 AliPHOSDigit *digit = gime->Digit(digitsList[iDigit]);
322 AliFatal(Form("Digit not found at the expected position %d!",iDigit));
325 amplList[iDigit] = (Short_t)digit->GetEnergy();
326 digiList[iDigit] = (Short_t)(digit->GetId());
327 //timeList[iDigit] = (Short_t)(digit->GetTime());
332 ec->SetPosition(xyz); //rec.point position in MARS
333 ec->SetE(rp->Energy()); //total particle energy
334 ec->SetClusterDisp(emcRP->GetDispersion()); //cluster dispersion
335 ec->SetPid (rp->GetPID()) ; //array of particle identification
336 ec->SetM02(emcRP->GetM2x()) ; //second moment M2x
337 ec->SetM20(emcRP->GetM2z()) ; //second moment M2z
338 ec->SetNExMax(emcRP->GetNExMax()); //number of local maxima
340 ec->SetEmcCpvDistance(-1); //not yet implemented
341 ec->SetClusterChi2(-1); //not yet implemented
342 ec->SetM11(-1) ; //not yet implemented
343 // TArrayS arrayAmpList(digitMult,amplList);
344 // TArrayS arrayTimeList(digitMult,timeList);
345 // TArrayS arrayIndexList(digitMult,digiList);
346 // ec->AddDigitAmplitude(arrayAmpList);
347 // ec->AddDigitTime(arrayTimeList);
348 // ec->AddDigitIndex(arrayIndexList);
349 //Distance to the nearest bad crystal
350 ec->SetDistanceToBadChannel(emcRP->GetDistanceToBadCrystal());
352 // add the track to the esd object
353 esd->AddCaloCluster(ec);
362 AliTracker* AliPHOSReconstructor::CreateTracker(AliRunLoader* runLoader) const
364 // creates the PHOS tracker
365 if (!runLoader) return NULL;
366 return new AliPHOSTracker(runLoader);