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 **************************************************************************/
16 /* $Id: AliITSUReconstructor.cxx 58442 2012-09-04 17:17:06Z masera $ */
18 ///////////////////////////////////////////////////////////////////////////////
20 // class for ITS reconstruction //
22 ///////////////////////////////////////////////////////////////////////////////
24 #include "Riostream.h"
25 #include "AliITSUReconstructor.h"
26 #include "AliITSURecoDet.h"
28 #include "AliRawReader.h"
29 #include "AliESDEvent.h"
31 #include "AliTracker.h"
32 #include "AliITSUTrackerCooked.h"
33 #include "AliITSUCATracker.h"
35 #include "AliITSUGeomTGeo.h"
36 #include "AliITSUSegmentationPix.h"
37 #include "AliITSUDigitPix.h"
38 #include "AliITSUClusterizer.h"
39 #include "AliITSUClusterPix.h"
40 #include "AliITSUVertexer.h"
43 ClassImp(AliITSUReconstructor)
45 //___________________________________________________________________________
46 AliITSUReconstructor::AliITSUReconstructor()
53 // Default constructor
57 //___________________________________________________________________________
58 AliITSUReconstructor::~AliITSUReconstructor()
62 if (!fGeom) return; // was not initialized
64 // same cluster finders and recpoint arrays might be attached to different layers
65 for (int i=fGeom->GetNLayers();i--;) {
66 TObject* clFinder = fClusterFinders.At(i);
68 while (fClusterFinders.Remove(clFinder)) {}
80 //______________________________________________________________________
81 void AliITSUReconstructor::Init()
83 // Initalize this reconstructor
84 // Note: fITS cannot be initialized here since it requires RecoParams (not available ar
85 // the moment of reconstructors initialization)
87 AliInfo("Initializing");
88 if (fGeom) AliFatal("was already done, something is wrong...");
90 fGeom = new AliITSUGeomTGeo(kTRUE,kTRUE);
91 AliITSUClusterPix::SetGeom(fGeom);
93 AliITSUClusterizer* clusPIX = 0;
94 fClusters = new TClonesArray*[fGeom->GetNLayers()];
96 for (int ilr=fGeom->GetNLayers();ilr--;) {
98 int tpDet = fGeom->GetLayerChipTypeID(ilr)/AliITSUGeomTGeo::kMaxSegmPerChipType;
99 if (tpDet == AliITSUGeomTGeo::kChipTypePix) {
100 if (!clusPIX) clusPIX = new AliITSUClusterizer();
101 fClusterFinders.AddAtAndExpand(clusPIX, ilr);
102 fClusters[ilr] = new TClonesArray(AliITSUClusterPix::Class());
104 // to expand the buffers to max.size
105 clusPIX->SetSegmentation((AliITSUSegmentationPix*)fGeom->GetSegmentation(ilr));
109 AliFatal(Form("ClusterFinder for detector type %d is not defined",tpDet));
116 //_____________________________________________________________________________
117 void AliITSUReconstructor::Reconstruct(TTree *digitsTree, TTree *clustersTree) const
119 // reconstruct clusters. If clustersTree is provided, write the tree
120 if (!digitsTree) return;
121 AliDebug(1,"ITSU Cluster finder (from digits tree) is initiated here \n");
123 // At the moment only pixel digits
124 TClonesArray *digArrPix = 0;
125 digitsTree->SetBranchAddress("ITSDigitsPix",&digArrPix);
127 // a new tree is created for each event: add each layer as separate branch
128 TBranch *lrBranch[fGeom->GetNLayers()];
130 for (int ilr=0;ilr<fGeom->GetNLayers();ilr++) {
132 if (clustersTree) { // do we write clusters tree?
133 int tp = fGeom->GetLayerChipTypeID(ilr)/AliITSUGeomTGeo::kMaxSegmPerChipType;
134 if (tp==AliITSUGeomTGeo::kChipTypePix) {
135 lrBranch[ilr] = clustersTree->Bronch(Form("ITSRecPoints%d",ilr),"TClonesArray",&fClusters[ilr]);
138 AliFatal(Form("Detector type %d is not defined",tp));
143 AliITSUClusterizer* clFinder = 0;
144 AliMagF* field = dynamic_cast<AliMagF*>(TGeoGlobalMagField::Instance()->GetField());
146 if (field == 0) AliError("Cannot get magnetic field from TGeoGlobalMagField");
147 else bz = field->SolenoidField();
148 const AliITSURecoParam* recPar = GetRecoParam();
150 for (int ilr=0;ilr<fGeom->GetNLayers();ilr++) {
152 fClusters[ilr]->Clear();
153 clFinder = (AliITSUClusterizer*)fClusterFinders[ilr];
154 clFinder->SetSegmentation((AliITSUSegmentationPix*)fGeom->GetSegmentation(ilr));
155 clFinder->SetLayerID(ilr);
156 clFinder->SetClusters(fClusters[ilr]);
157 clFinder->SetRecoParam(recPar); // RS: Do we need to set it for every event?
158 clFinder->SetAllowDiagonalClusterization(recPar->GetAllowDiagonalClusterization(ilr));
159 clFinder->PrepareLorentzAngleCorrection(bz);
161 int modF=fGeom->GetFirstChipIndex(ilr);
162 int modL=fGeom->GetLastChipIndex(ilr)+1;
163 for (int imod=modF;imod<modL;imod++) {
164 digitsTree->GetEntry(imod);
165 int ndig = digArrPix->GetEntries();
167 clFinder->SetVolID(imod);
168 clFinder->SetDigits(digArrPix);
169 clFinder->Clusterize();
172 // AliITSUClusterPix::SetSortMode( AliITSUClusterPix::SortModeIdTrkYZ());
173 // fClusters[ilr]->Sort();
174 AliDebug(1,Form(" -> Lr%d : %d Cluster",ilr,fClusters[ilr]->GetEntries()));
175 if (clustersTree) lrBranch[ilr]->Fill();
177 if (clustersTree) clustersTree->SetEntries();
181 //_____________________________________________________________________________
182 AliTracker* AliITSUReconstructor::CreateTracker() const
184 // create a ITS tracker
185 Int_t opt = GetRecoParam()->GetTracker();
186 Bool_t sa = GetRecoParam()->GetSAonly();
188 Info("CreateTracker","Creating the ITSU tracker %d in mode %d",opt,sa);
190 AliTracker *tracker=0;
193 tracker = new AliITSUTrackerGlo((AliITSUReconstructor*)this);
196 tracker = new AliITSUTrackerCooked((AliITSUReconstructor*)this);
197 ((AliITSUTrackerCooked*)tracker)->SetSAonly(sa);
200 tracker = new AliITSUCATracker((AliITSUReconstructor*)this);
201 ((AliITSUCATracker*)tracker)->SetSAonly(sa);
204 AliFatal("Undefined ITSU tracker type !");
211 //_____________________________________________________________________________
212 AliVertexer* AliITSUReconstructor::CreateVertexer() const
214 // create a ITS vertexer
216 AliInfo("Creating vertexer using tracklets with the first 3 ITS layers");
218 if (GetRecoParam()->GetEventSpecie() & AliRecoParam::kHighMult) {
219 return new AliITSUVertexer();
221 return new AliITSUVertexer(0.05,0.003,0.04,0.8,3);
225 //_____________________________________________________________________________
226 AliTrackleter* AliITSUReconstructor::CreateMultFinder() const
228 // create the SPD trackeleter for mult. reconstruction
229 // to be implemented for the upgrade
231 AliDebug(1,"ITSU MultFinder should be initiated here\n");
236 //_____________________________________________________________________________
237 AliTracker* AliITSUReconstructor::CreateTrackleter() const
239 // create the SPD trackeleter (for SPD PlaneEfficiency evaluation)
240 // to be implemented for the upgrade
242 AliDebug(1,"ITSU Trackleter should be initiated here\n");
248 //_____________________________________________________________________________
249 Int_t AliITSUReconstructor::LoadClusters(TTree* treeRP)
251 // read clusters from the tree, if it is provided
252 for (int ilr=fGeom->GetNLayers();ilr--;) {
253 if (!fClusters[ilr]) AliFatal(Form("Clusters array for layer %d is not defined",ilr));
254 TBranch* br = treeRP->GetBranch(Form("ITSRecPoints%d",ilr));
255 if (!br) AliFatal(Form("Provided cluster tree does not contain branch for layer %d",ilr));
256 br->SetAddress(&fClusters[ilr]);
258 treeRP->GetEntry(0); // we are still in 1 ev/tree mode...
264 //_____________________________________________________________________________
265 AliITSURecoDet* AliITSUReconstructor::GetITSInterface()
267 // Create reco oriented interface to geometry
268 if (fITS) return fITS;
270 fITS = new AliITSURecoDet(fGeom,"ITSURecoInterface");
271 int nLr = fITS->GetNLayersActive();
272 for (int ilr=nLr;ilr--;) fITS->GetLayerActive(ilr)->SetClusters(GetClusters(ilr));