]> git.uio.no Git - u/mrichter/AliRoot.git/blob - TPC/AliTPCReconstructor.cxx
e9fb3f19f3242e5389d14634603cd2d16df08ef7
[u/mrichter/AliRoot.git] / TPC / AliTPCReconstructor.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 // class for TPC reconstruction                                              //
21 //                                                                           //
22 ///////////////////////////////////////////////////////////////////////////////
23
24
25 #include "AliTPCReconstructor.h"
26 #include "AliRunLoader.h"
27 #include "AliRun.h"
28 #include "AliRawReader.h"
29 #include "AliTPCclustererMI.h"
30 #include "AliTPCtrackerMI.h"
31 #include "AliTPCpidESD.h"
32 #include "AliTPCParam.h"
33 #include "AliTPCParamSR.h"
34
35 ClassImp(AliTPCReconstructor)
36
37 Double_t AliTPCReconstructor::fgCtgRange = 1.05;
38 Double_t AliTPCReconstructor::fgMaxSnpTracker   = 0.95;   // max tangent in tracker - correspond to 3    
39 Double_t AliTPCReconstructor::fgMaxSnpTrack     = 0.999;  // tangent    
40 Int_t    AliTPCReconstructor::fgStreamLevel     = 0;      // stream (debug) level
41 //_____________________________________________________________________________
42 void AliTPCReconstructor::Reconstruct(AliRunLoader* runLoader) const
43 {
44 // reconstruct clusters
45
46   AliLoader* loader = runLoader->GetLoader("TPCLoader");
47   if (!loader) {
48     Error("Reconstruct", "TPC loader not found");
49     return;
50   }
51   loader->LoadRecPoints("recreate");
52   loader->LoadDigits("read");
53
54   AliTPCParam* param = GetTPCParam(runLoader);
55   if (!param) return;
56   AliTPCclustererMI clusterer(param);
57   Int_t nEvents = runLoader->GetNumberOfEvents();
58
59   for (Int_t iEvent = 0; iEvent < nEvents; iEvent++) {
60     runLoader->GetEvent(iEvent);
61
62     TTree* treeClusters = loader->TreeR();
63     if (!treeClusters) {
64       loader->MakeTree("R");
65       treeClusters = loader->TreeR();
66     }
67     TTree* treeDigits = loader->TreeD();
68     if (!treeDigits) {
69       Error("Reconstruct", "Can't get digits tree !");
70       return;
71     }
72
73     clusterer.SetInput(treeDigits);
74     clusterer.SetOutput(treeClusters);
75     clusterer.Digits2Clusters();
76          
77     loader->WriteRecPoints("OVERWRITE");
78   }
79
80   loader->UnloadRecPoints();
81   loader->UnloadDigits();
82 }
83
84 //_____________________________________________________________________________
85 void AliTPCReconstructor::Reconstruct(AliRunLoader* runLoader,
86                                       AliRawReader* rawReader) const
87 {
88 // reconstruct clusters from raw data
89
90   AliLoader* loader = runLoader->GetLoader("TPCLoader");
91   if (!loader) {
92     Error("Reconstruct", "TPC loader not found");
93     return;
94   }
95   loader->LoadRecPoints("recreate");
96
97   AliTPCParam* param = GetTPCParam(runLoader);
98   if (!param) {
99     AliWarning("Loading default TPC parameters !");
100     param = new AliTPCParamSR;
101   }
102   AliTPCclustererMI clusterer(param);
103
104   TString option = GetOption();
105   if (option.Contains("PedestalSubtraction"))
106     clusterer.SetPedSubtraction(kTRUE);
107   if (option.Contains("OldRCUFormat"))
108     clusterer.SetOldRCUFormat(kTRUE);
109  
110   Int_t iEvent = 0;
111   while (rawReader->NextEvent()) {
112     runLoader->GetEvent(iEvent++);
113
114     TTree* treeClusters = loader->TreeR();
115     if (!treeClusters) {
116       loader->MakeTree("R");
117       treeClusters = loader->TreeR();
118     }
119
120     clusterer.SetOutput(treeClusters);
121     clusterer.Digits2Clusters(rawReader);
122          
123     loader->WriteRecPoints("OVERWRITE");
124   }
125
126   loader->UnloadRecPoints();
127 }
128
129 //_____________________________________________________________________________
130 AliTracker* AliTPCReconstructor::CreateTracker(AliRunLoader* runLoader) const
131 {
132 // create a TPC tracker
133
134   AliTPCParam* param = GetTPCParam(runLoader);
135   if (!param) {
136     AliWarning("Loading default TPC parameters !");
137     param = new AliTPCParamSR;
138   }
139   param->ReadGeoMatrices();
140   return new AliTPCtrackerMI(param);
141 }
142
143 //_____________________________________________________________________________
144 void AliTPCReconstructor::FillESD(AliRunLoader* /*runLoader*/, 
145                                   AliESD* esd) const
146 {
147 // make PID
148
149   Double_t parTPC[] = {47., 0.10, 10.};
150   AliTPCpidESD tpcPID(parTPC);
151   tpcPID.MakePID(esd);
152 }
153
154
155 //_____________________________________________________________________________
156 AliTPCParam* AliTPCReconstructor::GetTPCParam(AliRunLoader* runLoader) const
157 {
158 // get the TPC parameters
159
160   TDirectory* saveDir = gDirectory;
161   runLoader->CdGAFile();
162
163   AliTPCParam* param = (AliTPCParam*) gDirectory->Get("75x40_100x60_150x60");
164   if (!param) Error("GetTPCParam", "no TPC parameters found");
165
166   saveDir->cd();
167   return param;
168 }