f80b98cb |
1 | //Author: Anders Strand Vestbo |
2 | //Last Modified: 28.6.01 |
3 | |
162e2e5a |
4 | #include <math.h> |
5 | #include <TTree.h> |
6 | #include <TFile.h> |
7 | |
8 | #include "GetGoodParticles.h" |
9 | #include "AliL3TrackArray.h" |
10 | #include "AliL3Logging.h" |
99e7186b |
11 | #include "AliL3HoughEval.h" |
4de874d1 |
12 | #include "AliL3HoughTransformer.h" |
99e7186b |
13 | #include "AliL3DigitData.h" |
4de874d1 |
14 | #include "AliL3HoughTrack.h" |
99e7186b |
15 | #include "AliL3Transform.h" |
16 | #include "AliL3Histogram.h" |
17 | #include "AliL3Defs.h" |
4de874d1 |
18 | |
19 | ClassImp(AliL3HoughEval) |
20 | |
4de874d1 |
21 | AliL3HoughEval::AliL3HoughEval() |
22 | { |
4fc9a6a4 |
23 | |
99e7186b |
24 | } |
4de874d1 |
25 | |
26 | AliL3HoughEval::AliL3HoughEval(AliL3HoughTransformer *transformer) |
27 | { |
4fc9a6a4 |
28 | |
4de874d1 |
29 | fHoughTransformer = transformer; |
30 | fTransform = new AliL3Transform(); |
99e7186b |
31 | |
32 | fSlice = fHoughTransformer->GetSlice(); |
33 | fPatch = fHoughTransformer->GetPatch(); |
34 | fNrows = NRows[fPatch][1] - NRows[fPatch][0] + 1; |
35 | fNEtaSegments = fHoughTransformer->GetNEtaSegments(); |
36 | fEtaMin = fHoughTransformer->GetEtaMin(); |
37 | fEtaMax = fHoughTransformer->GetEtaMax(); |
38 | fRemoveFoundTracks = kFALSE; |
f80b98cb |
39 | fNumOfPadsToLook = 1; |
99e7186b |
40 | fNumOfRowsToMiss = 1; |
41 | GenerateLUT(); |
4de874d1 |
42 | } |
43 | |
4de874d1 |
44 | AliL3HoughEval::~AliL3HoughEval() |
45 | { |
4de874d1 |
46 | if(fTransform) |
47 | delete fTransform; |
99e7186b |
48 | if(fRowPointers) |
49 | delete [] fRowPointers; |
4de874d1 |
50 | } |
51 | |
99e7186b |
52 | void AliL3HoughEval::GenerateLUT() |
4de874d1 |
53 | { |
99e7186b |
54 | //Generate a LUT, to limit the access to raw data |
55 | |
56 | fRowPointers = new AliL3DigitRowData*[fNrows]; |
57 | |
58 | AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fHoughTransformer->GetDataPointer(); |
59 | if(!tempPt) |
60 | printf("AliL3HoughEval::GenerateLUT : Zero data pointer\n"); |
f80b98cb |
61 | |
99e7186b |
62 | for(Int_t i=NRows[fPatch][0]; i<=NRows[fPatch][1]; i++) |
f80b98cb |
63 | { |
99e7186b |
64 | Int_t prow = i - NRows[fPatch][0]; |
65 | fRowPointers[prow] = tempPt; |
4fc9a6a4 |
66 | fHoughTransformer->UpdateDataPointer(tempPt); |
f80b98cb |
67 | } |
68 | |
99e7186b |
69 | } |
f80b98cb |
70 | |
99e7186b |
71 | Bool_t AliL3HoughEval::LookInsideRoad(AliL3HoughTrack *track,Int_t eta_index,Bool_t remove) |
72 | { |
73 | //Look at rawdata along the road specified by the track candidates. |
74 | //If track is good, return true, if not return false. |
162e2e5a |
75 | |
f80b98cb |
76 | Int_t sector,row; |
99e7186b |
77 | |
f80b98cb |
78 | Int_t nrow=0,npixs=0; |
79 | Float_t xyz[3]; |
99e7186b |
80 | |
81 | Double_t etaslice = (fEtaMax - fEtaMin)/fNEtaSegments; |
82 | for(Int_t padrow = NRows[fPatch][0]; padrow <= NRows[fPatch][1]; padrow++) |
f80b98cb |
83 | { |
99e7186b |
84 | Int_t prow = padrow - NRows[fPatch][0]; |
85 | |
f80b98cb |
86 | if(!track->GetCrossingPoint(padrow,xyz)) |
87 | { |
88 | printf("AliL3HoughEval::LookInsideRoad : Track does not cross line!!\n"); |
89 | continue; |
90 | } |
91 | |
99e7186b |
92 | fTransform->Slice2Sector(fSlice,padrow,sector,row); |
f80b98cb |
93 | fTransform->Local2Raw(xyz,sector,row); |
94 | npixs=0; |
99e7186b |
95 | |
96 | //Get the timebins for this pad |
97 | AliL3DigitRowData *tempPt = fRowPointers[prow]; |
98 | if(!tempPt) |
f80b98cb |
99 | { |
99e7186b |
100 | printf("AliL3HoughEval::LookInsideRoad : Zero data pointer\n"); |
101 | continue; |
f80b98cb |
102 | } |
103 | |
99e7186b |
104 | //Look at both sides of the pad: |
162e2e5a |
105 | for(Int_t p=(Int_t)rint(xyz[1])-fNumOfPadsToLook; p<=(Int_t)rint(xyz[1])+fNumOfPadsToLook; p++) |
99e7186b |
106 | { |
107 | AliL3DigitData *digPt = tempPt->fDigitData; |
108 | for(UInt_t j=0; j<tempPt->fNDigit; j++) |
109 | { |
110 | UChar_t pad = digPt[j].fPad; |
162e2e5a |
111 | |
99e7186b |
112 | if(pad < p) continue; |
113 | if(pad > p) break; |
114 | UShort_t time = digPt[j].fTime; |
162e2e5a |
115 | Double_t eta = fTransform->GetEta(padrow,pad,time); |
99e7186b |
116 | Int_t pixel_index = (Int_t)(eta/etaslice); |
117 | if(pixel_index > eta_index) continue; |
118 | if(pixel_index != eta_index) break; |
119 | if(remove) |
120 | digPt[j].fCharge = 0; //Delete the track from image |
121 | npixs++; |
122 | } |
123 | } |
124 | |
162e2e5a |
125 | if(npixs > 1)//At least 2 digits on this padrow |
f80b98cb |
126 | { |
127 | nrow++; |
128 | } |
129 | } |
99e7186b |
130 | if(remove) |
131 | return kTRUE; |
132 | |
133 | if(nrow >= fNrows - fNumOfRowsToMiss)//this was a good track |
f80b98cb |
134 | { |
135 | track->SetEtaIndex(eta_index); |
99e7186b |
136 | if(fRemoveFoundTracks) |
137 | LookInsideRoad(track,eta_index,kTRUE); |
f80b98cb |
138 | return kTRUE; |
139 | } |
140 | else |
141 | return kFALSE; |
142 | } |
143 | |
99e7186b |
144 | void AliL3HoughEval::DisplayEtaSlice(Int_t eta_index,AliL3Histogram *hist) |
f80b98cb |
145 | { |
99e7186b |
146 | //Display the current raw data inside the slice |
4de874d1 |
147 | |
99e7186b |
148 | if(!hist) |
4de874d1 |
149 | { |
99e7186b |
150 | printf("AliL3HoughEval::DisplayEtaSlice : No input histogram!\n"); |
4de874d1 |
151 | return; |
152 | } |
153 | |
99e7186b |
154 | Double_t etaslice = (fEtaMax - fEtaMin)/fNEtaSegments; |
155 | for(Int_t padrow = NRows[fPatch][0]; padrow <= NRows[fPatch][1]; padrow++) |
f80b98cb |
156 | { |
99e7186b |
157 | Int_t prow = padrow - NRows[fPatch][0]; |
158 | |
159 | AliL3DigitRowData *tempPt = fRowPointers[prow]; |
160 | if(!tempPt) |
f80b98cb |
161 | { |
99e7186b |
162 | printf("AliL3HoughEval::DisplayEtaSlice : Zero data pointer\n"); |
f80b98cb |
163 | continue; |
164 | } |
165 | |
99e7186b |
166 | AliL3DigitData *digPt = tempPt->fDigitData; |
167 | for(UInt_t j=0; j<tempPt->fNDigit; j++) |
f80b98cb |
168 | { |
99e7186b |
169 | UChar_t pad = digPt[j].fPad; |
170 | UChar_t charge = digPt[j].fCharge; |
171 | UShort_t time = digPt[j].fTime; |
4fc9a6a4 |
172 | if(charge < fHoughTransformer->GetThreshold()) continue; |
99e7186b |
173 | Float_t xyz[3]; |
4de874d1 |
174 | Int_t sector,row; |
99e7186b |
175 | fTransform->Slice2Sector(fSlice,padrow,sector,row); |
176 | fTransform->Raw2Local(xyz,sector,row,pad,time); |
177 | Double_t eta = fTransform->GetEta(xyz); |
178 | Int_t pixel_index = (Int_t)(eta/etaslice); |
179 | if(pixel_index != eta_index) continue; |
180 | hist->Fill(xyz[0],xyz[1],charge); |
f80b98cb |
181 | } |
182 | } |
f80b98cb |
183 | |
4de874d1 |
184 | } |
162e2e5a |
185 | |
186 | void AliL3HoughEval::CompareMC(AliL3TrackArray *tracks,Char_t *trackfile) |
187 | { |
188 | |
189 | struct GoodTrack goodtracks[15000]; |
190 | Int_t nt=0; |
191 | ifstream in(trackfile); |
192 | if(in) |
193 | { |
194 | printf("Reading good tracks from file %s\n",trackfile); |
195 | while (in>>goodtracks[nt].label>>goodtracks[nt].code>> |
196 | goodtracks[nt].px>>goodtracks[nt].py>>goodtracks[nt].pz>> |
197 | goodtracks[nt].pt>>goodtracks[nt].eta>>goodtracks[nt].nhits) |
198 | { |
199 | nt++; |
200 | if (nt==15000) |
201 | { |
202 | cerr<<"Too many good tracks"<<endl; |
203 | break; |
204 | } |
205 | } |
206 | if (!in.eof()) |
207 | { |
208 | LOG(AliL3Log::kError,"AliL3HoughEval::CompareMC","Input file") |
209 | <<"Error in file reading"<<ENDLOG; |
210 | return; |
211 | } |
212 | } |
213 | else |
214 | { |
215 | LOG(AliL3Log::kError,"AliL3HoughEval::CompareMC","Input") |
216 | <<"No input trackfile "<<trackfile<<ENDLOG; |
217 | } |
218 | |
219 | Int_t *particles = new Int_t[fNEtaSegments]; |
220 | Int_t *ftracks = new Int_t[fNEtaSegments]; |
221 | for(Int_t i=0; i<fNEtaSegments; i++) |
222 | { |
223 | particles[i]=0; |
224 | ftracks[i]=0; |
225 | } |
226 | |
227 | Double_t etaslice = (fEtaMax - fEtaMin)/fNEtaSegments; |
228 | for(Int_t i=0; i<tracks->GetNTracks(); i++) |
229 | { |
230 | AliL3HoughTrack *tr = (AliL3HoughTrack*)tracks->GetCheckedTrack(i); |
231 | if(!tr) continue; |
232 | Int_t trackindex = tr->GetEtaIndex(); |
233 | if(trackindex <0 || trackindex >= fNEtaSegments) continue; |
234 | ftracks[trackindex]++; |
235 | } |
236 | for(Int_t i=0; i<nt; i++) |
237 | { |
238 | if(goodtracks[i].nhits < 150) continue; |
239 | if(goodtracks[i].pt < 0.5) continue; |
240 | Int_t particleindex = (Int_t)(goodtracks[i].eta/etaslice); |
241 | if(particleindex < 0 || particleindex >= fNEtaSegments) continue; |
242 | particles[particleindex]++; |
243 | } |
244 | |
245 | Double_t found=0; |
246 | Double_t good =0; |
247 | for(Int_t i=0; i<fNEtaSegments; i++) |
248 | { |
249 | printf("Slice %d : Found tracks %d, good tracks %d\n",i,ftracks[i],particles[i]); |
250 | found += ftracks[i]; |
251 | good += particles[i]; |
252 | } |
253 | printf("And the total efficiency was: %f\n",found/good); |
254 | |
255 | delete [] particles; |
256 | delete [] ftracks; |
257 | |
258 | } |
259 | |