]>
Commit | Line | Data |
---|---|---|
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 | /////////////////////////////////////////////////////////////////////////////// | |
17 | // | |
18 | // Data container for relative ITS-TPC alignment analysis | |
19 | // Holds an array of AliRelAlignerKalman objects | |
20 | // and takes care of merging when processing data in parallel | |
21 | // | |
22 | // Origin: Mikolaj Krzewicki, Nikhef, Mikolaj.Krzewicki@cern.ch | |
23 | // | |
24 | ////////////////////////////////////////////////////////////////////////////// | |
25 | ||
26 | #include <TNamed.h> | |
27 | #include <TGraphErrors.h> | |
28 | #include <TAxis.h> | |
29 | #include <TTree.h> | |
30 | #include <TMath.h> | |
31 | #include <TObjArray.h> | |
32 | #include <TCollection.h> | |
33 | #include <AliESDEvent.h> | |
34 | #include <AliRelAlignerKalman.h> | |
35 | #include "AliRelAlignerKalmanArray.h" | |
36 | #include "TList.h" | |
37 | #include "TBrowser.h" | |
38 | ||
39 | ClassImp(AliRelAlignerKalmanArray) | |
40 | ||
41 | //______________________________________________________________________________ | |
42 | AliRelAlignerKalmanArray::AliRelAlignerKalmanArray(): | |
43 | TNamed("alignerArray","array of aligners"), | |
44 | fT0(0), | |
45 | fTimebinWidth(0), | |
46 | fSize(0), | |
47 | fOutRejSigmaOnMerge(10.), | |
48 | fOutRejSigmaOnSmooth(1.), | |
49 | fAlignerTemplate(), | |
50 | fPArray(NULL), | |
51 | fListOfGraphs(NULL) | |
52 | { | |
53 | //ctor | |
54 | } | |
55 | ||
56 | //______________________________________________________________________________ | |
57 | AliRelAlignerKalmanArray::AliRelAlignerKalmanArray(Int_t t0, Int_t tend, Int_t slotwidth): | |
58 | TNamed("alignerArray","array of aligners"), | |
59 | fT0(t0), | |
60 | fTimebinWidth(slotwidth), | |
61 | fSize(0), | |
62 | fOutRejSigmaOnMerge(10.), | |
63 | fOutRejSigmaOnSmooth(1.), | |
64 | fAlignerTemplate(), | |
65 | fPArray(NULL), | |
66 | fListOfGraphs(new TList) | |
67 | { | |
68 | //ctor | |
69 | if (slotwidth==0) fSize = 1; | |
70 | else fSize = (tend-t0)/slotwidth; | |
71 | fPArray = new AliRelAlignerKalman* [fSize]; | |
72 | if (fPArray) for (Int_t i=0;i<fSize;i++){fPArray[i]=NULL;}//fill with zeros | |
73 | else fSize=0; | |
74 | fListOfGraphs->SetName("graphs"); | |
75 | fListOfGraphs->SetOwner(kTRUE); | |
76 | } | |
77 | ||
78 | //______________________________________________________________________________ | |
79 | AliRelAlignerKalmanArray::AliRelAlignerKalmanArray( const AliRelAlignerKalmanArray& in): | |
80 | TNamed(in.GetName(), in.GetTitle()), | |
81 | fT0(in.fT0), | |
82 | fTimebinWidth(in.fTimebinWidth), | |
83 | fSize(in.fSize), | |
84 | fOutRejSigmaOnMerge(in.fOutRejSigmaOnMerge), | |
85 | fOutRejSigmaOnSmooth(in.fOutRejSigmaOnSmooth), | |
86 | fAlignerTemplate(in.fAlignerTemplate), | |
87 | fPArray(NULL), | |
88 | fListOfGraphs(new TList) | |
89 | { | |
90 | //copy ctor | |
91 | fPArray = new AliRelAlignerKalman* [fSize]; | |
92 | if (!fPArray) {fSize=0;return;} //if fail | |
93 | for (Int_t i=0;i<fSize;i++) | |
94 | { | |
95 | if (in.fPArray[i]) fPArray[i]=new AliRelAlignerKalman(*(in.fPArray[i])); | |
96 | else fPArray[i]=NULL; | |
97 | } | |
98 | fListOfGraphs->SetName("graphs"); | |
99 | fListOfGraphs->SetOwner(kTRUE); | |
100 | } | |
101 | ||
102 | //______________________________________________________________________________ | |
103 | AliRelAlignerKalmanArray::~AliRelAlignerKalmanArray() | |
104 | { | |
105 | //dtor | |
106 | ClearContents(); | |
107 | delete [] fPArray; | |
108 | delete fListOfGraphs; | |
109 | } | |
110 | ||
111 | //______________________________________________________________________________ | |
112 | void AliRelAlignerKalmanArray::SetupArray(Int_t t0, Int_t tend, Int_t slotwidth) | |
113 | { | |
114 | //setup array - clears old content | |
115 | ClearContents(); | |
116 | fT0=t0; | |
117 | fTimebinWidth=slotwidth; | |
118 | Int_t tmpsize; | |
119 | if (slotwidth==0) tmpsize = 1; | |
120 | else tmpsize = (tend-t0)/slotwidth; | |
121 | if (tmpsize!=fSize) | |
122 | { | |
123 | delete [] fPArray; | |
124 | fPArray=new AliRelAlignerKalman* [tmpsize]; | |
125 | if (fPArray) fSize=tmpsize; | |
126 | else fSize=0; | |
127 | } | |
128 | for (Int_t i=0;i<fSize;i++){fPArray[i]=NULL;}//fill with zeros | |
129 | } | |
130 | ||
131 | //______________________________________________________________________________ | |
132 | AliRelAlignerKalman* AliRelAlignerKalmanArray::GetAlignerTemplate() | |
133 | { | |
134 | //get aligner template | |
135 | return &fAlignerTemplate; | |
136 | } | |
137 | ||
138 | //______________________________________________________________________________ | |
139 | Long64_t AliRelAlignerKalmanArray::Merge( TCollection* list ) | |
140 | { | |
141 | //Merge all the arrays | |
142 | //tlihe merge is vertical, meaning matching entries in tree are merged | |
143 | ||
144 | AliRelAlignerKalmanArray *arrayFromList; | |
145 | if (!list) return 0; | |
146 | TIter next(list); | |
147 | while ( (arrayFromList = dynamic_cast<AliRelAlignerKalmanArray*>(next())) ) | |
148 | { | |
149 | if (fT0 != arrayFromList->fT0) continue; | |
150 | if (fTimebinWidth != arrayFromList->fTimebinWidth) continue; | |
151 | if (fSize != arrayFromList->fSize) continue; | |
152 | ||
153 | for (Int_t i=0; i<GetSize(); i++) | |
154 | { | |
155 | AliRelAlignerKalman* a1 = fPArray[i]; | |
156 | AliRelAlignerKalman* a2 = arrayFromList->fPArray[i]; | |
157 | if (a1 && a2) | |
158 | { | |
159 | a1->SetRejectOutliers(kFALSE); | |
160 | a1->SetOutRejSigma(fOutRejSigmaOnMerge); a1->Merge(a2); | |
161 | } | |
162 | else | |
163 | if (!a1 && a2) fPArray[i] = new AliRelAlignerKalman(*a2); | |
164 | } | |
165 | } | |
166 | fListOfGraphs->Delete(); | |
167 | return 0; | |
168 | } | |
169 | ||
170 | //______________________________________________________________________________ | |
171 | Int_t AliRelAlignerKalmanArray::Timebin( UInt_t timestamp ) const | |
172 | { | |
173 | //calculate binnumber for given timestamp | |
174 | if (fTimebinWidth==0) return 0; | |
175 | Int_t slot = (timestamp-fT0)/fTimebinWidth; //it's all integers! | |
176 | return slot; | |
177 | } | |
178 | ||
179 | //______________________________________________________________________________ | |
180 | AliRelAlignerKalman* AliRelAlignerKalmanArray::GetAligner(UInt_t ts) | |
181 | { | |
182 | //get the aligner for specified timestamp | |
183 | Int_t tb = Timebin(ts); | |
184 | if (tb<0) return NULL; | |
185 | if (tb>=fSize) return NULL; | |
186 | if (!fPArray) return NULL; | |
187 | if (!fPArray[tb]) | |
188 | { | |
189 | fPArray[tb] = new AliRelAlignerKalman( fAlignerTemplate ); | |
190 | fPArray[tb]->SetTimeStamp(fT0+tb*fTimebinWidth); | |
191 | } | |
192 | return fPArray[tb]; | |
193 | } | |
194 | ||
195 | //______________________________________________________________________________ | |
196 | AliRelAlignerKalman* AliRelAlignerKalmanArray::GetAligner(AliESDEvent* event) | |
197 | { | |
198 | //get the aligner for this event and set relevant info | |
199 | AliRelAlignerKalman* a = GetAligner(event->GetTimeStamp()); | |
200 | if (a) a->SetRunNumber(event->GetRunNumber()); | |
201 | if (a) a->SetMagField(event->GetMagneticField()); | |
202 | return a; | |
203 | } | |
204 | ||
205 | //______________________________________________________________________________ | |
206 | AliRelAlignerKalmanArray& AliRelAlignerKalmanArray::operator=(const AliRelAlignerKalmanArray& in) | |
207 | { | |
208 | //assignment operator | |
209 | if(&in == this) return *this; | |
210 | TNamed::operator=(in); | |
211 | ||
212 | if (fSize!=in.fSize) | |
213 | { | |
214 | //if sizes different, delete curent and make a new one with new size | |
215 | ClearContents(); | |
216 | delete [] fPArray; | |
217 | fPArray = new AliRelAlignerKalman* [in.fSize]; | |
218 | } | |
219 | ||
220 | fOutRejSigmaOnMerge=in.fOutRejSigmaOnMerge; | |
221 | fOutRejSigmaOnSmooth=in.fOutRejSigmaOnSmooth; | |
222 | ||
223 | if (!fPArray) fSize=0; | |
224 | else fSize = in.fSize; | |
225 | fTimebinWidth = in.fTimebinWidth; | |
226 | fT0 = in.fT0; | |
227 | for (Int_t i=0;i<fSize;i++) | |
228 | { | |
229 | if (in.fPArray[i]) fPArray[i] = new AliRelAlignerKalman(*(in.fPArray[i])); | |
230 | else fPArray[i]=NULL; | |
231 | } | |
232 | return *this; | |
233 | } | |
234 | ||
235 | //______________________________________________________________________________ | |
236 | void AliRelAlignerKalmanArray::ClearContents() | |
237 | { | |
238 | //clear contents of array | |
239 | for (Int_t i=0;i<fSize;i++) | |
240 | { | |
241 | if (fPArray[i]) delete fPArray[i]; | |
242 | } | |
243 | } | |
244 | ||
245 | //______________________________________________________________________________ | |
246 | AliRelAlignerKalman* AliRelAlignerKalmanArray::At( Int_t i ) const | |
247 | { | |
248 | //mimic TObjArray::At( Int_t i ) | |
249 | if (i>=fSize) {printf("AliRelAlignerKalmanArray::At index %i out of bounds, size=%i\n",i,fSize); return NULL;} | |
250 | if (i<0) return NULL; | |
251 | return fPArray[i]; | |
252 | } | |
253 | ||
254 | //______________________________________________________________________________ | |
255 | AliRelAlignerKalman* AliRelAlignerKalmanArray::operator[](Int_t i) const | |
256 | { | |
257 | //mimic TObjArray::operator[](Int_t) | |
258 | if (i>=fSize) {printf("AliRelAlignerKalmanArray::operator[] index %i out of bounds, size=%i\n",i,fSize); return NULL;} | |
259 | if (i<0) return NULL; | |
260 | return fPArray[i]; | |
261 | } | |
262 | ||
263 | //______________________________________________________________________________ | |
264 | AliRelAlignerKalman*& AliRelAlignerKalmanArray::operator[](Int_t i) | |
265 | { | |
266 | //mimic TObjArray::operator[](Int_t) can be used as lvalue | |
267 | if (i>=fSize||i<0) {Error("operator[]", "index %i out of bounds, size=%i\n",i,fSize); return fPArray[0];} | |
268 | return fPArray[i]; | |
269 | } | |
270 | ||
271 | //______________________________________________________________________________ | |
272 | AliRelAlignerKalman* AliRelAlignerKalmanArray::Last() const | |
273 | { | |
274 | //return aligner in last filled slot | |
275 | if (fSize==0) return NULL; | |
276 | return fPArray[fSize-1]; | |
277 | } | |
278 | ||
279 | //______________________________________________________________________________ | |
280 | void AliRelAlignerKalmanArray::Print(Option_t* option) const | |
281 | { | |
282 | // print some info | |
283 | TString optionStr(option); | |
284 | printf( "%s: t0: %i, tend: %i, width: %i, size: %i, entries: %i\n", | |
285 | GetName(), | |
286 | fT0, (fT0+(fSize)*fTimebinWidth), fTimebinWidth, | |
287 | fSize, GetEntries() ); | |
288 | if (optionStr.Contains("a")) | |
289 | for (Int_t i=0; i<fSize; i++) | |
290 | { | |
291 | AliRelAlignerKalman* al = fPArray[i]; | |
292 | if (!al) continue; | |
293 | al->Print(); | |
294 | } | |
295 | } | |
296 | ||
297 | //______________________________________________________________________________ | |
298 | Int_t AliRelAlignerKalmanArray::GetEntries() const | |
299 | { | |
300 | //get number of filled slots | |
301 | if (!fPArray) return 0; | |
302 | Int_t entries=0; | |
303 | for (Int_t i=0; i<fSize; i++) | |
304 | { | |
305 | if (fPArray[i]) entries++; | |
306 | } | |
307 | return entries; | |
308 | } | |
309 | ||
310 | //______________________________________________________________________________ | |
311 | void AliRelAlignerKalmanArray::FillTree( TTree* tree ) const | |
312 | { | |
313 | AliRelAlignerKalman* al = NULL; | |
314 | tree->Branch("aligner","AliRelAlignerKalman",&al); | |
315 | //fill a tree with filled slots | |
316 | for (Int_t i=0; i<fSize; i++) | |
317 | { | |
318 | al = fPArray[i]; | |
319 | if (al) tree->Fill(); | |
320 | } | |
321 | } | |
322 | ||
323 | //______________________________________________________________________________ | |
324 | TGraphErrors* AliRelAlignerKalmanArray::MakeGraph(Int_t iparam) const | |
325 | { | |
326 | //return a graph for the iparam-th parameter | |
327 | if (iparam>8 || iparam<0) | |
328 | { | |
329 | printf("wrong parameter number. must be from 0-8"); | |
330 | return NULL; | |
331 | } | |
332 | ||
333 | Int_t n=GetEntries(); | |
334 | TVectorD vx(n); | |
335 | TVectorD vy(n); | |
336 | TVectorD vex(n); | |
337 | TVectorD vey(n); | |
338 | Int_t entry=0; | |
339 | for (Int_t i=0; i<fSize; i++) | |
340 | { | |
341 | AliRelAlignerKalman* al = fPArray[i]; | |
342 | if (!al) continue; | |
343 | vx(entry) = al->GetTimeStamp(); | |
344 | vy(entry) = al->GetStateArr()[iparam]; | |
345 | TMatrixDSym* cm = al->GetStateCov(); | |
346 | vey(entry) = TMath::Sqrt((*cm)(iparam,iparam)); | |
347 | entry++; | |
348 | } | |
349 | ||
350 | TGraphErrors* graph = new TGraphErrors(vx,vy,vex,vey); | |
351 | TString graphtitle; | |
352 | TString graphtitley; | |
353 | switch (iparam) | |
354 | { | |
355 | case 0: | |
356 | graphtitle="rotation \\psi"; | |
357 | graphtitley="\\psi [rad]"; | |
358 | break; | |
359 | case 1: | |
360 | graphtitle="rotation \\theta"; | |
361 | graphtitley="\\theta [rad]"; | |
362 | break; | |
363 | case 2: | |
364 | graphtitle="rotation \\phi"; | |
365 | graphtitley="\\phi [rad]"; | |
366 | break; | |
367 | case 3: | |
368 | graphtitle="shift x"; | |
369 | graphtitley="x [cm]"; | |
370 | break; | |
371 | case 4: | |
372 | graphtitle="shift y"; | |
373 | graphtitley="y [cm]"; | |
374 | break; | |
375 | case 5: | |
376 | graphtitle="shift z"; | |
377 | graphtitley="z [cm]"; | |
378 | break; | |
379 | case 6: | |
380 | graphtitle="TPC vd correction"; | |
381 | graphtitley="correction factor []"; | |
382 | break; | |
383 | case 7: | |
384 | graphtitle="TPC t0 correction"; | |
385 | graphtitley="t0 correction [\\micros]"; | |
386 | break; | |
387 | case 8: | |
388 | graphtitle="TPC dv/dy"; | |
389 | graphtitley="dv/dy [cm/\\micros/m]"; | |
390 | break; | |
391 | } | |
392 | graph->SetName(graphtitle); | |
393 | graph->SetTitle(graphtitle); | |
394 | TAxis* xas = graph->GetXaxis(); | |
395 | TAxis* yas = graph->GetYaxis(); | |
396 | xas->SetTitle("time"); | |
397 | xas->SetTimeDisplay(1); | |
398 | yas->SetTitle(graphtitley); | |
399 | return graph; | |
400 | } | |
401 | ||
402 | //______________________________________________________________________________ | |
403 | AliRelAlignerKalmanArray* AliRelAlignerKalmanArray::MakeSmoothArray() const | |
404 | { | |
405 | //return a smoothed version of the data | |
406 | AliRelAlignerKalmanArray* outputarr = new AliRelAlignerKalmanArray(fT0, | |
407 | fT0+fSize*fTimebinWidth,fTimebinWidth); | |
408 | ||
409 | AliRelAlignerKalman* tmpaligner = outputarr->GetAlignerTemplate(); | |
410 | tmpaligner->SetOutRejSigma(fOutRejSigmaOnSmooth); | |
411 | //copy the first filled slot | |
412 | Int_t n=0; | |
413 | while (!fPArray[n]) {n++;} | |
414 | if (n==fSize) return NULL; | |
415 | *tmpaligner = *fPArray[n]; | |
416 | if (fPArray[n]->GetNUpdates()>10) | |
417 | (*outputarr)[n] = new AliRelAlignerKalman(*(fPArray[n])); | |
418 | //for the rest of slots use merge | |
419 | for (Int_t i=n+1;i<fSize;i++) | |
420 | { | |
421 | if (!fPArray[i]) continue; | |
422 | PropagateToTime(tmpaligner, fPArray[i]->GetTimeStamp()); | |
423 | if (tmpaligner->Merge(fPArray[i])) | |
424 | (*outputarr)[i] = new AliRelAlignerKalman(*tmpaligner); | |
425 | else | |
426 | (*outputarr)[i] = new AliRelAlignerKalman(*(fPArray[i])); | |
427 | } | |
428 | ||
429 | tmpaligner->Reset(); //clean up | |
430 | ||
431 | return outputarr; | |
432 | } | |
433 | ||
434 | //______________________________________________________________________________ | |
435 | void AliRelAlignerKalmanArray::PropagateToTime(AliRelAlignerKalman* al, UInt_t timestamp ) const | |
436 | { | |
437 | //propagate an aligner in time | |
438 | TMatrixDSym* cov = al->GetStateCov(); | |
439 | TMatrixDSym corr(TMatrixDSym::kZero,*cov); | |
440 | UInt_t dt = (timestamp>al->GetTimeStamp())? | |
441 | timestamp-al->GetTimeStamp():al->GetTimeStamp()-timestamp; | |
442 | //the propagation matrix | |
443 | //to be added to the covariance matrix of kalman filter | |
444 | corr(6,6) = dt*1e-8; //vdrift | |
445 | (*cov) += corr; | |
446 | } | |
447 | ||
448 | //______________________________________________________________________________ | |
449 | void AliRelAlignerKalmanArray::Browse(TBrowser* b ) | |
450 | { | |
451 | if (!b) return; | |
452 | if (!fListOfGraphs) | |
453 | { | |
454 | fListOfGraphs=new TList(); | |
455 | fListOfGraphs->SetName("graphs"); | |
456 | fListOfGraphs->SetOwner(kTRUE); | |
457 | } | |
458 | for (Int_t i=0;i<9;i++) | |
459 | { | |
460 | TGraphErrors* gr = dynamic_cast<TGraphErrors*>(fListOfGraphs->At(i)); | |
461 | if (gr) continue; | |
462 | fListOfGraphs->AddAt(MakeGraph(i),i); | |
463 | } | |
464 | b->Add(fListOfGraphs); | |
465 | } | |
466 |