1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 //-----------------------------------------------------------------------------
19 /// \class AliMUONAlignmentV5
20 /// Alignment class fro the ALICE DiMuon spectrometer
22 /// MUON specific alignment class which interface to AliMillepede.
23 /// For each track ProcessTrack calculates the local and global derivatives
24 /// at each hit and fill the corresponding local equations. Provide methods for
25 /// fixing or constraining detection elements for best results.
27 /// \author Bruce Becker, Javier Castillo
28 //-----------------------------------------------------------------------------
30 #include "AliMUONAlignment.h"
31 #include "AliMUONTrack.h"
32 #include "AliMUONRawCluster.h"
33 #include "AliMUONTrackParam.h"
34 #include "AliMUONHitForRec.h"
35 #include "AliMUONGeometryTransformer.h"
36 #include "AliMUONGeometryModuleTransformer.h"
37 #include "AliMUONGeometryDetElement.h"
38 #include "AliMUONGeometryBuilder.h"
39 #include "AliMUONConstants.h"
40 #include "AliMillepede.h"
42 #include "AliMpExMap.h"
48 ClassImp(AliMUONAlignment)
49 Int_t AliMUONAlignment::fgNDetElem = 4*2+4*2+18*2+26*2+26*2;
50 Int_t AliMUONAlignment::fgNDetElemCh[10] = {4,4,4,4,18,18,26,26,26,26};
51 Int_t AliMUONAlignment::fgSNDetElemCh[10] = {4,8,12,16,34,52,78,104,130,156};
52 Int_t AliMUONAlignment::fgNParCh = 3;
53 Int_t AliMUONAlignment::fgNCh = 10;
54 Int_t AliMUONAlignment::fgNSt = 5;
56 AliMUONAlignment::AliMUONAlignment()
68 fNGlobal(fgNDetElem*fgNParCh),
78 /// Default constructor, setting define alignment parameters
82 fDoF[0] = kTRUE; fDoF[1] = kTRUE; fDoF[2] = kTRUE;
83 fAllowVar[0] = 0.05; fAllowVar[1] = 0.05; fAllowVar[2] = 0.001;
85 AliInfo(Form("fAllowVar[0]: %f\t fAllowVar[1]: %f\t fPhi: %f\t fgNDetElem: %i\t fNGlobal: %i\t fNLocal: %i",fAllowVar[0],fAllowVar[1],fPhi,fgNDetElem,fNGlobal,fNLocal));
87 fMillepede = new AliMillepede();
89 Init(fNGlobal, fNLocal, fNStdDev);
92 AliInfo("Parameters initialized to zero");
96 AliMUONAlignment::~AliMUONAlignment() {
100 void AliMUONAlignment::Init(Int_t nGlobal, /* number of global paramers */
101 Int_t nLocal, /* number of local parameters */
102 Int_t nStdDev /* std dev cut */ )
104 /// Initialization of AliMillepede. Fix parameters, define constraints ...
105 fMillepede->InitMille(nGlobal,nLocal,nStdDev,fResCut,fResCutInitial);
107 Bool_t bStOnOff[5] = {kFALSE,kFALSE,kTRUE,kTRUE,kTRUE};
109 AllowVariations(bStOnOff);
111 // Fix parameters or add constraints here
112 for (Int_t iSt=0; iSt<5; iSt++)
113 if (!bStOnOff[iSt]) FixStation(iSt+1);
115 Bool_t bVarXYT[3] = {kFALSE,kTRUE,kFALSE};
116 Bool_t bDetTLBR[4] = {kFALSE,kTRUE,kFALSE,kTRUE};
118 AddConstraints(bStOnOff,bVarXYT,bDetTLBR);
119 bVarXYT[0] = kFALSE; bVarXYT[1] = kFALSE; bVarXYT[2] = kTRUE;
120 bDetTLBR[0] = kFALSE; bDetTLBR[1] = kTRUE; bDetTLBR[2] = kFALSE; bDetTLBR[3] = kFALSE;
121 AddConstraints(bStOnOff,bVarXYT,bDetTLBR);
122 bVarXYT[0] = kFALSE; bVarXYT[1] = kFALSE; bVarXYT[2] = kTRUE;
123 AddConstraints(bStOnOff,bVarXYT);
126 if (fStartFac>1) fMillepede->SetIterations(fStartFac);
129 void AliMUONAlignment::FixStation(Int_t iSt){
130 /// Fix all detection elements of station iSt
131 Int_t iDetElemFirst = (iSt>1) ? fgSNDetElemCh[2*(iSt-1)-1] : 0;
132 Int_t iDetElemLast = fgSNDetElemCh[2*(iSt)-1];
133 for (Int_t i = iDetElemFirst; i < iDetElemLast; i++){
134 FixParameter(i*fgNParCh+0, 0.0);
135 FixParameter(i*fgNParCh+1, 0.0);
136 FixParameter(i*fgNParCh+2, 0.0);
140 void AliMUONAlignment::SetNonLinear(Bool_t *lStOnOff,Bool_t *lVarXYT){
141 /// Set non linear parameter flag selected stations and degrees of freedom
142 for (Int_t i = 0; i < fgNDetElem; i++){
144 for (iCh=1; iCh<=fgNCh; iCh++){
145 if (i<fgSNDetElemCh[iCh-1]) break;
147 Int_t iSt = lStOnOff[(iCh-1)/2] ? (iCh+1)/2 : 0;
149 if (lVarXYT[0]) { // X constraint
150 SetNonLinear(i*fgNParCh+0);
152 if (lVarXYT[1]) { // Y constraint
153 SetNonLinear(i*fgNParCh+1);
155 if (lVarXYT[2]) { // T constraint
156 SetNonLinear(i*fgNParCh+2);
162 void AliMUONAlignment::AddConstraints(Bool_t *lStOnOff,Bool_t *lVarXYT){
163 /// Add constraint equations for selected stations and degrees of freedom
164 for (Int_t i = 0; i < fgNDetElem; i++){
166 for (iCh=1; iCh<=fgNCh; iCh++){
167 if (i<fgSNDetElemCh[iCh-1]) break;
169 Int_t iSt = lStOnOff[(iCh-1)/2] ? (iCh+1)/2 : 0;
171 if (lVarXYT[0]) { // X constraint
172 fConstraintX[i*fgNParCh+0]=1.0;
174 if (lVarXYT[1]) { // Y constraint
175 fConstraintY[i*fgNParCh+1]=1.0;
177 if (lVarXYT[2]) { // T constraint
178 fConstraintP[i*fgNParCh+2]=1.0;
182 if (lVarXYT[0]) { // X constraint
183 AddConstraint(fConstraintX,0.0);
185 if (lVarXYT[1]) { // Y constraint
186 AddConstraint(fConstraintY,0.0);
188 if (lVarXYT[2]) { // T constraint
189 AddConstraint(fConstraintP,0.0);
193 void AliMUONAlignment::AddConstraints(Bool_t *lStOnOff,Bool_t *lVarXYT, Bool_t *lDetTLBR){
194 /// Add constraint equations for selected stations, degrees of freedom detector half
195 for (Int_t i = 0; i < fgNDetElem; i++){
197 for (iCh=1; iCh<=fgNCh; iCh++){
198 if (i<fgSNDetElemCh[iCh-1]) break;
200 Int_t iSt = lStOnOff[(iCh-1)/2] ? (iCh+1)/2 : 0;
202 if (lVarXYT[0]) { // X constraint
203 if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintXT,0); // Top half
204 if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintXL,0); // Left half
205 if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintXB,0); // Bottom half
206 if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintXR,0); // Right half
208 if (lVarXYT[1]) { // X constraint
209 if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintYT,1); // Top half
210 if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintYL,1); // Left half
211 if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintYB,1); // Bottom half
212 if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintYR,1); // Right half
214 if (lVarXYT[2]) { // X constraint
215 if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintPT,2); // Top half
216 if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintPL,2); // Left half
217 if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintPB,2); // Bottom half
218 if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintPR,2); // Right half
222 if (lVarXYT[0]) { // X constraint
223 if (lDetTLBR[0]) AddConstraint(fConstraintXT,0.0); // Top half
224 if (lDetTLBR[1]) AddConstraint(fConstraintXL,0.0); // Left half
225 if (lDetTLBR[2]) AddConstraint(fConstraintXB,0.0); // Bottom half
226 if (lDetTLBR[3]) AddConstraint(fConstraintXR,0.0); // Right half
228 if (lVarXYT[1]) { // Y constraint
229 if (lDetTLBR[0]) AddConstraint(fConstraintYT,0.0); // Top half
230 if (lDetTLBR[1]) AddConstraint(fConstraintYL,0.0); // Left half
231 if (lDetTLBR[2]) AddConstraint(fConstraintYB,0.0); // Bottom half
232 if (lDetTLBR[3]) AddConstraint(fConstraintYR,0.0); // Right half
234 if (lVarXYT[2]) { // T constraint
235 if (lDetTLBR[0]) AddConstraint(fConstraintPT,0.0); // Top half
236 if (lDetTLBR[1]) AddConstraint(fConstraintPL,0.0); // Left half
237 if (lDetTLBR[2]) AddConstraint(fConstraintPB,0.0); // Bottom half
238 if (lDetTLBR[3]) AddConstraint(fConstraintPR,0.0); // Right half
242 void AliMUONAlignment::ConstrainT(Int_t lDetElem, Int_t lCh, Double_t *lConstraintT, Int_t iVar){
243 /// Set constrain equation for top half of spectrometer
244 Int_t lDetElemNumber = (lCh==1) ? lDetElem : lDetElem-fgSNDetElemCh[lCh-2];
245 if (lCh>=1 && lCh<=4){
246 if (lDetElemNumber==0 || lDetElemNumber==1){ // From track crossings
247 lConstraintT[lDetElem*fgNParCh+iVar]=1.0;
250 if (lCh>=5 && lCh<=6){
251 if (lDetElemNumber>=0&&lDetElemNumber<=9){
252 lConstraintT[lDetElem*fgNParCh+iVar]=1.0;
255 if (lCh>=7 && lCh<=10){
256 if (lDetElemNumber>=0&&lDetElemNumber<=13){
257 lConstraintT[lDetElem*fgNParCh+iVar]=1.0;
262 void AliMUONAlignment::ConstrainL(Int_t lDetElem, Int_t lCh, Double_t *lConstraintL, Int_t iVar){
263 /// Set constrain equation for left half of spectrometer
264 Int_t lDetElemNumber = (lCh==1) ? lDetElem : lDetElem-fgSNDetElemCh[lCh-2];
265 if (lCh>=1 && lCh<=4){
266 if (lDetElemNumber==1 || lDetElemNumber==2){ // From track crossings
267 lConstraintL[lDetElem*fgNParCh+iVar]=1.0;
270 if (lCh>=5 && lCh<=6){
271 if (lDetElemNumber>=5&&lDetElemNumber<=13){
272 lConstraintL[lDetElem*fgNParCh+iVar]=1.0;
275 if (lCh>=7 && lCh<=10){
276 if (lDetElemNumber>=7&&lDetElemNumber<=19){
277 lConstraintL[lDetElem*fgNParCh+iVar]=1.0;
282 void AliMUONAlignment::ConstrainB(Int_t lDetElem, Int_t lCh, Double_t *lConstraintB, Int_t iVar){
283 /// Set constrain equation for bottom half of spectrometer
284 Int_t lDetElemNumber = (lCh==1) ? lDetElem : lDetElem-fgSNDetElemCh[lCh-2];
285 if (lCh>=1 && lCh<=4){
286 if (lDetElemNumber==2 && lDetElemNumber==3){ // From track crossings
287 lConstraintB[lDetElem*fgNParCh+iVar]=1.0;
290 if (lCh>=5 && lCh<=6){
291 if ((lDetElemNumber>=9&&lDetElemNumber<=17) ||
292 (lDetElemNumber==0)){
293 lConstraintB[lDetElem*fgNParCh+iVar]=1.0;
296 if (lCh>=7 && lCh<=10){
297 if ((lDetElemNumber>=13&&lDetElemNumber<=25) ||
298 (lDetElemNumber==0)){
299 lConstraintB[lDetElem*fgNParCh+iVar]=1.0;
304 void AliMUONAlignment::ConstrainR(Int_t lDetElem, Int_t lCh, Double_t *lConstraintR, Int_t iVar){
305 /// Set constrain equation for right half of spectrometer
306 Int_t lDetElemNumber = (lCh==1) ? lDetElem : lDetElem-fgSNDetElemCh[lCh-2];
307 if (lCh>=1 && lCh<=4){
308 if (lDetElemNumber==0 && lDetElemNumber==3){ // From track crossings
309 lConstraintR[lDetElem*fgNParCh+iVar]=1.0;
312 if (lCh>=5 && lCh<=6){
313 if ((lDetElemNumber>=0&&lDetElemNumber<=4) ||
314 (lDetElemNumber>=14&&lDetElemNumber<=17)){
315 lConstraintR[lDetElem*fgNParCh+iVar]=1.0;
318 if (lCh>=7 && lCh<=10){
319 if ((lDetElemNumber>=0&&lDetElemNumber<=6) ||
320 (lDetElemNumber>=20&&lDetElemNumber<=25)){
321 lConstraintR[lDetElem*fgNParCh+iVar]=1.0;
326 void AliMUONAlignment::ResetConstraints(){
327 /// Reset all constraint equations
328 for (Int_t i = 0; i < fgNDetElem; i++){
329 fConstraintX[i*fgNParCh+0]=0.0;
330 fConstraintX[i*fgNParCh+1]=0.0;
331 fConstraintX[i*fgNParCh+2]=0.0;
332 fConstraintY[i*fgNParCh+0]=0.0;
333 fConstraintY[i*fgNParCh+1]=0.0;
334 fConstraintY[i*fgNParCh+2]=0.0;
335 fConstraintP[i*fgNParCh+0]=0.0;
336 fConstraintP[i*fgNParCh+1]=0.0;
337 fConstraintP[i*fgNParCh+2]=0.0;
338 fConstraintXT[i*fgNParCh+0]=0.0;
339 fConstraintXT[i*fgNParCh+1]=0.0;
340 fConstraintXT[i*fgNParCh+2]=0.0;
341 fConstraintYT[i*fgNParCh+0]=0.0;
342 fConstraintYT[i*fgNParCh+1]=0.0;
343 fConstraintYT[i*fgNParCh+2]=0.0;
344 fConstraintPT[i*fgNParCh+0]=0.0;
345 fConstraintPT[i*fgNParCh+1]=0.0;
346 fConstraintPT[i*fgNParCh+2]=0.0;
347 fConstraintXL[i*fgNParCh+0]=0.0;
348 fConstraintXL[i*fgNParCh+1]=0.0;
349 fConstraintXL[i*fgNParCh+2]=0.0;
350 fConstraintYL[i*fgNParCh+0]=0.0;
351 fConstraintYL[i*fgNParCh+1]=0.0;
352 fConstraintYL[i*fgNParCh+2]=0.0;
353 fConstraintPL[i*fgNParCh+0]=0.0;
354 fConstraintPL[i*fgNParCh+1]=0.0;
355 fConstraintPL[i*fgNParCh+2]=0.0;
356 fConstraintXB[i*fgNParCh+0]=0.0;
357 fConstraintXB[i*fgNParCh+1]=0.0;
358 fConstraintXB[i*fgNParCh+2]=0.0;
359 fConstraintYB[i*fgNParCh+0]=0.0;
360 fConstraintYB[i*fgNParCh+1]=0.0;
361 fConstraintYB[i*fgNParCh+2]=0.0;
362 fConstraintPB[i*fgNParCh+0]=0.0;
363 fConstraintPB[i*fgNParCh+1]=0.0;
364 fConstraintPB[i*fgNParCh+2]=0.0;
365 fConstraintXR[i*fgNParCh+0]=0.0;
366 fConstraintXR[i*fgNParCh+1]=0.0;
367 fConstraintXR[i*fgNParCh+2]=0.0;
368 fConstraintYR[i*fgNParCh+0]=0.0;
369 fConstraintYR[i*fgNParCh+1]=0.0;
370 fConstraintYR[i*fgNParCh+2]=0.0;
371 fConstraintPR[i*fgNParCh+0]=0.0;
372 fConstraintPR[i*fgNParCh+1]=0.0;
373 fConstraintPR[i*fgNParCh+2]=0.0;
377 void AliMUONAlignment::AddConstraint(Double_t *par, Double_t value) {
378 /// Constrain equation defined by par to value
379 fMillepede->SetGlobalConstraint(par, value);
380 AliInfo("Adding constraint");
383 void AliMUONAlignment::InitGlobalParameters(Double_t *par) {
384 /// Initialize global parameters with par array
385 fMillepede->SetGlobalParameters(par);
386 AliInfo("Init Global Parameters");
389 void AliMUONAlignment::FixParameter(Int_t iPar, Double_t value) {
390 /// Parameter iPar is encourage to vary in [-value;value].
391 /// If value == 0, parameter is fixed
392 fMillepede->SetParSigma(iPar, value);
393 if (value==0) AliInfo(Form("Parameter %i Fixed", iPar));
396 void AliMUONAlignment::ResetLocalEquation()
398 /// Reset the derivative vectors
399 for(int i=0; i<fNLocal; i++) {
400 fLocalDerivatives[i] = 0.0;
402 for(int i=0; i<fNGlobal; i++) {
403 fGlobalDerivatives[i] = 0.0;
407 void AliMUONAlignment::AllowVariations(Bool_t *bStOnOff) {
408 /// Set allowed variation for selected stations based on fDoF and fAllowVar
409 for (Int_t iSt=1; iSt<=5; iSt++) {
410 if (bStOnOff[iSt-1]) {
411 Int_t iDetElemFirst = (iSt>1) ? fgSNDetElemCh[2*(iSt-1)-1] : 0;
412 Int_t iDetElemLast = fgSNDetElemCh[2*(iSt)-1];
413 for (int i=0; i<fgNParCh; i++) {
414 AliDebug(1,Form("fDoF[%d]= %d",i,fDoF[i]));
416 for (Int_t j=iDetElemFirst; j<iDetElemLast; j++){
417 FixParameter(j*fgNParCh+i, fAllowVar[i]);
425 void AliMUONAlignment::SetNonLinear(Int_t iPar /* set non linear flag */ ) {
426 /// Set nonlinear flag for parameter iPar
427 fMillepede->SetNonLinear(iPar);
428 AliInfo(Form("Parameter %i set to non linear", iPar));
431 void AliMUONAlignment::LocalEquationX() {
432 /// Define local equation for current track and hit in x coor. measurement
433 // set local derivatives
434 SetLocalDerivative(0, fCosPhi);
435 SetLocalDerivative(1, fCosPhi * (fTrackPos[2] - fTrackPos0[2]));
436 SetLocalDerivative(2, fSinPhi);
437 SetLocalDerivative(3, fSinPhi * (fTrackPos[2] - fTrackPos0[2]));
439 // set global derivatives
440 SetGlobalDerivative(fDetElemNumber*fgNParCh+0, -1.);
441 SetGlobalDerivative(fDetElemNumber*fgNParCh+1, 0.);
443 SetGlobalDerivative(fDetElemNumber*fgNParCh+2,
444 -fSinPhi*(fTrackPos[0]-fDetElemPos[0])
445 +fCosPhi*(fTrackPos[1]-fDetElemPos[1]));
448 SetGlobalDerivative(fDetElemNumber*fgNParCh+2,
449 -fSinPhi*(fTrackPos0[0]+fTrackSlope0[0]*
450 (fTrackPos[2]-fTrackPos0[2])-fDetElemPos[0])
451 +fCosPhi*(fTrackPos0[1]+fTrackSlope0[1]*
452 (fTrackPos[2]-fTrackPos0[2])-fDetElemPos[1]));
455 fMillepede->SetLocalEquation(fGlobalDerivatives, fLocalDerivatives, fMeas[0], fSigma[0]);
458 void AliMUONAlignment::LocalEquationY() {
459 /// Define local equation for current track and hit in y coor. measurement
460 // set local derivatives
461 SetLocalDerivative(0,-fSinPhi);
462 SetLocalDerivative(1,-fSinPhi * (fTrackPos[2] - fTrackPos0[2]));
463 SetLocalDerivative(2, fCosPhi);
464 SetLocalDerivative(3, fCosPhi * (fTrackPos[2] - fTrackPos0[2]));
466 // set global derivatives
467 SetGlobalDerivative(fDetElemNumber*fgNParCh+0, 0.);
468 SetGlobalDerivative(fDetElemNumber*fgNParCh+1, -1.);
470 SetGlobalDerivative(fDetElemNumber*fgNParCh+2,
471 -fCosPhi*(fTrackPos[0]-fDetElemPos[0])
472 -fSinPhi*(fTrackPos[1]-fDetElemPos[1]));
475 SetGlobalDerivative(fDetElemNumber*fgNParCh+2,
476 -fCosPhi*(fTrackPos0[0]+fTrackSlope0[0]*
477 (fTrackPos[2]-fTrackPos0[2])-fDetElemPos[0])
478 -fSinPhi*(fTrackPos0[1]+fTrackSlope0[1]*
479 (fTrackPos[2]-fTrackPos0[2])-fDetElemPos[1]));
482 fMillepede->SetLocalEquation(fGlobalDerivatives, fLocalDerivatives, fMeas[1], fSigma[1]);
485 void AliMUONAlignment::FillRecPointData() {
486 /// Get information of current hit
487 fClustPos[0] = fRecHit->GetNonBendingCoor();
488 fClustPos[1] = fRecHit->GetBendingCoor();
489 fClustPos[2] = fRecHit->GetZ();
490 fTransform->Global2Local(fDetElemId,fClustPos[0],fClustPos[1],fClustPos[2],
491 fClustPosLoc[0],fClustPosLoc[1],fClustPosLoc[2]);
494 void AliMUONAlignment::FillTrackParamData() {
495 /// Get information of current track at current hit
496 fTrackPos[0] = fTrackParam->GetNonBendingCoor();
497 fTrackPos[1] = fTrackParam->GetBendingCoor();
498 fTrackPos[2] = fTrackParam->GetZ();
499 fTrackSlope[0] = fTrackParam->GetNonBendingSlope();
500 fTrackSlope[1] = fTrackParam->GetBendingSlope();
501 fTransform->Global2Local(fDetElemId,fTrackPos[0],fTrackPos[1],fTrackPos[2],
502 fTrackPosLoc[0],fTrackPosLoc[1],fTrackPosLoc[2]);
505 void AliMUONAlignment::FillDetElemData() {
506 /// Get information of current detection element
507 Double_t lDetElemLocX = 0.;
508 Double_t lDetElemLocY = 0.;
509 Double_t lDetElemLocZ = 0.;
510 fDetElemId = fRecHit->GetDetElemId();
511 fDetElemNumber = fDetElemId%100;
512 for (int iCh=0; iCh<fDetElemId/100-1; iCh++){
513 fDetElemNumber += fgNDetElemCh[iCh];
515 fTransform->Local2Global(fDetElemId,lDetElemLocX,lDetElemLocY,lDetElemLocZ,
516 fDetElemPos[0],fDetElemPos[1],fDetElemPos[2]);
517 if (fDetElemId/100 < 5){
527 void AliMUONAlignment::ProcessTrack(AliMUONTrack * track) {
528 /// Process track; Loop over hits and set local equations
530 // get tclones arrays.
531 fTrackParamAtHit = fTrack->GetTrackParamAtHit();
532 fHitForRecAtHit = fTrack->GetHitForRecAtHit();
534 // get size of arrays
535 Int_t nTrackParam = fTrackParamAtHit->GetEntries();
536 Int_t nHitForRec = fHitForRecAtHit->GetEntries();
537 AliInfo(Form("Number of track param entries : %i ", nTrackParam));
538 AliInfo(Form("Number of hit for rec entries : %i ", nHitForRec));
540 for(Int_t iHit=0; iHit<nHitForRec; iHit++) {
541 fRecHit = (AliMUONHitForRec *) fHitForRecAtHit->At(iHit);
542 fTrackParam = (AliMUONTrackParam *) fTrack->GetTrackParamAtHit()->At(iHit);
543 if (!fRecHit || !fTrackParam) continue;
544 // fill local variables for this position --> one measurement
547 FillTrackParamData();
548 // if (fDetElemId<500) continue;
549 fTrackPos0[0] = fTrackPos[0];
550 fTrackPos0[1] = fTrackPos[1];
551 fTrackPos0[2] = fTrackPos[2];
552 fTrackSlope0[0] = fTrackSlope[0];
553 fTrackSlope0[1] = fTrackSlope[1];
557 for(Int_t iHit=0; iHit<nHitForRec; iHit++) {
558 // and get new pointers
559 fRecHit = (AliMUONHitForRec *) fHitForRecAtHit->At(iHit);
560 fTrackParam = (AliMUONTrackParam *) fTrack->GetTrackParamAtHit()->At(iHit);
561 if (!fRecHit || !fTrackParam) continue;
562 // fill local variables for this position --> one measurement
565 FillTrackParamData();
566 // if (fDetElemId<500) continue;
567 AliDebug(1,Form("cluster: %i", iHit));
568 AliDebug(1,Form("x: %f\t y: %f\t z: %f\t DetElemID: %i\t ", fClustPos[0], fClustPos[1], fClustPos[2], fDetElemId));
569 AliDebug(1,Form("fDetElemPos[0]: %f\t fDetElemPos[1]: %f\t fDetElemPos[2]: %f\t DetElemID: %i\t ", fDetElemPos[0],fDetElemPos[1],fDetElemPos[2], fDetElemId));
571 AliDebug(1,Form("Track Parameter: %i", iHit));
572 AliDebug(1,Form("x: %f\t y: %f\t z: %f\t slopex: %f\t slopey: %f", fTrackPos[0], fTrackPos[1], fTrackPos[2], fTrackSlope[0], fTrackSlope[1]));
573 AliDebug(1,Form("x0: %f\t y0: %f\t z0: %f\t slopex0: %f\t slopey0: %f", fTrackPos0[0], fTrackPos0[1], fTrackPos0[2], fTrackSlope0[0], fTrackSlope0[1]));
575 fCosPhi = TMath::Cos(fPhi);
576 fSinPhi = TMath::Sin(fPhi);
578 fMeas[0] = fTrackPos[0] - fClustPos[0];
579 fMeas[1] = fTrackPos[1] - fClustPos[1];
582 fMeas[0] = - fClustPos[0];
583 fMeas[1] = - fClustPos[1];
585 AliDebug(1,Form("fMeas[0]: %f\t fMeas[1]: %f\t fSigma[0]: %f\t fSigma[1]: %f", fMeas[0], fMeas[1], fSigma[0], fSigma[1]));
586 // Set local equations
592 void AliMUONAlignment::LocalFit(Int_t iTrack, Double_t *lTrackParam, Int_t lSingleFit) {
593 /// Call local fit for this tracks
594 Int_t iRes = fMillepede->LocalFit(iTrack,lTrackParam,lSingleFit);
595 if (iRes && !lSingleFit) {
596 fMillepede->SetNLocalEquations(fMillepede->GetNLocalEquations()+1);
600 void AliMUONAlignment::GlobalFit(Double_t *parameters,Double_t *errors,Double_t *pulls) {
601 /// Call global fit; Global parameters are stored in parameters
602 fMillepede->GlobalFit(parameters,errors,pulls);
603 AliInfo("Done fitting global parameters!");
604 for (int iGlob=0; iGlob<fgNDetElem; iGlob++){
605 printf("%d\t %f\t %f\t %f \n",iGlob,parameters[iGlob*fgNParCh+0],parameters[iGlob*fgNParCh+1],parameters[iGlob*fgNParCh+2]);
609 Double_t AliMUONAlignment::GetParError(Int_t iPar) {
610 /// Get error of parameter iPar
611 Double_t lErr = fMillepede->GetParError(iPar);
615 void AliMUONAlignment::PrintGlobalParameters() {
616 /// Print global parameters
617 fMillepede->PrintGlobalParameters();
620 //_________________________________________________________________________
621 TGeoCombiTrans AliMUONAlignment::ReAlign(const TGeoCombiTrans & transform, double *detElemMisAlignment) const
623 /// Realign given transformation by given misalignment and return the misaligned transformation
625 Double_t cartMisAlig[3] = {0,0,0};
626 Double_t angMisAlig[3] = {0,0,0};
627 const Double_t *trans = transform.GetTranslation();
629 // check if the rotation we obtain is not NULL
630 if (transform.GetRotation()) {
631 rot = transform.GetRotation();
634 rot = new TGeoRotation("rot");
635 } // default constructor.
637 cartMisAlig[0] = -detElemMisAlignment[0];
638 cartMisAlig[1] = -detElemMisAlignment[1];
639 angMisAlig[2] = -detElemMisAlignment[2]*180./TMath::Pi();
641 TGeoTranslation newTrans(cartMisAlig[0] + trans[0], cartMisAlig[1] + trans[1], cartMisAlig[2] + trans[2]);
643 rot->RotateX(angMisAlig[0]);
644 rot->RotateY(angMisAlig[1]);
645 rot->RotateZ(angMisAlig[2]);
647 return TGeoCombiTrans(newTrans, *rot);
650 //______________________________________________________________________
651 AliMUONGeometryTransformer *
652 AliMUONAlignment::ReAlign(const AliMUONGeometryTransformer * transformer,
653 double *misAlignments, Bool_t verbose)
656 /////////////////////////////////////////////////////////////////////
657 // Takes the internal geometry module transformers, copies them
658 // and gets the Detection Elements from them.
659 // Takes misalignment parameters and applies these
660 // to the local transform of the Detection Element
661 // Obtains the global transform by multiplying the module transformer
662 // transformation with the local transformation
663 // Applies the global transform to a new detection element
664 // Adds the new detection element to a new module transformer
665 // Adds the new module transformer to a new geometry transformer
666 // Returns the new geometry transformer
668 Double_t lDetElemMisAlignment[3] = {0.,0.,0.};
670 AliMUONGeometryTransformer *newGeometryTransformer =
671 new AliMUONGeometryTransformer(kTRUE);
672 for (Int_t iMt = 0; iMt < transformer->GetNofModuleTransformers(); iMt++)
673 { // module transformers
675 const AliMUONGeometryModuleTransformer *kModuleTransformer =
676 transformer->GetModuleTransformer(iMt, true);
678 AliMUONGeometryModuleTransformer *newModuleTransformer =
679 new AliMUONGeometryModuleTransformer(iMt);
680 newGeometryTransformer->AddModuleTransformer(newModuleTransformer);
682 TGeoCombiTrans moduleTransform =
683 TGeoCombiTrans(*kModuleTransformer->GetTransformation());
684 TGeoCombiTrans *newModuleTransform = new TGeoCombiTrans(moduleTransform);
685 // same module transform as the previous one
686 // no mis align object created
687 newModuleTransformer->SetTransformation(moduleTransform);
689 AliMpExMap *detElements = kModuleTransformer->GetDetElementStore();
693 ("%i DEs in old GeometryStore %i",
694 detElements->GetSize(), iMt));
695 Int_t iBase = !iMt ? 0 : fgSNDetElemCh[iMt-1];
696 for (Int_t iDe = 0; iDe < detElements->GetSize(); iDe++)
697 { // detection elements.
698 AliMUONGeometryDetElement *detElement =
699 (AliMUONGeometryDetElement *) detElements->GetObject(iDe);
701 AliFatal("Detection element not found.");
703 /// make a new detection element
704 AliMUONGeometryDetElement *newDetElement =
705 new AliMUONGeometryDetElement(detElement->GetId(),
706 detElement->GetVolumePath());
707 for (int i=0; i<fgNParCh; i++) {
708 AliInfo(Form("iMt %i, iBase %i, iDe %i, iPar %i",iMt, iBase, iDe, (iBase+iDe)*fgNParCh+i));
709 lDetElemMisAlignment[i] =
710 (iMt<fgNCh) ? misAlignments[(iBase+iDe)*fgNParCh+i] : 0.;
712 // local transformation of this detection element.
713 TGeoCombiTrans localTransform
714 = TGeoCombiTrans(*detElement->GetLocalTransformation());
715 TGeoCombiTrans newLocalTransform = ReAlign(localTransform,lDetElemMisAlignment);
716 newDetElement->SetLocalTransformation(newLocalTransform);
718 // global transformation
719 TGeoHMatrix newGlobalTransform =
720 AliMUONGeometryBuilder::Multiply(*newModuleTransform,
722 newDetElement->SetGlobalTransformation(newGlobalTransform);
724 // add this det element to module
725 newModuleTransformer->GetDetElementStore()->Add(newDetElement->GetId(),
727 // Get delta transformation:
728 // Tdelta = Tnew * Told.inverse
729 TGeoHMatrix deltaTransform
730 = AliMUONGeometryBuilder::Multiply(
732 detElement->GetGlobalTransformation()->Inverse());
734 // Create mis alignment matrix
735 newGeometryTransformer
736 ->AddMisAlignDetElement(detElement->GetId(), deltaTransform);
739 AliInfo(Form("Added module transformer %i to the transformer", iMt));
740 newGeometryTransformer->AddModuleTransformer(newModuleTransformer);
742 return newGeometryTransformer;