Line data Source code
1 : /**************************************************************************
2 : * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 : * *
4 : * Author: The ALICE Off-line Project. *
5 : * Contributors are mentioned in the code where appropriate. *
6 : * *
7 : * Permission to use, copy, modify and distribute this software and its *
8 : * documentation strictly for non-commercial purposes is hereby granted *
9 : * without fee, provided that the above copyright notice appears in all *
10 : * copies and that both the copyright notice and this permission notice *
11 : * appear in the supporting documentation. The authors make no claims *
12 : * about the suitability of this software for any purpose. It is *
13 : * provided "as is" without express or implied warranty. *
14 : **************************************************************************/
15 :
16 : //-----------------------------------------------------------------------------
17 : /// \class AliMUONSurveyObj
18 : /// Base class for the survey processing of the ALICE DiMuon spectrometer
19 : ///
20 : /// This base object provides methods to process the survey+photogrammetry
21 : /// data of the Chambers (frames) and Detection Elements of the DiMuon
22 : /// Spectrometer and calculate their misalignments.
23 : ///
24 : /// \author Javier Castillo
25 : //-----------------------------------------------------------------------------
26 :
27 : #include <fstream>
28 :
29 : #include "TMath.h"
30 : #include "TVector3.h"
31 : #include "TGeoMatrix.h"
32 : #include "TFitter.h"
33 : #include "TMinuit.h"
34 : #include "TString.h"
35 : #include "TH2.h"
36 : #include "TF2.h"
37 : #include "TGraph2DErrors.h"
38 : #include "TArrayD.h"
39 :
40 : #include "AliLog.h"
41 : #include "AliSurveyPoint.h"
42 :
43 : #include "AliMUONSurveyObj.h"
44 : #include "AliMUONSurveyUtil.h"
45 :
46 : void SurveyFcn(int &npar, double *g, double &f, double *par, int iflag);
47 :
48 : /// \cond CLASSIMP
49 18 : ClassImp(AliMUONSurveyObj)
50 : /// \endcond
51 :
52 : AliMUONSurveyObj::AliMUONSurveyObj()
53 0 : : TObject()
54 0 : , fSTargets(0x0)
55 0 : , fGBTargets(0x0)
56 0 : , fLBTargets(0x0)
57 0 : , fLocalTrf(0x0)
58 0 : , fAlignTrf(0x0)
59 0 : , fBaseTrf(0x0)
60 0 : , fOwnerLocalTrf(kFALSE)
61 0 : , fOwnerAlignTrf(kTRUE)
62 0 : , fOwnerBaseTrf(kFALSE)
63 0 : , fUseCM(kTRUE)
64 0 : , fPlane(0x0)
65 0 : , fFitter(0x0)
66 0 : , fXMin(-400.)
67 0 : , fXMax(400.)
68 0 : , fYMin(-400.)
69 0 : , fYMax(400.)
70 0 : , fZMin(-2000.)
71 0 : , fZMax(2000.)
72 0 : {
73 : /// Default constructor
74 :
75 0 : fSTargets = new TObjArray();
76 0 : fSTargets->SetOwner(kFALSE);
77 0 : fGBTargets = new TObjArray();
78 0 : fGBTargets->SetOwner(kFALSE);
79 0 : fLBTargets = new TObjArray();
80 0 : fLBTargets->SetOwner(kFALSE);
81 :
82 0 : fAlignTrf = new TGeoCombiTrans();
83 :
84 0 : fFitter = new TFitter(100);
85 0 : }
86 :
87 0 : AliMUONSurveyObj::~AliMUONSurveyObj() {
88 : /// Destructor
89 0 : if(fSTargets) {
90 0 : fSTargets->Delete();
91 0 : fSTargets = 0x0;
92 0 : }
93 0 : if(fGBTargets) {
94 0 : fGBTargets->Delete();
95 0 : fGBTargets = 0x0;
96 0 : }
97 0 : if(fLBTargets) {
98 0 : fLBTargets->Delete();
99 0 : fLBTargets = 0x0;
100 0 : }
101 0 : if (fPlane) {
102 0 : fPlane->Delete();
103 0 : fPlane = 0x0;
104 0 : }
105 0 : if(fOwnerLocalTrf && fLocalTrf) {
106 0 : fLocalTrf->Delete();
107 0 : fLocalTrf = 0x0;
108 0 : }
109 0 : if(fOwnerAlignTrf && fAlignTrf) {
110 0 : fAlignTrf->Delete();
111 0 : fAlignTrf = 0x0;
112 0 : }
113 0 : if(fOwnerBaseTrf && fBaseTrf) {
114 0 : fBaseTrf->Delete();
115 0 : fBaseTrf = 0x0;
116 0 : }
117 0 : if (fFitter){
118 0 : fFitter->Delete();
119 0 : fFitter = 0x0;
120 0 : }
121 0 : }
122 :
123 : void AliMUONSurveyObj::AddStickerTarget(AliSurveyPoint *stPoint){
124 : /// Add sticker target
125 0 : if (fUseCM) {
126 0 : fSTargets->Add(ConvertPointUnits(stPoint,0.1));
127 0 : } else {
128 0 : fSTargets->Add(stPoint);
129 : }
130 0 : }
131 :
132 : void AliMUONSurveyObj::AddGButtonTarget(AliSurveyPoint *btPoint){
133 : /// Add global button target
134 0 : if (fUseCM) {
135 0 : fGBTargets->Add(ConvertPointUnits(btPoint,0.1));
136 0 : } else {
137 0 : fGBTargets->Add(btPoint);
138 : }
139 0 : }
140 :
141 : void AliMUONSurveyObj::AddLButtonTarget(AliSurveyPoint *btPoint){
142 : /// Add local button target target; AliSurveyPoint
143 0 : if (fUseCM) {
144 0 : fLBTargets->Add(ConvertPointUnits(btPoint,0.1));
145 0 : } else {
146 0 : fLBTargets->Add(btPoint);
147 : }
148 0 : }
149 :
150 : void AliMUONSurveyObj::AddLButtonTarget(TVector3 *btVector){
151 : /// Add local button target target; TVector3
152 0 : fLBTargets->Add(btVector);
153 0 : }
154 :
155 : Int_t AliMUONSurveyObj::AddStickerTargets(TObjArray *pArray, TString stBaseName, Int_t lTargetMax){
156 : /// Add a maximum of lTargetMax sticker targets with stBaseName from pArray
157 0 : if (!pArray) {
158 0 : AliError(Form("Survey points array is empty %p!",pArray));
159 0 : return 0;
160 : }
161 0 : if (stBaseName.IsNull()) {
162 0 : AliError(Form("Need base name for sticker targets %s!",stBaseName.Data()));
163 0 : return 0;
164 : }
165 :
166 : Int_t stIndex = 0;
167 : AliSurveyPoint *pointSST = 0x0;
168 :
169 0 : TString stNumber;
170 :
171 0 : for (int iPoint=0; iPoint<lTargetMax; iPoint++) {
172 0 : TString stFullName(stBaseName);
173 0 : stNumber = Form("%d",iPoint+1);
174 0 : if(lTargetMax>9&&iPoint+1<10) {
175 0 : stFullName+="0";
176 : }
177 0 : stFullName+=stNumber;
178 :
179 0 : pointSST = (AliSurveyPoint *)pArray->FindObject(stFullName.Data());
180 :
181 0 : if(pointSST) {
182 0 : AddStickerTarget(pointSST);
183 0 : AliInfo(Form("Added survey sticker target %s at index %d",pointSST->GetName(),stIndex));
184 0 : stIndex++;
185 0 : }
186 0 : }
187 :
188 0 : AliInfo(Form("Found %d sticker targets with base name %s",fSTargets->GetEntries(),stBaseName.Data()));
189 : return stIndex;
190 0 : }
191 :
192 : Int_t AliMUONSurveyObj::AddGButtonTargets(TObjArray *pArray, TString btBaseName, Int_t lTargetMax){
193 : /// Add a maximum of lTargetMax global button targets with stBaseName from pArray
194 0 : printf("%s \n",btBaseName.Data());
195 0 : if (!pArray) {
196 0 : AliError(Form("Survey points array is empty %p!",pArray));
197 0 : return 0;
198 : }
199 0 : if (btBaseName.IsNull()) {
200 0 : AliError(Form("Need base name for button targets %s!",btBaseName.Data()));
201 0 : return 0;
202 : }
203 :
204 : Int_t btIndex = 0;
205 : AliSurveyPoint *pointSBT = 0x0;
206 :
207 0 : TString btNumber;
208 :
209 0 : for (int iPoint=0; iPoint<lTargetMax; iPoint++) {
210 0 : TString btFullName(btBaseName);
211 0 : btNumber = Form("%d",iPoint+1);
212 0 : if(lTargetMax>9&&iPoint+1<10) {
213 0 : btFullName+="0";
214 : }
215 0 : btFullName+=btNumber;
216 0 : printf("%s \n",btFullName.Data());
217 0 : pointSBT = (AliSurveyPoint *)pArray->FindObject(btFullName.Data());
218 :
219 0 : if(pointSBT) {
220 0 : AddGButtonTarget(pointSBT);
221 0 : AliInfo(Form("Added survey button target %s at index %d",pointSBT->GetName(),btIndex));
222 0 : btIndex++;
223 0 : }
224 0 : }
225 :
226 0 : AliInfo(Form("Found %d button targets with base name %s",fGBTargets->GetEntries(),btBaseName.Data()));
227 : return btIndex;
228 0 : }
229 :
230 : Int_t AliMUONSurveyObj::AddLButtonTargets(TObjArray *pArray, TString btBaseName, Int_t lTargetMax){
231 : /// Add a maximum of lTargetMax local button targets with stBaseName from pArray
232 0 : printf("%s \n",btBaseName.Data());
233 0 : if (!pArray) {
234 0 : AliError(Form("Local points array is empty %p!",pArray));
235 0 : return 0;
236 : }
237 0 : if (btBaseName.IsNull()) {
238 0 : AliError(Form("Need base name for button targets %s!",btBaseName.Data()));
239 0 : return 0;
240 : }
241 :
242 : Int_t btIndex = 0;
243 : AliSurveyPoint *pointSBT = 0x0;
244 :
245 0 : TString btNumber;
246 :
247 0 : for (int iPoint=0; iPoint<lTargetMax; iPoint++) {
248 0 : TString btFullName(btBaseName);
249 0 : btNumber = Form("%d",iPoint+1);
250 0 : if(lTargetMax>9&&iPoint+1<10) {
251 0 : btFullName+="0";
252 : }
253 0 : btFullName+=btNumber;
254 0 : printf("%s \n",btFullName.Data());
255 0 : pointSBT = (AliSurveyPoint *)pArray->FindObject(btFullName.Data());
256 :
257 0 : if(pointSBT) {
258 0 : AddLButtonTarget(pointSBT);
259 0 : AliInfo(Form("Added local button target %s at index %d",pointSBT->GetName(),btIndex));
260 0 : btIndex++;
261 0 : }
262 0 : }
263 :
264 0 : AliInfo(Form("Found %d local button targets with base name %s",fLBTargets->GetEntries(),btBaseName.Data()));
265 : return btIndex;
266 0 : }
267 :
268 : Int_t AliMUONSurveyObj::GetNStickerTargets() {
269 : /// return number of sticker targets
270 0 : return fSTargets->GetEntriesFast();
271 : }
272 :
273 : AliSurveyPoint* AliMUONSurveyObj::GetStickerTarget(Int_t stIndex){
274 : /// return sticker target at stIndex
275 0 : if (stIndex<0||stIndex>=fSTargets->GetEntriesFast()) {
276 0 : AliError(Form("No sticker target at index %d",stIndex));
277 0 : return 0x0;
278 : }
279 : else {
280 0 : return (AliSurveyPoint*)fSTargets->At(stIndex);
281 : }
282 0 : }
283 :
284 : Int_t AliMUONSurveyObj::GetNGButtonTargets() {
285 : /// return number of global button targets
286 0 : return fGBTargets->GetEntriesFast();
287 : }
288 :
289 : AliSurveyPoint* AliMUONSurveyObj::GetGButtonTarget(Int_t btIndex){
290 : /// return global button target at btIndex
291 0 : if (btIndex<0||btIndex>=fGBTargets->GetEntriesFast()) {
292 0 : AliError(Form("No surveyed button target at index %d",btIndex));
293 0 : return 0x0;
294 : }
295 : else {
296 0 : return (AliSurveyPoint*)fGBTargets->At(btIndex);
297 : }
298 0 : }
299 :
300 : Int_t AliMUONSurveyObj::GetNLButtonTargets() {
301 : /// return number of local button targets
302 0 : return fGBTargets->GetEntriesFast();
303 : }
304 :
305 : AliSurveyPoint* AliMUONSurveyObj::GetLButtonTarget(Int_t btIndex){
306 : /// return local button target at btIndex
307 0 : if (btIndex<0||btIndex>=fLBTargets->GetEntriesFast()) {
308 0 : AliError(Form("No surveyed button target at index %d",btIndex));
309 0 : return 0x0;
310 : }
311 : else {
312 0 : if(fLBTargets->At(btIndex)->IsA()==TVector3::Class()){
313 0 : TVector3 *lBT = (TVector3*)fLBTargets->At(btIndex);
314 0 : TString str("B");
315 0 : return (new AliSurveyPoint(TString("local"),(float)lBT->X(),(float)lBT->Y(),(float)lBT->Z(),(float)0.,(float)0.,(float)0.,'B',kTRUE));
316 0 : } else if(fLBTargets->At(btIndex)->IsA()==AliSurveyPoint::Class()) {
317 0 : AliSurveyPoint *lBT = (AliSurveyPoint*)fLBTargets->At(btIndex);
318 : return lBT;
319 : } else {
320 0 : AliError(Form("Unexpected class %s ! Valid classes are TVector3 or AliSurveyPoint",fLBTargets->At(btIndex)->ClassName()));
321 0 : return 0;
322 : }
323 : }
324 0 : }
325 :
326 : void AliMUONSurveyObj::SetPlane(TString pName, Double_t xMin, Double_t xMax, Double_t yMin, Double_t yMax){
327 : /// Set the plane function for the plane fitting
328 0 : if(fPlane) {
329 0 : fPlane->Delete();
330 0 : fPlane = 0x0;
331 0 : }
332 0 : fPlane = new TF2(pName,this,&AliMUONSurveyObj::EqPlane,xMin,xMax,yMin,yMax,3,"AliMUONSurveyObj","EqPlane");
333 0 : }
334 :
335 : void AliMUONSurveyObj::SetPlaneParameters(Double_t p0, Double_t p1, Double_t p2) {
336 : /// Set the parameters of plane function for the plane fitting
337 0 : if (!fPlane) {
338 0 : AliError("Must use SetPlane before SetPlaneParameters!!!");
339 0 : }
340 : else {
341 0 : fPlane->SetParameter(0,p0);
342 0 : fPlane->SetParameter(1,p1);
343 0 : fPlane->SetParameter(2,p2);
344 : }
345 0 : }
346 :
347 : void AliMUONSurveyObj::DrawSTargets() {
348 : /// Draw a graph of the sticker targets
349 0 : TGraph2DErrors *gST = new TGraph2DErrors(3);
350 : AliSurveyPoint *pST = 0x0;
351 0 : for (Int_t iPoint=0; iPoint<GetNStickerTargets(); iPoint++) {
352 0 : pST = GetStickerTarget(iPoint);
353 : // pST->PrintPoint();
354 0 : gST->SetPoint(iPoint,pST->GetX(),pST->GetY(),pST->GetZ());
355 0 : gST->SetPointError(iPoint,pST->GetPrecisionX(),pST->GetPrecisionY(),pST->GetPrecisionZ());
356 : }
357 0 : gST->DrawClone("P0");
358 :
359 0 : delete gST;
360 0 : }
361 :
362 : Double_t AliMUONSurveyObj::FitPlane() {
363 : /// Fit plane to sticker targets
364 0 : if (!fPlane) {
365 0 : AliError("Must use SetPlane before FitPlane!!!");
366 0 : return 0.;
367 : }
368 0 : if (fSTargets->GetEntriesFast()<3) {
369 0 : AliError("Not enough sticker targets (%d) for plane fitting!!!");
370 0 : return 0.;
371 : }
372 :
373 0 : Double_t pl[3] = {0};
374 0 : Double_t pg[3] = {0};
375 :
376 0 : TGraph2DErrors *gST = new TGraph2DErrors(3);
377 : AliSurveyPoint *pST = 0x0;
378 0 : for (Int_t iPoint=0; iPoint<GetNStickerTargets(); iPoint++) {
379 0 : pST = GetStickerTarget(iPoint);
380 : // pST->PrintPoint();
381 0 : pg[0] = pST->GetX();
382 0 : pg[1] = pST->GetY();
383 0 : pg[2] = pST->GetZ();
384 0 : fBaseTrf->MasterToLocal(pg,pl);
385 0 : gST->SetPoint(iPoint,pl[0],pl[1],pl[2]);
386 0 : printf("%d %f %f %f\n",iPoint,pl[0],pl[1],pl[2]);
387 0 : gST->SetPointError(iPoint,pST->GetPrecisionX(),pST->GetPrecisionY(),pST->GetPrecisionZ());
388 : }
389 0 : gST->Fit(fPlane);
390 :
391 0 : delete gST;
392 :
393 0 : return fPlane->GetChisquare();
394 0 : }
395 :
396 : Double_t AliMUONSurveyObj::SurveyChi2(Double_t *par){
397 : /// Returns the chisquare between local2global transform of local button targets and their surveyed position
398 0 : TGeoTranslation transTemp;
399 0 : TGeoRotation rotTemp;
400 0 : TGeoCombiTrans trfTemp;
401 :
402 : Double_t lChi2=0.;
403 :
404 0 : trfTemp.SetTranslation(transTemp);
405 0 : trfTemp.SetRotation(rotTemp);
406 0 : trfTemp.Clear();
407 0 : trfTemp.RotateZ(TMath::RadToDeg()*par[5]);
408 0 : trfTemp.RotateY(TMath::RadToDeg()*par[4]);
409 0 : trfTemp.RotateX(TMath::RadToDeg()*par[3]);
410 0 : trfTemp.SetTranslation(par[0],par[1],par[2]);
411 :
412 0 : TGeoHMatrix matGlo = (*fBaseTrf)*trfTemp;
413 0 : TGeoCombiTrans trfGlo(matGlo);
414 :
415 0 : Double_t pl[3] = {0};
416 0 : Double_t pg[3] = {0};
417 :
418 0 : for(Int_t iPoint=0; iPoint<fGBTargets->GetEntries(); iPoint++){
419 0 : AliSurveyPoint *gBT = (AliSurveyPoint*)fGBTargets->At(iPoint);
420 0 : if(fLBTargets->At(iPoint)->IsA()==TVector3::Class()){
421 0 : TVector3 *lBT = (TVector3*)fLBTargets->At(iPoint);
422 0 : pl[0]=lBT->X();
423 0 : pl[1]=lBT->Y();
424 0 : pl[2]=lBT->Z();
425 0 : } else if(fLBTargets->At(iPoint)->IsA()==AliSurveyPoint::Class()) {
426 0 : AliSurveyPoint *lBT = (AliSurveyPoint*)fLBTargets->At(iPoint);
427 0 : pl[0]=lBT->GetX();
428 0 : pl[1]=lBT->GetY();
429 0 : pl[2]=lBT->GetZ();
430 : } else {
431 0 : AliError(Form("Unexpected class %s ! Valid classes are TVector3 or AliSurveyPoint",fLBTargets->At(iPoint)->ClassName()));
432 0 : return 0;
433 : }
434 :
435 0 : trfGlo.LocalToMaster(pl, pg);
436 : // printf("%d %f %f %f\n",iPoint,pg[0],pg[1],pg[2]);
437 0 : if(fLBTargets->At(iPoint)->IsA()==TVector3::Class()){
438 0 : lChi2 += (pg[0]-gBT->GetX())*(pg[0]-gBT->GetX())/(gBT->GetPrecisionX()*gBT->GetPrecisionX());
439 0 : lChi2 += (pg[1]-gBT->GetY())*(pg[1]-gBT->GetY())/(gBT->GetPrecisionY()*gBT->GetPrecisionY());
440 0 : lChi2 += (pg[2]-gBT->GetZ())*(pg[2]-gBT->GetZ())/(gBT->GetPrecisionZ()*gBT->GetPrecisionZ());
441 0 : } else if(fLBTargets->At(iPoint)->IsA()==AliSurveyPoint::Class()) {
442 0 : AliSurveyPoint *lBT = (AliSurveyPoint*)fLBTargets->At(iPoint);
443 0 : lChi2 += (pg[0]-gBT->GetX())*(pg[0]-gBT->GetX())/(gBT->GetPrecisionX()*gBT->GetPrecisionX()+lBT->GetPrecisionX()*lBT->GetPrecisionX());
444 0 : lChi2 += (pg[1]-gBT->GetY())*(pg[1]-gBT->GetY())/(gBT->GetPrecisionY()*gBT->GetPrecisionY()+lBT->GetPrecisionY()*lBT->GetPrecisionY());
445 0 : lChi2 += (pg[2]-gBT->GetZ())*(pg[2]-gBT->GetZ())/(gBT->GetPrecisionZ()*gBT->GetPrecisionZ()+lBT->GetPrecisionZ()*lBT->GetPrecisionZ());
446 : } else {
447 0 : AliError(Form("Unexpected class %s ! Valid classes are TVector3 or AliSurveyPoint",fLBTargets->At(iPoint)->ClassName()));
448 0 : return 0;
449 : }
450 0 : }
451 :
452 0 : return lChi2;
453 0 : }
454 :
455 : //_____________________________________________________________________________
456 : void SurveyFcn(int &npar, double *g, double &f, double *par, int iflag) {
457 : ///
458 : /// Standard function as needed by Minuit-like minimization procedures.
459 : /// For the set of parameters par calculates and returns chi-squared.
460 : ///
461 :
462 : // smuggle a C++ object into a C function
463 0 : AliMUONSurveyObj *aSurveyObj = (AliMUONSurveyObj*) gMinuit->GetObjectFit();
464 :
465 0 : f = aSurveyObj->SurveyChi2(par);
466 : if (iflag==3) {}
467 0 : if (npar) {}
468 : if (g) {} // no warnings about unused stuff...
469 :
470 0 : }
471 :
472 : //_____________________________________________________________________________
473 : Int_t AliMUONSurveyObj::SurveyToAlign(TGeoCombiTrans &quadTransf, Double_t *parErr, Double_t psi, Double_t tht, Double_t epsi, Double_t etht) {
474 : /// Main function to obtain the misalignments from the surveyed position of the button targets;
475 0 : if (fGBTargets->GetEntries()!=fLBTargets->GetEntries()){
476 0 : AliError(Form("Different number of button targets: %d survey points and %d local coord!",
477 : fGBTargets->GetEntries(),fLBTargets->GetEntries()));
478 0 : return 0;
479 : }
480 :
481 0 : TFitter fitter(100);
482 0 : gMinuit->SetObjectFit(this);
483 0 : fitter.SetFCN(SurveyFcn);
484 0 : fitter.SetParameter(0,"dx",0,0.1,-20,20);
485 0 : fitter.SetParameter(1,"dy",0,0.1,-20,20);
486 0 : fitter.SetParameter(2,"dz",0,0.1,0,0);
487 0 : fitter.SetParameter(3,"rx",psi,0.0001,-0.05,0.05);
488 0 : fitter.SetParameter(4,"ry",tht,0.0001,-0.05,0.05);
489 : // fitter.SetParameter(3,"rx",psi,0.0001,psi-5*epsi,psi+5*epsi);
490 : // fitter.SetParameter(4,"ry",tht,0.0001,tht-5*etht,tht+5*etht);
491 0 : fitter.SetParameter(5,"rz",0,0.0001,0,0);
492 :
493 0 : if(psi) fitter.FixParameter(3);
494 0 : if(tht) fitter.FixParameter(4);
495 :
496 0 : double arglist[100];
497 0 : arglist[0] = 2;
498 0 : fitter.ExecuteCommand("SET PRINT", arglist, 1);
499 0 : fitter.ExecuteCommand("SET ERR", arglist, 1);
500 0 : arglist[0]=0;
501 : //fitter.ExecuteCommand("SIMPLEX", arglist, 1);
502 : // fitter.ExecuteCommand("MINIMIZE", arglist, 1);
503 0 : fitter.ExecuteCommand("MIGRAD", arglist, 1);
504 0 : fitter.ExecuteCommand("IMPROVE", arglist, 1);
505 : // fitter.ExecuteCommand("MINOS", arglist, 1);
506 : // fitter.ExecuteCommand("CALL 3", arglist,0);
507 :
508 0 : for (int j=0; j<3; j++) printf("%10.3f ",fitter.GetParameter(j));
509 0 : for (int j=3; j<6; j++) printf("%10.3f ",1000*fitter.GetParameter(j));
510 0 : printf("\n");
511 0 : for (int j=0; j<3; j++) printf("%10.3f ",fitter.GetParError(j));
512 0 : for (int j=3; j<6; j++) printf("%10.3f ",1000*fitter.GetParError(j));
513 0 : printf("\n");
514 :
515 0 : quadTransf.Clear();
516 0 : quadTransf.RotateZ(TMath::RadToDeg()*fitter.GetParameter(5));
517 0 : quadTransf.RotateY(TMath::RadToDeg()*fitter.GetParameter(4));
518 0 : quadTransf.RotateX(TMath::RadToDeg()*fitter.GetParameter(3));
519 0 : quadTransf.SetTranslation(fitter.GetParameter(0),fitter.GetParameter(1),fitter.GetParameter(2));
520 :
521 0 : for(Int_t iPar=0; iPar<6; iPar++){
522 0 : parErr[iPar] = fitter.GetParError(iPar);
523 : }
524 0 : if(epsi) parErr[3] = epsi;
525 0 : if(etht) parErr[4] = etht;
526 :
527 : return 1;
528 :
529 0 : }
530 :
531 : //_____________________________________________________________________________
532 : Int_t AliMUONSurveyObj::SurveyToAlign(Double_t psi, Double_t tht, Double_t epsi, Double_t etht) {
533 : /// Main function to obtain the misalignments from the surveyed position of the button targets;
534 0 : if (fGBTargets->GetEntries()!=fLBTargets->GetEntries()){
535 0 : AliError(Form("Different number of button targets: %d survey points and %d local coord!",
536 : fGBTargets->GetEntries(),fLBTargets->GetEntries()));
537 0 : return 0;
538 : }
539 :
540 : // TFitter fitter(100);
541 0 : gMinuit->SetObjectFit(this);
542 0 : fFitter->SetFCN(SurveyFcn);
543 0 : fFitter->SetParameter(0,"dx",0,0.1,-20,20);
544 0 : fFitter->SetParameter(1,"dy",0,0.1,-20,20);
545 0 : fFitter->SetParameter(2,"dz",0,0.1,0,0);
546 0 : if(psi)
547 0 : fFitter->SetParameter(3,"rx",psi,epsi,psi-5*epsi,psi+5*epsi);
548 : else
549 0 : fFitter->SetParameter(3,"rx",psi,0.0001,-0.05,0.05);
550 0 : if(tht)
551 0 : fFitter->SetParameter(4,"ry",tht,etht,tht-5*etht,tht+5*etht);
552 : else
553 0 : fFitter->SetParameter(4,"ry",tht,0.0001,-0.05,0.05);
554 : // fFitter->SetParameter(3,"rx",psi,0.0001,psi-5*epsi,psi+5*epsi);
555 : // fFitter->SetParameter(4,"ry",tht,0.0001,tht-5*etht,tht+5*etht);
556 0 : fFitter->SetParameter(5,"rz",0,0.0001,0,0);
557 :
558 0 : if(psi) fFitter->FixParameter(3);
559 0 : if(tht) fFitter->FixParameter(4);
560 :
561 0 : double arglist[100];
562 0 : arglist[0] = 2;
563 0 : fFitter->ExecuteCommand("SET PRINT", arglist, 1);
564 0 : fFitter->ExecuteCommand("SET ERR", arglist, 1);
565 0 : arglist[0]=0;
566 : //fFitter->ExecuteCommand("SIMPLEX", arglist, 1);
567 : // fFitter->ExecuteCommand("MINIMIZE", arglist, 1);
568 0 : fFitter->ExecuteCommand("MIGRAD", arglist, 1);
569 0 : fFitter->ExecuteCommand("IMPROVE", arglist, 1);
570 : // fFitter->ExecuteCommand("MINOS", arglist, 1);
571 : // fFitter->ExecuteCommand("CALL 3", arglist,0);
572 :
573 0 : for (int j=0; j<3; j++) printf("%10.3f ",fFitter->GetParameter(j));
574 0 : for (int j=3; j<6; j++) printf("%10.3f ",1000*fFitter->GetParameter(j));
575 0 : printf("\n");
576 0 : for (int j=0; j<3; j++) printf("%10.3f ",fFitter->GetParError(j));
577 0 : for (int j=3; j<6; j++) printf("%10.3f ",1000*fFitter->GetParError(j));
578 0 : printf("\n");
579 :
580 0 : fAlignTrf->Clear();
581 0 : fAlignTrf->RotateZ(TMath::RadToDeg()*fFitter->GetParameter(5));
582 0 : fAlignTrf->RotateY(TMath::RadToDeg()*fFitter->GetParameter(4));
583 0 : fAlignTrf->RotateX(TMath::RadToDeg()*fFitter->GetParameter(3));
584 0 : fAlignTrf->SetTranslation(fFitter->GetParameter(0),fFitter->GetParameter(1),fFitter->GetParameter(2));
585 :
586 0 : if(epsi) fFitter->ReleaseParameter(3); // To get error
587 0 : if(etht) fFitter->ReleaseParameter(4); // To get error
588 :
589 0 : TGeoCombiTrans lGlobalTrf = TGeoCombiTrans((*fBaseTrf)*(*fAlignTrf));
590 : AliSurveyPoint *pointGBT;
591 : AliSurveyPoint *pointLBT;
592 0 : Double_t pl[3] = {0};
593 0 : Double_t pg[3] = {0};
594 0 : for (int iPoint=0; iPoint<GetNGButtonTargets(); iPoint++){
595 0 : pointGBT=GetGButtonTarget(iPoint);
596 0 : pointLBT=GetLButtonTarget(iPoint);
597 0 : pl[0] = pointLBT->GetX();
598 0 : pl[1] = pointLBT->GetY();
599 0 : pl[2] = pointLBT->GetZ();
600 0 : lGlobalTrf.LocalToMaster(pl,pg);
601 0 : printf("Point %d local: %.3f %.3f %.3f\n",iPoint,pl[0],pl[1],pl[2]);
602 0 : printf("Point %d global: %.3f %.3f %.3f\n",iPoint,pg[0],pg[1],pg[2]);
603 0 : printf("Point %d survey: %.3f %.3f %.3f\n",iPoint,pointGBT->GetX(),pointGBT->GetY(),pointGBT->GetZ());
604 : }
605 :
606 : return 1;
607 :
608 0 : }
609 :
610 : Double_t AliMUONSurveyObj::EvalFunction(const TF2 *lFunction, Int_t iP1, Int_t iP2, const Char_t *lCoord) {
611 : /// Evaluate the given function at the given points for the given coordinate
612 0 : if (!lFunction) {
613 0 : AliError("No function given!!!");
614 0 : return 0;
615 : }
616 0 : AliSurveyPoint *gP1 = GetGButtonTarget(iP1);
617 0 : AliSurveyPoint *gP2 = GetGButtonTarget(iP2);
618 :
619 0 : if(!gP1||!gP2){
620 0 : AliError("Missing global button target!!!");
621 0 : return 0;
622 : }
623 :
624 : // AliInfo(Form("Function %s parameters %f %f %f %f %f %f",lFunction->GetName(),lFunction->GetParameter(0),lFunction->GetParameter(1),lFunction->GetParameter(2),lFunction->GetParameter(3),lFunction->GetParameter(4),lFunction->GetParameter(5)));
625 0 : Double_t pl1[3] = {0};
626 0 : Double_t pl2[3] = {0};
627 0 : Double_t pg1[3] = {0};
628 0 : Double_t pg2[3] = {0};
629 :
630 0 : pg1[0] = gP1->GetX();
631 0 : pg1[1] = gP1->GetY();
632 0 : pg1[2] = gP1->GetZ();
633 0 : pg2[0] = gP2->GetX();
634 0 : pg2[1] = gP2->GetY();
635 0 : pg2[2] = gP2->GetZ();
636 :
637 0 : fBaseTrf->MasterToLocal(pg1,pl1);
638 0 : fBaseTrf->MasterToLocal(pg2,pl2);
639 :
640 : Double_t lVal = 0.;
641 0 : switch (lCoord[0]) {
642 : case 'X':
643 : {
644 0 : lVal = lFunction->Eval(pl1[0],pl2[0]);
645 : // lVal = lFunction->Eval(gP1->GetX(),gP2->GetX());
646 : // AliInfo(Form("case X, lVal = %f",lVal));
647 0 : return lVal;
648 : }
649 : case 'Y':
650 : {
651 0 : lVal = lFunction->Eval(pl1[1],pl2[1]);
652 : // lVal = lFunction->Eval(gP1->GetY(),gP2->GetY());
653 : // AliInfo(Form("case Y, lVal = %f",lVal));
654 0 : return lVal;
655 : }
656 : case 'Z':
657 : {
658 0 : lVal = lFunction->Eval(pl1[2],pl2[2]);
659 : // lVal = lFunction->Eval(gP1->GetZ(),gP2->GetZ());
660 : // AliInfo(Form("case Z, lVal = %f",lVal));
661 0 : return lVal;
662 : }
663 : default:
664 : {
665 0 : AliError(Form("Coordinate %s is not valid, options are X Y Z",lCoord));
666 0 : return 0;
667 : }
668 : }
669 0 : }
670 :
671 : void AliMUONSurveyObj::CalculateTranslation(TF2 *xFunc, TF2 *yFunc, TF2 *zFunc, Int_t iP1, Int_t iP2, Double_t *lCenTemp) {
672 : /// Calculate the center translation using analytic functions
673 0 : lCenTemp[0] = EvalFunction(xFunc,iP1,iP2,"X");
674 0 : lCenTemp[1] = EvalFunction(yFunc,iP1,iP2,"Y");
675 0 : lCenTemp[2] = EvalFunction(zFunc,iP1,iP2,"Z");
676 :
677 0 : }
678 :
679 : Double_t AliMUONSurveyObj::CalculateGlobalDiff(TGeoCombiTrans &lTransf, Int_t nPoints, TArrayD &lDiff){
680 : /// By hand computation of distance between local2global transform of target position and its surveyed position
681 0 : if (nPoints > GetNGButtonTargets()) {
682 0 : nPoints = GetNGButtonTargets();
683 0 : }
684 :
685 0 : for(Int_t iVal=0; iVal<nPoints*(3+1)+1; iVal++){
686 0 : lDiff[iVal] = 0.;
687 : }
688 :
689 0 : Double_t pl[3] = {0};
690 0 : Double_t pg[3] = {0};
691 0 : Double_t pml[3] = {0};
692 0 : Double_t pmg[3] = {0};
693 : AliSurveyPoint *gBT = 0x0;
694 0 : for(Int_t iPoint=0; iPoint<nPoints; iPoint++){
695 0 : gBT = GetGButtonTarget(iPoint);
696 0 : if(!gBT||!fLBTargets->At(iPoint)){
697 0 : AliError(Form("The local or global target %d is missing!",iPoint));
698 0 : lDiff[nPoints*(3+1)] = 1.e7;
699 0 : return lDiff[nPoints*(3+1)];
700 : }
701 0 : if(fLBTargets->At(iPoint)->IsA()==TVector3::Class()){
702 0 : TVector3 *lBT = (TVector3*)fLBTargets->At(iPoint);
703 0 : pl[0]=lBT->X();
704 0 : pl[1]=lBT->Y();
705 0 : pl[2]=lBT->Z();
706 0 : } else if(fLBTargets->At(iPoint)->IsA()==AliSurveyPoint::Class()) {
707 0 : AliSurveyPoint *lBT = (AliSurveyPoint*)fLBTargets->At(iPoint);
708 0 : pl[0]=lBT->GetX();
709 0 : pl[1]=lBT->GetY();
710 0 : pl[2]=lBT->GetZ();
711 : } else {
712 0 : AliError(Form("Unexpected class %s ! Valid classes are TVector3 or AliSurveyPoint",fLBTargets->At(iPoint)->ClassName()));
713 0 : return 0;
714 : }
715 :
716 0 : lTransf.LocalToMaster(pl,pg);
717 0 : pmg[0] = gBT->GetX();
718 0 : pmg[1] = gBT->GetY();
719 0 : pmg[2] = gBT->GetZ();
720 0 : fBaseTrf->MasterToLocal(pmg,pml);
721 : // printf("l %d %f %f %f\n",iPoint,pl[0],pl[1],pl[2]);
722 : // printf("g %d %f %f %f\n",iPoint,pg[0],pg[1],pg[2]);
723 : // printf("ml %d %f %f %f\n",iPoint,pml[0],pml[1],pml[2]);
724 : // printf("mg %d %f %f %f\n",iPoint,gBT->GetX(),gBT->GetY(),gBT->GetZ());
725 0 : lDiff[iPoint*(3+1)+0] = (pml[0]-pg[0]);
726 0 : lDiff[iPoint*(3+1)+1] = (pml[1]-pg[1]);
727 0 : lDiff[iPoint*(3+1)+2] = (pml[2]-pg[2]);
728 :
729 0 : lDiff[iPoint*(3+1)+3] = TMath::Sqrt(lDiff[iPoint*(3+1)+0]*lDiff[iPoint*(3+1)+0]+
730 0 : lDiff[iPoint*(3+1)+1]*lDiff[iPoint*(3+1)+1]+
731 0 : lDiff[iPoint*(3+1)+2]*lDiff[iPoint*(3+1)+2]);
732 :
733 0 : lDiff[nPoints*(3+1)] += lDiff[iPoint*(3+1)+3]*lDiff[iPoint*(3+1)+3];
734 : }
735 :
736 0 : lDiff[nPoints*(3+1)] = TMath::Sqrt(lDiff[nPoints*(3+1)]);
737 0 : return lDiff[nPoints*(3+1)];
738 0 : }
739 :
740 : Int_t AliMUONSurveyObj::CalculateBestTransf(Int_t iP1, Int_t iP2, Double_t *lXYZ, Double_t *lPTP) {
741 : /// By hand calculation of the best local to global transform using 2 button targets
742 0 : Double_t lPsi = lPTP[0];
743 0 : Double_t lTht = lPTP[1];
744 :
745 : Double_t pl1[3] = {0};
746 : Double_t pl2[3] = {0};
747 :
748 0 : if(!fLBTargets->At(iP1)||!fLBTargets->At(iP2)){
749 0 : AliError(Form("Local target %d or %d is missing!",iP1,iP2));
750 0 : return 0;
751 : }
752 :
753 0 : if(fLBTargets->At(iP1)->IsA()==TVector3::Class()){
754 0 : TVector3 *lBT1 = (TVector3*)fLBTargets->At(iP1);
755 0 : pl1[0]=lBT1->X();
756 0 : pl1[1]=lBT1->Y();
757 0 : pl1[2]=lBT1->Z();
758 0 : } else if(fLBTargets->At(iP1)->IsA()==AliSurveyPoint::Class()) {
759 0 : AliSurveyPoint *lBT1 = (AliSurveyPoint*)fLBTargets->At(iP1);
760 0 : pl1[0]=lBT1->GetX();
761 0 : pl1[1]=lBT1->GetY();
762 0 : pl1[2]=lBT1->GetZ();
763 : } else {
764 0 : AliError(Form("Unexpected class %s ! Valid classes are TVector3 or AliSurveyPoint",fLBTargets->At(iP1)->ClassName()));
765 0 : return 0;
766 : }
767 0 : if(fLBTargets->At(iP2)->IsA()==TVector3::Class()){
768 0 : TVector3 *lBT2 = (TVector3*)fLBTargets->At(iP2);
769 0 : pl2[0]=lBT2->X();
770 0 : pl2[1]=lBT2->Y();
771 0 : pl2[2]=lBT2->Z();
772 0 : } else if(fLBTargets->At(iP2)->IsA()==AliSurveyPoint::Class()) {
773 0 : AliSurveyPoint *lBT2 = (AliSurveyPoint*)fLBTargets->At(iP2);
774 0 : pl2[0]=lBT2->GetX();
775 0 : pl2[1]=lBT2->GetY();
776 0 : pl2[2]=lBT2->GetZ();
777 : } else {
778 0 : AliError(Form("Unexpected class %s ! Valid classes are TVector3 or AliSurveyPoint",fLBTargets->At(iP2)->ClassName()));
779 0 : return 0;
780 : }
781 :
782 :
783 0 : AliMUONSurveyUtil *surveyUtil = AliMUONSurveyUtil::Instance();
784 :
785 : // Xcenter functions
786 : const char *fxcName = "fXcn00";
787 0 : TF2 **fXc = new TF2*[2];
788 : fxcName = "fXcn";
789 0 : fXc[0] = new TF2(fxcName,surveyUtil,&AliMUONSurveyUtil::XnCenter,fXMin,fXMax,fYMin,fYMax,7,"AliMUONSurveyUtil","XnCenter");
790 : fxcName = "fXcp";
791 0 : fXc[1] = new TF2(fxcName,surveyUtil,&AliMUONSurveyUtil::XpCenter,fXMin,fXMax,fYMin,fYMax,7,"AliMUONSurveyUtil","XpCenter");
792 :
793 : // Ycenter functions
794 : const char *fycName = "fYcn00";
795 0 : TF2 **fYc = new TF2*[2];
796 : fycName = "fYcn";
797 0 : fYc[0] = new TF2(fycName,surveyUtil,&AliMUONSurveyUtil::YnCenter,fYMin,fYMax,fYMin,fYMax,8,"AliMUONSurveyUtil","YnCenter");
798 : fycName = "fYcp";
799 0 : fYc[1] = new TF2(fycName,surveyUtil,&AliMUONSurveyUtil::YpCenter,fYMin,fYMax,fYMin,fYMax,8,"AliMUONSurveyUtil","YpCenter");
800 :
801 : // Zcenter functions
802 : const char *fzcName = "fZcn00";
803 0 : TF2 **fZc = new TF2*[2];
804 : fzcName = "fZcn";
805 0 : fZc[0] = new TF2(fzcName,surveyUtil,&AliMUONSurveyUtil::ZnCenter,fZMin,fZMax,fZMin,fZMax,8,"AliMUONSurveyUtil","ZnCenter");
806 : fzcName = "fZcp";
807 0 : fZc[1] = new TF2(fzcName,surveyUtil,&AliMUONSurveyUtil::ZpCenter,fZMin,fZMax,fZMin,fZMax,8,"AliMUONSurveyUtil","ZpCenter");
808 :
809 : // Phi rotation using xglobal coords functions
810 : const char *fphixName = "fPhiXnn00";
811 0 : TF2 ***fPhiX = new TF2**[2];
812 0 : for (Int_t iX =0; iX<2; iX++) {
813 0 : fPhiX[iX] = new TF2*[2];
814 : }
815 : fphixName = "fPhiXnn";
816 0 : fPhiX[0][0] = new TF2(fphixName,surveyUtil,&AliMUONSurveyUtil::PhiXnn,fXMin,fXMax,fXMin,fXMax,7,"AliMUONSurveyUtil","PhiXnn");
817 : fphixName = "fPhiXnp";
818 0 : fPhiX[0][1] = new TF2(fphixName,surveyUtil,&AliMUONSurveyUtil::PhiXnp,fXMin,fXMax,fXMin,fXMax,7,"AliMUONSurveyUtil","PhiXnp");
819 : fphixName = "fPhiXpn";
820 0 : fPhiX[1][0] = new TF2(fphixName,surveyUtil,&AliMUONSurveyUtil::PhiXpn,fXMin,fXMax,fXMin,fXMax,7,"AliMUONSurveyUtil","PhiXpn");
821 : fphixName = "fPhiXpp";
822 0 : fPhiX[1][1] = new TF2(fphixName,surveyUtil,&AliMUONSurveyUtil::PhiXpp,fXMin,fXMax,fXMin,fXMax,7,"AliMUONSurveyUtil","PhiXpp");
823 :
824 : // Phi rotation using yglobal coords functions
825 : const char *fphiyName = "fPhiYnn00";
826 0 : TF2 ***fPhiY = new TF2**[2];
827 0 : for (Int_t iY =0; iY<2; iY++) {
828 0 : fPhiY[iY] = new TF2*[2];
829 : }
830 : fphiyName = "fPhiYnn";
831 0 : fPhiY[0][0] = new TF2(fphiyName,surveyUtil,&AliMUONSurveyUtil::PhiYnn,fYMin,fYMax,fYMin,fYMax,8,"AliMUONSurveyUtil","PhiYnn");
832 : fphiyName = "fPhiYnp";
833 0 : fPhiY[0][1] = new TF2(fphiyName,surveyUtil,&AliMUONSurveyUtil::PhiYnp,fYMin,fYMax,fYMin,fYMax,8,"AliMUONSurveyUtil","PhiYnp");
834 : fphiyName = "fPhiYpn";
835 0 : fPhiY[1][0] = new TF2(fphiyName,surveyUtil,&AliMUONSurveyUtil::PhiYpn,fYMin,fYMax,fYMin,fYMax,8,"AliMUONSurveyUtil","PhiYpn");
836 : fphiyName = "fPhiYpp";
837 0 : fPhiY[1][1] = new TF2(fphiyName,surveyUtil,&AliMUONSurveyUtil::PhiYpp,fYMin,fYMax,fYMin,fYMax,8,"AliMUONSurveyUtil","PhiYpp");
838 :
839 :
840 : // Set Parameters of functions
841 0 : for(Int_t iS=0; iS<2; iS++){
842 0 : fXc[iS]->SetParameters(pl1[0],pl1[1],pl1[2],pl2[0],pl2[1],pl2[2],lTht);
843 0 : fYc[iS]->SetParameters(pl1[0],pl1[1],pl1[2],pl2[0],pl2[1],pl2[2],lPsi,lTht);
844 0 : fZc[iS]->SetParameters(pl1[0],pl1[1],pl1[2],pl2[0],pl2[1],pl2[2],lPsi,lTht);
845 : // fXc[iS]->SetParameters(lBT1->X(),lBT1->Y(),lBT1->Z(),lBT2->X(),lBT2->Y(),lBT2->Z(),lTht);
846 : // fYc[iS]->SetParameters(lBT1->X(),lBT1->Y(),lBT1->Z(),lBT2->X(),lBT2->Y(),lBT2->Z(),lPsi,lTht);
847 : // fZc[iS]->SetParameters(lBT1->X(),lBT1->Y(),lBT1->Z(),lBT2->X(),lBT2->Y(),lBT2->Z(),lPsi,lTht);
848 0 : for(Int_t jS=0; jS<2; jS++){
849 0 : fPhiX[iS][jS]->SetParameters(pl1[0],pl1[1],pl1[2],pl2[0],pl2[1],pl2[2],lTht);
850 0 : fPhiY[iS][jS]->SetParameters(pl1[0],pl1[1],pl1[2],pl2[0],pl2[1],pl2[2],lPsi,lTht);
851 : // fPhiX[iS][jS]->SetParameters(lBT1->X(),lBT1->Y(),lBT1->Z(),lBT2->X(),lBT2->Y(),lBT2->Z(),lTht);
852 : // fPhiY[iS][jS]->SetParameters(lBT1->X(),lBT1->Y(),lBT1->Z(),lBT2->X(),lBT2->Y(),lBT2->Z(),lPsi,lTht);
853 : }
854 : }
855 :
856 0 : Double_t lCenTemp[3];
857 : Double_t lRotTemp[3];
858 :
859 0 : TGeoCombiTrans trfTemp;
860 :
861 0 : Int_t nPoints = GetNGButtonTargets();
862 :
863 0 : TArrayD lDiffTemp(nPoints*(3+1)+1);
864 0 : TArrayD lDiffMin(nPoints*(3+1)+1);
865 :
866 0 : for(Int_t i=0; i<nPoints*(3+1)+1; i++){
867 0 : lDiffMin[i]=1000000.;
868 0 : lDiffTemp[i]=0.;
869 : }
870 :
871 : //
872 : // Calculate Detection Element Center from button targets
873 : //
874 :
875 : // Trying 2x*2y*2z*(2phi+2phi) possibilities
876 0 : for(Int_t iX=0; iX<2; iX++){
877 0 : for(Int_t iY=0; iY<2; iY++){
878 0 : for(Int_t iZ=0; iZ<2; iZ++){
879 0 : CalculateTranslation(fXc[iX],fYc[iY],fZc[iZ],iP1,iP2,lCenTemp);
880 :
881 : lRotTemp[0] = lPsi;
882 : lRotTemp[1] = lTht;
883 0 : for(Int_t iP=0; iP<2; iP++){
884 0 : lRotTemp[2] = EvalFunction(fPhiX[iX][iP],iP1,iP2,"X");
885 :
886 0 : trfTemp.Clear();
887 0 : trfTemp.RotateZ(TMath::RadToDeg()*lRotTemp[2]);
888 0 : trfTemp.RotateY(TMath::RadToDeg()*lRotTemp[1]);
889 0 : trfTemp.RotateX(TMath::RadToDeg()*lRotTemp[0]);
890 0 : trfTemp.SetTranslation(lCenTemp);
891 :
892 0 : if(CalculateGlobalDiff(trfTemp,nPoints,lDiffTemp)<lDiffMin[nPoints*(3+1)]){
893 0 : printf("Diffs");
894 0 : for(Int_t i=0; i<nPoints*(3+1)+1; i++){
895 0 : printf(" %f",lDiffTemp[i]);
896 : }
897 0 : printf("\n");
898 0 : printf(" : mycenX%dY%dZ%d(%f,%f,%f); rotx%d(%f,%f,%f)\n",iX,iY,iZ,lCenTemp[0],lCenTemp[1],lCenTemp[2],iP,lRotTemp[0],lRotTemp[1],lRotTemp[2]);
899 0 : printf("Transformation improved ...\n");
900 0 : for (int i=0; i<3; i++) {
901 0 : lXYZ[i] = lCenTemp[i];
902 : }
903 0 : lPTP[2] = lRotTemp[2];
904 0 : for(Int_t i=0; i<nPoints*(3+1)+1; i++){
905 0 : lDiffMin[i]=lDiffTemp[i];
906 : }
907 0 : }
908 : }
909 0 : for(Int_t iP=0; iP<2; iP++){
910 0 : lRotTemp[2] = EvalFunction(fPhiY[iY][iP],iP1,iP2,"Y");
911 :
912 0 : trfTemp.Clear();
913 0 : trfTemp.RotateZ(TMath::RadToDeg()*lRotTemp[2]);
914 0 : trfTemp.RotateY(TMath::RadToDeg()*lRotTemp[1]);
915 0 : trfTemp.RotateX(TMath::RadToDeg()*lRotTemp[0]);
916 0 : trfTemp.SetTranslation(lCenTemp);
917 :
918 0 : if(CalculateGlobalDiff(trfTemp,nPoints,lDiffTemp)<lDiffMin[nPoints*(3+1)]){
919 0 : printf("Diffs");
920 0 : for(Int_t i=0; i<nPoints*(3+1)+1; i++){
921 0 : printf(" %f",lDiffTemp[i]);
922 : }
923 0 : printf("\n");
924 0 : printf(" : mycenX%dY%dZ%d(%f,%f,%f); roty%d(%f,%f,%f)\n",iX,iY,iZ,lCenTemp[0],lCenTemp[1],lCenTemp[2],iP,lRotTemp[0],lRotTemp[1],lRotTemp[2]);
925 0 : printf("Transformation improved ...\n");
926 0 : for (int i=0; i<3; i++) {
927 0 : lXYZ[i] = lCenTemp[i];
928 : }
929 0 : lPTP[2] = lRotTemp[2];
930 0 : for(Int_t i=0; i<nPoints*(3+1)+1; i++){
931 0 : lDiffMin[i]=lDiffTemp[i];
932 : }
933 0 : }
934 : }
935 : }
936 : }
937 : }
938 :
939 0 : for (Int_t i=0; i<2; i++) {
940 0 : delete fXc[i];
941 0 : delete fYc[i];
942 0 : delete fZc[i];
943 0 : for (Int_t j=0; j<2; j++) {
944 0 : delete fPhiX[i][j];
945 0 : delete fPhiY[i][j];
946 : }
947 0 : delete[] fPhiX[i];
948 0 : delete[] fPhiY[i];
949 : }
950 0 : delete[] fXc;
951 0 : delete[] fYc;
952 0 : delete[] fZc;
953 0 : delete[] fPhiX;
954 0 : delete[] fPhiY;
955 :
956 0 : if (lDiffMin[nPoints*(3+1)]>20) return 0;
957 :
958 0 : return 1;
959 0 : }
960 :
961 : void AliMUONSurveyObj::CalculateMeanTransf(Double_t *lXYZ, Double_t *lPTP) {
962 : /// By hand calculation of the mean (for nPairs of targets) of the best local to global transform using 2 button targets
963 : Double_t xce=0.;
964 : Double_t yce=0.;
965 : Double_t zce=0.;
966 : Double_t phi=0.;
967 :
968 : Int_t nPairs = 0;
969 0 : Int_t nPoints = GetNGButtonTargets();
970 : // Loop over all possible pairs of button tragets
971 0 : for(Int_t iP1=0; iP1<nPoints; iP1++){
972 0 : for(Int_t iP2=iP1+1; iP2<nPoints; iP2++){
973 0 : printf("%d and %d\n",iP1,iP2);
974 :
975 0 : if(CalculateBestTransf(iP1,iP2,lXYZ,lPTP)) {
976 0 : nPairs++;
977 :
978 0 : xce+=lXYZ[0];
979 0 : yce+=lXYZ[1];
980 0 : zce+=lXYZ[2];
981 0 : phi+=lPTP[2];
982 0 : }
983 : }
984 : }
985 :
986 0 : if (!nPairs) return;
987 :
988 0 : lXYZ[0]=xce/nPairs;
989 0 : lXYZ[1]=yce/nPairs;
990 0 : lXYZ[2]=zce/nPairs;
991 0 : lPTP[2]=phi/nPairs;
992 0 : }
993 :
994 : void AliMUONSurveyObj::PrintLocalTrf() {
995 : /// Print the local transformation
996 0 : Double_t lRotTemp[3];
997 0 : AliMUONSurveyUtil::MatrixToAngles(fLocalTrf->GetRotationMatrix(),lRotTemp);
998 0 : printf("(%.3f %.3f %.3f), (%.6f %.6f %.6f)\n",fLocalTrf->GetTranslation()[0],fLocalTrf->GetTranslation()[1],fLocalTrf->GetTranslation()[2],lRotTemp[0],lRotTemp[1],lRotTemp[2]);
999 0 : }
1000 :
1001 : void AliMUONSurveyObj::PrintAlignTrf() {
1002 : /// Print the alignment transformation
1003 0 : Double_t lRotTemp[3];
1004 0 : AliMUONSurveyUtil::MatrixToAngles(fAlignTrf->GetRotationMatrix(),lRotTemp);
1005 0 : printf("(%.3f %.3f %.3f), (%.6f %.6f %.6f)\n",fAlignTrf->GetTranslation()[0],fAlignTrf->GetTranslation()[1],fAlignTrf->GetTranslation()[2],lRotTemp[0],lRotTemp[1],lRotTemp[2]);
1006 0 : }
1007 :
1008 : void AliMUONSurveyObj::FillSTHistograms(TString baseNameC, TH2 *hSTc, TString baseNameA, TH2 *hSTa) {
1009 : /// Fill sticker target histograms for monitoring
1010 0 : if(baseNameC.IsNull()||!hSTc){
1011 0 : AliError("Need base name for points on side C and/or a histogram for them!");
1012 0 : return;
1013 : }
1014 : AliSurveyPoint *pointST = 0x0;
1015 0 : for (Int_t iPoint=0; iPoint<GetNStickerTargets(); iPoint++) {
1016 0 : pointST = GetStickerTarget(iPoint);
1017 0 : if (!pointST) continue;
1018 0 : if (pointST->GetPointName().Contains(baseNameC)){
1019 0 : hSTc->Fill(pointST->GetX(),pointST->GetY(),-pointST->GetZ());
1020 0 : } else if ((!baseNameA.IsNull()) &&
1021 0 : (pointST->GetPointName().Contains(baseNameA))) {
1022 0 : if (!hSTa){
1023 0 : AliError("Base name for points on side A provided but no histogram for them!");
1024 0 : continue;
1025 : }
1026 0 : hSTa->Fill(pointST->GetX(),pointST->GetY(),-pointST->GetZ());
1027 0 : }
1028 : }
1029 0 : }
1030 :
1031 : Double_t AliMUONSurveyObj::GetAlignResX() {
1032 : /// Returns the uncertainty of the x translation parameter
1033 0 : if(!fFitter) {
1034 0 : AliError("There is no fitter for this object! X resolution will be 0.");
1035 0 : return 0.;
1036 : }
1037 0 : return fFitter->GetParError(0);
1038 0 : }
1039 :
1040 : Double_t AliMUONSurveyObj::GetAlignResY() {
1041 : /// Returns the uncertainty of the y translation parameter
1042 0 : if(!fFitter) {
1043 0 : AliError("There is no fitter for this object! Y resolution will be 0.");
1044 0 : return 0.;
1045 : }
1046 0 : return fFitter->GetParError(1);
1047 0 : }
1048 :
1049 : AliSurveyPoint* AliMUONSurveyObj::ConvertPointUnits(AliSurveyPoint *stPoint, Float_t lFactor) {
1050 : /// Return the AliSurveyPoint with new units. Default is from mm -> cm
1051 0 : return new AliSurveyPoint(stPoint->GetPointName(),
1052 0 : lFactor*stPoint->GetX(),lFactor*stPoint->GetY(),lFactor*stPoint->GetZ(),
1053 0 : lFactor*stPoint->GetPrecisionX(),lFactor*stPoint->GetPrecisionY(),lFactor*stPoint->GetPrecisionZ(),
1054 0 : stPoint->GetType(), stPoint->GetTarget());
1055 0 : }
|