Line data Source code
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 : // An AliTRDalignment object contains the alignment data (3 shifts and 3 //
20 : // tilts) for all the alignable volumes of the TRD, i.e. for 18 supermodules //
21 : // and 540 chambers. The class provides simple tools for reading and writing //
22 : // these data in different formats, and for generating fake data that can be //
23 : // used to simulate misalignment. //
24 : // The six alignment variables have the following meaning: //
25 : // shift in rphi //
26 : // shift in z //
27 : // shift in r //
28 : // tilt around rphi //
29 : // tilt around z //
30 : // tilt around r //
31 : // The shifts are in cm and the tilts are in degrees. //
32 : // The currently supported formats are: //
33 : // - ascii //
34 : // - root file containing a TClonesArray of alignment objects //
35 : // - offline conditions database //
36 : // - OCDB-like root file //
37 : // - geometry file (like misaligned_geometry.root) //
38 : // //
39 : // Some examples of usage (in an aliroot session): //
40 : // AliTRDalignment a,b,c,d,e; //
41 : // double xsm[]={0,0,0,-70,0,0}; //
42 : // double xch[]={0,0,-50,0,0,0}; //
43 : // a.SetSm(4,xsm); //
44 : // a.SetCh(120,xch); //
45 : // a.WriteAscii("kuku.dat"); //
46 : // TGeoManager::Import("geometry.root"); a.WriteRoot("kuku.root"); //
47 : // TGeoManager::Import("geometry.root"); a.WriteDB("kukudb.root",0,0); //
48 : // TGeoManager::Import("geometry.root"); //
49 : // a.WriteDB("local://$ALICE_ROOT/OCDB", "TRD/Align/Data", 0,0); //
50 : // TGeoManager::Import("geometry.root"); a.WriteGeo("kukugeometry.root"); //
51 : // //
52 : // b.ReadAscii("kuku.dat"); //
53 : // TGeoManager::Import("geometry.root"); c.ReadRoot("kuku.root"); //
54 : // TGeoManager::Import("geometry.root"); d.ReadDB("kukudb.root"); //
55 : // TGeoManager::Import("kukugeometry.root"); e.ReadCurrentGeo(); //
56 : // //
57 : // e.PrintSm(4); //
58 : // e.PrintCh(120); //
59 : // a.PrintRMS(); //
60 : // b.PrintRMS(); //
61 : // e.PrintRMS(); //
62 : // //
63 : // //
64 : // D.Miskowiec, November 2006 //
65 : // //
66 : ///////////////////////////////////////////////////////////////////////////////
67 :
68 : #include <iostream>
69 : #include <fstream>
70 :
71 : #include "TMath.h"
72 : #include "TFile.h"
73 : #include "TGeoManager.h"
74 : #include "TGeoPhysicalNode.h"
75 : #include "TClonesArray.h"
76 : #include "TString.h"
77 : #include "TFitter.h"
78 : #include "TMinuit.h"
79 :
80 : #include "AliLog.h"
81 : #include "AliAlignObj.h"
82 : #include "AliAlignObjParams.h"
83 : #include "AliCDBManager.h"
84 : #include "AliCDBStorage.h"
85 : #include "AliCDBMetaData.h"
86 : #include "AliCDBEntry.h"
87 : #include "AliSurveyObj.h"
88 : #include "AliSurveyPoint.h"
89 :
90 : #include "AliTRDalignment.h"
91 :
92 : void trdAlignmentFcn(Int_t &npar, Double_t *gin, Double_t &f, Double_t *x, Int_t iflag);
93 :
94 : using std::ostream;
95 : using std::fstream;
96 48 : ClassImp(AliTRDalignment)
97 :
98 : //_____________________________________________________________________________
99 : AliTRDalignment::AliTRDalignment()
100 0 : :TObject()
101 0 : ,fComment()
102 0 : ,fRan(0)
103 0 : {
104 : //
105 : // constructor
106 : //
107 :
108 0 : SetZero();
109 :
110 0 : for (int i=0; i<18; i++) for (int j=0; j<2; j++) for (int k=0; k<2; k++) for (int l=0; l<2; l++) {
111 0 : fSurveyX[i][j][k][l] = 0.0;
112 0 : fSurveyY[i][j][k][l] = 0.0;
113 0 : fSurveyZ[i][j][k][l] = 0.0;
114 0 : fSurveyEX[i][j][k][l] = 0.0;
115 0 : fSurveyEY[i][j][k][l] = 0.0;
116 0 : fSurveyEZ[i][j][k][l] = 0.0;
117 : }
118 :
119 : // Initialize the nominal positions of the survey points
120 : // in the local frame of supermodule (where y is the long side,
121 : // z corresponds to the radius in lab, and x to the phi in lab).
122 : // Four survey marks are on each z-side of the supermodule.
123 : // A B
124 : // ----o-----------o---- x |
125 : // \ / |
126 : // \ / |
127 : // \ / |
128 : // \ / |
129 : // ---o-----o--- -------------->
130 : // C D y
131 : //
132 : // For the purpose of this explanation lets define the origin such that
133 : // the supermodule occupies 0 < x < 77.9 cm. Then the coordinates (x,y)
134 : // are (in cm)
135 : // A (76.2,-30.25)
136 : // B (76.2,+30.25)
137 : // C ( 2.2,-22.5 )
138 : // D ( 2.2,+22.5 )
139 : //
140 :
141 0 : double x[2] = {22.5,30.25}; // lab phi, or tracking-y
142 0 : double y[2] = {353.0, -353.0}; // lab z; inc. 2 cm survey target offset
143 0 : double z[2] = {-(77.9/2.0-2.0),77.9/2.0-1.5}; // lab r, or better tracking-x
144 :
145 0 : for (int j=0; j<2; j++) for (int k=0; k<2; k++) for (int l=0; l<2; l++) {
146 0 : fSurveyX0[j][k][l] = -TMath::Power(-1,l) * x[k];
147 0 : fSurveyY0[j][k][l] = y[j];
148 0 : fSurveyZ0[j][k][l] = z[k];
149 : }
150 :
151 0 : for (int i=0; i<1000; i++) {
152 0 : fIbuffer[i] = 0;
153 0 : fDbuffer[i] = 0.0;
154 : }
155 :
156 0 : }
157 :
158 : //_____________________________________________________________________________
159 : AliTRDalignment::AliTRDalignment(const AliTRDalignment& source)
160 0 : :TObject(source)
161 0 : ,fComment(source.fComment)
162 0 : ,fRan(source.fRan)
163 0 : {
164 : //
165 : // copy constructor
166 : //
167 :
168 0 : for (int i=0; i<18; i++) SetSm(i,source.fSm[i]);
169 0 : for (int i=0; i<540; i++) SetCh(i,source.fCh[i]);
170 0 : for (int i=0; i<18; i++) for (int j=0; j<2; j++) for (int k=0; k<2; k++) for (int l=0; l<2; l++) {
171 0 : fSurveyX[i][j][k][l] = source.fSurveyX[i][j][k][l];
172 0 : fSurveyY[i][j][k][l] = source.fSurveyY[i][j][k][l];
173 0 : fSurveyZ[i][j][k][l] = source.fSurveyZ[i][j][k][l];
174 0 : fSurveyEX[i][j][k][l] = source.fSurveyEX[i][j][k][l];
175 0 : fSurveyEY[i][j][k][l] = source.fSurveyEY[i][j][k][l];
176 0 : fSurveyEZ[i][j][k][l] = source.fSurveyEZ[i][j][k][l];
177 : }
178 0 : for (int j=0; j<2; j++) for (int k=0; k<2; k++) for (int l=0; l<2; l++) {
179 0 : fSurveyX0[j][k][l] = source.fSurveyX0[j][k][l];
180 0 : fSurveyY0[j][k][l] = source.fSurveyY0[j][k][l];
181 0 : fSurveyZ0[j][k][l] = source.fSurveyZ0[j][k][l];
182 : }
183 0 : for (int i=0; i<1000; i++) {
184 0 : fIbuffer[i] = 0;
185 0 : fDbuffer[i] = 0.0;
186 : }
187 :
188 0 : }
189 :
190 : //_____________________________________________________________________________
191 : AliTRDalignment& AliTRDalignment::operator=(const AliTRDalignment &source)
192 : {
193 : //
194 : // assignment operator
195 : //
196 :
197 0 : if (this != &source) {
198 0 : for (int i = 0; i < 18; i++) SetSm(i,source.fSm[i]);
199 0 : for (int i = 0; i < 540; i++) SetCh(i,source.fCh[i]);
200 0 : for (int i=0; i<18; i++) for (int j=0; j<2; j++) for (int k=0; k<2; k++) for (int l=0; l<2; l++) {
201 0 : fSurveyX[i][j][k][l] = source.fSurveyX[i][j][k][l];
202 0 : fSurveyY[i][j][k][l] = source.fSurveyY[i][j][k][l];
203 0 : fSurveyZ[i][j][k][l] = source.fSurveyZ[i][j][k][l];
204 0 : fSurveyEX[i][j][k][l] = source.fSurveyEX[i][j][k][l];
205 0 : fSurveyEY[i][j][k][l] = source.fSurveyEY[i][j][k][l];
206 0 : fSurveyEZ[i][j][k][l] = source.fSurveyEZ[i][j][k][l];
207 : }
208 0 : for (int j=0; j<2; j++) for (int k=0; k<2; k++) for (int l=0; l<2; l++) {
209 0 : fSurveyX0[j][k][l] = source.fSurveyX0[j][k][l];
210 0 : fSurveyY0[j][k][l] = source.fSurveyY0[j][k][l];
211 0 : fSurveyZ0[j][k][l] = source.fSurveyZ0[j][k][l];
212 : }
213 0 : fComment = source.fComment;
214 0 : }
215 :
216 0 : return *this;
217 :
218 : }
219 :
220 : //_____________________________________________________________________________
221 : AliTRDalignment& AliTRDalignment::operator*=(double fac)
222 : {
223 : //
224 : // multiplication operator
225 : //
226 :
227 0 : for (int i = 0; i < 18; i++) for (int j = 0; j < 6; j++) this->fSm[i][j] *= fac;
228 0 : for (int i = 0; i < 540; i++) for (int j = 0; j < 6; j++) this->fCh[i][j] *= fac;
229 :
230 0 : return *this;
231 :
232 : }
233 :
234 : //_____________________________________________________________________________
235 : AliTRDalignment& AliTRDalignment::operator+=(const AliTRDalignment &source)
236 : {
237 : //
238 : // addition operator
239 : //
240 :
241 0 : for (int i = 0; i < 18; i++) for (int j = 0; j < 6; j++) this->fSm[i][j] += source.fSm[i][j];
242 0 : for (int i = 0; i < 540; i++) for (int j = 0; j < 6; j++) this->fCh[i][j] += source.fCh[i][j];
243 :
244 0 : return *this;
245 :
246 : }
247 :
248 : //_____________________________________________________________________________
249 : AliTRDalignment& AliTRDalignment::operator-=(const AliTRDalignment &source)
250 : {
251 : //
252 : // subtraction operator
253 : //
254 :
255 0 : for (int i = 0; i < 18; i++) for (int j = 0; j < 6; j++) fSm[i][j] -= source.fSm[i][j];
256 0 : for (int i = 0; i < 540; i++) for (int j = 0; j < 6; j++) fCh[i][j] -= source.fCh[i][j];
257 :
258 0 : return *this;
259 :
260 : }
261 :
262 : //_____________________________________________________________________________
263 : Bool_t AliTRDalignment::operator==(const AliTRDalignment &source) const
264 : {
265 : //
266 : // comparison operator
267 : //
268 :
269 : Bool_t areEqual = 1;
270 :
271 0 : for (int i = 0; i < 18; i++) for (int j = 0; j < 6; j++) areEqual &= (fSm[i][j] == source.fSm[i][j]);
272 0 : for (int i = 0; i < 540; i++) for (int j = 0; j < 6; j++) areEqual &= (fCh[i][j] == source.fCh[i][j]);
273 :
274 0 : return areEqual;
275 :
276 : }
277 :
278 : //_____________________________________________________________________________
279 : void AliTRDalignment::SetSmZero()
280 : {
281 : //
282 : // reset to zero supermodule data
283 : //
284 :
285 0 : memset(&fSm[0][0],0,sizeof(fSm));
286 :
287 0 : }
288 :
289 : //_____________________________________________________________________________
290 : void AliTRDalignment::SetChZero()
291 : {
292 : //
293 : // reset to zero chamber data
294 : //
295 :
296 0 : memset(&fCh[0][0],0,sizeof(fCh));
297 :
298 0 : }
299 :
300 : //_____________________________________________________________________________
301 : void AliTRDalignment::SetSmRandom(double a[6])
302 : {
303 : //
304 : // generate random gaussian supermodule data with sigmas a
305 : //
306 :
307 0 : double x[6];
308 0 : double xmax[6]={999, 0.6, 999, 999, 999, 999};
309 :
310 0 : for (int i = 0; i < 18; i++) {
311 0 : for (int j = 0; j < 6; j++) {
312 0 : do {x[j] = fRan.Gaus(0,a[j]);} while (TMath::Abs(x[j]) > xmax[j]);
313 : }
314 0 : SetSm(i,x);
315 : //PrintSm(i);
316 : }
317 :
318 0 : }
319 :
320 : //_____________________________________________________________________________
321 : void AliTRDalignment::SetChRandom(double a[6])
322 : {
323 : //
324 : // generate random gaussian chamber data with sigmas a
325 : //
326 :
327 0 : double x[6];
328 :
329 0 : for (int i = 0; i < 540; i++) {
330 0 : fRan.Rannor(x[0],x[1]);
331 0 : fRan.Rannor(x[2],x[3]);
332 0 : fRan.Rannor(x[4],x[5]);
333 0 : for (int j = 0; j < 6; j++) x[j] *= a[j];
334 0 : SetCh(i,x);
335 : //PrintCh(i);
336 : }
337 :
338 0 : }
339 :
340 : //_____________________________________________________________________________
341 : void AliTRDalignment::SetSmFull()
342 : {
343 : //
344 : // generate random gaussian supermodule data similar to the misalignment
345 : // expected from the mechanical precision
346 : //
347 :
348 0 : double a[6];
349 :
350 0 : a[0] = 0.3; // phi
351 0 : a[1] = 0.3; // z
352 0 : a[2] = 0.3; // r
353 0 : a[3] = 0.4/1000.0 / TMath::Pi()*180.0; // phi
354 0 : a[4] = 2.0/1000.0 / TMath::Pi()*180.0; // z
355 0 : a[5] = 0.4/1000.0 / TMath::Pi()*180.0; // r
356 :
357 0 : SetSmRandom(a);
358 :
359 0 : }
360 :
361 : //_____________________________________________________________________________
362 : void AliTRDalignment::SetChFull()
363 : {
364 : //
365 : // generate random gaussian chamber data similar to the misalignment
366 : // expected from the mechanical precision
367 : //
368 :
369 0 : double a[6];
370 :
371 0 : a[0] = 0.1; // phi
372 0 : a[1] = 0.1; // z
373 0 : a[2] = 0.1; // r
374 0 : a[3] = 1.0/1000.0 / TMath::Pi()*180.0; // phi
375 0 : a[4] = 1.0/1000.0 / TMath::Pi()*180.0; // z
376 0 : a[5] = 0.7/1000.0 / TMath::Pi()*180.0; // r
377 :
378 0 : SetChRandom(a);
379 :
380 0 : }
381 :
382 : //_____________________________________________________________________________
383 : void AliTRDalignment::SetSmResidual()
384 : {
385 : //
386 : // generate random gaussian supermodule data similar to the misalignment
387 : // remaining after full calibration
388 : // I assume that it will be negligible
389 : //
390 :
391 0 : SetSmZero();
392 :
393 0 : }
394 :
395 : //_____________________________________________________________________________
396 : void AliTRDalignment::SetChResidual()
397 : {
398 : //
399 : // generate random gaussian chamber data similar to the misalignment
400 : // remaining after full calibration
401 : //
402 :
403 0 : double a[6];
404 :
405 0 : a[0] = 0.002; // phi
406 0 : a[1] = 0.003; // z
407 0 : a[2] = 0.007; // r
408 0 : a[3] = 0.3/1000.0 / TMath::Pi()*180.0; // phi
409 0 : a[4] = 0.3/1000.0 / TMath::Pi()*180.0; // z
410 0 : a[5] = 0.1/1000.0 / TMath::Pi()*180.0; // r
411 :
412 0 : SetChRandom(a);
413 :
414 0 : }
415 :
416 : //_____________________________________________________________________________
417 : void AliTRDalignment::PrintSm(int i, FILE * const fp) const
418 : {
419 : //
420 : // print the supermodule data
421 : //
422 :
423 0 : fprintf(fp,"%4d %11.4f %11.4f %11.4f %11.5f %11.5f %11.5f %6d %s\n"
424 0 : ,i,fSm[i][0],fSm[i][1],fSm[i][2],fSm[i][3],fSm[i][4],fSm[i][5]
425 0 : ,0,GetSmName(i));
426 :
427 0 : }
428 :
429 : //_____________________________________________________________________________
430 : void AliTRDalignment::PrintCh(int i, FILE * const fp) const
431 : {
432 : //
433 : // print the chamber data
434 : //
435 :
436 0 : fprintf(fp,"%4d %11.4f %11.4f %11.4f %11.5f %11.5f %11.5f %6d %s\n"
437 0 : ,i,fCh[i][0],fCh[i][1],fCh[i][2],fCh[i][3],fCh[i][4],fCh[i][5]
438 0 : ,GetVoi(i),GetChName(i));
439 :
440 0 : }
441 :
442 : //_____________________________________________________________________________
443 : void AliTRDalignment::ReadAscii(const char * const filename)
444 : {
445 : //
446 : // read the alignment data from ascii file
447 : //
448 :
449 0 : double x[6]; // alignment data
450 0 : int volid; // volume id
451 0 : std::string syna; // symbolic name
452 0 : int j; // dummy index
453 :
454 0 : fstream fi(filename,fstream::in);
455 0 : if (!fi) {
456 0 : AliError(Form("cannot open input file %s",filename));
457 0 : return;
458 : }
459 :
460 : // supermodules
461 :
462 0 : for (int i = 0; i < 18; i++) {
463 0 : fi>>j>>x[0]>>x[1]>>x[2]>>x[3]>>x[4]>>x[5]>>volid>>syna;
464 0 : if (j != i) AliError(Form("sm %d expected, %d found",i,j));
465 0 : if (volid != 0) AliError(Form("sm %d volume id %d expected, %d found",i,0,volid));
466 0 : std::string symnam = GetSmName(i);
467 0 : if (syna != symnam) AliError(Form("sm %d name %s expected, %s found",i,symnam.data(),syna.data()));
468 0 : SetSm(i,x);
469 0 : }
470 :
471 : // chambers
472 :
473 0 : for (int i = 0; i < 540; i++) {
474 0 : fi>>j>>x[0]>>x[1]>>x[2]>>x[3]>>x[4]>>x[5]>>volid>>syna;
475 0 : if (j != i) AliError(Form("ch %d expected, %d found",i,j));
476 0 : if (volid != GetVoi(i)) AliError(Form("ch %d volume id %d expected, %d found",i,GetVoi(i),volid));
477 0 : std::string symnam = GetChName(i);
478 0 : if (syna != symnam) AliError(Form("ch %d name %s expected, %s found",i,symnam.data(),syna.data()));
479 0 : SetCh(i,x);
480 0 : }
481 :
482 0 : fi.close();
483 :
484 0 : }
485 :
486 : //_____________________________________________________________________________
487 : void AliTRDalignment::ReadCurrentGeo()
488 : {
489 : //
490 : // use currently loaded geometry to determine misalignment by comparing
491 : // original and misaligned matrix of the last node
492 : // Now, original, does not mean "ideal". It is the matrix before the alignment.
493 : // So, if alignment was applied more than once, the numbers extracted will
494 : // represent just the last alignment. -- check this!
495 : //
496 :
497 : TGeoPNEntry *pne;
498 0 : TGeoHMatrix *ideSm[18]; // ideal
499 0 : TGeoHMatrix *misSm[18]; // misaligned
500 0 : for (int i = 0; i < 18; i++) if ((pne = gGeoManager->GetAlignableEntry(GetSmName(i)))) {
501 :
502 : // read misaligned and original matrices
503 :
504 0 : TGeoPhysicalNode *node = pne->GetPhysicalNode();
505 0 : if (!node) AliError(Form("physical node entry %s has no physical node",GetSmName(i)));
506 0 : if (!node) continue;
507 0 : misSm[i] = new TGeoHMatrix(*node->GetNode(node->GetLevel())->GetMatrix());
508 0 : ideSm[i] = new TGeoHMatrix(*node->GetOriginalMatrix());
509 :
510 : // calculate the local misalignment matrices as inverse misaligned times ideal
511 :
512 0 : TGeoHMatrix mat(ideSm[i]->Inverse());
513 0 : mat.Multiply(misSm[i]);
514 0 : double *tra = mat.GetTranslation();
515 0 : double *rot = mat.GetRotationMatrix();
516 0 : double pars[6];
517 0 : pars[0] = tra[0];
518 0 : pars[1] = tra[1];
519 0 : pars[2] = tra[2];
520 0 : if (TMath::Abs(rot[0])<1e-7 || TMath::Abs(rot[8])<1e-7) AliError("Failed to extract roll-pitch-yall angles!");
521 0 : double raddeg = TMath::RadToDeg();
522 0 : pars[3] = raddeg * TMath::ATan2(-rot[5],rot[8]);
523 0 : pars[4] = raddeg * TMath::ASin(rot[2]);
524 0 : pars[5] = raddeg * TMath::ATan2(-rot[1],rot[0]);
525 0 : SetSm(i,pars);
526 :
527 : // cleanup
528 :
529 0 : delete ideSm[i];
530 0 : delete misSm[i];
531 0 : }
532 :
533 0 : TGeoHMatrix *ideCh[540]; // ideal
534 0 : TGeoHMatrix *misCh[540]; // misaligned
535 0 : for (int i = 0; i < 540; i++) if ((pne = gGeoManager->GetAlignableEntry(GetChName(i)))) {
536 :
537 : // read misaligned and original matrices
538 :
539 0 : TGeoPhysicalNode *node = pne->GetPhysicalNode();
540 0 : if (!node) AliError(Form("physical node entry %s has no physical node",GetChName(i)));
541 0 : if (!node) continue;
542 0 : misCh[i] = new TGeoHMatrix(*node->GetNode(node->GetLevel())->GetMatrix());
543 0 : ideCh[i] = new TGeoHMatrix(*node->GetOriginalMatrix());
544 :
545 : // calculate the local misalignment matrices as inverse misaligned times ideal
546 :
547 0 : TGeoHMatrix mat(ideCh[i]->Inverse());
548 0 : mat.Multiply(misCh[i]);
549 0 : double *tra = mat.GetTranslation();
550 0 : double *rot = mat.GetRotationMatrix();
551 0 : double pars[6];
552 0 : pars[0] = tra[0];
553 0 : pars[1] = tra[1];
554 0 : pars[2] = tra[2];
555 0 : if(TMath::Abs(rot[0])<1e-7 || TMath::Abs(rot[8])<1e-7) {
556 0 : AliError("Failed to extract roll-pitch-yall angles!");
557 0 : return;
558 : }
559 0 : double raddeg = TMath::RadToDeg();
560 0 : pars[3] = raddeg * TMath::ATan2(-rot[5],rot[8]);
561 0 : pars[4] = raddeg * TMath::ASin(rot[2]);
562 0 : pars[5] = raddeg * TMath::ATan2(-rot[1],rot[0]);
563 0 : SetCh(i,pars);
564 :
565 : // cleanup
566 0 : delete ideCh[i];
567 0 : delete misCh[i];
568 0 : }
569 :
570 0 : return;
571 :
572 0 : }
573 :
574 : //_____________________________________________________________________________
575 : void AliTRDalignment::ReadRoot(const char * const filename)
576 : {
577 : //
578 : // read the alignment data from root file
579 : //
580 :
581 0 : TFile fi(filename,"READ");
582 :
583 0 : if (fi.IsOpen()) {
584 0 : TClonesArray *ar = (TClonesArray*) fi.Get("TRDAlignObjs");
585 0 : ArToNumbers(ar);
586 0 : fi.Close();
587 0 : }
588 0 : else AliError(Form("cannot open input file %s",filename));
589 :
590 : return;
591 :
592 0 : }
593 :
594 : //_____________________________________________________________________________
595 : void AliTRDalignment::ReadDB(const char * const filename)
596 : {
597 : //
598 : // read the alignment data from database file
599 : //
600 :
601 0 : TFile fi(filename,"READ");
602 :
603 0 : if (fi.IsOpen()) {
604 0 : AliCDBEntry *e = (AliCDBEntry *) fi.Get("AliCDBEntry");
605 0 : e->PrintMetaData();
606 0 : fComment.SetString(e->GetMetaData()->GetComment());
607 0 : TClonesArray *ar = (TClonesArray *) e->GetObject();
608 0 : ArToNumbers(ar);
609 0 : fi.Close();
610 0 : }
611 0 : else AliError(Form("cannot open input file %s",filename));
612 :
613 : return;
614 :
615 0 : }
616 :
617 : //_____________________________________________________________________________
618 : void AliTRDalignment::ReadDB(const char * const db, const char * const path,
619 : int run, int version, int subversion)
620 : {
621 : //
622 : // read the alignment data from database
623 : //
624 :
625 0 : AliCDBManager *cdb = AliCDBManager::Instance();
626 0 : AliCDBStorage *storLoc = cdb->GetStorage(db);
627 0 : AliCDBEntry *e = storLoc->Get(path,run,version,subversion);
628 0 : if (e) {
629 0 : e->PrintMetaData();
630 0 : fComment.SetString(e->GetMetaData()->GetComment());
631 0 : TClonesArray *ar = (TClonesArray *) e->GetObject();
632 0 : ArToNumbers(ar);
633 0 : }
634 0 : }
635 :
636 : //_____________________________________________________________________________
637 : Bool_t AliTRDalignment::DecodeSurveyPointName(TString pna, Int_t &sm, Int_t &iz,
638 : Int_t &ir, Int_t &iphi) {
639 : // decode the survey point name and extract the sm, z, r and phi indices
640 :
641 0 : if (pna(0,6)!="TRD_sm") {
642 0 : AliError(Form("unexpected point name: %s",pna.Data()));
643 0 : return kFALSE;
644 : }
645 0 : sm = atoi(pna(6,2).Data()); // supermodule number
646 0 : iz = -1;
647 0 : if (pna(8) == 'a') iz=0; // anticlockwise, positive z
648 0 : if (pna(8) == 'c') iz=1; // clockwise, negative z
649 0 : ir = -1;
650 0 : if (pna(9) == 'l') ir=0; // low radius
651 0 : if (pna(9) == 'h') ir=1; // high radius
652 0 : iphi = -1;
653 0 : if (pna(10) == '0') iphi = 0; // low phi within supermodule
654 0 : if (pna(10) == '1') iphi = 1; // high phi within supermodule
655 0 : if (sm>=0 && sm<18 && iz>=0 && iz<2 && ir>=0 && ir<2 && iphi>=0 && iphi<2) return kTRUE;
656 0 : AliError(Form("cannot decode point name: %s",pna.Data()));
657 0 : return kFALSE;
658 0 : }
659 :
660 : //_____________________________________________________________________________
661 : void AliTRDalignment::ReadSurveyReport(const char * const filename)
662 : {
663 : //
664 : // Read survey report and store the numbers in fSurveyX, fSurveyY, fSurveyZ,
665 : // and fSurveyE. Store the survey info in the fComment.
666 : // Each supermodule has 8 survey points. The point names look like
667 : // TRD_sm08ah0 and have the following meaning.
668 : //
669 : // sm00..17 mean supermodule 0 through 17, following the phi.
670 : // Supermodule 00 is between phi=0 and phi=20 degrees.
671 : //
672 : // a or c denotes the anticlockwise and clockwise end of the supermodule
673 : // in z. Clockwise end is where z is negative and where the muon arm sits.
674 : //
675 : // l or h denote low radius and high radius holes
676 : //
677 : // 0 or 1 denote the hole at smaller and at larger phi, respectively.
678 : //
679 :
680 : // read the survey file
681 :
682 0 : fstream in(filename,fstream::in);
683 0 : if (!in) {
684 0 : AliError(Form("cannot open input file %s",filename));
685 0 : return;
686 : }
687 :
688 : // loop through the lines of the file until the beginning of data
689 :
690 0 : TString title,date,subdetector,url,version,observations,system,units;
691 : while (1) {
692 0 : char pee=in.peek();
693 0 : if (pee==EOF) break;
694 0 : TString line;
695 0 : line.ReadLine(in);
696 0 : if (line.Contains("Title:")) title.ReadLine(in);
697 0 : if (line.Contains("Date:")) date.ReadLine(in);
698 0 : if (line.Contains("Subdetector:")) subdetector.ReadLine(in);
699 0 : if (line.Contains("URL:")) url.ReadLine(in);
700 0 : if (line.Contains("Version:")) version.ReadLine(in);
701 0 : if (line.Contains("Observations:")) observations.ReadLine(in);
702 0 : if (line.Contains("System:")) system.ReadLine(in);
703 0 : if (line.Contains("Units:")) units.ReadLine(in);
704 0 : if (line.Contains("Data:")) break;
705 0 : }
706 :
707 : // check what we found so far (watch out, they have \r at the end)
708 :
709 0 : std::cout<<"title .........."<<title<<std::endl;
710 0 : std::cout<<"date ..........."<<date<<std::endl;
711 0 : std::cout<<"subdetector ...."<<subdetector<<std::endl;
712 0 : std::cout<<"url ............"<<url<<std::endl;
713 0 : std::cout<<"version ........"<<version<<std::endl;
714 0 : std::cout<<"observations ..."<<observations<<std::endl;
715 0 : std::cout<<"system ........."<<system<<std::endl;
716 0 : std::cout<<"units .........."<<units<<std::endl;
717 :
718 0 : if (!subdetector.Contains("TRD")) {
719 0 : AliWarning(Form("Not a TRD survey file, subdetector = %s",subdetector.Data()));
720 0 : return;
721 : }
722 : double tocm = 0; // we want to have it in cm
723 0 : if (units.Contains("mm")) tocm = 0.1;
724 0 : else if (units.Contains("cm")) tocm = 1.0;
725 0 : else if (units.Contains("m")) tocm = 100.0;
726 0 : else if (units.Contains("pc")) tocm = 3.24078e-15;
727 : else {
728 0 : AliError(Form("unexpected units: %s",units.Data()));
729 0 : return;
730 : }
731 0 : if (!system.Contains("ALICEPH")) {
732 0 : AliError(Form("wrong system: %s, should be ALICEPH",system.Data()));
733 0 : return;
734 : }
735 :
736 : // scan the rest of the file which should contain list of surveyed points
737 : // for every point, decode the point name and store the numbers in the right
738 : // place in the arrays fSurveyX etc.
739 :
740 : while (1) {
741 0 : TString pna; // point name
742 0 : char type, target;
743 0 : double x,y,z,precision;
744 :
745 0 : in >> pna >> x >> y >> z >> type >> target >> precision;
746 0 : if (in.fail()) break;
747 0 : Int_t i,j,k,l;
748 0 : if (DecodeSurveyPointName(pna,i,j,k,l)) {
749 0 : fSurveyX[i][j][k][l] = tocm*x;
750 0 : fSurveyY[i][j][k][l] = tocm*y;
751 0 : fSurveyZ[i][j][k][l] = tocm*z;
752 0 : fSurveyEX[i][j][k][l] = precision/10; // "precision" is supposed to be in mm
753 0 : fSurveyEY[i][j][k][l] = precision/10; // "precision" is supposed to be in mm
754 0 : fSurveyEZ[i][j][k][l] = precision/10; // "precision" is supposed to be in mm
755 : // if, at some point, separate precision numbers for x,y,z show up in the
756 : // survey reports the function will fail here
757 0 : printf("decoded %s %02d %d %d %d %8.2f %8.2f %8.2f %6.2f %6.2f %6.2f\n",
758 0 : pna.Data(), i, j, k, l,
759 0 : fSurveyX[i][j][k][l], fSurveyY[i][j][k][l], fSurveyZ[i][j][k][l],
760 0 : fSurveyEX[i][j][k][l], fSurveyEY[i][j][k][l], fSurveyEZ[i][j][k][l]);
761 0 : } else AliError(Form("cannot decode point name: %s",pna.Data()));
762 0 : }
763 0 : in.close();
764 0 : TString info = "Survey "+title+" "+date+" "+url+" "+version+" "+observations;
765 0 : info.ReplaceAll("\r","");
766 0 : fComment.SetString(info.Data());
767 :
768 0 : }
769 :
770 : //_____________________________________________________________________________
771 : void AliTRDalignment::ReadSurveyReport(const AliSurveyObj * const so)
772 : {
773 : //
774 : // Read survey report and store the numbers in fSurveyX, fSurveyY, fSurveyZ,
775 : // and fSurveyE. Store the survey info in the fComment.
776 : // Each supermodule has 8 survey points. The point names look like
777 : // TRD_sm08ah0 and have the following meaning.
778 : //
779 : // sm00..17 mean supermodule 0 through 17, following the phi.
780 : // Supermodule 00 is between phi=0 and phi=20 degrees.
781 : //
782 : // a or c denotes the anticlockwise and clockwise end of the supermodule
783 : // in z. Clockwise end is where z is negative and where the muon arm sits.
784 : //
785 : // l or h denote low radius and high radius holes
786 : //
787 : // 0 or 1 denote the hole at smaller and at larger phi, respectively.
788 : //
789 :
790 : // read and process the data from the survey object
791 :
792 0 : Int_t size = so->GetEntries();
793 0 : printf("-> %d\n", size);
794 :
795 0 : TString title = so->GetReportTitle();
796 0 : TString date = so->GetReportDate();
797 0 : TString subdetector = so->GetDetector();
798 0 : TString url = so->GetURL();
799 0 : TString report = so->GetReportNumber();
800 0 : TString version = so->GetReportVersion();
801 0 : TString observations = so->GetObservations();
802 0 : TString system = so->GetCoordSys();
803 0 : TString units = so->GetUnits();
804 :
805 : // check what we found so far (watch out, they have \r at the end)
806 :
807 0 : std::cout<<"title .........."<<title<<std::endl;
808 0 : std::cout<<"date ..........."<<date<<std::endl;
809 0 : std::cout<<"subdetector ...."<<subdetector<<std::endl;
810 0 : std::cout<<"url ............"<<url<<std::endl;
811 0 : std::cout<<"version ........"<<version<<std::endl;
812 0 : std::cout<<"observations ..."<<observations<<std::endl;
813 0 : std::cout<<"system ........."<<system<<std::endl;
814 0 : std::cout<<"units .........."<<units<<std::endl;
815 :
816 0 : if (!subdetector.Contains("TRD")) {
817 0 : AliWarning(Form("Not a TRD survey file, subdetector = %s",subdetector.Data()));
818 0 : return;
819 : }
820 : double tocm = 0; // we want to have it in cm
821 0 : if (units.Contains("mm")) tocm = 0.1;
822 0 : else if (units.Contains("cm")) tocm = 1.0;
823 0 : else if (units.Contains("m")) tocm = 100.0;
824 0 : else if (units.Contains("pc")) tocm = 3.24078e-15;
825 : else {
826 0 : AliError(Form("unexpected units: %s",units.Data()));
827 0 : return;
828 : }
829 0 : if (!system.Contains("ALICEPH")) {
830 0 : AliError(Form("wrong system: %s, should be ALICEPH",system.Data()));
831 0 : return;
832 : }
833 :
834 : // for every survey point, decode the point name and store the numbers in
835 : // the right place in the arrays fSurveyX etc.
836 :
837 0 : TObjArray *points = so->GetData();
838 0 : for (int ip = 0; ip<points->GetEntries(); ++ip) {
839 0 : AliSurveyPoint *po = (AliSurveyPoint *) points->At(ip);
840 0 : TString pna = po->GetPointName();
841 0 : Int_t i,j,k,l;
842 0 : if (DecodeSurveyPointName(pna,i,j,k,l)) {
843 0 : fSurveyX[i][j][k][l] = tocm*po->GetX();
844 0 : fSurveyY[i][j][k][l] = tocm*po->GetY();
845 0 : fSurveyZ[i][j][k][l] = tocm*po->GetZ();
846 0 : fSurveyEX[i][j][k][l] = po->GetPrecisionX()/10; // "precision" is supposed to be in mm
847 0 : fSurveyEY[i][j][k][l] = po->GetPrecisionY()/10;
848 0 : fSurveyEZ[i][j][k][l] = po->GetPrecisionZ()/10;
849 0 : printf("decoded %s %02d %d %d %d %8.2f %8.2f %8.2f %6.2f %6.2f %6.2f\n",
850 0 : pna.Data(), i, j, k, l,
851 0 : fSurveyX[i][j][k][l], fSurveyY[i][j][k][l], fSurveyZ[i][j][k][l],
852 0 : fSurveyEX[i][j][k][l], fSurveyEY[i][j][k][l], fSurveyEZ[i][j][k][l]);
853 0 : } else AliError(Form("cannot decode point name: %s",pna.Data()));
854 0 : }
855 :
856 0 : TString info = "Survey "+title+" "+date+" "+url+" "+report+" "+version+" "+observations;
857 0 : info.ReplaceAll("\r","");
858 0 : fComment.SetString(info.Data());
859 0 : }
860 :
861 : //_____________________________________________________________________________
862 : double AliTRDalignment::SurveyChi2(int i, const double * const a) {
863 :
864 : //
865 : // Compare the survey results to the ideal positions of the survey marks
866 : // in the local frame of supermodule. When transforming, use the alignment
867 : // parameters a[6]. Return chi-squared.
868 : //
869 :
870 0 : if (!IsGeoLoaded()) return 0;
871 0 : printf("Survey of supermodule %d\n",i);
872 0 : AliAlignObjParams al(GetSmName(i),0,a[0],a[1],a[2],a[3],a[4],a[5],0);
873 :
874 0 : TGeoPNEntry *pne = gGeoManager->GetAlignableEntry(GetSmName(i));
875 0 : if (!pne) AliError(Form("no such physical node entry: %s",GetSmName(i)));
876 0 : TGeoPhysicalNode *node = pne->GetPhysicalNode();
877 0 : if (!node) {
878 0 : AliWarning(Form("physical node entry %s has no physical node; making a new one",GetSmName(i)));
879 0 : node = gGeoManager->MakeAlignablePN(pne);
880 0 : }
881 :
882 : // al.ApplyToGeometry();
883 : // node = pne->GetPhysicalNode(); // changed in the meantime
884 : // TGeoHMatrix *ma = node->GetMatrix();
885 :
886 : // a less destructive method (it does not modify geometry), gives the same result:
887 :
888 0 : TGeoHMatrix *ma = new TGeoHMatrix();
889 0 : al.GetLocalMatrix(*ma);
890 0 : ma->MultiplyLeft(node->GetMatrix()); // global trafo, modified by a[]
891 :
892 : double chi2=0;
893 0 : printf(" sm z r phi x (lab phi) y (lab z) z (lab r) all in cm\n");
894 0 : for (int j=0; j<2; j++) for (int k=0; k<2; k++) for (int l=0; l<2; l++) {
895 0 : if (fSurveyEX[i][j][k][l] == 0.0
896 0 : && fSurveyEY[i][j][k][l] == 0.0
897 0 : && fSurveyEZ[i][j][k][l] == 0.0) continue; // no data for this survey point
898 0 : double master[3] = {fSurveyX[i][j][k][l],fSurveyY[i][j][k][l],fSurveyZ[i][j][k][l]};
899 0 : double local[3];
900 0 : ma->MasterToLocal(master,local);
901 0 : double dx = local[0]-fSurveyX0[j][k][l];
902 0 : double dy = local[1]-fSurveyY0[j][k][l];
903 0 : double dz = local[2]-fSurveyZ0[j][k][l];
904 0 : chi2 += dx*dx/fSurveyEX[i][j][k][l]/fSurveyEX[i][j][k][l];
905 0 : chi2 += dy*dy/fSurveyEY[i][j][k][l]/fSurveyEY[i][j][k][l];
906 0 : chi2 += dz*dz/fSurveyEZ[i][j][k][l]/fSurveyEZ[i][j][k][l];
907 0 : printf("local survey %3d %3d %3d %3d %12.3f %12.3f %12.3f\n",i,j,k,l,local[0],local[1],local[2]);
908 0 : printf("local ideal %12.3f %12.3f %12.3f\n",fSurveyX0[j][k][l],
909 0 : fSurveyY0[j][k][l],fSurveyZ0[j][k][l]);
910 0 : printf("difference %12.3f %12.3f %12.3f\n",dx,dy,dz);
911 0 : }
912 0 : printf("chi2 = %.2f\n",chi2);
913 : return chi2;
914 0 : }
915 :
916 : //_____________________________________________________________________________
917 : void trdAlignmentFcn(int &npar, double *g, double &f, double *par, int iflag) {
918 :
919 : //
920 : // Standard function as needed by Minuit-like minimization procedures.
921 : // For the set of parameters par calculates and returns chi-squared.
922 : //
923 :
924 : // smuggle a C++ object into a C function
925 0 : AliTRDalignment *alignment = (AliTRDalignment*) gMinuit->GetObjectFit();
926 :
927 0 : f = alignment->SurveyChi2(par);
928 : if (iflag==3) {}
929 0 : if (npar) {}
930 : if (g) {} // no warnings about unused stuff...
931 :
932 0 : }
933 :
934 : //_____________________________________________________________________________
935 : void AliTRDalignment::SurveyToAlignment(int i, const char * const flag) {
936 :
937 : //
938 : // Find the supermodule alignment parameters needed to make the survey
939 : // results coincide with the ideal positions of the survey marks.
940 : // The string flag should look like "101000"; the six characters corresponds
941 : // to the six alignment parameters and 0/1 mean that the parameter should
942 : // be fixed/released in the fit.
943 :
944 0 : if (strlen(flag)!=6) {
945 0 : AliError(Form("unexpected flag: %s",flag));
946 0 : return;
947 : }
948 :
949 0 : printf("Finding alignment matrix for supermodule %d\n",i);
950 0 : fIbuffer[0] = i; // store the sm number in the buffer so minuit can see it
951 :
952 0 : TFitter fitter(100);
953 0 : gMinuit->SetObjectFit(this);
954 0 : fitter.SetFCN(trdAlignmentFcn);
955 0 : fitter.SetParameter(0,"dx",0,0.5,0,0);
956 0 : fitter.SetParameter(1,"dy",0,0.5,0,0);
957 0 : fitter.SetParameter(2,"dz",0,0.5,0,0);
958 0 : fitter.SetParameter(3,"rx",0,0.1,0,0);
959 0 : fitter.SetParameter(4,"ry",0,0.1,0,0);
960 0 : fitter.SetParameter(5,"rz",0,0.1,0,0);
961 :
962 0 : for (int j=0; j<6; j++) if (flag[j]=='0') fitter.FixParameter(j);
963 :
964 0 : double arglist[100];
965 0 : arglist[0] = 2;
966 0 : fitter.ExecuteCommand("SET PRINT", arglist, 1);
967 0 : fitter.ExecuteCommand("SET ERR", arglist, 1);
968 0 : arglist[0]=50;
969 : //fitter.ExecuteCommand("SIMPLEX", arglist, 1);
970 0 : fitter.ExecuteCommand("MINIMIZE", arglist, 1);
971 0 : fitter.ExecuteCommand("CALL 3", arglist,0);
972 0 : double a[6];
973 0 : for (int j=0; j<6; j++) a[j] = fitter.GetParameter(j);
974 0 : SetSm(i,a);
975 0 : for (int j=0; j<6; j++) printf("%10.3f ",fitter.GetParameter(j));
976 0 : printf("\n");
977 0 : for (int j=0; j<6; j++) printf("%10.3f ",fitter.GetParError(j));
978 0 : printf("\n");
979 :
980 0 : }
981 :
982 : //_____________________________________________________________________________
983 : void AliTRDalignment::ReadAny(const char * const filename)
984 : {
985 : //
986 : // read the alignment data from any kind of file
987 : //
988 :
989 0 : TString fist(filename);
990 0 : if (fist.EndsWith(".txt")) ReadAscii(filename);
991 0 : if (fist.EndsWith(".dat")) ReadAscii(filename);
992 0 : if (fist.EndsWith(".root")) {
993 0 : if (fist.Contains("Run")) ReadDB(filename);
994 0 : else ReadRoot(filename);
995 : }
996 :
997 0 : }
998 :
999 : //_____________________________________________________________________________
1000 : void AliTRDalignment::WriteAscii(const char * const filename) const
1001 : {
1002 : //
1003 : // store the alignment data on ascii file
1004 : //
1005 :
1006 0 : FILE *fp = fopen(filename, "w");
1007 0 : if (!fp) {
1008 0 : AliError(Form("cannot open output file %s",filename));
1009 0 : return;
1010 : }
1011 :
1012 0 : PrintSm(fp);
1013 0 : PrintCh(fp);
1014 :
1015 0 : fclose(fp);
1016 :
1017 0 : }
1018 :
1019 : //_____________________________________________________________________________
1020 : void AliTRDalignment::WriteRoot(const char * const filename)
1021 : {
1022 : //
1023 : // store the alignment data on root file
1024 : //
1025 :
1026 0 : TClonesArray *ar = new TClonesArray("AliAlignObjParams",10000);
1027 0 : NumbersToAr(ar);
1028 0 : TFile fo(filename,"RECREATE");
1029 0 : if (fo.IsOpen()) {
1030 0 : fo.cd();
1031 0 : fo.WriteObject(ar,"TRDAlignObjs","kSingleKey");
1032 0 : fo.Close();
1033 : }
1034 0 : else AliError(Form("cannot open output file %s",filename));
1035 :
1036 0 : delete ar;
1037 :
1038 0 : }
1039 :
1040 : //_____________________________________________________________________________
1041 : void AliTRDalignment::WriteDB(const char * const filename, int run0, int run1, int ver, int subver)
1042 : {
1043 : //
1044 : // dumping on a DB-like file
1045 : //
1046 :
1047 0 : TClonesArray *ar = new TClonesArray("AliAlignObjParams",10000);
1048 0 : NumbersToAr(ar);
1049 : const Char_t *path = "TRD/Align/Data";
1050 0 : AliCDBId id(path,run0,run1);
1051 0 : AliCDBMetaData *md = new AliCDBMetaData();
1052 0 : md->SetResponsible("Dariusz Miskowiec");
1053 0 : md->SetComment(fComment.GetString().Data());
1054 0 : AliCDBEntry *e = new AliCDBEntry(ar, id, md);
1055 0 : e->SetVersion(ver);
1056 0 : e->SetSubVersion(subver);
1057 0 : TFile fi(filename,"RECREATE");
1058 0 : if (fi.IsOpen()) {
1059 0 : e->Write();
1060 0 : fi.Close();
1061 : }
1062 0 : else AliError(Form("cannot open input file %s",filename));
1063 :
1064 0 : delete e;
1065 0 : delete md;
1066 0 : delete ar;
1067 :
1068 : return;
1069 :
1070 0 : }
1071 :
1072 : //_____________________________________________________________________________
1073 : void AliTRDalignment::WriteDB(char * const db, const char * const path, int run0, int run1)
1074 : {
1075 : //
1076 : // store the alignment data in database
1077 : //
1078 :
1079 0 : TClonesArray *ar = new TClonesArray("AliAlignObjParams",10000);
1080 0 : NumbersToAr(ar);
1081 0 : AliCDBManager *cdb = AliCDBManager::Instance();
1082 0 : AliCDBStorage *storLoc = cdb->GetStorage(db);
1083 0 : AliCDBMetaData *md = new AliCDBMetaData();
1084 0 : md->SetResponsible("Dariusz Miskowiec");
1085 0 : md->SetComment(fComment.GetString().Data());
1086 0 : AliCDBId id(path,run0,run1);
1087 0 : storLoc->Put(ar,id,md);
1088 0 : md->Delete();
1089 0 : delete ar;
1090 :
1091 0 : }
1092 :
1093 : //_____________________________________________________________________________
1094 : void AliTRDalignment::WriteGeo(char *filename)
1095 : {
1096 : //
1097 : // apply misalignment to current geometry and store the
1098 : // resulting geometry on a root file
1099 : //
1100 :
1101 0 : TClonesArray *ar = new TClonesArray("AliAlignObjParams",10000);
1102 0 : NumbersToAr(ar);
1103 0 : delete ar;
1104 0 : gGeoManager->Export(filename);
1105 :
1106 0 : }
1107 :
1108 : //_____________________________________________________________________________
1109 : double AliTRDalignment::GetSmRMS(int xyz) const
1110 : {
1111 : //
1112 : // rms fSm[][xyz]
1113 : //
1114 :
1115 : double s1 = 0.0;
1116 : double s2 = 0.0;
1117 0 : for (int i = 0; i < 18; i++) {
1118 0 : s1 += fSm[i][xyz];
1119 0 : s2 += fSm[i][xyz]*fSm[i][xyz];
1120 : }
1121 0 : double rms2 = s2/18.0 - s1*s1/18.0/18.0;
1122 :
1123 0 : return rms2>0 ? sqrt(rms2) : 0.0;
1124 :
1125 : }
1126 :
1127 : //_____________________________________________________________________________
1128 : double AliTRDalignment::GetChRMS(int xyz) const
1129 : {
1130 : //
1131 : // rms fCh[][xyz]
1132 : //
1133 :
1134 : double s1 =0.0;
1135 : double s2 =0.0;
1136 0 : for (int i = 0; i < 540; i++) {
1137 0 : s1 += fCh[i][xyz];
1138 0 : s2 += fCh[i][xyz]*fCh[i][xyz];
1139 : }
1140 0 : double rms2 = s2/540.0 - s1*s1/540.0/540.0;
1141 :
1142 0 : return rms2>0 ? sqrt(rms2) : 0.0;
1143 :
1144 : }
1145 :
1146 : //_____________________________________________________________________________
1147 : void AliTRDalignment::PrintSmRMS() const
1148 : {
1149 : //
1150 : // dump rms of fSm
1151 : //
1152 :
1153 0 : printf(" %11.4f %11.4f %11.4f %11.5f %11.5f %11.5f supermodule rms\n"
1154 0 : ,GetSmRMS(0),GetSmRMS(1),GetSmRMS(2),GetSmRMS(3),GetSmRMS(4),GetSmRMS(5));
1155 :
1156 0 : }
1157 :
1158 : //_____________________________________________________________________________
1159 : void AliTRDalignment::PrintChRMS() const
1160 : {
1161 : //
1162 : // dump rms of fCh
1163 : //
1164 :
1165 0 : printf(" %11.4f %11.4f %11.4f %11.5f %11.5f %11.5f chamber rms\n"
1166 0 : ,GetChRMS(0),GetChRMS(1),GetChRMS(2),GetChRMS(3),GetChRMS(4),GetChRMS(5));
1167 :
1168 0 : }
1169 :
1170 : //_____________________________________________________________________________
1171 : void AliTRDalignment::ArToNumbers(TClonesArray * const ar)
1172 : {
1173 : //
1174 : // for each of the alignment objects in array ar extract the six local
1175 : // alignment parameters; recognize by name to which supermodule or chamber
1176 : // the alignment object pertains; set the respective fSm or fCh
1177 : //
1178 :
1179 0 : ar->Sort();
1180 0 : if (!IsGeoLoaded()) return;
1181 0 : for (int i = 0; i < ar->GetEntries(); i++) {
1182 0 : AliAlignObj *aao = (AliAlignObj *) ar->At(i);
1183 0 : aao->ApplyToGeometry();
1184 : }
1185 0 : SetZero();
1186 0 : ReadCurrentGeo();
1187 :
1188 0 : }
1189 :
1190 : //_____________________________________________________________________________
1191 : void AliTRDalignment::NumbersToAr(TClonesArray * const ar)
1192 : {
1193 : //
1194 : // build array of AliAlignObj objects based on fSm and fCh data
1195 : // at the same time, apply misalignment to the currently loaded geometry
1196 : // it is important to apply misalignment of supermodules before creating
1197 : // alignment objects for chambers
1198 : //
1199 :
1200 0 : if (!IsGeoLoaded()) return;
1201 : TClonesArray &alobj = *ar;
1202 : int nobj = 0;
1203 0 : for (int i = 0; i < 18; i++) {
1204 0 : new(alobj[nobj]) AliAlignObjParams(GetSmName(i)
1205 : ,0
1206 0 : ,fSm[i][0],fSm[i][1],fSm[i][2]
1207 0 : ,fSm[i][3],fSm[i][4],fSm[i][5]
1208 : ,0);
1209 0 : ((AliAlignObj *) alobj[nobj])->ApplyToGeometry();
1210 0 : nobj++;
1211 : }
1212 :
1213 0 : for (int i = 0; i < 540; i++) {
1214 0 : if (gGeoManager->GetAlignableEntry(GetChName(i))) {
1215 0 : new(alobj[nobj]) AliAlignObjParams(GetChName(i)
1216 0 : ,GetVoi(i)
1217 0 : ,fCh[i][0],fCh[i][1],fCh[i][2]
1218 0 : ,fCh[i][3],fCh[i][4],fCh[i][5]
1219 : ,0);
1220 0 : ((AliAlignObj *) alobj[nobj])->ApplyToGeometry();
1221 0 : nobj++;
1222 0 : }
1223 : }
1224 0 : AliInfo("current geometry modified");
1225 :
1226 0 : }
1227 :
1228 : //_____________________________________________________________________________
1229 : int AliTRDalignment::IsGeoLoaded()
1230 : {
1231 : //
1232 : // check whether a geometry is loaded
1233 : // issue a warning if geometry is not ideal
1234 : //
1235 :
1236 0 : if (gGeoManager) {
1237 0 : if (gGeoManager->GetListOfPhysicalNodes()->GetEntries()) AliWarning("current geometry is not ideal");
1238 0 : return 1;
1239 : } else {
1240 0 : AliError("first load geometry by calling TGeoManager::Import(filename)");
1241 0 : return 0;
1242 : }
1243 :
1244 0 : }
1245 :
1246 : //_____________________________________________________________________________
|
