]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/AliITSAlignMille2Module.cxx
Optimisation
[u/mrichter/AliRoot.git] / ITS / AliITSAlignMille2Module.cxx
CommitLineData
6526a72c 1/************************************************************************** \r
2 * Copyright(c) 2007-2009, ALICE Experiment at CERN, All rights reserved. * \r
3 * * \r
4 * Author: The ALICE Off-line Project. * \r
5 * Contributors are mentioned in the code where appropriate. * \r
6 * * \r
7 * Permission to use, copy, modify and distribute this software and its * \r
8 * documentation strictly for non-commercial purposes is hereby granted * \r
9 * without fee, provided that the above copyright notice appears in all * \r
10 * copies and that both the copyright notice and this permission notice * \r
11 * appear in the supporting documentation. The authors make no claims * \r
12 * about the suitability of this software for any purpose. It is * \r
13 * provided "as is" without express or implied warranty. * \r
14 **************************************************************************/ \r
15 \r
16/* $Id$ */ \r
17//----------------------------------------------------------------------------- \r
18/// \class AliITSAlignMille2Module\r
19/// Alignment class for the ALICE ITS detector \r
20/// \r
21/// This class is used by AliITSAlignMille to build custom supermodules \r
22/// made of ITS sensitive modules. These supermodules are then aligned\r
23/// \r
24/// Custom supermodules must have VolumeID > 14335\r
25///\r
26/// \author M. Lunardon \r
27//----------------------------------------------------------------------------- \r
28 \r
29#include <TGeoManager.h> \r
30#include <TGeoMatrix.h> \r
31 \r
32#include "AliITSAlignMille2Module.h" \r
33#include "AliITSgeomTGeo.h" \r
34#include "AliGeomManager.h" \r
35#include "AliAlignObjParams.h" \r
36#include "AliLog.h" \r
37#include "AliITSAlignMille2.h"\r
38 \r
39/// \cond CLASSIMP \r
40ClassImp(AliITSAlignMille2Module) \r
41/// \endcond \r
42\r
43#define CORHW_\r
44\r
45AliAlignObjParams AliITSAlignMille2Module::fgTempAlignObj;\r
8102b2c9 46const Float_t AliITSAlignMille2Module::fgkDummyConstraint = 1e-2;//1.E3;\r
6526a72c 47 \r
48//-------------------------------------------------------------\r
49AliITSAlignMille2Module::AliITSAlignMille2Module() : \r
50 TNamed(), \r
51 fNSensVol(0), \r
52 fIndex(-1), \r
53 fDetType(-1),\r
54 fVolumeID(0),\r
55 fNParTot(0),\r
56 fNParFree(0),\r
57 fParOffs(0),\r
58 fNProcPoints(0),\r
59 fParVals(0),\r
60 fParErrs(0),\r
61 fParCstr(0),\r
62 fSensVolIndex(0),\r
63 fSensVolVolumeID(0),\r
64 fMatrix(NULL),\r
65 fSensVolMatrix(NULL),\r
66 fSensVolModifMatrix(NULL),\r
67 fParent(NULL),\r
68 fChildren(0)\r
69{ \r
70 /// void constructor \r
71 fMatrix = new TGeoHMatrix; \r
72 fSensVolMatrix = new TGeoHMatrix; \r
73 fSensVolModifMatrix = new TGeoHMatrix; \r
74 fSensVolIndex.Set(1);\r
75 fSensVolVolumeID.Set(1);\r
76 fSigmaFactor[0]=fSigmaFactor[1]=fSigmaFactor[2]=1.0;\r
77} \r
78\r
79//-------------------------------------------------------------\r
6be22b3f 80AliITSAlignMille2Module::AliITSAlignMille2Module(Int_t index, UShort_t volid, const char* symname,\r
81 const TGeoHMatrix *m, Int_t nsv, const UShort_t *volidsv) : \r
6526a72c 82 TNamed(), \r
83 fNSensVol(0), \r
84 fIndex(-1), \r
85 fDetType(-1), \r
86 fVolumeID(0),\r
87 fNParTot(0),\r
88 fNParFree(0),\r
89 fParOffs(kMaxParGeom),\r
90 fNProcPoints(0),\r
91 fParVals(0),\r
92 fParErrs(0),\r
93 fParCstr(0),\r
94 fSensVolIndex(0),\r
95 fSensVolVolumeID(0), \r
96 fMatrix(NULL),\r
97 fSensVolMatrix(NULL),\r
98 fSensVolModifMatrix(NULL),\r
99 fParent(NULL),\r
100 fChildren(0)\r
101{ \r
102 /// void constructor \r
103 fMatrix = new TGeoHMatrix; \r
104 fSensVolMatrix = new TGeoHMatrix; \r
105 fSensVolModifMatrix = new TGeoHMatrix; \r
106 fSigmaFactor[0]=fSigmaFactor[1]=fSigmaFactor[2]=1.0;\r
107 for (int i=kMaxParGeom;i--;) fParOffs[i] = -1;\r
108 if (Set(index,volid,symname,m,nsv,volidsv)) {\r
109 AliInfo("Error in AliITSAlignMille2Module::Set() - initializing void supermodule...");\r
110 }\r
111 AssignDetType();\r
112} \r
113\r
114//-------------------------------------------------------------\r
115AliITSAlignMille2Module::AliITSAlignMille2Module(UShort_t volid) : \r
116 TNamed(), \r
117 fNSensVol(0), \r
118 fIndex(-1), \r
119 fDetType(-1),\r
120 fVolumeID(0),\r
121 fNParTot(0),\r
122 fNParFree(0),\r
123 fParOffs(kMaxParGeom),\r
124 fNProcPoints(0),\r
125 fParVals(0),\r
126 fParErrs(0),\r
127 fParCstr(0), \r
128 fSensVolIndex(0),\r
129 fSensVolVolumeID(0),\r
130 fMatrix(NULL),\r
131 fSensVolMatrix(NULL),\r
132 fSensVolModifMatrix(NULL),\r
133 fParent(NULL),\r
134 fChildren(0)\r
135{ \r
136 /// simple constructor building a supermodule from a single sensitive volume \r
137 fMatrix = new TGeoHMatrix; \r
138 fSensVolMatrix = new TGeoHMatrix; \r
139 fSensVolModifMatrix = new TGeoHMatrix; \r
140 // temporary align object, just use the rotation...\r
141 fSensVolIndex.Set(1);\r
142 fSensVolVolumeID.Set(1);\r
143 fSigmaFactor[0]=fSigmaFactor[1]=fSigmaFactor[2]=1.0;\r
144 for (int i=kMaxParGeom;i--;) fParOffs[i] = -1;\r
145 //\r
146 fIndex = GetIndexFromVolumeID(volid); \r
147 if (fIndex>=0 && gGeoManager) { // good sensitive module and geometry loaded\r
148 SetName(AliGeomManager::SymName(volid));\r
149 fVolumeID = volid;\r
150 AddSensitiveVolume(volid);\r
151 SetSensorsProvided(kTRUE);\r
152 if (SensVolMatrix(volid, fMatrix))\r
153 AliInfo("Matrix not defined");\r
154 }\r
155 else {\r
156 AliInfo("Wrong VolumeID or Geometry not loaded - initializing void supermodule...");\r
157 }\r
158 AssignDetType();\r
159} \r
160\r
161\r
162//_____________________________________________________________________________\r
163AliITSAlignMille2Module::AliITSAlignMille2Module(const AliITSAlignMille2Module &m) :\r
164 TNamed(m),\r
165 fNSensVol(m.fNSensVol),\r
166 fIndex(m.fIndex), \r
167 fDetType(m.fDetType),\r
168 fVolumeID(m.fVolumeID),\r
169 fNParTot(m.fNParTot),\r
170 fNParFree(m.fNParFree),\r
171 fParOffs(m.fNParTot),\r
172 fNProcPoints(0),\r
173 fParVals(0),\r
174 fParErrs(0),\r
175 fParCstr(0), \r
176 fSensVolIndex(m.fSensVolIndex),\r
177 fSensVolVolumeID(m.fSensVolVolumeID),\r
178 fMatrix(new TGeoHMatrix(*m.GetMatrix())),\r
179 fSensVolMatrix(new TGeoHMatrix),\r
180 fSensVolModifMatrix(new TGeoHMatrix),\r
181 fParent(m.fParent),\r
182 fChildren(0)\r
183{\r
184 // Copy constructor\r
185 fSensVolIndex = m.fSensVolIndex;\r
186 fSensVolVolumeID = m.fSensVolVolumeID;\r
187 for (int i=m.fNParTot;i--;) fParOffs[i] = m.fParOffs[i];\r
188 for (int i=3;i--;) fSigmaFactor[i] = m.fSigmaFactor[i];\r
189 if (fNParTot) {\r
190 fParVals = new Float_t[fNParTot];\r
191 fParErrs = new Float_t[fNParTot];\r
192 fParCstr = new Float_t[fNParTot];\r
193 for (int i=fNParTot;i--;) {\r
194 fParVals[i] = m.fParVals[i];\r
195 fParErrs[i] = m.fParErrs[i];\r
196 fParCstr[i] = m.fParCstr[i];\r
197 }\r
198 }\r
199}\r
200\r
201//_____________________________________________________________________________\r
202AliITSAlignMille2Module& AliITSAlignMille2Module::operator=(const AliITSAlignMille2Module &m) \r
203{\r
204 // operator =\r
205 //\r
206 if(this==&m) return *this;\r
207 ((TNamed *)this)->operator=(m);\r
208 //\r
209 fNSensVol=m.fNSensVol;\r
210 fIndex=m.fIndex;\r
211 fDetType = m.fDetType;\r
212 fVolumeID=m.fVolumeID;\r
213 fNParTot = m.fNParTot;\r
214 fNParFree = m.fNParFree;\r
215 delete[] fParVals; fParVals = 0;\r
216 delete[] fParErrs; fParErrs = 0;\r
217 delete[] fParCstr; fParCstr = 0;\r
218 //\r
219 if (fNParTot) {\r
220 fParVals = new Float_t[fNParTot];\r
221 fParErrs = new Float_t[fNParTot];\r
8102b2c9 222 fParCstr = new Float_t[fNParTot];\r
6526a72c 223 for (int i=m.GetNParTot();i--;) {\r
224 fParVals[i] = m.fParVals[i];\r
225 fParErrs[i] = m.fParErrs[i];\r
226 fParCstr[i] = m.fParCstr[i];\r
227 }\r
228 }\r
229 //\r
230 fParOffs.Set(fNParTot);\r
231 for (int i=0;i<fNParTot;i++) fParOffs[i] = m.fParOffs[i];\r
232 for (int i=3;i--;) fSigmaFactor[i] = m.fSigmaFactor[i];\r
233 if (fMatrix) delete fMatrix;\r
234 fMatrix=new TGeoHMatrix(*m.GetMatrix());\r
235 fSensVolIndex = m.fSensVolIndex;\r
236 fSensVolVolumeID = m.fSensVolVolumeID;\r
237 fParent = m.fParent;\r
238 fChildren.Clear();\r
239 for (int i=0;i<m.GetNChildren();i++) fChildren.Add(m.GetChild(i));\r
240 return *this;\r
241}\r
242\r
243\r
244//-------------------------------------------------------------\r
245AliITSAlignMille2Module::~AliITSAlignMille2Module() { \r
246 /// Destructor \r
247 delete fMatrix; \r
248 delete fSensVolMatrix; \r
249 delete fSensVolModifMatrix; \r
250 delete[] fParVals;\r
251 delete[] fParErrs;\r
252 delete[] fParCstr;\r
253 fChildren.Clear();\r
254} \r
255\r
256//-------------------------------------------------------------\r
6be22b3f 257Int_t AliITSAlignMille2Module::Set(Int_t index, UShort_t volid, const char* symname, \r
258 const TGeoHMatrix *m, Int_t nsv, const UShort_t *volidsv) \r
6526a72c 259{\r
260 // initialize a custom supermodule\r
261 // index, volid, symname and matrix must be given\r
262 // if (volidsv) add nsv sensitive volumes to the supermodules\r
263 // return 0 if success\r
264\r
265 if (index<2198) {\r
266 AliInfo("Index must be >= 2198");\r
267 return -1;\r
268 }\r
269 if (volid<14336) {\r
270 AliInfo("VolumeID must be >= 14336");\r
271 return -2;\r
272 }\r
273 \r
274 if (!symname) return -3;\r
275 for (Int_t i=0; i<2198; i++) {\r
276 if (!strcmp(symname,AliITSgeomTGeo::GetSymName(i))) {\r
277 AliInfo("Symname already used by a Sensitive Volume");\r
278 return -3;\r
279 }\r
280 }\r
281 \r
282 if (!m) return -4;\r
283\r
284 // can initialize needed stuffs\r
285 fIndex = index;\r
286 fVolumeID = volid;\r
287 SetName(symname);\r
288 //\r
289 (*fMatrix) = (*m);\r
290 //\r
291 fSensVolIndex.Set(nsv);\r
292 fSensVolVolumeID.Set(nsv);\r
293 // add sensitive volumes\r
294 for (Int_t i=0; i<nsv; i++) AddSensitiveVolume(volidsv[i]);\r
295\r
296 return 0;\r
297}\r
298\r
299//-------------------------------------------------------------\r
300void AliITSAlignMille2Module::SetFreeDOF(Int_t dof,Double_t cstr)\r
301{\r
ef24eb3b 302 if (AliITSAlignMille2::IsZero(cstr)) fParCstr[dof] = 0; // fixed parameter\r
303 else if (cstr>0) fParCstr[dof] = fgkDummyConstraint+1.; // the parameter is free and unconstrained\r
304 else fParCstr[dof] = -cstr; // the parameter is free but constrained\r
6526a72c 305}\r
306\r
307//-------------------------------------------------------------\r
308Bool_t AliITSAlignMille2Module::IsSensor(UShort_t voluid) \r
309{\r
310 // Does this volid correspond to sensor ?\r
311 AliGeomManager::ELayerID layId = AliGeomManager::VolUIDToLayerSafe(voluid);\r
312 if (layId>0 && layId<7) {\r
313 Int_t mId = Int_t(voluid & 0x7ff);\r
314 if( mId>=0 && mId<AliGeomManager::LayerSize(layId)) return kTRUE;\r
315 }\r
316 return kFALSE;\r
317}\r
318\r
319//-------------------------------------------------------------\r
320Int_t AliITSAlignMille2Module::GetIndexFromVolumeID(UShort_t voluid) {\r
321 /// index from volume ID\r
322 AliGeomManager::ELayerID lay = AliGeomManager::VolUIDToLayer(voluid);\r
323 if (lay<1|| lay>6) return -1;\r
324 Int_t idx=Int_t(voluid)-2048*lay;\r
325 if (idx>=AliGeomManager::LayerSize(lay)) return -1;\r
326 for (Int_t ilay=1; ilay<lay; ilay++) \r
327 idx += AliGeomManager::LayerSize(ilay);\r
328 return idx;\r
329}\r
330\r
331//-------------------------------------------------------------\r
332void AliITSAlignMille2Module::AddSensitiveVolume(UShort_t voluid)\r
333{\r
334 /// add a sensitive volume to this supermodule\r
335 if (GetIndexFromVolumeID(voluid)<0) return; // bad volid\r
336 //\r
337 // in principle, the correct size of fSensVol... arrays was set outside but check anyway\r
8fd71c0a 338 if (fSensVolVolumeID.GetSize()<fNSensVol+1) {\r
6526a72c 339 fSensVolVolumeID.Set(fNSensVol+1);\r
340 fSensVolIndex.Set(fNSensVol+1);\r
341 }\r
342 //\r
343 fSensVolVolumeID[fNSensVol] = Short_t(voluid);\r
344 fSensVolIndex[fNSensVol] = GetIndexFromVolumeID(voluid);\r
345 fNSensVol++;\r
346}\r
347\r
348//-------------------------------------------------------------\r
349void AliITSAlignMille2Module::DelSensitiveVolume(Int_t at)\r
350{\r
351 // Suppress sensor at position "at"\r
352 // in fact we are swapping with the last valid one \r
353 int lastValid = --fNSensVol;\r
354 int tmpv = fSensVolIndex[at];\r
355 fSensVolIndex[at] = fSensVolIndex[lastValid];\r
356 tmpv = fSensVolVolumeID[at];\r
357 fSensVolVolumeID[at] = fSensVolVolumeID[lastValid];\r
358 fSensVolVolumeID[lastValid] = tmpv;\r
359 //\r
360}\r
361\r
362//-------------------------------------------------------------\r
363Bool_t AliITSAlignMille2Module::IsIn(UShort_t voluid) const \r
364{\r
365 /// check if voluid is defined\r
366 if (!voluid) return kFALSE; // only positive voluid are accepted\r
367 for (Int_t i=0; i<fNSensVol; i++) if (UShort_t(fSensVolVolumeID[i])==voluid) return kTRUE;\r
368 return kFALSE;\r
369}\r
370\r
371//-------------------------------------------------------------\r
372Bool_t AliITSAlignMille2Module::BelongsTo(AliITSAlignMille2Module* parent) const\r
373{\r
374 /// check if parent contains the sensors of this volume\r
375 if (fNSensVol<1 || fNSensVol>=parent->GetNSensitiveVolumes()) return kFALSE;\r
376 return parent->IsIn( fSensVolVolumeID[0] );\r
377}\r
378\r
379//-------------------------------------------------------------\r
b80c197e 380TGeoHMatrix *AliITSAlignMille2Module::GetSensitiveVolumeModifiedMatrix(UShort_t voluid, const Double_t *delta,Bool_t local)\r
6526a72c 381{\r
382 // modify the original TGeoHMatrix of the sensitive module 'voluid' according\r
383 // with a delta transform. applied to the supermodule matrix\r
384 // return NULL if error\r
385\r
386 if (!IsIn(voluid)) return NULL;\r
387 if (!gGeoManager) return NULL;\r
388\r
389 // prepare the TGeoHMatrix\r
390 Double_t tr[3],ang[3];\r
391 tr[0]=delta[0]; // in centimeter\r
392 tr[1]=delta[1]; \r
393 tr[2]=delta[2];\r
394 ang[0]=delta[3]; // psi (X) in deg\r
395 ang[1]=delta[4]; // theta (Y)\r
396 ang[2]=delta[5]; // phi (Z)\r
397 //\r
398 fgTempAlignObj.SetRotation(ang[0],ang[1],ang[2]);\r
399 fgTempAlignObj.SetTranslation(tr[0],tr[1],tr[2]);\r
400 AliDebug(3,Form("Delta angles: psi=%f theta=%f phi=%f",ang[0],ang[1],ang[2]));\r
401 TGeoHMatrix hm;\r
402 fgTempAlignObj.GetMatrix(hm);\r
403 //printf("\n0: delta matrix\n");hm.Print();\r
404\r
405 // 1) start setting fSensVolModif = fSensVol\r
406 if (SensVolMatrix(voluid, fSensVolModifMatrix)) return NULL;\r
407 //\r
408 if (local) {\r
409 // 2) set fSensVolModif = SensVolRel\r
410 fSensVolModifMatrix->MultiplyLeft( &fMatrix->Inverse() );\r
411 // 3) multiply left by delta\r
412 fSensVolModifMatrix->MultiplyLeft( &hm );\r
413 // 4) multiply left by fMatrix\r
414 fSensVolModifMatrix->MultiplyLeft( fMatrix );\r
415 }\r
416 else fSensVolModifMatrix->MultiplyLeft( &hm );\r
417 //\r
418 return fSensVolModifMatrix;\r
419}\r
420\r
421//-------------------------------------------------------------\r
b80c197e 422AliAlignObjParams *AliITSAlignMille2Module::GetSensitiveVolumeMisalignment(UShort_t voluid, const Double_t *deltalocal)\r
6526a72c 423{\r
424 // calculate misalignment of sens.vol. 'voluid' according with a displacement 'deltalocal'\r
425 // of the mother volume. The misalignment is returned as AliAlignObjParams object\r
426\r
427 if (!IsIn(voluid)) return NULL;\r
428 if (!gGeoManager) return NULL;\r
429 \r
430 // prepare the TGeoHMatrix\r
431 Double_t tr[3],ang[3];\r
432 tr[0]=deltalocal[0]; // in centimeter\r
433 tr[1]=deltalocal[1]; \r
434 tr[2]=deltalocal[2];\r
435 ang[0]=deltalocal[3]; // psi (X) in deg\r
436 ang[1]=deltalocal[4]; // theta (Y)\r
437 ang[2]=deltalocal[5]; // phi (Z)\r
438 //\r
439 fgTempAlignObj.SetRotation(ang[0],ang[1],ang[2]);\r
440 fgTempAlignObj.SetTranslation(tr[0],tr[1],tr[2]);\r
441 AliDebug(3,Form("Delta angles: psi=%f theta=%f phi=%f",ang[0],ang[1],ang[2]));\r
442 //\r
443 return GetSensitiveVolumeMisalignment(voluid,&fgTempAlignObj);\r
444}\r
445\r
446//-------------------------------------------------------------\r
b80c197e 447AliAlignObjParams *AliITSAlignMille2Module::GetSensitiveVolumeMisalignment(UShort_t voluid, const AliAlignObjParams *a)\r
6526a72c 448{\r
449 // return the misalignment of the sens. vol. 'voluid' corresponding with \r
450 // a misalignment 'a' in the mother volume\r
451 // return NULL if error\r
452\r
453 // Gsv = Gg * Gg-1 * Gsv -> Lsv,g = Gg-1 * Gsv\r
454 // G'sv = Gg * Dg * Lsv,g === Gsv * Dsv\r
455 // Gg * Dg * Gg-1 * Gsv = Gsv * Gsv-1 * Gg * Dg * Gg-1 * Gsv\r
456 //\r
457 // => Dsv = (Gsv-1 * Gg * Dg * Gg-1 * Gsv)\r
458 //\r
459\r
460 if (!IsIn(voluid)) return NULL;\r
461 if (!gGeoManager) return NULL;\r
462\r
463 //a->Print("");\r
464\r
465 // prepare the Delta matrix Dg\r
466 TGeoHMatrix dg;\r
467 a->GetMatrix(dg);\r
468 //dg.Print();\r
469\r
470 // 1) start setting fSensVolModif = Gsv\r
471 if (SensVolMatrix(voluid, fSensVolModifMatrix)) return NULL;\r
472 //printf("\n1: modif=orig del sensvol\n");fSensVolModifMatrix->Print();\r
473\r
474 // 2) set fSensVolModif = Gg-1 * Gsv\r
475 fSensVolModifMatrix->MultiplyLeft( &fMatrix->Inverse() );\r
476 //printf("\n2: modif=relative del sensvol\n");fSensVolModifMatrix->Print();\r
477 \r
478 // 3) set fSensVolModif = Dg * Gg-1 * Gsv\r
479 fSensVolModifMatrix->MultiplyLeft( &dg );\r
480 //printf("\n3: modif= delta*relative\n");fSensVolModifMatrix->Print();\r
481 \r
482 // 4) set fSensVolModif = Gg * Dg * Gg-1 * Gsv\r
483 fSensVolModifMatrix->MultiplyLeft( fMatrix );\r
484 //printf("\n4: modif=quasi,manca il Gsv-1...\n");fSensVolModifMatrix->Print();\r
485\r
486 // 5) set fSensVolModif = Gsv-1 * Gg * Dg * Gg-1 * Gsv\r
487 if (SensVolMatrix(voluid, &dg)) return NULL;\r
488 fSensVolModifMatrix->MultiplyLeft( &dg.Inverse() );\r
489 //printf("\n5: modif=finale\n");fSensVolModifMatrix->Print();\r
490 //\r
491 // >> RS\r
492 // 6) mo' fSensVolModif dovrebbe essere la Dsv(loc) t.c. G'sv = Gsv*Dsv(loc)\r
493 // per trasformarla in Dsv(loc rispetto al Gsv0, non modificato) dovrebbe essere:\r
494 // Dsv(loc) -> Dpre * Dsv(loc) * Dpre-1\r
495 //TGeoHMatrix dpre; // dpre = Gsv0-1*Gsv\r
496 //if (SensVolOrigGlobalMatrix(voluid, &dg)) return NULL;\r
497 //if (SensVolMatrix(voluid, &dpre)) return NULL;\r
498 //dpre.MultiplyLeft( &dg.Inverse() );\r
499 //fSensVolModifMatrix->Multiply( &dpre.Inverse() );\r
500 //fSensVolModifMatrix->MultiplyLeft( &dpre );\r
501 // direi che NON FUNZIONA!!!! \r
502\r
503 // << RS\r
504\r
505 // reset align object (may not be needed...)\r
506 fgTempAlignObj.SetVolUID(0);\r
507 fgTempAlignObj.SetSymName("");\r
508 fgTempAlignObj.SetTranslation(0,0,0);\r
509 fgTempAlignObj.SetRotation(0,0,0);\r
510 //\r
511 // >> RS\r
512#ifdef CORHW_\r
513 // correction for SPD y-shift\r
514 if (voluid>=2048 && voluid<4256) {\r
515 TGeoHMatrix deltay;\r
516 double dy[3]={0.,0.0081,0.};\r
517 deltay.SetTranslation(dy);\r
518 fSensVolModifMatrix->MultiplyLeft( &deltay );\r
519 fSensVolModifMatrix->Multiply( &deltay.Inverse() );\r
520 }\r
521#endif\r
522 // << RS\r
523 if (!fgTempAlignObj.SetMatrix(*fSensVolModifMatrix)) return NULL;\r
524 fgTempAlignObj.SetVolUID(voluid);\r
525 fgTempAlignObj.SetSymName(AliGeomManager::SymName(voluid));\r
526 //\r
527 return &fgTempAlignObj;\r
528}\r
529\r
530// >> RS\r
531//-------------------------------------------------------------\r
b80c197e 532AliAlignObjParams *AliITSAlignMille2Module::GetSensitiveVolumeTotalMisalignment(UShort_t voluid, const Double_t *deltalocal)\r
6526a72c 533{\r
534 // calculate misalignment of sens.vol. 'voluid' according with a displacement 'deltalocal'\r
535 // of the mother volume. The misalignment is returned as AliAlignObjParams object including\r
536 // the (evenctual) prealignment => no merging needed\r
537\r
538 if (!IsIn(voluid)) return NULL;\r
539 if (!gGeoManager) return NULL;\r
540 \r
541 // prepare the TGeoHMatrix\r
542 Double_t tr[3],ang[3];\r
543 tr[0]=deltalocal[0]; // in centimeter\r
544 tr[1]=deltalocal[1]; \r
545 tr[2]=deltalocal[2];\r
546 ang[0]=deltalocal[3]; // psi (X) in deg\r
547 ang[1]=deltalocal[4]; // theta (Y)\r
548 ang[2]=deltalocal[5]; // phi (Z)\r
549\r
550 // reset align object (may not be needed...)\r
551 fgTempAlignObj.SetVolUID(0);\r
552 fgTempAlignObj.SetSymName("");\r
553 fgTempAlignObj.SetRotation(ang[0],ang[1],ang[2]);\r
554 fgTempAlignObj.SetTranslation(tr[0],tr[1],tr[2]);\r
555 AliDebug(3,Form("Delta angles: psi=%f theta=%f phi=%f",ang[0],ang[1],ang[2]));\r
556\r
557 // Gsv = Gg * Gg-1 * Gsv -> Lsv,g = Gg-1 * Gsv\r
558 // G'sv = Gg * Dg * Lsv,g === DGsv * Gsv \r
559 //\r
560 // => Dsv = (G0sv-1 * Gg * Dg * Gg-1 * GMsv) //\r
561 //\r
562\r
563 // prepare the Delta matrix Dg\r
564 TGeoHMatrix dg;\r
565 fgTempAlignObj.GetMatrix(dg);\r
566 //dg.Print();\r
567\r
568 // 1) start setting fSensVolModif = Gsv\r
569 if (SensVolMatrix(voluid, fSensVolModifMatrix)) return NULL;\r
570 //printf("\n1: modif=orig del sensvol\n");fSensVolModifMatrix->Print();\r
571\r
572 // 2) set fSensVolModif = Gg-1 * Gsv\r
573 fSensVolModifMatrix->MultiplyLeft( &fMatrix->Inverse() );\r
574 //printf("\n2: modif=relative del sensvol\n");fSensVolModifMatrix->Print();\r
575 \r
576 // 3) set fSensVolModif = Dg * Gg-1 * Gsv\r
577 fSensVolModifMatrix->MultiplyLeft( &dg );\r
578 //printf("\n3: modif= delta*relative\n");fSensVolModifMatrix->Print();\r
579 \r
580 // 4) set fSensVolModif = Gg * Dg * Gg-1 * Gsv\r
581 fSensVolModifMatrix->MultiplyLeft( fMatrix );\r
582 //printf("\n4: modif=quasi,manca il Gsv-1...\n");fSensVolModifMatrix->Print();\r
583\r
584 // 5) set fSensVolModif = G0sv-1 * Gg * Dg * Gg-1 * Gsv\r
585 // qui usa l'orig anziche' la prealigned...\r
586 if (SensVolOrigGlobalMatrix(voluid, &dg)) return NULL;\r
587 fSensVolModifMatrix->MultiplyLeft( &dg.Inverse() );\r
588 //printf("\n5: modif=finale\n");fSensVolModifMatrix->Print();\r
589\r
590 // reset align object (may not be needed...)\r
591 fgTempAlignObj.SetVolUID(0);\r
592 fgTempAlignObj.SetSymName("");\r
593 fgTempAlignObj.SetTranslation(0,0,0);\r
594 fgTempAlignObj.SetRotation(0,0,0);\r
595\r
596#ifdef CORHW_\r
597 // correction for SPD y-shift\r
598 if (voluid>=2048 && voluid<4256) {\r
599 TGeoHMatrix deltay;\r
600 double dy[3]={0.,0.0081,0.};\r
601 deltay.SetTranslation(dy);\r
602 fSensVolModifMatrix->MultiplyLeft( &deltay );\r
603 fSensVolModifMatrix->Multiply( &deltay.Inverse() );\r
604 }\r
605#endif\r
606 if (!fgTempAlignObj.SetMatrix(*fSensVolModifMatrix)) return NULL;\r
607 fgTempAlignObj.SetVolUID(voluid);\r
608 fgTempAlignObj.SetSymName(AliGeomManager::SymName(voluid));\r
609\r
610 \r
611 //fgTempAlignObj.Print("");\r
612\r
613 return &fgTempAlignObj;\r
614}\r
615//-------------------------------------------------------------\r
616\r
617//-------------------------------------------------------------\r
b80c197e 618AliAlignObjParams *AliITSAlignMille2Module::GetSensitiveVolumeGlobalMisalignment(UShort_t voluid, const Double_t *deltalocal)\r
6526a72c 619{\r
620 // calculate misalignment of sens.vol. 'voluid' according with a displacement 'deltalocal'\r
621 // of the mother volume. The misalignment is returned as AliAlignObjParams object\r
622\r
623 if (!IsIn(voluid)) return NULL;\r
624 if (!gGeoManager) return NULL;\r
625 \r
626 // prepare the TGeoHMatrix\r
627 Double_t tr[3],ang[3];\r
628 tr[0]=deltalocal[0]; // in centimeter\r
629 tr[1]=deltalocal[1]; \r
630 tr[2]=deltalocal[2];\r
631 ang[0]=deltalocal[3]; // psi (X) in deg\r
632 ang[1]=deltalocal[4]; // theta (Y)\r
633 ang[2]=deltalocal[5]; // phi (Z)\r
634\r
635 // reset align object (may not be needed...)\r
636 fgTempAlignObj.SetTranslation(0,0,0);\r
637 fgTempAlignObj.SetRotation(0,0,0);\r
638\r
639 fgTempAlignObj.SetRotation(ang[0],ang[1],ang[2]);\r
640 fgTempAlignObj.SetTranslation(tr[0],tr[1],tr[2]);\r
641 AliDebug(3,Form("Delta angles: psi=%f theta=%f phi=%f",ang[0],ang[1],ang[2]));\r
642\r
643 // Gsv = Gg * Gg-1 * Gsv -> Lsv,g = Gg-1 * Gsv\r
644 // G'sv = Gg * Dg * Lsv,g === DGsv * Gsv \r
645 //\r
646 // => DGsv = (Gg * Dg * Gg-1)\r
647 //\r
648\r
649 // prepare the Delta matrix Dg\r
650 TGeoHMatrix dg;\r
651 fgTempAlignObj.GetMatrix(dg);\r
652 //dg.Print();\r
653\r
654 dg.MultiplyLeft( fMatrix );\r
655 dg.Multiply( &fMatrix->Inverse() );\r
656\r
657 // reset align object (may not be needed...)\r
658 fgTempAlignObj.SetTranslation(0,0,0);\r
659 fgTempAlignObj.SetRotation(0,0,0);\r
660\r
661 fgTempAlignObj.SetVolUID(voluid);\r
662 fgTempAlignObj.SetSymName(AliGeomManager::SymName(voluid));\r
663\r
664 if (!fgTempAlignObj.SetMatrix(dg)) return NULL;\r
665 \r
666 //fgTempAlignObj.Print("");\r
667\r
668 return &fgTempAlignObj;\r
669}\r
670// << RS\r
671\r
672//-------------------------------------------------------------\r
673TGeoHMatrix *AliITSAlignMille2Module::GetSensitiveVolumeMatrix(UShort_t voluid)\r
674{\r
675 // return TGeoHMatrix of the sens.vol. 'voluid' in the current geometry\r
676 if (SensVolMatrix(voluid,fSensVolMatrix)) return NULL;\r
677 return fSensVolMatrix;\r
678}\r
679\r
680//-------------------------------------------------------------\r
681TGeoHMatrix *AliITSAlignMille2Module::GetSensitiveVolumeOrigGlobalMatrix(UShort_t voluid)\r
682{\r
683 // return original ideal position (from AliGeomManager::GetOrigGlobalMatrix())\r
684 if (SensVolOrigGlobalMatrix(voluid,fSensVolMatrix)) return NULL;\r
685 return fSensVolMatrix;\r
686}\r
687//-------------------------------------------------------------\r
688Int_t AliITSAlignMille2Module::SensVolMatrix(UShort_t volid, TGeoHMatrix *m) \r
689{\r
690 // set matrix for sensitive modules (SPD corrected)\r
691 // return 0 if success\r
692 Double_t rot[9];\r
693 Int_t idx=GetIndexFromVolumeID(volid);\r
694 if (idx<0) return -1;\r
695 if (!AliITSgeomTGeo::GetRotation(idx,rot)) return -2;\r
696 m->SetRotation(rot);\r
697 Double_t oLoc[3]={0,0,0};\r
698 Double_t oGlo[3]={0,0,0};\r
699 if (!AliITSgeomTGeo::LocalToGlobal(idx,oLoc,oGlo)) return -3;\r
700 m->SetTranslation(oGlo);\r
701 return 0;\r
702}\r
703\r
704//-------------------------------------------------------------\r
705Int_t AliITSAlignMille2Module::SensVolOrigGlobalMatrix(UShort_t volid, TGeoHMatrix *m) \r
706{\r
707 // set original global matrix for sensitive modules (SPD corrected)\r
708 // return 0 if success\r
709 Int_t idx=GetIndexFromVolumeID(volid);\r
710 if (idx<0) return -1;\r
711 TGeoHMatrix mo;\r
712 if (!AliGeomManager::GetOrigGlobalMatrix(AliGeomManager::SymName(volid),mo)) return -1;\r
713 (*m)=mo;\r
714 //\r
715#ifdef CORHW_\r
716 // SPD y-shift by 81 mu\r
717 if (volid<5000) { \r
718 Double_t oLoc[3]={0.0,0.0081,0.0};\r
719 Double_t oGlo[3]={0,0,0};\r
720 m->LocalToMaster(oLoc,oGlo);\r
721 m->SetTranslation(oGlo);\r
722 }\r
723#endif\r
724 return 0;\r
725}\r
726\r
727//-------------------------------------------------------------\r
728UShort_t AliITSAlignMille2Module::GetVolumeIDFromSymname(const Char_t *symname) {\r
729 /// volume ID from symname\r
730 if (!symname) return 0;\r
731\r
732 for (UShort_t voluid=2000; voluid<13300; voluid++) {\r
733 Int_t modId;\r
734 AliGeomManager::ELayerID layerId = AliGeomManager::VolUIDToLayer(voluid,modId);\r
735 if (layerId>0 && layerId<7 && modId>=0 && modId<AliGeomManager::LayerSize(layerId)) {\r
736 if (!strcmp(symname,AliGeomManager::SymName(layerId,modId))) return voluid;\r
737 }\r
738 }\r
739\r
740 return 0;\r
741}\r
742\r
743//-------------------------------------------------------------\r
744UShort_t AliITSAlignMille2Module::GetVolumeIDFromIndex(Int_t index) {\r
745 /// volume ID from index\r
746 if (index<0 || index>2197) return 0;\r
747 return GetVolumeIDFromSymname(AliITSgeomTGeo::GetSymName(index));\r
748}\r
749\r
750//-------------------------------------------------------------\r
751void AliITSAlignMille2Module::Print(Option_t*) const \r
752{\r
b80c197e 753 // print data\r
6526a72c 754 //\r
755 const char* typeName[] = {"SPD","SDD","SSD"};\r
756 printf("*** ITS SuperModule for AliITSAlignMille ***\n");\r
757 printf("symname : %s (type: %s)\n",GetName(),fDetType<0 ? "N/A":typeName[fDetType]);\r
758 printf("parent : %s | %d children\n",fParent ? fParent->GetName() : "N/A",GetNChildren());\r
759 printf("volumeID : %4d | index : %4d | Geom.Params are %s\n",fVolumeID,fIndex,\r
760 GeomParamsGlobal() ? "Global":"Local");\r
761 printf("Factors : X=%.2f Y=%.2f Z=%.2f\n"\r
762 "DOF: %cTx:%5d| %cTy:%5d| %cTz:%5d| %cPsi:%5d| %cTheta:%5d| %cPhi:%5d|",\r
763 fSigmaFactor[0],fSigmaFactor[1],fSigmaFactor[2],\r
ef24eb3b 764 IsFreeDOF(kDOFTX) ? '+':'-',GetParOffset(kDOFTX),IsFreeDOF(kDOFTY) ? '+':'-',GetParOffset(kDOFTY),\r
765 IsFreeDOF(kDOFTZ) ? '+':'-',GetParOffset(kDOFTZ),IsFreeDOF(kDOFPS) ? '+':'-',GetParOffset(kDOFPS),\r
766 IsFreeDOF(kDOFTH) ? '+':'-',GetParOffset(kDOFTH),IsFreeDOF(kDOFPH) ? '+':'-',GetParOffset(kDOFPH));\r
767 if (IsSDD()) {\r
768 printf("%cT0:%5d| %cDVl:%5d| %cDVr:%5d|",IsFreeDOF(kDOFT0)?'+':'-',GetParOffset(kDOFT0),\r
769 IsFreeDOF(kDOFDVL)?'+':'-',GetParOffset(kDOFDVL),IsFreeDOF(kDOFDVR)?'+':'-',GetParOffset(kDOFDVR));\r
770 if (IsVDriftLRSame()) printf("(dVL=dVR)");\r
771 }\r
6526a72c 772 printf("\n");\r
773 fMatrix->Print();\r
774 printf("%4d Sensitive volumes | %6d Processed Points\n",fNSensVol,fNProcPoints);\r
775 for (Int_t i=0; i<fNSensVol; i++) printf(" voluid[%d] = %d\n",i,UShort_t(fSensVolVolumeID[i]));\r
776}\r
777\r
778//-------------------------------------------------------------\r
779Bool_t AliITSAlignMille2Module::IsAlignable() const\r
780{\r
b80c197e 781 // it it alignable?\r
6526a72c 782 TGeoManager* geoManager = AliGeomManager::GetGeometry();\r
783 if (!geoManager) {\r
784 AliInfo("Couldn't initialize geometry");\r
785 return kFALSE;\r
786 }\r
787 return geoManager->GetAlignableEntry(GetName())!=0;\r
788}\r
789\r
790//-------------------------------------------------------------\r
791void AliITSAlignMille2Module::GetLocalMatrix(TGeoHMatrix &mat) const\r
792{\r
793 // return the local matrix for transformation to its parent\r
794 mat = *fMatrix;\r
795 if (fParent) mat.MultiplyLeft( &fParent->GetMatrix()->Inverse() );\r
796}\r
797\r
798//-------------------------------------------------------------\r
799void AliITSAlignMille2Module::AssignDetType()\r
800{\r
b80c197e 801 // assign the detector type\r
6526a72c 802 TString tp = GetName();\r
803 if (tp.Contains("SPD",TString::kIgnoreCase)) fDetType = kSPD;\r
804 else if (tp.Contains("SDD",TString::kIgnoreCase)) fDetType = kSDD;\r
805 else if (tp.Contains("SSD",TString::kIgnoreCase)) fDetType = kSSD;\r
806 else fDetType = -1;\r
ef24eb3b 807 fNParTot = IsSDD() ? kMaxParTot:kMaxParGeom;\r
6526a72c 808 fNParFree = 0;\r
809 fParVals = new Float_t[fNParTot];\r
810 fParErrs = new Float_t[fNParTot]; \r
811 fParCstr = new Float_t[fNParTot]; \r
812 if (fParOffs.GetSize()<fNParTot) fParOffs.Set(fNParTot);\r
813 for (int i=fNParTot;i--;) {\r
814 fParVals[i] = fParErrs[i] = 0.; \r
815 fParCstr[i] = 0.;\r
816 fParOffs[i] = -1;\r
817 }\r
818}\r
819\r
820//-------------------------------------------------------------\r
821void AliITSAlignMille2Module::EvaluateDOF()\r
822{\r
b80c197e 823 // count d.o.f.\r
6526a72c 824 fNParFree = 0;\r
825 for (int i=fNParTot;i--;) if (IsFreeDOF(i)) fNParFree++;\r
826}\r
827\r
828//-------------------------------------------------------------\r
829void AliITSAlignMille2Module::GetSensVolGlobalParams(UShort_t volid,Double_t *t, Double_t *r)\r
830{\r
831 // return global parameters of the sensor volid\r
832 for (int i=3;i--;) t[i] = r[i] = 0.;\r
833 if (SensVolMatrix(volid,fSensVolMatrix)) return; \r
834 fgTempAlignObj.SetMatrix(*fSensVolMatrix);\r
835 fgTempAlignObj.GetPars(t,r);\r
836}\r
837\r
838//-------------------------------------------------------------\r
839void AliITSAlignMille2Module::GetSensVolLocalParams(UShort_t volid,Double_t *t, Double_t *r)\r
840{\r
841 // return parameters of the sensor volid in the current module\r
842 for (int i=3;i--;) t[i] = r[i] = 0.;\r
843 if (SensVolMatrix(volid,fSensVolMatrix)) return; \r
844 fSensVolMatrix->MultiplyLeft( &fMatrix->Inverse() );\r
845 fgTempAlignObj.SetMatrix(*fSensVolMatrix);\r
846 fgTempAlignObj.GetPars(t,r);\r
847}\r
848\r
849//-------------------------------------------------------------\r
b80c197e 850void AliITSAlignMille2Module::GetSensVolGlobalParams(UShort_t volid,const Double_t* loct, const Double_t* locr,Double_t *t, Double_t *r)\r
6526a72c 851{\r
852 // return global parameters of the sensor volid modified by the localDelta params\r
853 for (int i=3;i--;) t[i] = r[i] = 0.;\r
854 if (SensVolMatrix(volid,fSensVolMatrix)) return; \r
855 fgTempAlignObj.SetTranslation(loct[0],loct[1],loct[2]);\r
856 fgTempAlignObj.SetRotation(locr[0],locr[1],locr[2]);\r
857 //\r
858 fgTempAlignObj.GetMatrix(*fSensVolModifMatrix); // obtain local delta\r
859 fSensVolModifMatrix->MultiplyLeft( fSensVolMatrix ); // obtain global delta\r
860 fgTempAlignObj.SetMatrix(*fSensVolModifMatrix);\r
861 fgTempAlignObj.GetPars(t,r); // obtain global params\r
862}\r
863\r
864//-------------------------------------------------------------\r
b80c197e 865void AliITSAlignMille2Module::GetSensVolLocalParams(UShort_t volid,const Double_t* loct,const Double_t* locr,Double_t *t, Double_t *r)\r
6526a72c 866{\r
867 // return parameters of the sensor volid (modified by the localDelta params) in the current volume\r
868 for (int i=3;i--;) t[i] = r[i] = 0.;\r
869 if (SensVolMatrix(volid,fSensVolMatrix)) return; \r
870 fgTempAlignObj.SetTranslation(loct[0],loct[1],loct[2]);\r
871 fgTempAlignObj.SetRotation(locr[0],locr[1],locr[2]);\r
872 //\r
873 fgTempAlignObj.GetMatrix(*fSensVolModifMatrix); // obtain local delta\r
874 fSensVolModifMatrix->MultiplyLeft( fSensVolMatrix ); // obtain global delta\r
875 fSensVolModifMatrix->MultiplyLeft( &fMatrix->Inverse() ); // obtain delta in current volume\r
876 fgTempAlignObj.SetMatrix(*fSensVolModifMatrix);\r
877 fgTempAlignObj.GetPars(t,r); // obtain params\r
878}\r
879\r
880//-------------------------------------------------------------\r
881void AliITSAlignMille2Module::SetParVals(Double_t *vl,Int_t npar)\r
882{\r
b80c197e 883 // set parameters\r
6526a72c 884 for (int i=TMath::Min(npar,(Int_t)fNParTot);i--;) fParVals[i] = vl[i];\r
885}\r
886\r
887//-------------------------------------------------------------\r
888void AliITSAlignMille2Module::GetGeomParamsGlo(Double_t *pars)\r
889{\r
890 // recompute parameters from local to global frame\r
891 //\r
892 // is there anything to do?\r
893 if (GeomParamsGlobal()) {for (int i=kMaxParGeom;i--;) pars[i] = fParVals[i]; return;}\r
894 //\r
895 // IMPORTANT: It is assumed that the parents params are defined in a same way (local or global)\r
896 // as for the current module. Since in the mp2 the modules are stored from parents to children,\r
897 // it is safe to call this method in loop staring from the lowest level child, i.e. from the end\r
898 // of the modules array.\r
899 //\r
900 // DeltaGlobal = (ModifParents)*DeltaLocal*(ModifParents)^-1 \r
901 //\r
902 *fSensVolMatrix = *fMatrix; // current global matrix\r
903 AliITSAlignMille2Module* parent = GetParent();\r
904 while (parent) {\r
905 if (parent->GeomParamsGlobal()) {\r
906 AliError("Cannot convert params to Global when the parents are already Global\n");\r
907 for (int i=kMaxParGeom;i--;) pars[i] = 0;\r
908 return;\r
909 }\r
910 fSensVolMatrix->MultiplyLeft( &parent->GetMatrix()->Inverse() ); // Local Matrix\r
911 Float_t *parpar = parent->GetParVals();\r
912 fgTempAlignObj.SetTranslation(parpar[0],parpar[1],parpar[2]);\r
913 fgTempAlignObj.SetRotation(parpar[3],parpar[4],parpar[5]);\r
914 fgTempAlignObj.GetMatrix(*fSensVolModifMatrix);\r
915 fSensVolMatrix->MultiplyLeft(fSensVolModifMatrix);\r
916 fSensVolMatrix->MultiplyLeft(parent->GetMatrix()); // global matrix after parents modifications\r
917 parent = parent->GetParent();\r
918 }\r
919 //\r
920 fgTempAlignObj.SetTranslation(fParVals[0],fParVals[1],fParVals[2]);\r
921 fgTempAlignObj.SetRotation(fParVals[3],fParVals[4],fParVals[5]);\r
922 fgTempAlignObj.GetMatrix(*fSensVolModifMatrix); // local delta matrix\r
923 fSensVolModifMatrix->Multiply( &fSensVolMatrix->Inverse() );\r
924 fSensVolModifMatrix->MultiplyLeft( fSensVolMatrix );\r
925 fgTempAlignObj.SetMatrix( *fSensVolModifMatrix ); // global delta matrix\r
926 fgTempAlignObj.GetPars(pars,pars+3);\r
927 //\r
928}\r
929\r
930//-------------------------------------------------------------\r
931void AliITSAlignMille2Module::GetGeomParamsLoc(Double_t *pars)\r
932{\r
933 // recompute parameters from global to local frame\r
934 //\r
935 // is there anything to do?\r
936 if (!GeomParamsGlobal()) {for (int i=kMaxParGeom;i--;) pars[i] = fParVals[i]; return;}\r
937 //\r
938 // IMPORTANT: It is assumed that the parents params are defined in a same way (local or global)\r
939 // as for the current module. Since in the mp2 the modules are stored from parents to children,\r
940 // it is safe to call this method in loop staring from the lowest level child, i.e. from the end\r
941 // of the modules array.\r
942 //\r
943 // DeltaLocal = (DeltaParents*GlobalMat)^-1*DeltaGlobal*(DeltaParents*GlobalMat)\r
944 //\r
945 AliITSAlignMille2Module* parent = GetParent();\r
946 fgTempAlignObj.SetTranslation(0.,0.,0.);\r
947 fgTempAlignObj.SetRotation(0.,0.,0.);\r
948 fgTempAlignObj.GetMatrix(*fSensVolMatrix); // get no-shift matrix\r
949 //\r
950 while (parent) { // accumulate the product of parents global modifications\r
951 if (!parent->GeomParamsGlobal()) {\r
952 AliError("Cannot convert params to Local when the parents are already Local\n");\r
953 for (int i=kMaxParGeom;i--;) pars[i] = 0;\r
954 return;\r
955 }\r
956 Float_t *parpar = parent->GetParVals();\r
957 fgTempAlignObj.SetTranslation(parpar[0],parpar[1],parpar[2]);\r
958 fgTempAlignObj.SetRotation(parpar[3],parpar[4],parpar[5]);\r
959 fgTempAlignObj.GetMatrix(*fSensVolModifMatrix);\r
960 fSensVolMatrix->Multiply(fSensVolModifMatrix); \r
961 parent = parent->GetParent();\r
962 }\r
963 // global matrix after parents modifications\r
964 fSensVolMatrix->Multiply(fMatrix);\r
965 //\r
966 fgTempAlignObj.SetTranslation(fParVals[0],fParVals[1],fParVals[2]);\r
967 fgTempAlignObj.SetRotation(fParVals[3],fParVals[4],fParVals[5]);\r
968 fgTempAlignObj.GetMatrix(*fSensVolModifMatrix); // global delta matrix\r
969 fSensVolModifMatrix->MultiplyLeft( &fSensVolMatrix->Inverse() );\r
970 fSensVolModifMatrix->Multiply( fSensVolMatrix );\r
971 fgTempAlignObj.SetMatrix( *fSensVolModifMatrix ); // local delta matrix\r
972 fgTempAlignObj.GetPars(pars,pars+3);\r
973 //\r
974}\r
975\r
6be22b3f 976\r
977//-------------------------------------------------------------\r
978void AliITSAlignMille2Module::CalcDerivDPosDPar(Int_t sensVol,const Double_t* pl, Double_t *deriv)\r
979{\r
980 // calculate jacobian of the global position vs Parameters (dPos/dParam) \r
981 // for the point in the sensor sensVol\r
982 const double kDel = 0.01;\r
983 double pos0[3],pos1[3],pos2[3],pos3[3];\r
984 double delta[kMaxParGeom];\r
985 //\r
986 for (int ip=kMaxParGeom;ip--;) delta[ip] = 0;\r
987 //\r
988 for (int ip=kMaxParGeom;ip--;) {\r
989 //\r
990 delta[ip] -= kDel;\r
991 GetSensitiveVolumeModifiedMatrix(sensVol,delta,!GeomParamsGlobal())->LocalToMaster(pl,pos0); \r
992 delta[ip] += kDel/2;\r
993 GetSensitiveVolumeModifiedMatrix(sensVol,delta,!GeomParamsGlobal())->LocalToMaster(pl,pos1); \r
994 delta[ip] += kDel;\r
995 GetSensitiveVolumeModifiedMatrix(sensVol,delta,!GeomParamsGlobal())->LocalToMaster(pl,pos2); \r
996 delta[ip] += kDel/2;\r
997 GetSensitiveVolumeModifiedMatrix(sensVol,delta,!GeomParamsGlobal())->LocalToMaster(pl,pos3); \r
998 //\r
999 delta[ip] = 0;\r
1000 double *curd = deriv + ip*3;\r
1001 for (int i=3;i--;) curd[i] = (8.*(pos2[i]-pos1[i]) - (pos3[i]-pos0[i]))/6./kDel;\r
1002 }\r
1003 //\r
1004}\r
1005\r
6526a72c 1006//-------------------------------------------------------------\r
1007void AliITSAlignMille2Module::CalcDerivGloLoc(Int_t idx,Double_t *deriv)\r
1008{\r
1009 // calculate derivative of global params vs local param idx: deriv[j] = dParGlo[j]/dParLoc[idx]\r
1010 Double_t lpar[kMaxParGeom];\r
1011 for (int i=kMaxParGeom;i--;) lpar[i] = 0.;\r
1012 // using f(x+h),f(x-h),f(x+h/2),f(x-h2)...\r
1013 Double_t par1[kMaxParGeom]; // f(x-h)\r
1014 Double_t par2[kMaxParGeom]; // f(x-h/2)\r
1015 Double_t par3[kMaxParGeom]; // f(x+h/2)\r
1016 Double_t par4[kMaxParGeom]; // f(x+h)\r
1017 //\r
1018 const Double_t dpar = 1e-3;\r
1019 //\r
1020 // first values\r
1021 lpar[idx] -= dpar;\r
1022 GetGlobalParams(lpar,lpar+3, par1,par1+3);\r
1023 //\r
1024 // second values\r
1025 lpar[idx] += dpar/2;\r
1026 GetGlobalParams(lpar,lpar+3, par2,par2+3);\r
1027 //\r
1028 // third values\r
1029 lpar[idx] += dpar;\r
1030 GetGlobalParams(lpar,lpar+3, par3,par3+3);\r
1031 //\r
1032 // fourth values\r
1033 lpar[idx] += dpar/2;\r
1034 GetGlobalParams(lpar,lpar+3, par4,par4+3);\r
1035 //\r
1036 Double_t h2 = 1./(2.*dpar);\r
1037 for (int i=kMaxParGeom;i--;) {\r
1038 Double_t d0 = par4[i]-par1[i];\r
1039 Double_t d2 = 2.*(par3[i]-par2[i]);\r
1040 deriv[i] = h2*(4*d2 - d0)/3.;\r
1041 if (TMath::Abs(deriv[i]) < 1.0e-9) deriv[i] = 0.0;\r
1042 }\r
1043 //\r
1044}\r
1045\r
1046//-------------------------------------------------------------\r
1047void AliITSAlignMille2Module::CalcDerivLocGlo(Double_t *deriv)\r
1048{\r
1049 // calculate derivative of local params vs global params: deriv[i][j] = dParLoc[i]/dParGlo[j]\r
1050 Double_t gpar[kMaxParGeom];\r
1051 for (int i=kMaxParGeom;i--;) gpar[i] = 0.;\r
1052 // using f(x+h),f(x-h),f(x+h/2),f(x-h2)...\r
1053 Double_t par1[kMaxParGeom]; // f(x-h)\r
1054 Double_t par2[kMaxParGeom]; // f(x-h/2)\r
1055 Double_t par3[kMaxParGeom]; // f(x+h/2)\r
1056 Double_t par4[kMaxParGeom]; // f(x+h)\r
1057 //\r
1058 const Double_t dpar = 1e-3;\r
1059 //\r
1060 for (int ig=kMaxParGeom;ig--;) {\r
1061 // first values\r
1062 gpar[ig] -= dpar;\r
1063 GetLocalParams(gpar,gpar+3, par1,par1+3);\r
1064 //\r
1065 // second values\r
1066 gpar[ig] += dpar/2;\r
1067 GetLocalParams(gpar,gpar+3, par2,par2+3);\r
1068 //\r
1069 // third values\r
1070 gpar[ig] += dpar;\r
1071 GetLocalParams(gpar,gpar+3, par3,par3+3);\r
1072 //\r
1073 // fourth values\r
1074 gpar[ig] += dpar/2;\r
1075 GetLocalParams(gpar,gpar+3, par4,par4+3);\r
1076 //\r
1077 Double_t h2 = 1./(2.*dpar);\r
1078 for (int i=kMaxParGeom;i--;) {\r
1079 Double_t d0 = par4[i]-par1[i];\r
1080 Double_t d2 = 2.*(par3[i]-par2[i]);\r
1081 int idig = i*kMaxParGeom + ig;\r
1082 deriv[idig] = h2*(4*d2 - d0)/3.;\r
1083 if (TMath::Abs(deriv[idig]) < 1.0e-9) deriv[idig] = 0.0;\r
1084 }\r
1085 }\r
1086 //\r
1087}\r
1088\r
1089//________________________________________________________________________________________________________\r
1090void AliITSAlignMille2Module::CalcDerivGloLoc(Int_t sensVol,Int_t paridx,Double_t* derivative)\r
1091{\r
1092 /// calculate numerically the derivatives of global params vs local param paridx for sensor sensVol: dPglob/dPloc_paridx\r
1093 //\r
1094 Double_t lpar[kMaxParGeom];\r
1095 for (int i=kMaxParGeom;i--;) lpar[i] = 0.;\r
1096 // using f(x+h),f(x-h),f(x+h/2),f(x-h2)...\r
1097 Double_t par1[kMaxParGeom]; // f(x-h)\r
1098 Double_t par2[kMaxParGeom]; // f(x-h/2)\r
1099 Double_t par3[kMaxParGeom]; // f(x+h/2)\r
1100 Double_t par4[kMaxParGeom]; // f(x+h)\r
1101 //\r
1102 const Double_t dpar = 1e-3;\r
1103 //\r
1104 // first values\r
1105 lpar[paridx] -= dpar;\r
1106 GetSensVolGlobalParams(sensVol,lpar,lpar+3, par1,par1+3);\r
1107 //\r
1108 // second values\r
1109 lpar[paridx] += dpar/2;\r
1110 GetSensVolGlobalParams(sensVol,lpar,lpar+3, par2,par2+3);\r
1111 //\r
1112 // third values\r
1113 lpar[paridx] += dpar;\r
1114 GetSensVolGlobalParams(sensVol,lpar,lpar+3, par3,par3+3);\r
1115 //\r
1116 // fourth values\r
1117 lpar[paridx] += dpar/2;\r
1118 GetSensVolGlobalParams(sensVol,lpar,lpar+3, par4,par4+3);\r
1119 //\r
1120 Double_t h2 = 1./(2.*dpar);\r
1121 for (int i=kMaxParGeom;i--;) {\r
1122 Double_t d0 = par4[i]-par1[i];\r
1123 Double_t d2 = 2.*(par3[i]-par2[i]);\r
1124 derivative[i] = h2*(4*d2 - d0)/3.;\r
1125 if (TMath::Abs(derivative[i]) < 1.0e-9) derivative[i] = 0.0;\r
1126 }\r
1127 //\r
1128}\r
1129\r
1130//________________________________________________________________________________________________________\r
1131void AliITSAlignMille2Module::CalcDerivCurLoc(Int_t sensVol,Int_t paridx,Double_t* derivative) \r
1132{\r
1133 /// calculate numerically the derivatives of sensor params in the current volume vs sensor local param paridx\r
1134 //\r
1135 Double_t lpar[kMaxParGeom];\r
1136 for (int i=kMaxParGeom;i--;) lpar[i] = 0.;\r
1137 // using f(x+h),f(x-h),f(x+h/2),f(x-h2)...\r
1138 Double_t par1[kMaxParGeom]; // f(x-h)\r
1139 Double_t par2[kMaxParGeom]; // f(x-h/2)\r
1140 Double_t par3[kMaxParGeom]; // f(x+h/2)\r
1141 Double_t par4[kMaxParGeom]; // f(x+h)\r
1142 //\r
1143 const Double_t dpar = 1e-3;\r
1144 //\r
1145 // first values\r
1146 lpar[paridx] -= dpar;\r
1147 GetSensVolLocalParams(sensVol,lpar,lpar+3, par1,par1+3);\r
1148 //\r
1149 // second values\r
1150 lpar[paridx] += dpar/2;\r
1151 GetSensVolLocalParams(sensVol,lpar,lpar+3, par2,par2+3);\r
1152 //\r
1153 // third values\r
1154 lpar[paridx] += dpar;\r
1155 GetSensVolLocalParams(sensVol,lpar,lpar+3, par3,par3+3);\r
1156 //\r
1157 // fourth values\r
1158 lpar[paridx] += dpar/2;\r
1159 GetSensVolLocalParams(sensVol,lpar,lpar+3, par4,par4+3);\r
1160 //\r
1161 Double_t h2 = 1./(2.*dpar);\r
1162 for (int i=kMaxParGeom;i--;) {\r
1163 Double_t d0 = par4[i]-par1[i];\r
1164 Double_t d2 = 2.*(par3[i]-par2[i]);\r
1165 derivative[i] = h2*(4*d2 - d0)/3.;\r
1166 if (TMath::Abs(derivative[i]) < 1.0e-9) derivative[i] = 0.0;\r
1167 }\r
1168 //\r
1169}\r
1170\r
1171\r
1172//-------------------------------------------------------------\r
1173void AliITSAlignMille2Module::GetGlobalParams(Double_t *t, Double_t *r)\r
1174{\r
1175 // global parameters of the module\r
1176 fgTempAlignObj.SetMatrix( *fMatrix );\r
1177 fgTempAlignObj.GetPars(t,r);\r
1178}\r
1179\r
1180//-------------------------------------------------------------\r
1181void AliITSAlignMille2Module::GetGlobalParams(const Double_t* loct, const Double_t* locr, Double_t *t, Double_t *r)\r
1182{\r
1183 // global parameters of the module after the modification by local loct,locr\r
1184 fgTempAlignObj.SetTranslation(loct[0],loct[1],loct[2]);\r
1185 fgTempAlignObj.SetRotation(locr[0],locr[1],locr[2]);\r
1186 fgTempAlignObj.GetMatrix(*fSensVolModifMatrix); \r
1187 *fSensVolMatrix = *fMatrix;\r
1188 fSensVolMatrix->Multiply(fSensVolModifMatrix);\r
1189 fgTempAlignObj.SetMatrix(*fSensVolMatrix);\r
1190 fgTempAlignObj.GetPars(t,r);\r
1191}\r
1192\r
1193//-------------------------------------------------------------\r
1194void AliITSAlignMille2Module::GetLocalParams(const Double_t* glot, const Double_t* glor, Double_t *t, Double_t *r)\r
1195{\r
1196 // obtain local delta parameters from global delta params\r
1197 fgTempAlignObj.SetTranslation(glot[0],glot[1],glot[2]);\r
1198 fgTempAlignObj.SetRotation(glor[0],glor[1],glor[2]);\r
1199 fgTempAlignObj.GetMatrix(*fSensVolMatrix); \r
1200 fSensVolMatrix->Multiply( fMatrix );\r
1201 fSensVolMatrix->MultiplyLeft( &fMatrix->Inverse() );\r
1202 fgTempAlignObj.SetMatrix(*fSensVolMatrix);\r
1203 fgTempAlignObj.GetPars(t,r);\r
1204}\r