ATLAS Offline Software
Loading...
Searching...
No Matches
InDetAlignCog.cxx
Go to the documentation of this file.
1/*
2 Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
3*/
4
11
12
13#include "InDetAlignGenAlgs/InDetAlignCog.h"
14#include "InDetReadoutGeometry/SiDetectorElement.h"
15#include "Identifier/Identifier.h"
16#include "Identifier/IdentifierHash.h"
17#include "InDetIdentifier/PixelID.h"
18#include "InDetIdentifier/SCT_ID.h"
19#include "InDetIdentifier/TRT_ID.h"
20#include "TRT_ReadoutGeometry/TRT_BaseElement.h"
21
22#include "AthenaKernel/IOVTime.h"
23#include "AthenaBaseComps/AthCheckMacros.h"
24
25#include "GeoPrimitives/CLHEPtoEigenConverter.h"
26
27#include "StoreGate/ReadCondHandleKey.h"
28
29#include <limits>
30
31namespace {
32 inline
33 Amg::AngleAxis3D _rotationX(double angle) {
34 return Amg::AngleAxis3D(angle, Amg::Vector3D(1., 0.,0.));
35 }
36
37 inline
38 Amg::AngleAxis3D _rotationY(double angle) {
39 return Amg::AngleAxis3D(angle, Amg::Vector3D(0.,1.,0.));
40 }
41
42 inline
43 Amg::AngleAxis3D _rotationZ(double angle) {
44 return Amg::AngleAxis3D(angle, Amg::Vector3D(0.,0.,1.));
45 }
46
47 inline Amg::Transform3D makeAffine3d(double angle_x, double angle_y, double angle_z, const Amg::Vector3D &translation_vector) {
48 Amg::Translation3D temp_translation(translation_vector);
49 Amg::Transform3D affine = temp_translation * _rotationX(angle_x);
50 affine *= _rotationY(angle_y);
51 affine *= _rotationZ(angle_z);
52 return affine;
53 }
54 inline Amg::Transform3D makeAffine3d(const AmgVector(6) &trans) {
55 return makeAffine3d(trans[3],trans[4],trans[5],Amg::Vector3D(trans[0],trans[1],trans[2]));
56 }
57 typedef double DoubleArray6_t[6];
58 inline Amg::Transform3D makeAffine3d(const DoubleArray6_t &trans) {
59 return makeAffine3d(trans[3],trans[4],trans[5],Amg::Vector3D(trans[0],trans[1],trans[2]));
60 }
61
62
63 template <class T_Matrix, class T_Vector>
64 bool solve(T_Matrix &A, T_Vector &x, const T_Vector &b) {
65 Eigen::LDLT< T_Matrix > llt;
66 // Eigen::LLT< T_Matrix > llt;
67 llt.compute(A);
68 if (llt.info() != Eigen::Success) return false;
69 x=b;
70 llt.solveInPlace(x);
71 return true;
72 }
73
74
75}
76
77#include <iostream>
78#include <iomanip>
79
80static const double tenmu = 0.010;
81static const double onemrad = 0.001;
82
83InDetAlignCog::InDetAlignCog(const std::string& name, ISvcLocator* pSvcLocator)
84 : AthAlgorithm(name, pSvcLocator),
85 m_pixid(nullptr),
86 m_sctid(nullptr),
87 m_trtid(nullptr),
88 m_IDAlignDBTool("InDetAlignDBTool"),
89 m_TRTAlignDbTool("TRT_AlignDbSvc",name),
90 //m_TRTAlignDbTool("TRTAlignDbTool"),
91 m_firstEvent(true),
92 //m_useChi2(true),
93 m_sigXpixB(0.1),
94 m_sigYpixB(0.1),
95 m_sigZpixB(0.1),
96 m_sigXpixE(0.1),
97 m_sigYpixE(0.1),
98 m_sigZpixE(0.1),
99 m_sigXsctB(0.1),
100 m_sigYsctB(0.1),
101 m_sigZsctB(0.1),
102 m_sigXsctE(0.1),
103 m_sigYsctE(0.1),
104 m_sigZsctE(0.1),
105 m_sigXtrtB(0.1),
106 m_sigYtrtB(0.1),
107 m_sigZtrtB(0.1),
108 m_sigXtrtE(0.1),
109 m_sigYtrtE(0.1),
110 m_sigZtrtE(0.1),
111 m_doTX(true),
112 m_doTY(true),
113 m_doTZ(true),
114 m_doRX(true),
115 m_doRY(true),
116 m_doRZ(true),
117 m_doCoG(true),
118 m_doL1(false),
119 m_traX(0.0),
120 m_traY(0.0),
121 m_traZ(0.0),
122 m_rotX(0.0),
123 m_rotY(0.0),
124 m_rotZ(0.0),
125 m_counter{}
126{
127
128
129 declareProperty("Det",m_det=99);// Pixel, SCT and TRT
130 declareProperty("SiBec",m_Si_bec=99);// silicon barrel and end-caps
131 declareProperty("TRT_Bec",m_TRT_bec=99);// TRT barrel and end-caps
132 declareProperty("SiLayer",m_Si_layer=99);// all silicon layers
133 declareProperty("TRT_Layer",m_TRT_layer=99);// all TRT layers
134
135 declareProperty("PrintFullMatrix", m_fullMatrix=false);
136 declareProperty("PrintDB", m_printDB=false);
137 declareProperty("ErrorTranslation",m_errTrans=1e-3); // one micron
138 declareProperty("ErrorRotation",m_errRot=1e-6); // one microrad
139 declareProperty("InDetAlignDBTool",m_IDAlignDBTool);
140 declareProperty("TRTAlignDbTool",m_TRTAlignDbTool);
141
142 declareProperty("SiTextOutput",m_SiTxtOutput=false);
143 declareProperty("TRT_TextOutput",m_TRT_TxtOutput=false);
144 declareProperty("SiTextFile",m_sitxtfile="cog_si_file.txt");
145 declareProperty("TRT_TextFile",m_trt_txtfile="cog_trt_file.txt");
146
147 declareProperty("SQLiteTag",m_SQLiteTag="cog_tag");
148
149 declareProperty("UseChi2",m_useChi2=true);
150 declareProperty("SXpixBarrel",m_sigXpixB=0.1);
151 declareProperty("SYpixBarrel",m_sigYpixB=0.1);
152 declareProperty("SZpixBarrel",m_sigZpixB=0.1);
153 declareProperty("SXpixEndcap",m_sigXpixE=0.1);
154 declareProperty("SYpixEndcap",m_sigYpixE=0.1);
155 declareProperty("SZpixEndcap",m_sigZpixE=0.1);
156 declareProperty("SXsctBarrel",m_sigXsctB=0.1);
157 declareProperty("SYsctBarrel",m_sigYsctB=0.1);
158 declareProperty("SZsctBarrel",m_sigZsctB=0.1);
159 declareProperty("SXsctEndcap",m_sigXsctE=0.1);
160 declareProperty("SYsctEndcap",m_sigYsctE=0.1);
161 declareProperty("SZsctEndcap",m_sigZsctE=0.1);
162 declareProperty("SXtrtBarrel",m_sigXtrtB=0.1);
163 declareProperty("SYtrtBarrel",m_sigYtrtB=0.1);
164 declareProperty("SZtrtBarrel",m_sigZtrtB=0.1);
165 declareProperty("SXtrtEndcap",m_sigXtrtE=0.1);
166 declareProperty("SYtrtEndcap",m_sigYtrtE=0.1);
167 declareProperty("SZtrtEndcap",m_sigZtrtE=0.1);
168
169 declareProperty("DoTX",m_doTX=true);
170 declareProperty("DoTY",m_doTY=true);
171 declareProperty("DoTZ",m_doTZ=true);
172 declareProperty("DoRX",m_doRX=true);
173 declareProperty("DoRY",m_doRY=true);
174 declareProperty("DoRZ",m_doRZ=true);
175
176 declareProperty("DoCoG",m_doCoG=true);
177 declareProperty("DoL1",m_doL1=false);
178 declareProperty("TraX",m_traX=0.0);
179 declareProperty("TraY",m_traY=0.0);
180 declareProperty("TraZ",m_traZ=0.0);
181 declareProperty("RotX",m_rotX=0.0);
182 declareProperty("RotY",m_rotY=0.0);
183 declareProperty("RotZ",m_rotZ=0.0);
184
185}
186
187//===================================================
188// initialize
189//===================================================
190StatusCode InDetAlignCog::initialize(){
191
192 ATH_MSG_DEBUG( "initialize()" );
193
194 // Get DetectorStore service
195 ATH_CHECK( detStore().retrieve() );
196
197 // get Pixel helper
198 ATH_CHECK( detStore()->retrieve(m_pixid));
199
200 // get SCT helper
201 ATH_CHECK( detStore()->retrieve(m_sctid));
202
203 // get TRT helper
204 ATH_CHECK( detStore()->retrieve(m_trtid));
205
206 // Get InDetAlignDBTool
207 ATH_CHECK( m_IDAlignDBTool.retrieve() );
208 ATH_MSG_DEBUG ( "Retrieved tool " << m_IDAlignDBTool );
209
210 // Get TRTAlignDBTool
211 ATH_CHECK( m_TRTAlignDbTool.retrieve() );
212
213 // ReadCondHandleKey
214 ATH_CHECK(m_trtDetEleContKey.initialize(m_det==99 || m_det==3));
215 ATH_CHECK(m_pixelDetEleCollKey.initialize(m_det==99 || m_det==1 || m_det==12));
216 ATH_CHECK(m_SCTDetEleCollKey.initialize(m_det==99 || m_det==2 || m_det==12));
217
218 ATH_MSG_DEBUG ( "Retrieved tool " << m_TRTAlignDbTool );
219
220 return StatusCode::SUCCESS;
221}
222
223
224//===================================================
225// execute
226//===================================================
227StatusCode InDetAlignCog::execute() {
228 ATH_MSG_DEBUG( "execute()" );
229
230 const InDetDD::SiDetectorElementCollection* pixelElements(nullptr);
231 if (m_det==99 || m_det==1 || m_det==12) {
232 SG::ReadCondHandle<InDetDD::SiDetectorElementCollection> pixelDetEleHandle(m_pixelDetEleCollKey);
233 pixelElements = *pixelDetEleHandle;
234 if (not pixelDetEleHandle.isValid() or pixelElements==nullptr) {
235 ATH_MSG_FATAL(m_pixelDetEleCollKey.fullKey() << " is not available.");
236 return StatusCode::FAILURE;
237 }
238 }
239
240 const InDetDD::SiDetectorElementCollection* sctElements(nullptr);
241 if (m_det==99 || m_det==2 || m_det==12) {
242 SG::ReadCondHandle<InDetDD::SiDetectorElementCollection> sctDetEleHandle(m_SCTDetEleCollKey);
243 sctElements = *sctDetEleHandle;
244 if (not sctDetEleHandle.isValid() or sctElements==nullptr) {
245 ATH_MSG_FATAL(m_SCTDetEleCollKey.fullKey() << " is not available.");
246 return StatusCode::FAILURE;
247 }
248 }
249
250 const InDetDD::TRT_DetElementCollection* trtElements(nullptr);
251 if (m_det==99 || m_det==3) {
252 SG::ReadCondHandle<InDetDD::TRT_DetElementContainer> trtDetEleHandle(m_trtDetEleContKey);
253 trtElements = trtDetEleHandle->getElements();
254 if (not trtDetEleHandle.isValid() or trtElements==nullptr) {
255 ATH_MSG_FATAL(m_trtDetEleContKey.fullKey() << " is not available.");
256 return StatusCode::FAILURE;
257 }
258 }
259
260 if (m_firstEvent) {
261 m_firstEvent = false;
262 m_counter = 0;
263
264 if(m_printDB){
265 m_IDAlignDBTool->printDB(2);
266 // m_TRTAlignDbTool->printCondObjects(); // displacements are printed in microns...
267 }
268
269 // create vectors and matrices Vi,Ai,Mi :
270 Params_t params;
271
272 // initialize transforms to identity
273 for( int i=0; i<6; i++){ m_cog[i]=0.0;}
274 for( int i=0; i<6; i++){ m_resglob[i]=0.0;}
275 m_CoG.setIdentity();
276 m_ResGlob.setIdentity();
277
278 // first loop to calculate cog
279 //StatusCode sc;
280 if(m_det==99 || m_det==1 || m_det==12) ATH_CHECK( getSiElements(pixelElements,false,params) );
281 if(m_det==99 || m_det==2 || m_det==12) ATH_CHECK( getSiElements(sctElements,false,params) );
282 if(m_det==99 || m_det==3) ATH_CHECK( getTRT_Elements(trtElements,false, params) );
283 //if(sc.isFailure())
284 // ATH_MSG_ERROR( "Problem getting elements from managers" );
285 if( !m_useChi2 ) {
286
287
288 // normalization of m_cog
289 if (m_counter !=0.){
290 for( int i=0; i<6; i++){m_cog[i]/=(double) m_counter;}
291 } else {
292 ATH_MSG_WARNING("m_counter is zero, giving undefined behaviour");
293 for( int i=0; i<6; i++){m_cog[i]=std::numeric_limits<double>::infinity();}
294 }
295
296
297 // convert to HepGeom::Transform3D:
298 m_CoG = makeAffine3d(m_cog);
299
300 // scaling by -1 translation and rotation angles
301 scaleTransform(m_CoG,-1.);
302
303 // show results
304 ATH_MSG_DEBUG( "+++++++++++++++++++++++++++++++++++++++++++++" );
305 ATH_MSG_DEBUG( "Processed " << m_counter << " elements" );
306 ATH_MSG_DEBUG( "Center-of-gravity : " << printTransform(m_CoG) );
307 ATH_MSG_DEBUG( "+++++++++++++++++++++++++++++++++++++++++++++" );
308
309 } else {
310
311 ATH_MSG_INFO( "+++++++++++++++++++++++++++++++++++++++++++++" );
312 ATH_MSG_INFO( "Processed " << m_counter << " elements" );
313 ATH_MSG_INFO( "+++++++++++++++++++++++++++++++++++++++++++++" );
314
315 // solve for the CoG correction:
316 // solve M1 * A1 = V1 for A1
317 bool success = solve(params.m_M1,params.m_A1, params.m_V1 );
318 if( success ) {
319 // (params.m_A1) = (params.m_M1)*(params.m_V1);
320 m_CoG = makeAffine3d(params.m_A1);
321
322 ATH_MSG_INFO( "+++++++++++++++++++++++++++++++++++++++++++++" );
323 ATH_MSG_INFO( "Inversion of matrix M1 successful." );
324 ATH_MSG_INFO( "Residual global transformation : " << printTransform(m_CoG) );
325 ATH_MSG_INFO( "+++++++++++++++++++++++++++++++++++++++++++++" );
326 } else {
327 ATH_MSG_ERROR( "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
328 ATH_MSG_ERROR( "Inversion of matrix M1 failed!!!");
329 ATH_MSG_ERROR( "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
330 }
331 }
332
333 /*
334 Below is cross-check: calculate residual global transformation
335 after cog substraction
336 */
337
338 // re-initialize normalization factor
339 m_counter=0;
340
341 // second loop to compute residual transform after substracting cog
342 if(m_det==99 || m_det==1 || m_det==12) ATH_CHECK( getSiElements(pixelElements,true, params) );
343 if(m_det==99 || m_det==2 || m_det==12) ATH_CHECK( getSiElements(sctElements,true, params) );
344 if(m_det==99 || m_det==3) ATH_CHECK( getTRT_Elements(trtElements,true, params) );
345 // if(sc.isFailure())
346 // ATH_MSG_ERROR( "Problem getting elements from managers" );
347 if( !m_useChi2 ) {
348 if (m_counter==0) throw std::logic_error("No Si-elements.");
349 // normalization of m_resglob
350 for( int i=0; i<6; i++){m_resglob[i]/=(double) m_counter;}
351
352 // convert to HepGeom::Transform3D:
353 m_ResGlob = makeAffine3d(m_resglob);
354
355
356 // show results
357 ATH_MSG_DEBUG( "+++++++++++++++++++++++++++++++++++++++++++++" );
358 ATH_MSG_DEBUG( "Processed " << m_counter << " elements" );
359 ATH_MSG_DEBUG( "Residual global transformation : " << printTransform(m_ResGlob) );
360 ATH_MSG_DEBUG( "+++++++++++++++++++++++++++++++++++++++++++++" );
361
362 } else {
363
364 ATH_MSG_INFO( "+++++++++++++++++++++++++++++++++++++++++++++" );
365 ATH_MSG_INFO( "Processed " << m_counter << " elements" );
366 ATH_MSG_INFO( "+++++++++++++++++++++++++++++++++++++++++++++" );
367
368 // solve for the CoG correction:
369 // solve M2 * A2 = V2 for A2
370 bool success = solve(params.m_M2,params.m_A2, params.m_V2 );
371 if( success ) {
372 m_ResGlob = makeAffine3d(params.m_A2);
373 ATH_MSG_INFO( "+++++++++++++++++++++++++++++++++++++++++++++" );
374 ATH_MSG_INFO( "Inversion of matrix M2 successful." );
375 ATH_MSG_INFO( "Residual global transformation : " << printTransform(m_ResGlob));
376 ATH_MSG_INFO( "+++++++++++++++++++++++++++++++++++++++++++++" );
377 } else {
378 ATH_MSG_ERROR( "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" );
379 ATH_MSG_ERROR( "Inversion of matrix M2 failed!!!" );
380 ATH_MSG_ERROR( "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" );
381 }
382 }
383
384
385
386 // check if remaining global transformation is small enough
387 // so that ID can be finally shifted
388 if(testIdentity(m_ResGlob,m_errRot,m_errTrans) || !m_doCoG)
389 {
390 // Decide which CoG DoF's are to be actually corrected:
391 if( m_doCoG ) {
392 ATH_MSG_INFO( "<<< CoG correction will be applied. >>>" );
393 enableCoG(m_CoG, m_doTX, m_doTY, m_doTZ, m_doRX, m_doRY, m_doRZ);
394 } else {
395 ATH_MSG_INFO( "<<< CoG correction will be skipped. >>>" );
396 m_CoG.setIdentity();
397 }
398 // Decide if an arbitrary L1 is to be done:
399 if( m_doL1 ) {
400 ATH_MSG_INFO( "<<< Predefined L1 transformation will be added. >>>" );
401 addL1();
402 }
403 // Update the CoG:
404 if(shiftIDbyCog().isFailure()){
405 ATH_MSG_ERROR( "Problem updating ID level 1 constants !!!" );
406 return StatusCode::FAILURE;
407 }
408 }
409
410 }
411 return StatusCode::SUCCESS;
412}
413
414//===================================================
415// finalize
416//===================================================
417StatusCode InDetAlignCog::finalize() {
418 ATH_MSG_DEBUG ( "finalize()" );
419
420 // output new si constants to root file
421 if(m_IDAlignDBTool->outputObjs().isFailure()){
422 ATH_MSG_ERROR( "Write of silicon alignableTransforms fails" );
423 return StatusCode::FAILURE;
424 }
425
426 // Save a text file with Si new constants
427 if(m_SiTxtOutput) m_IDAlignDBTool->writeFile(false,m_sitxtfile);
428
429 // Fill the SQLite file (silicon)
430 m_IDAlignDBTool->fillDB(m_SQLiteTag,
431 IOVTime::MINRUN,
432 IOVTime::MINEVENT,
433 IOVTime::MAXRUN,
434 IOVTime::MAXEVENT);
435
436 // output new TRT constants to root file
437 if(m_TRTAlignDbTool->streamOutAlignObjects().isFailure()){
438 ATH_MSG_ERROR( "Write of TRT alignableTransforms fails" );
439 return StatusCode::FAILURE;
440 }
441
442 // Save a text file with TRT new constants
443 if(m_TRT_TxtOutput)
444 if(m_TRTAlignDbTool->writeAlignTextFile(m_trt_txtfile).isFailure())
445 ATH_MSG_ERROR( "Write of TRT new constants txt file fails" );
446
447 // Fill the SQLite file (TRT)
448 if(m_TRTAlignDbTool->registerAlignObjects(m_SQLiteTag+"_TRT",
449 IOVTime::MINRUN,
450 IOVTime::MINEVENT,
451 IOVTime::MAXRUN,
452 IOVTime::MAXEVENT).isFailure())
453 ATH_MSG_ERROR( "Write of TRT new constants db file fails" );
454 return StatusCode::SUCCESS;
455}
456
457//===================================================
458// getSiElements
459//===================================================
460StatusCode InDetAlignCog::getSiElements(const InDetDD::SiDetectorElementCollection *elements,
461 bool cog_already_calculated,
462 InDetAlignCog::Params_t &params
463 ){
464
465 for (const InDetDD::SiDetectorElement *element: *elements) {
466 // @TODO can element be null ?
467 if (element) {
468 Identifier id = element->identify();
469 int det = 0;
470 int bec = 0;
471 int layer_disk = 0;
472 int phi_module = 0;
473 int eta_module = 0;
474 int side = 0;
475
476 if(element->isPixel()) {
477 det = 1;
478 bec = m_pixid->barrel_ec(id);
479 layer_disk = m_pixid->layer_disk(id);
480 phi_module = m_pixid->phi_module(id);
481 eta_module = m_pixid->eta_module(id);
482 side = 0;
483 } else { // SCT
484 det = 2;
485 bec = m_sctid->barrel_ec(id);
486 layer_disk = m_sctid->layer_disk(id);
487 phi_module = m_sctid->phi_module(id);
488 eta_module = m_sctid->eta_module(id);
489 side = m_sctid->side(id);
490 }
491
493 double sigmas[6]={0.1,0.1,0.1,1.,1.,1.};
494 if( det==1 ) { // Pixel
495 if( bec==0 ) { // barrel
496 sigmas[0]=m_sigXpixB; sigmas[1]=m_sigYpixB; sigmas[2]=m_sigZpixB;
497 } else { // EC's
498 sigmas[0]=m_sigXpixE; sigmas[1]=m_sigYpixE; sigmas[2]=m_sigZpixE;
499 }
500 } else { // SCT
501 if( bec==0 ) { // barrel
502 sigmas[0]=m_sigXsctB; sigmas[1]=m_sigYsctB; sigmas[2]=m_sigZsctB;
503 } else { // EC's
504 sigmas[0]=m_sigXsctE; sigmas[1]=m_sigYsctE; sigmas[2]=m_sigZsctE;
505 }
506 }
507
508 ATH_MSG_DEBUG( cog_already_calculated << " " << m_counter );
509
510 // perform selections
511 if((m_det==1 || m_det==2) && det!=m_det) continue;
512 if(m_Si_bec!=99 && bec!=m_Si_bec) continue;
513 if(m_Si_layer!=99 && layer_disk!=m_Si_layer) continue;
514
515 if(element->isPixel())
516 ATH_MSG_VERBOSE(std::setw(4)
517 << m_counter << " Module " << m_pixid->show_to_string(id));
518 else
519 ATH_MSG_VERBOSE(std::setw(4)
520 << m_counter << " Module " << m_sctid->show_to_string(id));
521
522 m_counter++;
523
524 // default Local-to-global transform (i.e. without misalignments)
525 const Amg::Transform3D defTransform(element->defModuleTransform());
526 ATH_MSG_VERBOSE( "defTransform " << 2 << " " << det << " " << bec << " " << layer_disk << " "
527 << phi_module << " " << eta_module << " " << side << " "
528 << printTransform(defTransform) );
529
530 // Local-to-global module transform
531 const Amg::Transform3D transform(element->moduleTransform());
532 ATH_MSG_VERBOSE("transform " << 2 << " " << det << " " << bec << " " << layer_disk << " "
533 << phi_module << " " << eta_module << " " << side << " "
534 << printTransform(transform) );
535
536 // corrections in local frame
537 const Amg::Transform3D localDelta = defTransform.inverse()*transform;
538
539 ATH_MSG_VERBOSE("localDelta " << 2 << " " << det << " " << bec << " " << layer_disk << " "
540 << phi_module << " " << eta_module << " " << side << " "
541 << printTransform(localDelta));
542
543 prepareDerivative(defTransform, m_useChi2);
544
545 if(!cog_already_calculated){
546 m_useChi2 ? accumulateChi2(localDelta,params.m_M1,params.m_V1,sigmas) : accumulate(localDelta,m_cog);
547 }
548 else{
549 // coglocalDelta : local transform corresponding to global cog transform
550 // recall that m_CoG already scaled by -1 after first pass
551 if((m_det==1 || m_det==2) && det!=m_det) continue;
552 if(m_Si_bec!=99 && bec!=m_Si_bec) continue;
553 if(m_Si_layer!=99 && layer_disk!=m_Si_layer) continue;
554 Amg::Transform3D coglocalDelta = defTransform.inverse() * m_CoG * defTransform;
555
556 // newlocalDelta : localDelta (alignment corrections) + coglocalDelta
557 Amg::Transform3D newlocalDelta = sumTransforms(localDelta, coglocalDelta);
558
559 m_useChi2 ? accumulateChi2(newlocalDelta,params.m_M2,params.m_V2,sigmas) : accumulate(newlocalDelta,m_resglob);
560 }
561 }
562 }
563 return StatusCode::SUCCESS;
564}
565
566//===================================================
567// getTRT_Elements
568//
569// From TRTAlignDbTool.h
570//-----------------------
571// TRT alignment transforms are built up like this:
572//
573// aligned = level-2 * level-1 * nominal
574//
575// The 'level-2' correction is the global transform.
576// There is a separate transform for each 'bec'.
577// However, note that the barrel is bec=-1: the +1
578// is ignored in the geometry.
579//
580// The 'level-1' transform is a transform per module.
581//===================================================
582StatusCode InDetAlignCog::getTRT_Elements(const InDetDD::TRT_DetElementCollection *elements,bool cog_already_calculated, InDetAlignCog::Params_t &params){
583 ATH_MSG_DEBUG( "in getTRT_Elements " );
584
585
586 for (const InDetDD::TRT_BaseElement *element: *elements) {
587 // @TODO can element be null ?
588 if (element) {
589 Identifier id = element->identify();
590
591 int bec = m_trtid->barrel_ec(id);
592 int phi_module = m_trtid->phi_module(id);
593 int layer_wheel = m_trtid->layer_or_wheel(id);
594 int straw_layer = m_trtid->straw_layer(id);
595
596 // perform selections
597 if(m_TRT_bec!=99 && bec!=m_TRT_bec) continue;
598 if(m_TRT_layer!=99 && layer_wheel!=m_TRT_layer) continue;
599
600 // Skip A Side of
601 if(bec == 1) continue;
602
603 ATH_MSG_VERBOSE(std::setw(4)
604 << m_counter << " Module " << m_trtid->show_to_string(id) );
605
606 m_counter++;
607
609 double sigmas[6]={0.1,0.1,0.1,1.,1.,1.};
610 if( bec==0 ) { // barrel
611 sigmas[0]=m_sigXtrtB; sigmas[1]=m_sigYtrtB; sigmas[2]=m_sigZtrtB;
612 } else { // EC's
613 sigmas[0]=m_sigXtrtE; sigmas[1]=m_sigYtrtE; sigmas[2]=m_sigZtrtE;
614 }
615
616
617 // Get Default Transform (of module in barrel, layer in endcap)
618 const Amg::Transform3D defTransform = element->defTransform();
619 ATH_MSG_VERBOSE("defTransform 2 3" << " " << bec << " " << phi_module << " "
620 << layer_wheel << " " << straw_layer << " " << printTransform(defTransform));
621
622 // retrieve level 1 transform (module level)
623 const Amg::Transform3D t1 = getL1Transform(bec);
624
625 // const HepGeom::Transform3D t1 = m_TRTAlignDbTool->getTrans(id);
626 ATH_MSG_VERBOSE( "t1 2 3" << " " << bec << " " << phi_module << " "
627 << layer_wheel << " " << straw_layer << " " << printTransform(t1));
628
629 // retrieve level 2 transform (barrel/endcaps)
630 const Amg::Transform3D t2 = m_TRTAlignDbTool->getAlignmentTransform(id,2);
631 // const HepGeom::Transform3D t2 = m_TRTAlignDbTool->getGlobalTrans(bec);
632 ATH_MSG_VERBOSE( "t2 2 3" << " " << bec << " " << phi_module << " "
633 << layer_wheel << " " << straw_layer << " " << printTransform(t2));
634
635 // build up globalDelta, i.e, the combination in GLOB of t2 and t1
636 const Amg::Transform3D globalDelta = (t2)*(t1);
637 ATH_MSG_VERBOSE("globalDelta 2 3" << " " << bec << " " << phi_module << " "
638 << layer_wheel << " " << straw_layer << " " << printTransform(globalDelta));
639
640 // equivalent 'local' transform (at the module level)
641 const Amg::Transform3D localDelta = defTransform.inverse() * globalDelta * defTransform;
642
643 ATH_MSG_VERBOSE("localDelta 2 3" << " " << bec << " " << phi_module << " "
644 << layer_wheel << " " << straw_layer << " " << printTransform(localDelta));
645
646 prepareDerivative(defTransform, m_useChi2);
647
648 if(!cog_already_calculated){
649 m_useChi2 ? accumulateChi2(localDelta,params.m_M1,params.m_V1,sigmas) : accumulate(localDelta,m_cog);
650 }
651 else{
652 // coglocalDelta : local transform corresponding to global cog transform
653 // recall that m_CoG already scaled by -1 after first pass
654 if(m_TRT_bec!=99 && bec!=m_TRT_bec) continue;
655 if(m_TRT_layer!=99 && layer_wheel!=m_TRT_layer) continue;
656
657 Amg::Transform3D coglocalDelta = defTransform.inverse() * m_CoG * defTransform;
658 Amg::Transform3D newlocalDelta = sumTransforms(localDelta, coglocalDelta);
659
660 m_useChi2 ? accumulateChi2(newlocalDelta,params.m_M2,params.m_V2,sigmas) : accumulate(newlocalDelta,m_resglob);
661 }
662 }
663 }
664 return StatusCode::SUCCESS;
665}
666
667//=====================================================================
668// get the correct L1 id for the COG
669//=====================================================================
670const Amg::Transform3D InDetAlignCog::getL1Transform(int bec){
671
672 Identifier thisL1Identifier = m_trtid->module_id(bec,0,0);
673
674 const Amg::Transform3D t1 = m_TRTAlignDbTool->getAlignmentTransform(thisL1Identifier,1);
675 return t1;
676}
677
678//=====================================================================
679// accumulate
680//=====================================================================
681void InDetAlignCog::accumulate(const Amg::Transform3D &trans,
682 double* val){
683
684 // ATH_MSG_DEBUG("in accumulate for transform " << printransform(trans));
685
686
687 double alpha[6], beta[6], gamma[6];
688 for(int i=0; i<6; i++){
689 alpha[i]=0;
690 beta[i]=0;
691 gamma[i]=0;
692 }
693
694 // extract parameters from input transform
695 double x,y,z,a,b,g;
696
697 const int x_i=0;
698 const int y_i=1;
699 const int z_i=2;
700
701 const Amg::Vector3D translation_part(trans.translation());
702 x = translation_part[x_i];
703 y = translation_part[y_i];
704 z = translation_part[z_i];
705
706 const Amg::Transform3D &transAMG = trans ;
707
708 m_IDAlignDBTool->extractAlphaBetaGamma(transAMG,a,b,g);
709
710 const Amg::Vector3D glob_x_trans = m_glob_x.translation();
711 const Amg::Vector3D glob_y_trans = m_glob_y.translation();
712 const Amg::Vector3D glob_z_trans = m_glob_z.translation();
713
714 const Amg::Vector3D grot_x_trans = m_grot_x.translation();
715 const Amg::Vector3D grot_y_trans = m_grot_y.translation();
716 const Amg::Vector3D grot_z_trans = m_grot_z.translation();
717
718 val[0] += glob_x_trans[x_i]/tenmu*x + glob_y_trans[x_i]/tenmu*y + glob_z_trans[x_i]/tenmu*z;
719 val[0] += grot_x_trans[x_i]/onemrad*a + grot_y_trans[x_i]/onemrad*b + grot_z_trans[x_i]/onemrad*g;
720
721 val[1] += glob_x_trans[y_i]/tenmu*x + glob_y_trans[y_i]/tenmu*y + glob_z_trans[y_i]/tenmu*z;
722 val[1] += grot_x_trans[y_i]/onemrad*a + grot_y_trans[y_i]/onemrad*b + grot_z_trans[y_i]/onemrad*g;
723
724 val[2] += glob_x_trans[z_i]/tenmu*x + glob_y_trans[z_i]/tenmu*y + glob_z_trans[z_i]/tenmu*z;
725 val[2] += grot_x_trans[z_i]/onemrad*a + grot_y_trans[z_i]/onemrad*b + grot_z_trans[z_i]/onemrad*g;
726
727 // extract alpha, beta and gamma from each of the 6 global transformation
728 m_IDAlignDBTool->extractAlphaBetaGamma(m_glob_x, alpha[0], beta[0], gamma[0]);
729 m_IDAlignDBTool->extractAlphaBetaGamma(m_glob_y, alpha[1], beta[1], gamma[1]);
730 m_IDAlignDBTool->extractAlphaBetaGamma(m_glob_z, alpha[2], beta[2], gamma[2]);
731 m_IDAlignDBTool->extractAlphaBetaGamma(m_grot_x, alpha[3], beta[3], gamma[3]);
732 m_IDAlignDBTool->extractAlphaBetaGamma(m_grot_y, alpha[4], beta[4], gamma[4]);
733 m_IDAlignDBTool->extractAlphaBetaGamma(m_grot_z, alpha[5], beta[5], gamma[5]);
734
735 val[3] += alpha[0]/tenmu*x + alpha[1]/tenmu*y + alpha[2]/tenmu*z;
736 val[3] += alpha[3]/onemrad*a + alpha[4]/onemrad*b + alpha[5]/onemrad*g;
737
738 val[4] += beta[0]/tenmu*x + beta[1]/tenmu*y + beta[2]/tenmu*z;
739 val[4] += beta[3]/onemrad*a + beta[4]/onemrad*b + beta[5]/onemrad*g;
740
741 val[5] += gamma[0]/tenmu*x + gamma[1]/tenmu*y + gamma[2]/tenmu*z;
742 val[5] += gamma[3]/onemrad*a + gamma[4]/onemrad*b + gamma[5]/onemrad*g;
743
744}
745
746//=====================================================================
747// accumulate Chi2 derivatives
748//=====================================================================
749void InDetAlignCog::accumulateChi2(const Amg::Transform3D &trans, AmgSymMatrix(6)& M, AmgVector(6)& V, const double* sigmas){
750
751 ATH_MSG_DEBUG("in accumulateChi2 for transform " << printTransform(trans));
752
753
755 double delta[6];
756
757 Amg::Vector3D translation_part(trans.translation());
758 delta[0] = translation_part[0];
759 delta[1] = translation_part[1];
760 delta[2] = translation_part[2];
761 m_IDAlignDBTool->extractAlphaBetaGamma(trans ,delta[3],delta[4],delta[5]);
762
764
765 AmgMatrix(6,6) jacobian;
766 // jacobian.setZero();
767
768 // extract the translations
769 const Amg::Vector3D glob_x_trans = m_glob_x.translation();
770 const Amg::Vector3D glob_y_trans = m_glob_y.translation();
771 const Amg::Vector3D glob_z_trans = m_glob_z.translation();
772
773 const Amg::Vector3D grot_x_trans = m_grot_x.translation();
774 const Amg::Vector3D grot_y_trans = m_grot_y.translation();
775 const Amg::Vector3D grot_z_trans = m_grot_z.translation();
776
777 // extract alpha, beta and gamma from each of the 6 global transformation
778 double alpha[6], beta[6], gamma[6];
779 m_IDAlignDBTool->extractAlphaBetaGamma(m_glob_x, alpha[0], beta[0], gamma[0]);
780 m_IDAlignDBTool->extractAlphaBetaGamma(m_glob_y, alpha[1], beta[1], gamma[1]);
781 m_IDAlignDBTool->extractAlphaBetaGamma(m_glob_z, alpha[2], beta[2], gamma[2]);
782 m_IDAlignDBTool->extractAlphaBetaGamma(m_grot_x, alpha[3], beta[3], gamma[3]);
783 m_IDAlignDBTool->extractAlphaBetaGamma(m_grot_y, alpha[4], beta[4], gamma[4]);
784 m_IDAlignDBTool->extractAlphaBetaGamma(m_grot_z, alpha[5], beta[5], gamma[5]);
785
786 // local translations:
787 jacobian(0,0) = glob_x_trans[0]/tenmu;
788 jacobian(0,1) = glob_y_trans[0]/tenmu;
789 jacobian(0,2) = glob_z_trans[0]/tenmu;
790 jacobian(0,3) = grot_x_trans[0]/onemrad;
791 jacobian(0,4) = grot_y_trans[0]/onemrad;
792 jacobian(0,5) = grot_z_trans[0]/onemrad;
793 jacobian(1,0) = glob_x_trans[1]/tenmu;
794 jacobian(1,1) = glob_y_trans[1]/tenmu;
795 jacobian(1,2) = glob_z_trans[1]/tenmu;
796 jacobian(1,3) = grot_x_trans[1]/onemrad;
797 jacobian(1,4) = grot_y_trans[1]/onemrad;
798 jacobian(1,5) = grot_z_trans[1]/onemrad;
799 jacobian(2,0) = glob_x_trans[2]/tenmu;
800 jacobian(2,1) = glob_y_trans[2]/tenmu;
801 jacobian(2,2) = glob_z_trans[2]/tenmu;
802 jacobian(2,3) = grot_x_trans[2]/onemrad;
803 jacobian(2,4) = grot_y_trans[2]/onemrad;
804 jacobian(2,5) = grot_z_trans[2]/onemrad;
805
806 // local rotations:
807 for (unsigned int i=3; i<6; ++i) {
808 for (unsigned int j=0; j<3; ++j) {
809 jacobian(i,j)=0;
810 }
811 }
812 jacobian(3,3) = alpha[3]/onemrad;
813 jacobian(3,4) = alpha[4]/onemrad;
814 jacobian(3,5) = alpha[5]/onemrad;
815 jacobian(4,3) = beta[3]/onemrad;
816 jacobian(4,4) = beta[4]/onemrad;
817 jacobian(4,5) = beta[5]/onemrad;
818 jacobian(5,3) = gamma[3]/onemrad;
819 jacobian(5,4) = gamma[4]/onemrad;
820 jacobian(5,5) = gamma[5]/onemrad;
821
822 // and the other six are null by construction.
823
825
826 for (int i=0; i<6; i++) {
827 for (int k=0; k<3; k++) {
828 V[i] -= 2./(sigmas[k]*sigmas[k])*jacobian(k,i)*delta[k];
829 for (int j=0; j<=i; j++) {
830 M(i,j) += 2./(sigmas[k]*sigmas[k])*jacobian(k,i)*jacobian(k,j);
831 }
832 }
833 }
834
835}
836
837//===================================================
838// prepareDerivative
839//===================================================
840void InDetAlignCog::prepareDerivative(const Amg::Transform3D &trans, const bool dLdG){
841
842 ATH_MSG_DEBUG("in prepareDerivative for transform " << printTransform(trans));
843
844 Amg::Transform3D id = Amg::Transform3D::Identity();
845
846 const Amg::Transform3D epsilon_x = Amg::Translation3D(tenmu,0,0) * id ;
847 const Amg::Transform3D epsilon_y = Amg::Translation3D(0,tenmu,0) * id ;
848 const Amg::Transform3D epsilon_z = Amg::Translation3D(0,0,tenmu) * id ;
849
850 const Amg::Transform3D epsilon_a = id * _rotationX(onemrad);
851 const Amg::Transform3D epsilon_b = id * _rotationY(onemrad);
852 const Amg::Transform3D epsilon_g = id * _rotationZ(onemrad);
853
854 Amg::Transform3D inv_trans = trans.inverse();
855 if( !dLdG ) {
857 m_glob_x = trans * epsilon_x * inv_trans;
858 m_glob_y = trans * epsilon_y * inv_trans;
859 m_glob_z = trans * epsilon_z * inv_trans;
860 m_grot_x = trans * epsilon_a * inv_trans;
861 m_grot_y = trans * epsilon_b * inv_trans;
862 m_grot_z = trans * epsilon_g * inv_trans;
863 } else {
865 m_glob_x = inv_trans * epsilon_x * trans;
866 m_glob_y = inv_trans * epsilon_y * trans;
867 m_glob_z = inv_trans * epsilon_z * trans;
868 m_grot_x = inv_trans * epsilon_a * trans;
869 m_grot_y = inv_trans * epsilon_b * trans;
870 m_grot_z = inv_trans * epsilon_g * trans;
871 }
872
873}
874
875//===================================================
876// shiftIDbyCog
877//===================================================
878StatusCode InDetAlignCog::shiftIDbyCog(){
879 ATH_MSG_DEBUG("in ShiftIDbyCog with det = " << m_det );
880
881 const Amg::Transform3D &cogAMG = m_CoG;
882
883 // update level 1 silicon constants
884 if(m_det==99 || m_det==12 || m_det==1 || m_det==2) {
885
886 ATH_MSG_DEBUG("Will try to update Pixel and SCT level 1 constants...");
887
888 const int level = 1;
889 if(!(m_IDAlignDBTool->tweakTrans(m_pixid->wafer_id(0,0,0,0),level,cogAMG)) || // whole Pixel
890 !(m_IDAlignDBTool->tweakTrans(m_sctid->wafer_id(-2,0,0,0,0),level,cogAMG)) || // SCT EC-C
891 !(m_IDAlignDBTool->tweakTrans(m_sctid->wafer_id(0,0,0,0,0),level,cogAMG)) || // SCT barrel
892 !(m_IDAlignDBTool->tweakTrans(m_sctid->wafer_id(2,0,0,0,0),level,cogAMG))){ // SCT EC-A
893 ATH_MSG_ERROR( "Could not update level 1 silicon constants !!" );
894 return StatusCode::FAILURE;
895 }
896 else{
897 ATH_MSG_DEBUG( "Successful update of level 1 silicon constants " );
898 }
899 }
900
901 // update level 2 TRT constants
902 if(m_det==99 || m_det==3){
903 ATH_MSG_DEBUG( "Will try to update TRT level 1 constants...");
904
905 if(//StatusCode::SUCCESS!=m_TRTAlignDbTool->tweakTrans2(m_trtid->module_id(-1,0,0),cogAMG) || // TRT barrel
906 StatusCode::SUCCESS!=m_TRTAlignDbTool->tweakAlignTransform(m_trtid->module_id(-1,0,0),cogAMG,1) || // TRT barrel
907 //StatusCode::SUCCESS!=m_TRTAlignDbTool->tweakTrans2(m_trtid->module_id(2,0,0),cogAMG) || // TRT pos. end-cap
908 StatusCode::SUCCESS!=m_TRTAlignDbTool->tweakAlignTransform(m_trtid->module_id(2,0,0),cogAMG,1) || // TRT pos. end-cap
909 //StatusCode::SUCCESS!=m_TRTAlignDbTool->tweakTrans2(m_trtid->module_id(-2,0,0),cogAMG)){ // TRT neg. end-cap
910 StatusCode::SUCCESS!=m_TRTAlignDbTool->tweakAlignTransform(m_trtid->module_id(-2,0,0),cogAMG,1)){ // TRT neg. end-cap
911 ATH_MSG_ERROR("Could not update level 1 TRT constants !!" );
912 return StatusCode::FAILURE;
913 }
914 else{
915 ATH_MSG_DEBUG( "Successful update of level 1 TRT constants " );
916 }
917 }
918 return StatusCode::SUCCESS;
919}
920
921//=====================================================================
922// enableCog
923//=====================================================================
924void InDetAlignCog::enableCoG(Amg::Transform3D & trans,
925 bool dotx, bool doty, bool dotz, bool dorx, bool dory, bool dorz){
926 ATH_MSG_DEBUG("in enableCoG with decisions " << dotx << doty << dotz << dorx << dory << dorz );
927
928
929 Amg::Vector3D vec = trans.translation();
930 if( !dotx ) vec[0]=0.0;
931 if( !doty ) vec[1]=0.0;
932 if( !dotz ) vec[2]=0.0;
933
934 double a,b,g;
935 m_IDAlignDBTool->extractAlphaBetaGamma(trans,a,b,g);
936 if( !dorx ) a=0.0;
937 if( !dory ) b=0.0;
938 if( !dorz ) g=0.0;
939
940 trans = makeAffine3d( a, b, g, vec);
941}
942
943//=====================================================================
944// add L1 transformation
945//=====================================================================
946void InDetAlignCog::addL1(){
947 ATH_MSG_DEBUG("in addL1... " );
948
949
950 Amg::Transform3D L1transform(makeAffine3d(m_rotX, m_rotY, m_rotZ, Amg::Vector3D(m_traX, m_traY, m_traZ) ));
951
952 ATH_MSG_INFO( "+++++++++++++++++++++++++++++++++++++++++++++" );
953 ATH_MSG_INFO( "An L1 transformation will be added:" );
954 ATH_MSG_INFO( printTransform(L1transform) );
955 ATH_MSG_INFO( "+++++++++++++++++++++++++++++++++++++++++++++" );
956
957 // add the two contributions:
958 Amg::Transform3D trans = sumTransforms( m_CoG, L1transform );
959
960 // substitute the original m_CoG:
961 m_CoG = trans;
962}
963
964//=====================================================================
965// normalizeCog
966//=====================================================================
967StatusCode InDetAlignCog::normalizeTransform(Amg::Transform3D & trans,
968 const int norm){
969 ATH_MSG_DEBUG("in normalizeTransform with factor " << norm );
970
971 if(norm==0){
972 ATH_MSG_ERROR( "norm factor is null !!!" );
973 return StatusCode::FAILURE;
974 }
975
976 Amg::Vector3D vec = trans.translation();
977 vec /= (double) norm;
978
979 double a,b,g;
980 m_IDAlignDBTool->extractAlphaBetaGamma(trans,a,b,g);
981 a /= (double) norm;
982 b /= (double) norm;
983 g /= (double) norm;
984
985 trans = makeAffine3d( a, b, g, vec);
986
987 return StatusCode::SUCCESS;
988}
989
990//=====================================================================
991// scaleTransform
992//=====================================================================
993void InDetAlignCog::scaleTransform(Amg::Transform3D & trans,
994 const float scale){
995
996 Amg::Vector3D vec = trans.translation();
997 vec *= scale;
998
999 double a,b,g;
1000 m_IDAlignDBTool->extractAlphaBetaGamma(trans,a,b,g);
1001 a *= scale;
1002 b *= scale;
1003 g *= scale;
1004
1005 trans = makeAffine3d(a,b,g,vec);
1006}
1007
1008//=====================================================================
1009// sumTransform
1010//=====================================================================
1011Amg::Transform3D InDetAlignCog::sumTransforms(const Amg::Transform3D & trans1,
1012 const Amg::Transform3D & trans2) const {
1013 double a1,b1,g1;
1014 double a2,b2,g2;
1015
1016 Amg::Vector3D vec1 = trans1.translation();
1017 Amg::Vector3D vec2 = trans2.translation();
1018
1019 m_IDAlignDBTool->extractAlphaBetaGamma(trans1,a1,b1,g1);
1020 m_IDAlignDBTool->extractAlphaBetaGamma(trans2,a2,b2,g2);
1021
1022 Amg::Vector3D sumvec = vec1+vec2;
1023 return makeAffine3d(a1+a2,b1+b2,g1+g2,sumvec);
1024}
1025
1026//===================================================
1027// printTransform
1028//===================================================
1029std::string InDetAlignCog::printTransform(const Amg::Transform3D &trans) const{
1030 std::ostringstream ostr;
1031
1032 AmgMatrix(3,3) rotation( trans.rotation());
1033 Amg::Vector3D translation( trans.translation());
1034 if(m_fullMatrix) {
1035 ostr.setf(std::ios::fixed);
1036 ostr.setf(std::ios::showpoint);
1037
1038 ostr << std::endl << std::endl;
1039 ostr << std::setw(10) << std::setprecision(6) << rotation(0,0)
1040 << std::setw(12) << std::setprecision(6) << rotation(0,1)
1041 << std::setw(12) << std::setprecision(6) << rotation(0,2)
1042 << " Tx = " << std::setw(10) << std::setprecision(6) << translation[0] << std::endl;
1043
1044 ostr << std::setw(10) << std::setprecision(6) << rotation(1,0)
1045 << std::setw(12) << std::setprecision(6) << rotation(1,1)
1046 << std::setw(12) << std::setprecision(6) << rotation(1,2)
1047 << " Ty = " << std::setw(10) << std::setprecision(6) << translation[1] << std::endl;
1048
1049 ostr << std::setw(10) << std::setprecision(6) << rotation(2,0)
1050 << std::setw(12) << std::setprecision(6) << rotation(2,1)
1051 << std::setw(12) << std::setprecision(6) << rotation(2,2)
1052 << " Tz = " << std::setw(10) << std::setprecision(6) << translation[2];
1053
1054 ostr << std::endl;
1055 }
1056 else{
1057 ostr << "(" << translation[0] << "," << translation[1] << "," << translation[2] << ")mm ";
1058 double alpha=0, beta=0, gamma=0;
1059 m_IDAlignDBTool->extractAlphaBetaGamma(trans, alpha, beta, gamma);
1060 ostr << "( A=" << 1000*alpha << ", B=" << 1000*beta << ", G=" << 1000*gamma << ")mrad";
1061 }
1062 return ostr.str();
1063}
1064
1065
1066
1067
1068//=====================================================================
1069// testIdentity
1070//=====================================================================
1071bool InDetAlignCog::testIdentity(const Amg::Transform3D &transform,
1072 double errRot,
1073 double errTrans) const {
1074 ATH_MSG_DEBUG ("in testIdentity for transform " << printTransform(transform));
1075 ATH_MSG_DEBUG ("errTrans=" << m_errTrans << " errRot=" << m_errRot );
1076
1077 bool pass = true;
1078
1079 const Amg::Transform3D & t1 = transform;
1080 const Amg::Transform3D t2(Amg::Transform3D::Identity()); // identity
1081
1082 // Rotation/Scale
1083 for (int i=0; i<3; i++){
1084 for (int j=0; j<3; j++){
1085 double diff = std::abs(t1(i,j) - t2(i,j));
1086 if (diff > errRot) pass = false;
1087 }
1088 }
1089 // Translation
1090 for (int i=0; i<3; i++){
1091 double diff = std::abs(t1(i,3) - t2(i,3));
1092 if (diff > errTrans) pass = false;
1093 }
1094
1095 if (pass) {
1096 ATH_MSG_DEBUG( "Remaining global transform within tolerances. Ok." );
1097 }
1098 else {
1099 ATH_MSG_WARNING( "Remaining global transform outside tolerances."
1100 " No level 1 updates will be done" );
1101 }
1102 return pass;
1103}
1104
1105
1106
ATH_CHECK
#define ATH_CHECK
Evaluate an expression and check for errors.
Definition AthCheckMacros.h:40
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition AthMsgStreamMacros.h:33
ATH_MSG_FATAL
#define ATH_MSG_FATAL(x)
Definition AthMsgStreamMacros.h:34
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition AthMsgStreamMacros.h:31
ATH_MSG_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition AthMsgStreamMacros.h:28
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition AthMsgStreamMacros.h:32
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition AthMsgStreamMacros.h:29
vec
std::vector< size_t > vec
Definition CombinationsGeneratorTest.cxx:9
AmgSymMatrix
#define AmgSymMatrix(dim)
Definition EventPrimitives.h:50
AmgVector
#define AmgVector(rows)
Definition EventPrimitives.h:55
AmgMatrix
#define AmgMatrix(rows, cols)
Definition EventPrimitives.h:49
IOVTime.h
Basic time unit for IOVSvc.
onemrad
static const double onemrad
Definition InDetAlignCog.cxx:81
tenmu
static const double tenmu
Definition InDetAlignCog.cxx:80
InDetAlignCog.h
Algorithm for ID cog calculation.
a
static Double_t a
Definition LArPhysWaveHECTool.cxx:38
PixelID.h
This is an Identifier helper class for the Pixel subdetector.
diff
void diff(const Jet &rJet1, const Jet &rJet2, std::map< std::string, double > varDiff)
Difference between jets - Non-Class function required by trigger.
Definition Jet.cxx:631
SCT_ID.h
This is an Identifier helper class for the SCT subdetector.
angle
double angle(const GeoTrf::Vector2D &a, const GeoTrf::Vector2D &b)
Definition TRTDetectorFactory_Full.cxx:71
TRT_BaseElement.h
TRT_ID.h
This is an Identifier helper class for the TRT subdetector.
y
#define y
x
#define x
z
#define z
AthAlgorithm::AthAlgorithm
AthAlgorithm(const std::string &name, ISvcLocator *pSvcLocator)
Constructor with parameters:
Definition AthAlgorithm.cxx:24
AthCommonDataStore< AthCommonMsg< Algorithm > >::declareProperty
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)
Definition AthCommonDataStore.h:145
AthCommonDataStore< AthCommonMsg< Algorithm > >::detStore
const ServiceHandle< StoreGateSvc > & detStore() const
Definition AthCommonDataStore.h:95
IOVTime::MAXRUN
static constexpr uint32_t MAXRUN
Definition IOVTime.h:48
IOVTime::MINEVENT
static constexpr uint32_t MINEVENT
Definition IOVTime.h:50
IOVTime::MAXEVENT
static constexpr uint32_t MAXEVENT
Definition IOVTime.h:51
IOVTime::MINRUN
static constexpr uint32_t MINRUN
Definition IOVTime.h:44
InDetAlignCog::Params_t
Definition InDetAlignCog.h:49
InDetAlignCog::m_det
int m_det
Pixel=1, SCT=2, Pixel+SCT=12, TRT=3, all (silicon and TRT)=99.
Definition InDetAlignCog.h:109
InDetAlignCog::m_pixid
const PixelID * m_pixid
Definition InDetAlignCog.h:93
InDetAlignCog::m_TRT_layer
int m_TRT_layer
a particular TRT layer or all (TRT)=99
Definition InDetAlignCog.h:113
InDetAlignCog::m_SiTxtOutput
bool m_SiTxtOutput
output Si constants to txt file ?
Definition InDetAlignCog.h:118
InDetAlignCog::getTRT_Elements
StatusCode getTRT_Elements(const InDetDD::TRT_DetElementCollection *, const bool, InDetAlignCog::Params_t &params)
Definition InDetAlignCog.cxx:582
InDetAlignCog::m_pixelDetEleCollKey
SG::ReadCondHandleKey< InDetDD::SiDetectorElementCollection > m_pixelDetEleCollKey
Definition InDetAlignCog.h:104
InDetAlignCog::m_glob_x
Amg::Transform3D m_glob_x
Definition InDetAlignCog.h:169
InDetAlignCog::m_trtid
const TRT_ID * m_trtid
Definition InDetAlignCog.h:95
InDetAlignCog::m_errTrans
double m_errTrans
acceptable value for residual global translation
Definition InDetAlignCog.h:122
InDetAlignCog::m_resglob
double m_resglob[6]
Definition InDetAlignCog.h:181
InDetAlignCog::finalize
StatusCode finalize()
Definition InDetAlignCog.cxx:417
InDetAlignCog::m_traZ
double m_traZ
translation/rotation values (CLHEP::mm, CLHEP::rad) for the arbitrary transformation
Definition InDetAlignCog.h:156
InDetAlignCog::m_SCTDetEleCollKey
SG::ReadCondHandleKey< InDetDD::SiDetectorElementCollection > m_SCTDetEleCollKey
Definition InDetAlignCog.h:105
InDetAlignCog::getL1Transform
const Amg::Transform3D getL1Transform(int bec)
Definition InDetAlignCog.cxx:670
InDetAlignCog::m_trtDetEleContKey
SG::ReadCondHandleKey< InDetDD::TRT_DetElementContainer > m_trtDetEleContKey
Definition InDetAlignCog.h:106
InDetAlignCog::m_TRTAlignDbTool
ServiceHandle< ITRT_AlignDbSvc > m_TRTAlignDbTool
Definition InDetAlignCog.h:102
InDetAlignCog::m_grot_z
Amg::Transform3D m_grot_z
Definition InDetAlignCog.h:174
InDetAlignCog::execute
StatusCode execute()
Definition InDetAlignCog.cxx:227
InDetAlignCog::m_glob_z
Amg::Transform3D m_glob_z
Definition InDetAlignCog.h:171
InDetAlignCog::m_sigXsctB
double m_sigXsctB
assumed error for SCT barrel local X matchnig in the Xi2 method
Definition InDetAlignCog.h:130
InDetAlignCog::m_TRT_TxtOutput
bool m_TRT_TxtOutput
output TRT constants to txt file ?
Definition InDetAlignCog.h:119
InDetAlignCog::scaleTransform
void scaleTransform(Amg::Transform3D &, const float)
Definition InDetAlignCog.cxx:993
InDetAlignCog::normalizeTransform
StatusCode normalizeTransform(Amg::Transform3D &, const int)
Definition InDetAlignCog.cxx:967
InDetAlignCog::m_CoG
Amg::Transform3D m_CoG
Definition InDetAlignCog.h:178
InDetAlignCog::m_ResGlob
Amg::Transform3D m_ResGlob
Definition InDetAlignCog.h:182
InDetAlignCog::m_glob_y
Amg::Transform3D m_glob_y
Definition InDetAlignCog.h:170
InDetAlignCog::printTransform
std::string printTransform(const Amg::Transform3D &) const
Definition InDetAlignCog.cxx:1029
InDetAlignCog::m_errRot
double m_errRot
acceptable value for residual global rotation angles
Definition InDetAlignCog.h:123
InDetAlignCog::m_sigXpixB
double m_sigXpixB
assumed error for Pixel barrel local X matchnig in the Xi2 method
Definition InDetAlignCog.h:124
InDetAlignCog::prepareDerivative
void prepareDerivative(const Amg::Transform3D &, const bool=false)
Definition InDetAlignCog.cxx:840
InDetAlignCog::sumTransforms
Amg::Transform3D sumTransforms(const Amg::Transform3D &, const Amg::Transform3D &) const
Definition InDetAlignCog.cxx:1011
InDetAlignCog::getSiElements
StatusCode getSiElements(const InDetDD::SiDetectorElementCollection *, const bool, InDetAlignCog::Params_t &params)
Definition InDetAlignCog.cxx:460
InDetAlignCog::m_doL1
bool m_doL1
enable/disable introducing the arbitrary L1 correction to the output objects
Definition InDetAlignCog.h:151
InDetAlignCog::testIdentity
bool testIdentity(const Amg::Transform3D &, double, double) const
Definition InDetAlignCog.cxx:1071
InDetAlignCog::shiftIDbyCog
StatusCode shiftIDbyCog()
Definition InDetAlignCog.cxx:878
InDetAlignCog::m_TRT_bec
int m_TRT_bec
Barrel=-1, Endcaps=+-2, all (TRT)=99.
Definition InDetAlignCog.h:111
InDetAlignCog::m_doCoG
bool m_doCoG
enable/disable introducing the CoG correction to the output objects
Definition InDetAlignCog.h:150
InDetAlignCog::accumulate
void accumulate(const Amg::Transform3D &, double *)
Definition InDetAlignCog.cxx:681
InDetAlignCog::initialize
StatusCode initialize()
Definition InDetAlignCog.cxx:190
InDetAlignCog::m_sigXtrtB
double m_sigXtrtB
assumed error for TRT barrel local X matchnig in the Xi2 method
Definition InDetAlignCog.h:136
InDetAlignCog::m_SQLiteTag
std::string m_SQLiteTag
SQLite tag name.
Definition InDetAlignCog.h:166
InDetAlignCog::accumulateChi2
void accumulateChi2(const Amg::Transform3D &, AmgSymMatrix(6)&, AmgVector(6)&, const double *)
Definition InDetAlignCog.cxx:749
InDetAlignCog::m_cog
double m_cog[6]
Definition InDetAlignCog.h:177
InDetAlignCog::m_IDAlignDBTool
ToolHandle< IInDetAlignDBTool > m_IDAlignDBTool
Definition InDetAlignCog.h:98
InDetAlignCog::m_grot_x
Amg::Transform3D m_grot_x
Definition InDetAlignCog.h:172
InDetAlignCog::m_sitxtfile
std::string m_sitxtfile
text file with dump of Si alignment constants after cog shift
Definition InDetAlignCog.h:164
InDetAlignCog::m_doTZ
bool m_doTZ
enable/disable writing of indivitual DoF's to the db
Definition InDetAlignCog.h:145
InDetAlignCog::m_trt_txtfile
std::string m_trt_txtfile
text file with dump of TRT alignment constants after cog shift
Definition InDetAlignCog.h:165
InDetAlignCog::m_Si_bec
int m_Si_bec
Barrel=1, Endcaps=+-2, all (silicon)=99.
Definition InDetAlignCog.h:110
InDetAlignCog::enableCoG
void enableCoG(Amg::Transform3D &, bool, bool, bool, bool, bool, bool)
Definition InDetAlignCog.cxx:924
InDetAlignCog::InDetAlignCog
InDetAlignCog(const std::string &name, ISvcLocator *pSvcLocator)
Definition InDetAlignCog.cxx:83
InDetAlignCog::m_counter
int m_counter
normalization factor
Definition InDetAlignCog.h:162
InDetAlignCog::m_grot_y
Amg::Transform3D m_grot_y
Definition InDetAlignCog.h:173
InDetAlignCog::m_sigXpixE
double m_sigXpixE
assumed error for Pixel endcap local X matchnig in the Xi2 method
Definition InDetAlignCog.h:127
InDetAlignCog::m_sctid
const SCT_ID * m_sctid
Definition InDetAlignCog.h:94
InDetAlignCog::m_Si_layer
int m_Si_layer
a particular silicon layer or all (silicon)=99
Definition InDetAlignCog.h:112
InDetDD::SiDetectorElementCollection
Class to hold the SiDetectorElement objects to be put in the detector store.
Definition SiDetectorElementCollection.h:27
InDetDD::SiDetectorElement
Class to hold geometrical description of a silicon detector element.
Definition SiDetectorElement.h:109
InDetDD::TRT_BaseElement
Virtual base class of TRT readout elements.
Definition TRT_BaseElement.h:52
InDetDD::TRT_DetElementCollection
Class to hold collection of TRT detector elements.
Definition TRT_DetElementCollection.h:26
Amg::AngleAxis3D
Eigen::AngleAxisd AngleAxis3D
Definition GeoPrimitives.h:45
Amg::Transform3D
Eigen::Affine3d Transform3D
Definition GeoPrimitives.h:46
Amg::Vector3D
Eigen::Matrix< double, 3, 1 > Vector3D
Definition GeoPrimitives.h:47
Amg::Translation3D
Eigen::Translation< double, 3 > Translation3D
Definition GeoPrimitives.h:44
A
hold the test vectors and ease the comparison
Generated on for ATLAS Offline Software by doxygen 1.14.0