]> git.uio.no Git - u/mrichter/AliRoot.git/blame - MUON/MUONSurveyCh1.C
Using send pressure sensor by default
[u/mrichter/AliRoot.git] / MUON / MUONSurveyCh1.C
CommitLineData
ba8b0266 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id:$ */
17
18/// \ingroup macros
52ce8802 19/// \file MUONSurveyCh1.C
ba8b0266 20/// \brief Macro to process survey and photogrammetry data of chamber 1
21///
22/// Macro loads the survey data from .txt file using AliSurveyObj.
23/// Macro uses AliMUONSurvey... classes.
24/// The transformations of the detection elements or chambers can be obtained
25/// in three ways:
26/// - A: Fit of local to global transformation using the fixed button targets.
27/// - B: Fit a plane to the sticker targets -> psi, theta
28/// Use above psi and theta and fit remaining 4 parameters using the fixed
29/// button targets
30/// - C: Fit a plane to the sticker targets -> psi, theta
31/// Use above psi and theta to calculate xc, yc, zc and phi by solving
32/// the equations from a local to global transformation of the fixed button
33/// targets
34///
35/// Methods A and B are prefered to C, and B is better if sticker targets are
36/// available and lie on a plane!
37/// For slats only methods B and C can be used.
38/// Various histograms are filled and printed for monitoring.
39/// MisAlignment object is then created.
40///
41/// \author Javier Castillo
42
43#if !defined(__CINT__) || defined(__MAKECINT__)
44
45#include "AliMUONGeometryTransformer.h"
46#include "AliMUONGeometryModuleTransformer.h"
47#include "AliMUONGeometryDetElement.h"
48#include "AliMUONGeometryMisAligner.h"
49#include "AliMUONSurveyChamber.h"
50#include "AliMUONSurveyDetElem.h"
51#include "AliMUONSurveyUtil.h"
52
53#include "AliSurveyObj.h"
54#include "AliSurveyPoint.h"
55#include "AliGeomManager.h"
56#include "AliCDBManager.h"
57#include "AliCDBMetaData.h"
58#include "AliCDBId.h"
59#include "AliAlignObjMatrix.h"
60#include "AliAlignObj.h"
61
62#include <TROOT.h>
63#include <TGeoManager.h>
64#include <TClonesArray.h>
65#include <TObjArray.h>
66#include <TArrayD.h>
67#include <TMath.h>
68#include <TString.h>
69#include <TFitter.h>
70#include <TH2.h>
71#include <TF2.h>
72#include <TGraphErrors.h>
73#include <TCanvas.h>
74#include <TPad.h>
75#include <TLine.h>
76#include <TPostScript.h>
77#include <TPaveLabel.h>
78#include <TStyle.h>
79#include <TFile.h>
80#include <TMatrixDSym.h>
81
82#include <fstream>
83
84#endif
85
86void MUONSurveyCh1() {
87
88 TString str;
89 TString title;
90
91 Bool_t bMonitor = kTRUE;
92 Bool_t bOCDB = kTRUE;
93 Bool_t saveps = kFALSE;
94 const int cWidth = (int)(700*(29./21.));
95 const int cHeight = 700;
96 const int filetype = 112; // landscape
97
98 Int_t chId = 0;
99 Int_t nChs = 2;
100 Int_t nDetElems = 4;
101 Int_t nDetElemsI = 2;
102 // Int_t nDetElemsO = 2;
103 Int_t iDetElemToDetElemId[4] = {100,101,102,103};
104 Int_t iDetElemPseudoIdToDetElem[4] = {0,1,2,3};
105 Int_t iDetElemsOfChamber[2][2] = {{0,3},
106 {1,2}};
107
108 TObjArray *myChamberArray = new TObjArray();
109 myChamberArray->Add(new AliMUONSurveyChamber(chId));
110 myChamberArray->Add(new AliMUONSurveyChamber(chId));
111
112 AliMUONSurveyChamber *myChamberI = 0x0;
113 AliMUONSurveyChamber *myChamberO = 0x0;
114 AliMUONSurveyChamber *myChamber = 0x0;
115 AliMUONSurveyDetElem *myDetElem = 0x0;
116
117 myChamberI = (AliMUONSurveyChamber*)myChamberArray->At(0);
118 myChamberI->GetSurveyObj()->FillFromLocalFile("~/phgwrk/Reports/AliceSt1_TC1_1171d.txt");
119 myChamberO = (AliMUONSurveyChamber*)myChamberArray->At(1);
120 myChamberO->GetSurveyObj()->FillFromLocalFile("~/phgwrk/Reports/AliceSt1_TC1_1171d.txt");
121
122 myChamber = myChamberI;
123 myChamber->PrintSurveyReport();
124
125 // Chamber & DetElems button targets local coordinates
126 AliSurveyObj *lSO = new AliSurveyObj();
127 lSO->FillFromLocalFile("$ALICE_ROOT/MUON/data/MUONTrackingTargetsLocal.txt");
128
129 // Set survey targets for chambers
130 myChamberI->AddGButtonTargets(TString("T4_520"),3);
131 myChamberI->AddStickerTargets(TString("T4_620"),6);
132 myChamberI->AddLButtonTargets(lSO->GetData(),Form("%d_520",chId+1),3);
133 myChamberI->AddGButtonTargets(TString("T4_510"),3);
134 myChamberI->AddStickerTargets(TString("T4_610"),6);
135 myChamberI->AddLButtonTargets(lSO->GetData(),Form("%d_510",chId+1),3);
136 myChamberO->AddGButtonTargets(TString("T4_520"),3);
137 myChamberO->AddStickerTargets(TString("T4_620"),6);
138 myChamberO->AddLButtonTargets(lSO->GetData(),Form("%d_520",chId+1),3);
139 myChamberO->AddGButtonTargets(TString("T4_510"),3);
140 myChamberO->AddStickerTargets(TString("T4_610"),6);
141 myChamberO->AddLButtonTargets(lSO->GetData(),Form("%d_510",chId+1),3);
142
143 // Set survey targets for detection elements
144 for (int iCh =0; iCh<=1; iCh++) {
145 myChamber = (AliMUONSurveyChamber*)myChamberArray->At(iCh);
146 for (int iDetElemI=0; iDetElemI<nDetElemsI; iDetElemI++){
147 Int_t iDetElem = iDetElemsOfChamber[iCh][iDetElemI];
148 myChamber->AddSurveyDetElem(iDetElemToDetElemId[iDetElem]);
149 TString baseName;
150 myDetElem = myChamber->GetDetElem(iDetElemI);
151 if (myDetElem) {
152 baseName = Form("T4_53%d",iDetElem+1);
153 myDetElem->AddGButtonTargets(baseName,3);
154 baseName = Form("%d_53%d",chId+1,iDetElem+1);
155 myDetElem->AddLButtonTargets(lSO->GetData(),baseName,3);
156 }
157 }
158 }
159 printf("All targets added! \n");
160
161 Double_t **lCenDetElem = new Double_t*[nDetElems+nChs];
162 Double_t **lRotDetElem = new Double_t*[nDetElems+nChs];
163 Double_t **lDiffCenDetElem0 = new Double_t*[nDetElems+nChs];
164 Double_t **lDiffRotDetElem0 = new Double_t*[nDetElems+nChs];
165 Double_t **lDiffThCenDetElem0 = new Double_t*[nDetElems+nChs];
166 Double_t **lDiffThRotDetElem0 = new Double_t*[nDetElems+nChs];
167 Double_t **lDeltaDiffCenDetElem0 = new Double_t*[nDetElems+nChs];
168 Double_t **lDeltaDiffRotDetElem0 = new Double_t*[nDetElems+nChs];
169
170 for (int iDetElem=0; iDetElem<nDetElems+nChs; iDetElem++) {
171 lCenDetElem[iDetElem] = new Double_t[3];
172 lRotDetElem[iDetElem] = new Double_t[3];
173 lDiffCenDetElem0[iDetElem] = new Double_t[3];
174 lDiffRotDetElem0[iDetElem] = new Double_t[3];
175 lDiffThCenDetElem0[iDetElem] = new Double_t[3];
176 lDiffThRotDetElem0[iDetElem] = new Double_t[3];
177 lDeltaDiffCenDetElem0[iDetElem] = new Double_t[3];
178 lDeltaDiffRotDetElem0[iDetElem] = new Double_t[3];
179 }
180
181 TGeoCombiTrans dtrfDetElem[nDetElems+nChs];
182 TGeoCombiTrans localTrfDetElem[nDetElems+nChs];
183 TGeoCombiTrans localTrfThDetElem[nDetElems+nChs];
184
185 // Import TGeo geometry
186 char* geoFilename = "geometry.root";
187 if ( ! AliGeomManager::GetGeometry() ) {
188 AliGeomManager::LoadGeometry(geoFilename);
189 if (! AliGeomManager::GetGeometry() ) {
190 printf("MUONSurveyCh%d: getting geometry from file %s failed\n", chId+1,geoFilename);
191 return;
192 }
193 cout << "geometry imported" << endl;
194 }
195
196 AliMUONGeometryTransformer *transform = new AliMUONGeometryTransformer();
197 transform->LoadGeometryData();
198
199 TGeoCombiTrans trfThChamber;
200 TGeoCombiTrans trfThDetElem;
201
202 for (int iCh =0; iCh<=1; iCh++) {
203 myChamber = (AliMUONSurveyChamber*)myChamberArray->At(iCh);
204
205 trfThChamber = TGeoCombiTrans(*transform->GetModuleTransformerByDEId(iDetElemToDetElemId[iCh*nDetElemsI])->GetTransformation());
ba8b0266 206 trfThChamber.Print();
207 myChamber->SetLocalTransformation(new TGeoCombiTrans(trfThChamber),kTRUE);
208
209 // Set Chamber plane function
210 cout << "Setting plane for Chamber" << iCh+1 << " ..." << endl;
211 myChamber->SetPlane(Form("fChamber%d",iCh+1));
212 myChamber->SetPlaneParameters(0.,0.,0.);
213
214 // Fit a plane to sticker targets
215 Double_t lCChi2 = myChamber->FitPlane();
216 cout << "... chi2 = " << lCChi2 << " ..." << endl;
217
218 // Get best transformation from fitting method
219 // (local to global transformation)
220 cout << "Trying fitting method for chamber " << iCh << endl;
221 myChamber->SurveyToAlign();
222 // myChamber->SurveyToAlign(myChamber->GetPlane()->GetParameter(0),myChamber->GetPlane()->GetParameter(1),myChamber->GetPlane()->GetParError(0),myChamber->GetPlane()->GetParError(1));
223
224 for (int iDetElemI=0; iDetElemI<nDetElemsI; iDetElemI++){
225 myDetElem = myChamber->GetDetElem(iDetElemI);
226 Int_t iDetElem = iDetElemsOfChamber[iCh][iDetElemI];
227
228 trfThDetElem.Clear();
229 trfThDetElem = TGeoCombiTrans(*transform->GetDetElement(iDetElemToDetElemId[iDetElem])->GetLocalTransformation());
ba8b0266 230 trfThDetElem.Print();
231
232 myDetElem->SetLocalTransformation(new TGeoCombiTrans(trfThDetElem),kTRUE);
233
234 for (int iCor=0; iCor<3; iCor++){
235 lCenDetElem[iDetElem][iCor] = 0;
236 lRotDetElem[iDetElem][iCor] = 0;
237 }
238
239 if (bMonitor){
240 // MONITOR: Draw graph with sticker targets for plane fitting
241 myDetElem->DrawSTargets();
242 gPad->Update();
243 }
244
245 // Get DetElem transformation.
246 // Psi and Theta are obtained by fitting a plane to the sticker targets.
247 // Then Xc, Yc, Zc and Phi are obtained by solving the equations to the ref.
248 // syst. transformation of the button targets
249
250 // Set DetElem plane function
251 cout << "Setting plane for DetElem" << iDetElem+1 << " ..." << endl;
252 myDetElem->SetPlane(Form("fDetElem%d",iDetElem+1));
f7d36813 253 myDetElem->SetPlaneParameters(0.,0.,3.);
254
ba8b0266 255 // Fit a plane to sticker targets
256 Double_t lChi2 = myDetElem->FitPlane();
257 cout << "... chi2 = " << lChi2 << " ..." << endl;
258
259 lRotDetElem[iDetElem][0]=TMath::ATan(myDetElem->GetPlane()->GetParameter(0));
260 lRotDetElem[iDetElem][1]=TMath::ATan(myDetElem->GetPlane()->GetParameter(1));
261
262 // Calculate Mean transformation using previous plane fit
263 // and pairs of button targets
264 myDetElem->CalculateMeanTransf(lCenDetElem[iDetElem],lRotDetElem[iDetElem]);
265
266 cout << "DetElem" << iDetElem+1 << " : mycen(" << lCenDetElem[iDetElem][0] << "," << lCenDetElem[iDetElem][1] << "," << lCenDetElem[iDetElem][2] << "); rot(" << lRotDetElem[iDetElem][0] << "," << lRotDetElem[iDetElem][1] << "," << lRotDetElem[iDetElem][2] << ")" << endl;
267
268
269 // Get best transformation from fitting method
270 // (local to global transformation)
271 cout << "Trying fitting method for DetElemId " << iDetElemToDetElemId[iDetElem] << endl;
272 myDetElem->SurveyToAlign();
273 // myDetElem->SurveyToAlign(lRotDetElem[iDetElem][0],lRotDetElem[iDetElem][1],myDetElem->GetPlane()->GetParError(0),myDetElem->GetPlane()->GetParError(1));
274 }
275 }
276
277 // Print found transformation of Detection Element (plane fit + loop)
278 for (int iCh =0; iCh<=1; iCh++) {
279 for (int iDetElemI=0; iDetElemI<nDetElemsI; iDetElemI++){
280 Int_t iDetElem = iDetElemsOfChamber[iCh][iDetElemI];
281 cout << "DetElem" << iDetElem+1 << " : mycen(" << lCenDetElem[iDetElem][0] << "," << lCenDetElem[iDetElem][1] << "," << lCenDetElem[iDetElem][2] << "); rot(" << lRotDetElem[iDetElem][0] << "," << lRotDetElem[iDetElem][1] << "," << lRotDetElem[iDetElem][2] << ")" << endl;
282 }
283 }
284
285 // Print Theoretical transformation of Detection Element
286 for (int iCh =0; iCh<=1; iCh++) {
287 myChamber = (AliMUONSurveyChamber*)myChamberArray->At(iCh);
288 for (int iDetElemI=0; iDetElemI<nDetElemsI; iDetElemI++){
289 myChamber->GetDetElem(iDetElemI)->PrintLocalTrf();
290 }
291 }
292
293 // Print found delta transformation of Detection Element
294 for (int iCh =0; iCh<=1; iCh++) {
295 myChamber = (AliMUONSurveyChamber*)myChamberArray->At(iCh);
296 for (int iDetElemI=0; iDetElemI<nDetElemsI; iDetElemI++){
297 myChamber->GetDetElem(iDetElemI)->PrintAlignTrf();
298 }
299 }
300
301 //
302 // Compare transformations to expected ones
303 //
304 Int_t iDetElemToPos[18] = {0, 1, 2, 3};
305
306 TGraphErrors *gDeltaDiffCenXDetElem0 = new TGraphErrors(nDetElems);
307 TGraphErrors *gDeltaDiffCenYDetElem0 = new TGraphErrors(nDetElems);
308 TGraphErrors *gDeltaDiffCenZDetElem0 = new TGraphErrors(nDetElems);
309 TGraphErrors *gDeltaDiffPsiDetElem0 = new TGraphErrors(nDetElems);
310 TGraphErrors *gDeltaDiffThtDetElem0 = new TGraphErrors(nDetElems);
311 TGraphErrors *gDeltaDiffPhiDetElem0 = new TGraphErrors(nDetElems);
312
313 for (int iCh =0; iCh<=1; iCh++) {
314 myChamber = (AliMUONSurveyChamber*)myChamberArray->At(iCh);
315 for (int iDetElemI=0; iDetElemI<nDetElemsI; iDetElemI++){
316 myChamber->GetDetElem(iDetElemI)->GetAlignTrf()->Print();
317 }
318 }
319
320 for (int iCh =0; iCh<=1; iCh++) {
321 myChamber = (AliMUONSurveyChamber*)myChamberArray->At(iCh);
322 // Store delta transformations to use for CDB entry creation
323 dtrfDetElem[nDetElems+iCh].Clear();
324 dtrfDetElem[nDetElems+iCh] = *(myChamber->GetAlignTrf());
325
326 // Those are for checks and visualizations
327 localTrfDetElem[nDetElems+iCh].Clear();
328 localTrfDetElem[nDetElems+iCh] = (*(myChamber->GetLocalTrf())*(*(myChamber->GetAlignTrf())));
329 localTrfDetElem[nDetElems+iCh].Print();
330 lDiffCenDetElem0[nDetElems+iCh] = (Double_t*)localTrfDetElem[nDetElems+iCh].GetTranslation();
331 AliMUONSurveyUtil::MatrixToAngles(localTrfDetElem[nDetElems+iCh].GetRotationMatrix(),lDiffRotDetElem0[nDetElems+iCh]);
332
333 localTrfThDetElem[nDetElems+iCh].Clear();
334 localTrfThDetElem[nDetElems+iCh] = (*(myChamber->GetLocalTrf()));
335 localTrfThDetElem[nDetElems+iCh].Print();
336 lDiffThCenDetElem0[nDetElems+iCh] = (Double_t*)localTrfThDetElem[nDetElems+iCh].GetTranslation();
337 AliMUONSurveyUtil::MatrixToAngles(localTrfThDetElem[nDetElems+iCh].GetRotationMatrix(),lDiffThRotDetElem0[nDetElems+iCh]);
338
339 for (int iCor=0; iCor<3; iCor++){
340 lDeltaDiffCenDetElem0[nDetElems+iCh][iCor] = lDiffCenDetElem0[nDetElems+iCh][iCor]-lDiffThCenDetElem0[nDetElems+iCh][iCor];
341 lDeltaDiffRotDetElem0[nDetElems+iCh][iCor] = lDiffRotDetElem0[nDetElems+iCh][iCor]-lDiffThRotDetElem0[nDetElems+iCh][iCor];
342 if (lDeltaDiffRotDetElem0[nDetElems+iCh][iCor]>TMath::Pi()) lDeltaDiffRotDetElem0[nDetElems+iCh][iCor]-=TMath::TwoPi();
343 if (lDeltaDiffRotDetElem0[nDetElems+iCh][iCor]<-TMath::Pi()) lDeltaDiffRotDetElem0[nDetElems+iCh][iCor]+=TMath::TwoPi();
344 }
345
346 for (int iDetElemI=0; iDetElemI<nDetElemsI; iDetElemI++){
347 Int_t iDetElem = iDetElemsOfChamber[iCh][iDetElemI];
348 myDetElem = myChamber->GetDetElem(iDetElemI);
349 // Store delta transformations to use for CDB entry creation
350 dtrfDetElem[iDetElem].Clear();
351 dtrfDetElem[iDetElem] = *(myDetElem->GetAlignTrf());
352
353 // Those are for checks and visualizations
354 localTrfDetElem[iDetElem].Clear();
355 localTrfDetElem[iDetElem] = (*(myDetElem->GetLocalTrf())*(*(myDetElem->GetAlignTrf())));
356 // localTrfDetElem[iDetElem] = (*(myDetElem->GetBaseTrf())*(*(myDetElem->GetAlignTrf())));
357 localTrfDetElem[iDetElem].Print();
358 lDiffCenDetElem0[iDetElem] = (Double_t*)localTrfDetElem[iDetElem].GetTranslation();
359 AliMUONSurveyUtil::MatrixToAngles(localTrfDetElem[iDetElem].GetRotationMatrix(),lDiffRotDetElem0[iDetElem]);
360// lDiffCenDetElem0[iDetElem] = lCenDetElem[iDetElem];
361// lDiffRotDetElem0[iDetElem] = lRotDetElem[iDetElem];
362
363 localTrfThDetElem[iDetElem].Clear();
364 localTrfThDetElem[iDetElem] = (*(myDetElem->GetLocalTrf()));
365 // localTrfThDetElem[iDetElem] = (*(myDetElem->GetBaseTrf()));
366 localTrfThDetElem[iDetElem].Print();
367 lDiffThCenDetElem0[iDetElem] = (Double_t*)localTrfThDetElem[iDetElem].GetTranslation();
368 AliMUONSurveyUtil::MatrixToAngles(localTrfThDetElem[iDetElem].GetRotationMatrix(),lDiffThRotDetElem0[iDetElem]);
369
370 for (int iCor=0; iCor<3; iCor++){
371 lDeltaDiffCenDetElem0[iDetElem][iCor] = lDiffCenDetElem0[iDetElem][iCor]-lDiffThCenDetElem0[iDetElem][iCor];
372 lDeltaDiffRotDetElem0[iDetElem][iCor] = lDiffRotDetElem0[iDetElem][iCor]-lDiffThRotDetElem0[iDetElem][iCor];
373 if (lDeltaDiffRotDetElem0[iDetElem][iCor]>TMath::Pi()) lDeltaDiffRotDetElem0[iDetElem][iCor]-=TMath::TwoPi();
374 if (lDeltaDiffRotDetElem0[iDetElem][iCor]<-TMath::Pi()) lDeltaDiffRotDetElem0[iDetElem][iCor]+=TMath::TwoPi();
375 }
376
f7d36813 377 gDeltaDiffCenXDetElem0->SetPoint(iDetElem,1e1*lDeltaDiffCenDetElem0[iDetElem][0],iDetElemToPos[iDetElem]+1);
378 gDeltaDiffCenYDetElem0->SetPoint(iDetElem,1e1*lDeltaDiffCenDetElem0[iDetElem][1],iDetElemToPos[iDetElem]+1);
379 gDeltaDiffCenZDetElem0->SetPoint(iDetElem,1e1*lDeltaDiffCenDetElem0[iDetElem][2],iDetElemToPos[iDetElem]+1);
ba8b0266 380 gDeltaDiffPsiDetElem0->SetPoint(iDetElem,1e3*lDeltaDiffRotDetElem0[iDetElem][0],iDetElemToPos[iDetElem]+1);
381 gDeltaDiffThtDetElem0->SetPoint(iDetElem,1e3*lDeltaDiffRotDetElem0[iDetElem][1],iDetElemToPos[iDetElem]+1);
382 gDeltaDiffPhiDetElem0->SetPoint(iDetElem,1e3*lDeltaDiffRotDetElem0[iDetElem][2],iDetElemToPos[iDetElem]+1);
f7d36813 383 gDeltaDiffCenXDetElem0->SetPointError(iDetElem,1e1*myDetElem->GetFitter()->GetParError(0),0.);
384 gDeltaDiffCenYDetElem0->SetPointError(iDetElem,1e1*myDetElem->GetFitter()->GetParError(1),0.);
385 gDeltaDiffCenZDetElem0->SetPointError(iDetElem,1e1*myDetElem->GetFitter()->GetParError(2),0.);
ba8b0266 386 gDeltaDiffPsiDetElem0->SetPointError(iDetElem,1e3*myDetElem->GetFitter()->GetParError(3),0.);
387 gDeltaDiffThtDetElem0->SetPointError(iDetElem,1e3*myDetElem->GetFitter()->GetParError(4),0.);
388 gDeltaDiffPhiDetElem0->SetPointError(iDetElem,1e3*myDetElem->GetFitter()->GetParError(5),0.);
389 }
390 }
391
392 // Apply the found misalignments to the geometry
393 AliMUONGeometryTransformer *newTransform = AliMUONSurveyUtil::ReAlign(transform,chId,nDetElems,iDetElemPseudoIdToDetElem,dtrfDetElem,true);
394 newTransform->WriteTransformations(Form("transform2ReAlignSurveyCh%d.dat",chId+1));
395
396 if(bOCDB){
397 // Generate realigned data in local cdb
398 const TClonesArray* array = newTransform->GetMisAlignmentData();
399
400 // Set the alignment resolution in the align objects for this chamber
f7d36813 401 Double_t chResX = myChamberI->GetAlignResX();
402 Double_t chResY = myChamberI->GetAlignResY();
403 Double_t deResX = (myChamberI->GetMeanDetElemAlignResX()+myChamberO->GetMeanDetElemAlignResX())/2.;
404 Double_t deResY = (myChamberI->GetMeanDetElemAlignResY()+myChamberO->GetMeanDetElemAlignResY())/2.;
405
ba8b0266 406 printf("Chamber alignment resolutions: resX=%f , resY=%f\n",chResX,chResY);
407 printf("Detection Elements alignment resolutions: resX=%f , resY=%f\n",deResX,deResY);
408 chResX = TMath::Sqrt(0.1*0.1+chResX*chResX);
409 chResY = TMath::Sqrt(0.1*0.1+chResY*chResY);
410 AliMUONSurveyUtil::SetAlignmentResolution(array,chId,chResX,chResY,deResX,deResY);
411
412 // CDB manager
413 AliCDBManager* cdbManager = AliCDBManager::Instance();
414 cdbManager->SetDefaultStorage(Form("local://ReAlignSurveyCh%dCDB",chId+1));
415
416 AliCDBMetaData* cdbData = new AliCDBMetaData();
417 cdbData->SetResponsible("Dimuon Offline project");
418 cdbData->SetComment("MUON alignment objects with survey realignment");
419 AliCDBId id("MUON/Align/Data", 0, AliCDBRunRange::Infinity());
420 cdbManager->Put(const_cast<TClonesArray*>(array), id, cdbData);
421 }
422
423 for(Int_t iCor=0; iCor<3; iCor++){
424 for(Int_t iDetElem=0; iDetElem<nDetElems; iDetElem++){
425 cout << lDeltaDiffCenDetElem0[iDetElem][iCor] << " " << lDiffCenDetElem0[iDetElem][iCor] << " " << lDiffThCenDetElem0[iDetElem][iCor] << endl;
426 }
427 cout << endl;
428 }
429 cout << endl;
430 for(Int_t iCor=0; iCor<3; iCor++){
431 for(Int_t iDetElem=0; iDetElem<nDetElems; iDetElem++){
432 cout << lDeltaDiffRotDetElem0[iDetElem][iCor] << " " << lDiffRotDetElem0[iDetElem][iCor] << " " << lDiffThRotDetElem0[iDetElem][iCor] << endl;
433 }
434 cout << endl;
435 }
f7d36813 436
437 // Plot the found misalignements (Units here are mm and mrad!)
ba8b0266 438 TH1F *myDetElemDeltaDiffCenX = new TH1F("myDetElemDeltaDiffCenX","myDetElemDeltaDiffCenX",100,-10,10);
439 myDetElemDeltaDiffCenX->SetMaximum(nDetElems+1);
440 myDetElemDeltaDiffCenX->SetMinimum(0);
441 TH1F *myDetElemDeltaDiffCenY = new TH1F("myDetElemDeltaDiffCenY","myDetElemDeltaDiffCenY",100,-10,10);
442 myDetElemDeltaDiffCenY->SetMaximum(nDetElems+1);
443 myDetElemDeltaDiffCenY->SetMinimum(0);
444 TH1F *myDetElemDeltaDiffCenZ = new TH1F("myDetElemDeltaDiffCenZ","myDetElemDeltaDiffCenZ",100,-20,20);
445 myDetElemDeltaDiffCenZ->SetMaximum(nDetElems+1);
446 myDetElemDeltaDiffCenZ->SetMinimum(0);
447
f7d36813 448 TH1F *myDetElemDeltaDiffRotX = new TH1F("myDetElemDeltaDiffRotX","myDetElemDeltaDiffRotX",100,-10,10);
ba8b0266 449 myDetElemDeltaDiffRotX->SetMaximum(nDetElems+1);
450 myDetElemDeltaDiffRotX->SetMinimum(0);
f7d36813 451 TH1F *myDetElemDeltaDiffRotY = new TH1F("myDetElemDeltaDiffRotY","myDetElemDeltaDiffRotY",100,-10,10);
ba8b0266 452 myDetElemDeltaDiffRotY->SetMaximum(nDetElems+1);
453 myDetElemDeltaDiffRotY->SetMinimum(0);
f7d36813 454 TH1F *myDetElemDeltaDiffRotZ = new TH1F("myDetElemDeltaDiffRotZ","myDetElemDeltaDiffRotZ",100,-6,6);
ba8b0266 455 myDetElemDeltaDiffRotZ->SetMaximum(nDetElems+1);
456 myDetElemDeltaDiffRotZ->SetMinimum(0);
457
458 //
459 // ******** Starting plots
460 //
461 TCanvas *canvas;
462 TPad *pad;
463 TPaveLabel *theTitle;
464 gStyle->SetPalette(1);
465
466 TPostScript *ps = 0;
467
468 if( saveps ){
469 ps = new TPostScript(Form("surveyChamber%d",chId+1),filetype);
470 ps->NewPage();
471 }
472
473 // Observed misalignments
474 str = Form("Chamber %d",chId+1);
475 TCanvas *cvn2 = new TCanvas("cvn2",str,cWidth,cHeight);
476 canvas = cvn2;
477 canvas->Range(0,0,21,29);
478
479 title = Form(" MisAlignments Chamber %d - L2G - In Single Frame ",chId+1);
480 TPaveLabel *theTitle2 = new TPaveLabel(3,27.0,18,28.5,title,"br");
481 theTitle = theTitle2;
482 theTitle->SetFillColor(18);
483 theTitle->SetTextFont(32);
484 theTitle->SetTextSize(0.4);
485 theTitle->SetTextColor(1);
486 theTitle->Draw();
487
488 TPad *pad2 = new TPad("pad2","pad2",0.01,0.01,0.98,0.91,0);
489 pad = pad2;
490 pad->Draw();
491 TLine *ch0Line = new TLine(0,1,0,2);
492 TLine *ch1Line = new TLine(0,1,0,2);
493 ch1Line->SetLineStyle(2);
494 pad->Divide(3,2);
495
496 pad->cd(1);
497 myDetElemDeltaDiffCenX->Draw();
498 myDetElemDeltaDiffCenX->SetXTitle("#Delta[xc_{i}^{m}-xc_{i}^{th}] (mm)");
499 myDetElemDeltaDiffCenX->SetYTitle("DetElem arbitrary ordering");
500 gDeltaDiffCenXDetElem0->SetMarkerStyle(20);
501 gDeltaDiffCenXDetElem0->Draw("P");
f7d36813 502 ch0Line->DrawLine(1e1*lDeltaDiffCenDetElem0[nDetElems+0][0],0.5,1e1*lDeltaDiffCenDetElem0[nDetElems+0][0],1.5);
503 ch0Line->DrawLine(1e1*lDeltaDiffCenDetElem0[nDetElems+0][0],3.5,1e1*lDeltaDiffCenDetElem0[nDetElems+0][0],4.5);
504 ch1Line->DrawLine(1e1*lDeltaDiffCenDetElem0[nDetElems+1][0],1.5,1e1*lDeltaDiffCenDetElem0[nDetElems+1][0],3.5);
ba8b0266 505
506 pad->cd(2);
507 myDetElemDeltaDiffCenY->Draw();
508 myDetElemDeltaDiffCenY->SetXTitle("#Delta[yc_{i}^{m}-yc_{i}^{th}] (mm)");
509 myDetElemDeltaDiffCenY->SetYTitle("DetElem arbitrary ordering");
510 gDeltaDiffCenYDetElem0->SetMarkerStyle(20);
511 gDeltaDiffCenYDetElem0->Draw("P");
f7d36813 512 ch0Line->DrawLine(1e1*lDeltaDiffCenDetElem0[nDetElems+0][1],0.5,1e1*lDeltaDiffCenDetElem0[nDetElems+0][1],1.5);
513 ch0Line->DrawLine(1e1*lDeltaDiffCenDetElem0[nDetElems+0][1],3.5,1e1*lDeltaDiffCenDetElem0[nDetElems+0][1],4.5);
514 ch1Line->DrawLine(1e1*lDeltaDiffCenDetElem0[nDetElems+1][1],1.5,1e1*lDeltaDiffCenDetElem0[nDetElems+1][1],3.5);
ba8b0266 515
516 pad->cd(3);
517 myDetElemDeltaDiffCenZ->Draw();
518 myDetElemDeltaDiffCenZ->SetXTitle("#Delta[zc_{i}^{m}-zc_{i}^{th}] (mm)");
519 myDetElemDeltaDiffCenZ->SetYTitle("DetElem arbitrary ordering");
520 gDeltaDiffCenZDetElem0->SetMarkerStyle(20);
521 gDeltaDiffCenZDetElem0->Draw("P");
f7d36813 522 ch0Line->DrawLine(1e1*lDeltaDiffCenDetElem0[nDetElems+0][2],0.5,1e1*lDeltaDiffCenDetElem0[nDetElems+0][2],1.5);
523 ch0Line->DrawLine(1e1*lDeltaDiffCenDetElem0[nDetElems+0][2],3.5,1e1*lDeltaDiffCenDetElem0[nDetElems+0][2],4.5);
524 ch1Line->DrawLine(1e1*lDeltaDiffCenDetElem0[nDetElems+1][2],1.5,1e1*lDeltaDiffCenDetElem0[nDetElems+1][2],3.5);
ba8b0266 525
526 pad->cd(4);
527 myDetElemDeltaDiffRotX->Draw();
528 myDetElemDeltaDiffRotX->SetXTitle("#Delta[#psi_{i}^{m}-#psi_{i}^{th}] (mrad)");
529 myDetElemDeltaDiffRotX->SetYTitle("DetElem arbitrary ordering");
530 gDeltaDiffPsiDetElem0->SetMarkerStyle(20);
531 gDeltaDiffPsiDetElem0->Draw("P");
532 ch0Line->DrawLine(1e3*lDeltaDiffRotDetElem0[nDetElems+0][0],0.5,1e3*lDeltaDiffRotDetElem0[nDetElems+0][0],1.5);
533 ch0Line->DrawLine(1e3*lDeltaDiffRotDetElem0[nDetElems+0][0],3.5,1e3*lDeltaDiffRotDetElem0[nDetElems+0][0],4.5);
534 ch1Line->DrawLine(1e3*lDeltaDiffRotDetElem0[nDetElems+1][0],1.5,1e3*lDeltaDiffRotDetElem0[nDetElems+1][0],3.5);
535
536 pad->cd(5);
537 myDetElemDeltaDiffRotY->Draw();
538 myDetElemDeltaDiffRotY->SetXTitle("#Delta[#theta_{i}^{m}-#theta_{i}^{th}] (mrad)");
539 myDetElemDeltaDiffRotY->SetYTitle("DetElem arbitrary ordering");
540 gDeltaDiffThtDetElem0->SetMarkerStyle(20);
541 gDeltaDiffThtDetElem0->Draw("P");
542 ch0Line->DrawLine(1e3*lDeltaDiffRotDetElem0[nDetElems+0][1],0.5,1e3*lDeltaDiffRotDetElem0[nDetElems+0][1],1.5);
543 ch0Line->DrawLine(1e3*lDeltaDiffRotDetElem0[nDetElems+0][1],3.5,1e3*lDeltaDiffRotDetElem0[nDetElems+0][1],4.5);
544 ch1Line->DrawLine(1e3*lDeltaDiffRotDetElem0[nDetElems+1][1],1.5,1e3*lDeltaDiffRotDetElem0[nDetElems+1][1],3.5);
545
546 pad->cd(6);
547 myDetElemDeltaDiffRotZ->Draw();
548 myDetElemDeltaDiffRotZ->SetXTitle("#Delta[#phi_{i}^{m}-#phi_{i}^{th}] (mrad)");
549 myDetElemDeltaDiffRotZ->SetYTitle("DetElem arbitrary ordering");
550 gDeltaDiffPhiDetElem0->SetMarkerStyle(20);
551 gDeltaDiffPhiDetElem0->Draw("P");
552 ch0Line->DrawLine(1e3*lDeltaDiffRotDetElem0[nDetElems+0][2],0.5,1e3*lDeltaDiffRotDetElem0[nDetElems+0][2],1.5);
553 ch0Line->DrawLine(1e3*lDeltaDiffRotDetElem0[nDetElems+0][2],3.5,1e3*lDeltaDiffRotDetElem0[nDetElems+0][2],4.5);
554 ch1Line->DrawLine(1e3*lDeltaDiffRotDetElem0[nDetElems+1][2],1.5,1e3*lDeltaDiffRotDetElem0[nDetElems+1][2],3.5);
555
556 pad->Update();
557
558 if(bMonitor){
559 // MONITOR: Histograms for monitoring
f7d36813 560 TH2F *hCPSTa = new TH2F("hCPSTa","hCPSTa",100,-200,200,100,-200,200);
561 TH2F *hCPSTc = new TH2F("hCPSTc","hCPSTc",100,-200,200,100,-200,200);
562 TH2F *hSSTa = new TH2F("hSSTa","hSSTa",100,-200,200,100,-200,200);
563 TH2F *hSSTc = new TH2F("hSSTc","hSSTc",100,-200,200,100,-200,200);
564 TH2F *hSSTap = new TH2F("hSSTap","hSSTap",800,-200,200,800,-200,200);
565 TH2F *hSSTcp = new TH2F("hSSTcp","hSSTcp",800,-200,200,800,-200,200);
ba8b0266 566
567 // MONITOR: Fill histograms with chambers and slats sticker target positions
568 for (int iCh =0; iCh<=1; iCh++) {
569 myChamber = (AliMUONSurveyChamber*)myChamberArray->At(iCh);
570 cout << "Filling histograms of sticker target for chamber" << iCh+1 << " ..." << endl;
571 myChamber->FillCPSTHistograms(Form("T4_6%d",2-iCh),hCPSTc);
f7d36813 572 myChamber->FillDESTHistograms(TString("T4_46"),hSSTc);
ba8b0266 573
574 for (int iDetElemI=0; iDetElemI<nDetElemsI; iDetElemI++){
575 // Int_t iDetElem = iDetElemsOfChamber[iCh][iDetElemI];
576 myDetElem = myChamber->GetDetElem(iDetElemI);
577 // MONITOR: Fill slat plane for fit monitor.
578 Double_t pl[3] = {0};
579 Double_t pg[3] = {0};
580 AliSurveyPoint *pointSBT0 = myDetElem->GetLButtonTarget(0);
581 AliSurveyPoint *pointSBT1 = myDetElem->GetLButtonTarget(1);
582 if(pointSBT0&&pointSBT1) {
583 if (pointSBT0->GetX()>pointSBT1->GetX()){
584 pointSBT0=myDetElem->GetLButtonTarget(1);
585 pointSBT1=myDetElem->GetLButtonTarget(0);
586 }
587 Double_t lX = pointSBT0->GetX();
588 while(lX<pointSBT1->GetX()){
f7d36813 589 Double_t lY = pointSBT0->GetY()-20;
590 while(lY<pointSBT0->GetY()+20){
ba8b0266 591 pl[0] = lX; pl[1] = lY; pl[2] = 0.;
592 (TGeoCombiTrans((*(myDetElem->GetBaseTrf()))*(*(myDetElem->GetAlignTrf())))).LocalToMaster(pl,pg);
593 if(myDetElem->GetGButtonTarget(0)->GetPointName().Contains("A")){
f7d36813 594 if (hSSTap->GetBinContent(hSSTap->FindBin(pg[0],pg[1]))==0)
595 hSSTap->Fill(pg[0],pg[1],-pg[2]);
ba8b0266 596 }
597 else {
f7d36813 598 if (hSSTcp->GetBinContent(hSSTcp->FindBin(pg[0],pg[1]))==0)
599 hSSTcp->Fill(pg[0],pg[1],-pg[2]);
ba8b0266 600 }
601 lY+=hSSTap->GetYaxis()->GetBinWidth(1);
602 }
603 lX+=hSSTap->GetXaxis()->GetBinWidth(1);
604 }
605 }
606 }
607 }
608
609 if( saveps ){
610 ps->NewPage();
611 }
612
613 // View from side A
614 str = Form("Chamber %d - side A",chId+1);
615 TCanvas *cvn0 = new TCanvas("cvn0",str,cWidth,cHeight);
616 canvas = cvn0;
617 canvas->Range(0,0,21,29);
618
619 title = Form(" Deformations of chamber %d - View from side A ",chId+1);
620 TPaveLabel *theTitle0 = new TPaveLabel(3,27.0,18,28.5,title,"br");
621 theTitle = theTitle0;
622 theTitle->SetFillColor(18);
623 theTitle->SetTextFont(32);
624 theTitle->SetTextSize(0.4);
625 theTitle->SetTextColor(1);
626 theTitle->Draw();
627
628 TPad *pad0 = new TPad("pad0","pad0",0.01,0.01,0.98,0.91,0);
629 pad = pad0;
630 pad->Draw();
631 pad->Divide(2,2);
632
633 pad->cd(1);
f7d36813 634 hCPSTa->SetMinimum(950);
635 hCPSTa->SetMaximum(975);
ba8b0266 636 hCPSTa->Draw("lego2z");
637
638 pad->cd(2);
f7d36813 639 hSSTa->SetMinimum(950);
640 hSSTa->SetMaximum(980);
ba8b0266 641 hSSTa->Draw("lego2z");
642
643 pad->cd(3);
644
645 pad->cd(4);
f7d36813 646 hSSTap->SetMinimum(950);
647 hSSTap->SetMaximum(980);
ba8b0266 648 hSSTap->Draw("surf2z");
649
650 pad->Update();
651 if(saveps){
652 ps->NewPage();
653 }
654
655 // Inv Mass, Multiplicity
656 str = Form("Chamber %d - side C",chId+1);
657 TCanvas *cvn1 = new TCanvas("cvn1",str,cWidth,cHeight);
658 canvas = cvn1;
659 canvas->Range(0,0,21,29);
660
661 title = Form(" Deformations of chamber %d - View from side C ",chId+1);
662 TPaveLabel *theTitle1 = new TPaveLabel(3,27.0,18,28.5,title,"br");
663 theTitle = theTitle1;
664 theTitle->SetFillColor(18);
665 theTitle->SetTextFont(32);
666 theTitle->SetTextSize(0.4);
667 theTitle->SetTextColor(1);
668 theTitle->Draw();
669
670 TPad *pad1 = new TPad("pad1","pad1",0.01,0.01,0.98,0.91,0);
671 pad = pad1;
672 pad->Draw();
673 pad->Divide(2,2);
674
f7d36813 675 Double_t lMin, lMax;
676
ba8b0266 677 pad->cd(1);
f7d36813 678 lMin = hCPSTc->GetMinimum(0);
679 hCPSTc->SetMinimum(TMath::Floor(lMin));
680 lMax = hCPSTc->GetMaximum();
681 hCPSTc->SetMaximum(TMath::Ceil(lMax));
ba8b0266 682 hCPSTc->Draw("lego2z");
683
684 pad->cd(2);
f7d36813 685// lMin = hSSTc->GetMinimum(0);
686// hSSTc->SetMinimum(TMath::Floor(lMin));
687// lMax = hSSTc->GetMaximum();
688// hSSTc->SetMaximum(TMath::Ceil(lMax));
689// hSSTc->Draw("lego2z");
690
ba8b0266 691
692 pad->cd(3);
693
694 pad->cd(4);
f7d36813 695 lMin = hSSTcp->GetMinimum(0);
696 hSSTcp->SetMinimum(TMath::Floor(lMin));
697 lMax = hSSTcp->GetMaximum();
698 hSSTcp->SetMaximum(TMath::Ceil(lMax));
ba8b0266 699 hSSTcp->Draw("surf2z");
700 }
701
702 pad->Update();
703 if( saveps ){
704 ps->Close();
705 }
706
707 TFile *hFile = new TFile(Form("spCH%d_L2G_ISF.root",chId+1),"RECREATE");
708 gDeltaDiffCenXDetElem0->Write("gDeltaDiffCenXDetElem0");
709 gDeltaDiffCenYDetElem0->Write("gDeltaDiffCenYDetElem0");
710 gDeltaDiffCenZDetElem0->Write("gDeltaDiffCenZDetElem0");
711 gDeltaDiffPsiDetElem0->Write("gDeltaDiffPsiDetElem0");
712 gDeltaDiffThtDetElem0->Write("gDeltaDiffThtDetElem0");
713 gDeltaDiffPhiDetElem0->Write("gDeltaDiffPhiDetElem0");
714 myDetElemDeltaDiffCenX->Write();
715 myDetElemDeltaDiffCenY->Write();
716 myDetElemDeltaDiffCenZ->Write();
717 myDetElemDeltaDiffRotX->Write();
718 myDetElemDeltaDiffRotY->Write();
719 myDetElemDeltaDiffRotZ->Write();
720 hFile->Close();
721 hFile->Delete();
722}