]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/AliITSResidualsAnalysis.cxx
Added a protection against bad (wrt pedestal values) channels, and removed an inadequ...
[u/mrichter/AliRoot.git] / ITS / AliITSResidualsAnalysis.cxx
CommitLineData
146f76de 1/**************************************************************************
2 * Copyright(c) 1998-2008, 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// Implementation of the tracks residuals analysis class
18// It provides an access to the track space points
19// written along the esd tracks. The class enables
20// the user to plug any track fitter (deriving from
21// AliTrackFitter class) and minimization fo the
22// track residual sums (deriving from the AliTrackResiduals).
23//-----------------------------------------------------------------
24
25#include <Riostream.h>
146f76de 26#include <TFile.h>
27#include <TMath.h>
28#include <TArrayI.h>
29#include <TClonesArray.h>
30#include <TNtuple.h>
31#include <TTree.h>
32#include <TF1.h>
33#include <TH1F.h>
34#include <TH2F.h>
35#include <TCanvas.h>
36#include <TGraphErrors.h>
37#include "TGeoMatrix.h"
38#include "TGeoManager.h"
39#include "TGeoPhysicalNode.h"
bad3394b 40#include "TMatrixDSym.h"
146f76de 41#include "TMatrixDSymEigen.h"
bad3394b 42#include "TMatrixD.h"
146f76de 43#include "TString.h"
44
45#include "AliAlignmentTracks.h"
46#include "AliTrackPointArray.h"
47#include "AliAlignObjParams.h"
48#include "AliTrackResiduals.h"
49#include "AliTrackFitter.h"
50#include "AliTrackFitterKalman.h"
51#include "AliTrackFitterRieman.h"
52#include "AliTrackResiduals.h"
53#include "AliTrackResidualsChi2.h"
54#include "AliTrackResidualsFast.h"
55#include "AliLog.h"
bad3394b 56#include "AliITSgeomTGeo.h"
146f76de 57
58#include "AliITSResidualsAnalysis.h"
59
146f76de 60ClassImp(AliITSResidualsAnalysis)
61
62//____________________________________________________________________________
63 AliITSResidualsAnalysis::AliITSResidualsAnalysis():
64 AliAlignmentTracks(),
65 fnHist(0),
66 fnPhi(0),
67 fnZ(0),
68 fvolidsToBin(0),
69 fLastVolVolid(0),
70 fCoordToBinTable(0),
71 fVolResHistRPHI(0),
72 fResHistZ(0),
bad3394b 73 fResHistX(0),
74 fResHistXLocsddL(0),
75 fResHistXLocsddR(0),
76 fHistCoordGlobY(0),
146f76de 77 fPullHistRPHI(0),
78 fPullHistZ(0),
79 fTrackDirPhi(0),
80 fTrackDirLambda(0),
81 fTrackDirLambda2(0),
82 fTrackDirAlpha(0),
83 fTrackDirPhiAll(0),
84 fTrackDirLambdaAll(0),
85 fTrackDirLambda2All(0),
86 fTrackDirAlphaAll(0),
87 fTrackDir(0),
88 fTrackDirAll(0),
89 fTrackDir2All(0),
90 fTrackDirXZAll(0),
91 fResHistGlob(0),
92 fhistCorrVol(0),
93 fVolNTracks(0),
94 fhEmpty(0),
95 fhistVolNptsUsed(0),
96 fhistVolUsed(0),
97 fSigmaVolZ(0),
98 fsingleLayer(0),
99 fWriteHist(0),
100 fpTrackVolIDs(0),
101 fVolVolids(0),
102 fVolUsed(0),
103 fRealignObjFileIsOpen(kFALSE),
104 fClonesArray(0),
105 fAliTrackPoints("AliTrackPoints.root"),
106 fGeom("geometry.root")
146f76de 107{
108
109 //
110 // Defaults
111 //
112
113
114}
115
116//____________________________________________________________________________
bad3394b 117AliITSResidualsAnalysis::AliITSResidualsAnalysis(const TString aliTrackPoints,
146f76de 118 const TString geom):
119 AliAlignmentTracks(),
120 fnHist(0),
121 fnPhi(0),
122 fnZ(0),
123 fvolidsToBin(0),
124 fLastVolVolid(0),
125 fCoordToBinTable(0),
126 fVolResHistRPHI(0),
127 fResHistZ(0),
bad3394b 128 fResHistX(0),
129 fResHistXLocsddL(0),
130 fResHistXLocsddR(0),
131 fHistCoordGlobY(0),
146f76de 132 fPullHistRPHI(0),
133 fPullHistZ(0),
134 fTrackDirPhi(0),
135 fTrackDirLambda(0),
136 fTrackDirLambda2(0),
137 fTrackDirAlpha(0),
138 fTrackDirPhiAll(0),
139 fTrackDirLambdaAll(0),
140 fTrackDirLambda2All(0),
141 fTrackDirAlphaAll(0),
142 fTrackDir(0),
143 fTrackDirAll(0),
144 fTrackDir2All(0),
145 fTrackDirXZAll(0),
146 fResHistGlob(0),
147 fhistCorrVol(0),
148 fVolNTracks(0),
149 fhEmpty(0),
150 fhistVolNptsUsed(0),
151 fhistVolUsed(0),
152 fSigmaVolZ(0),
153 fsingleLayer(0),
154 fWriteHist(0),
155 fpTrackVolIDs(0),
156 fVolVolids(0),
157 fVolUsed(0),
158 fRealignObjFileIsOpen(kFALSE),
159 fClonesArray(0),
160 fAliTrackPoints(aliTrackPoints),
161 fGeom(geom)
162{
163 //
bad3394b 164 // Standard Constructor (alitrackpoints)
146f76de 165 //
166
167
168}
169
170//____________________________________________________________________________
171AliITSResidualsAnalysis::AliITSResidualsAnalysis(const TArrayI *volIDs):
172 AliAlignmentTracks(),
173 fnHist(0),
174 fnPhi(0),
175 fnZ(0),
176 fvolidsToBin(0),
177 fLastVolVolid(0),
178 fCoordToBinTable(0),
179 fVolResHistRPHI(0),
180 fResHistZ(0),
bad3394b 181 fResHistX(0),
182 fResHistXLocsddL(0),
183 fResHistXLocsddR(0),
184 fHistCoordGlobY(0),
146f76de 185 fPullHistRPHI(0),
186 fPullHistZ(0),
187 fTrackDirPhi(0),
188 fTrackDirLambda(0),
189 fTrackDirLambda2(0),
190 fTrackDirAlpha(0),
191 fTrackDirPhiAll(0),
192 fTrackDirLambdaAll(0),
193 fTrackDirLambda2All(0),
194 fTrackDirAlphaAll(0),
195 fTrackDir(0),
196 fTrackDirAll(0),
197 fTrackDir2All(0),
198 fTrackDirXZAll(0),
199 fResHistGlob(0),
200 fhistCorrVol(0),
201 fVolNTracks(0),
202 fhEmpty(0),
203 fhistVolNptsUsed(0),
204 fhistVolUsed(0),
205 fSigmaVolZ(0),
206 fsingleLayer(0),
207 fWriteHist(0),
208 fpTrackVolIDs(0),
209 fVolVolids(0),
210 fVolUsed(0),
211 fRealignObjFileIsOpen(kFALSE),
212 fClonesArray(0),
213 fAliTrackPoints("AliTrackPoints.root"),
214 fGeom("geometry.root")
215
216{
217 //
218 // Original Constructor
219 //
bad3394b 220 if(!gGeoManager) AliGeomManager::LoadGeometry(GetFileNameGeometry());
146f76de 221 InitHistograms(volIDs);
222
223}
224
146f76de 225//____________________________________________________________________________
226AliITSResidualsAnalysis::~AliITSResidualsAnalysis()
227{
228 //
229 // Destructor
230 //
231
232 if(fvolidsToBin) delete[] fvolidsToBin;
233 if(fLastVolVolid) delete[] fLastVolVolid;
234 if(fCoordToBinTable) delete[] fCoordToBinTable;
bad3394b 235 if(fHistCoordGlobY) delete[] fHistCoordGlobY;
146f76de 236 if(fVolResHistRPHI) delete fVolResHistRPHI;
bad3394b 237 if(fResHistZ) delete fResHistZ;
238 if(fResHistX){
239 for(Int_t i=0; i<fnHist; i++) delete fResHistX[i];
240 delete [] fResHistX;
241 }
242 if(fResHistXLocsddL){
243 for(Int_t i=0; i<fnHist; i++) delete fResHistXLocsddL[i];
244 delete [] fResHistXLocsddL;
245 }
246 if(fResHistXLocsddR){
247 for(Int_t i=0; i<fnHist; i++) delete fResHistXLocsddR[i];
248 delete [] fResHistXLocsddR;
249 }
146f76de 250 if(fPullHistRPHI) delete fPullHistRPHI;
251 if(fPullHistZ) delete fPullHistZ;
252 if(fTrackDirPhi) delete fTrackDirPhi;
253 if(fTrackDirLambda) delete fTrackDirLambda;
254 if(fTrackDirLambda2) delete fTrackDirLambda2;
255 if(fTrackDirAlpha) delete fTrackDirAlpha;
256 if(fTrackDirPhiAll) delete fTrackDirPhiAll;
257 if(fTrackDirLambdaAll) delete fTrackDirLambdaAll;
258 if(fTrackDirLambda2All) delete fTrackDirLambda2All;
259 if(fTrackDirAlphaAll) delete fTrackDirAlphaAll;
260 if(fTrackDir) delete fTrackDir;
261 if(fTrackDirAll) delete fTrackDirAll;
262 if(fTrackDir2All) delete fTrackDir2All;
263 if(fTrackDirXZAll) delete fTrackDirXZAll;
264 if(fResHistGlob) delete fResHistGlob;
265 if(fhistCorrVol) delete fhistCorrVol;
266 if(fVolNTracks) delete fVolNTracks;
267 if(fhEmpty) delete fhEmpty;
268 if(fhistVolNptsUsed) delete fhistVolNptsUsed;
269 if(fhistVolUsed) delete fhistVolUsed;
270 if(fSigmaVolZ) delete fSigmaVolZ;
271 if(fpTrackVolIDs) delete fpTrackVolIDs;
272 if(fVolVolids) delete fVolVolids;
273 if(fVolUsed) delete fVolUsed;
274 if(fClonesArray) delete fClonesArray;
275
276 SetFileNameTrackPoints("");
277 SetFileNameGeometry("");
278
279}
280
281//____________________________________________________________________________
282void AliITSResidualsAnalysis::InitHistograms(const TArrayI *volIDs)
283{
284 //
285 // Method that sets and creates the required hisstrograms
286 // with the correct binning (it dos not fill them)
287 //
288
bad3394b 289
290 if(!gGeoManager) AliGeomManager::LoadGeometry(GetFileNameGeometry());
291
146f76de 292 TString histnameRPHI="HistRPHI_volID_",aux;
293 TString histnameZ="HistZ_volID_";
bad3394b 294 TString histnameX="HistX_volID_";
146f76de 295 TString histnameGlob="HistGlob_volID_";
296 TString histnameCorrVol="HistCorrVol_volID";
297 TString histnamePullRPHI="HistPullRPHI_volID_";
298 TString histnamePullZ="HistPullZ_volID_";
299
300 TString histnameDirPhi="HistTrackDirPhi_volID_";
301 TString histnameDirLambda="HistTrackDirLambda_volID_";
302 TString histnameDirLambda2="HistTrackDirLambda2_volID_";
303 TString histnameDirAlpha="HistTrackDirAlpha_volID_";
304 TString histnameDir="HistTrackDir_volID_";
305
146f76de 306 fnHist=volIDs->GetSize();
307 fVolResHistRPHI=new TH1F*[fnHist];
308 fResHistGlob=new TH1F*[fnHist];
309 fResHistZ=new TH1F*[fnHist];
bad3394b 310 fResHistX=new TH1F*[fnHist];
311 fResHistXLocsddL=new TH1F*[fnHist];
312 fResHistXLocsddR=new TH1F*[fnHist];
313 fHistCoordGlobY=new TH1F*[fnHist];
146f76de 314 fPullHistRPHI=new TH1F*[fnHist];
315 fPullHistZ=new TH1F*[fnHist];
316 fhistCorrVol=new TH2F*[fnHist];
146f76de 317
318 fTrackDirPhi=new TH1F*[fnHist];
319 fTrackDirLambda=new TH1F*[fnHist];
320 fTrackDirLambda2=new TH1F*[fnHist];
321 fTrackDirAlpha=new TH1F*[fnHist];
146f76de 322
323 fTrackDirPhiAll=new TH1F("fTrackDirPhiAll","fTrackDirPhiAll",100,-180,180);
324 fTrackDirLambdaAll=new TH1F("fTrackDirLambdaAll","fTrackDirLambdaAll",100,-180,180);
325 fTrackDirLambda2All=new TH1F("fTrackDirLambda2All","fTrackDirLambda2All",100,0,180);
326 fTrackDirAlphaAll=new TH1F("fTrackDirAlphaAll","fTrackDirAlphaAll",100,-180,180);
327
328 fTrackDirAll=new TH2F("fTrackDirAll","Hist with trakcs directions",100,-180,180,100,-180,180);
329 fTrackDir2All=new TH2F("fTrackDir2All","Hist with trakcs directions",100,-180,180,100,-180,180);
bad3394b 330 fTrackDirXZAll=new TH2F("fTrackDirXZAll","Hist with trakcs directions from XZ ",100,-3.,3.,100,-3.,3.);
146f76de 331
332 fTrackDir=new TH2F*[fnHist];
333
bad3394b 334 Bool_t binning;
335 Float_t **binningZPhi;
336 Float_t *binningz;
337 Float_t *binningphi;
146f76de 338
bad3394b 339 binningZPhi=CheckSingleLayer(volIDs);
340 fvolidsToBin=new Int_t*[fnPhi*fnZ];
341 binningphi=binningZPhi[0];
342 binningz=binningZPhi[1];
343 binning=SetBinning(volIDs,binningphi,binningz);
344
345 if(binning){ //ONLY FOR A SINGLE LAYER!
146f76de 346 fVolNTracks=new TH2F("fVolNTracks","Hist with N tracks passing through a given module==(r,phi) zone",fnPhi,binningphi,fnZ,binningz);
347 fhistVolNptsUsed=new TH2F("fhistVolNptsUsed","Hist with N points used for given module==(r,phi) ",fnPhi,binningphi,fnZ,binningz);
348 fhistVolUsed=new TH2F("fhistVolUsed","Hist with N modules used for a given module==(r,phi) zone",fnPhi,binningphi,fnZ,binningz);
349 fSigmaVolZ=new TH2F("fSigmaVolZ","Hist with Sigma of Residual distribution for each module",fnPhi,binningphi,fnZ,binningz);
350 fhEmpty=new TH2F("fhEmpty","Hist for getting binning",fnPhi,binningphi,fnZ,binningz);
351 fVolNTracks->SetXTitle("Volume #phi");
352 fVolNTracks->SetYTitle("Volume z ");
353 fhistVolNptsUsed->SetXTitle("Volume #phi");
354 fhistVolNptsUsed->SetYTitle("Volume z ");
355 fhistVolUsed->SetXTitle("Volume #phi");
356 fhistVolUsed->SetYTitle("Volume z ");
357 fSigmaVolZ->SetXTitle("Volume #phi");
358 fSigmaVolZ->SetYTitle("Volume z ");
bad3394b 359 } else{
146f76de 360 fVolNTracks=new TH2F("fVolNTracks","Hist with N tracks passing through a given module==(r,phi) zone",50,-3.2,3.2,100,-80,80);
361 fhistVolNptsUsed=new TH2F("fhistVolNptsUsed","Hist with N points used for given module==(r,phi) ",50,-3.2,3.2,100,-80,80);
362 fhistVolUsed=new TH2F("fhistVolUsed","Hist with N modules used for a given module==(r,phi) zone",50,-3.2,3.2,100,-80,80);
363 fSigmaVolZ=new TH2F("fSigmaVolZ","Hist with Sigma of Residual distribution for each module",50,-3.2,3.2,100,-80,80);
364 fhEmpty=new TH2F("fhEmpty","Hist for getting binning",50,-3.2,3.2,100,-80,80);
365 fVolNTracks->SetXTitle("Volume #phi");
366 fVolNTracks->SetYTitle("Volume z ");
367 fhistVolNptsUsed->SetXTitle("Volume #phi");
368 fhistVolNptsUsed->SetYTitle("Volume z ");
369 fhistVolUsed->SetXTitle("Volume #phi");
370 fhistVolUsed->SetYTitle("Volume z ");
371 fSigmaVolZ->SetXTitle("Volume #phi");
372 fSigmaVolZ->SetYTitle("Volume z ");
373 }
bad3394b 374
146f76de 375 fpTrackVolIDs=new TArrayI(fnHist);
376 fVolUsed=new TArrayI*[fnHist];
377 fVolVolids=new TArrayI*[fnHist];
378 fLastVolVolid=new Int_t[fnHist];
379
380 for (Int_t nhist=0;nhist<fnHist;nhist++){
381 fpTrackVolIDs->AddAt(volIDs->At(nhist),nhist);
382 aux=histnameRPHI;
383 aux+=volIDs->At(nhist);
5fd5f7a6 384 fVolResHistRPHI[nhist]=new TH1F("histname","histname",20000,-5.0,5.0);
146f76de 385 fVolResHistRPHI[nhist]->SetName(aux.Data());
386 fVolResHistRPHI[nhist]->SetTitle(aux.Data());
387
388 aux=histnameZ;
389 aux+=volIDs->At(nhist);
5fd5f7a6 390 fResHistZ[nhist]=new TH1F("histname","histname",20000,-5.0,5.0);
146f76de 391 fResHistZ[nhist]->SetName(aux.Data());
392 fResHistZ[nhist]->SetTitle(aux.Data());
393
bad3394b 394 aux=histnameX;
395 aux+=volIDs->At(nhist);
5fd5f7a6 396 fResHistX[nhist]=new TH1F("histname","histname",20000,-5.0,5.0);
bad3394b 397 fResHistX[nhist]->SetName(aux.Data());
398 fResHistX[nhist]->SetTitle(aux.Data());
399
400 aux=histnameX;
401 aux+=volIDs->At(nhist);
402 aux.Append("LocalSDDLeft");
5fd5f7a6 403 fResHistXLocsddL[nhist]=new TH1F("histname","histname",20000,-5.0,5.0);
bad3394b 404 fResHistXLocsddL[nhist]->SetName(aux.Data());
405 fResHistXLocsddL[nhist]->SetTitle(aux.Data());
406 aux=histnameX;
407
408 aux+=volIDs->At(nhist);
409 aux.Append("LocalSDDRight");
5fd5f7a6 410 fResHistXLocsddR[nhist]=new TH1F("histname","histname",20000,-5.0,5.0);
bad3394b 411 fResHistXLocsddR[nhist]->SetName(aux.Data());
412 fResHistXLocsddR[nhist]->SetTitle(aux.Data());
413
414 aux="fHistCoordGlobY";
415 fHistCoordGlobY[nhist]=new TH1F("histname","histname",24000,-30.,30.);
416 fHistCoordGlobY[nhist]->SetName(aux.Data());
417 fHistCoordGlobY[nhist]->SetTitle(aux.Data());
418
146f76de 419 aux=histnamePullRPHI;
420 aux+=volIDs->At(nhist);
421 fPullHistRPHI[nhist]=new TH1F("histname","histname",100,-7.,7.);
422 fPullHistRPHI[nhist]->SetName(aux.Data());
423 fPullHistRPHI[nhist]->SetTitle(aux.Data());
424
425 aux=histnamePullZ;
426 aux+=volIDs->At(nhist);
427 fPullHistZ[nhist]=new TH1F("histname","histname",100,-7.,7.);
428 fPullHistZ[nhist]->SetName(aux.Data());
429 fPullHistZ[nhist]->SetTitle(aux.Data());
430
431 aux=histnameDirPhi;
432 aux+=volIDs->At(nhist);
433 fTrackDirPhi[nhist]=new TH1F("histname","histname",100,-180,180);
434 fTrackDirPhi[nhist]->SetName(aux.Data());
435 fTrackDirPhi[nhist]->SetTitle(aux.Data());
436
437 aux=histnameDirLambda;
438 aux+=volIDs->At(nhist);
439 fTrackDirLambda[nhist]=new TH1F("histname","histname",100,0,180);
440 fTrackDirLambda[nhist]->SetName(aux.Data());
441 fTrackDirLambda[nhist]->SetTitle(aux.Data());
442
443 aux=histnameDirLambda2;
444 aux+=volIDs->At(nhist);
445 fTrackDirLambda2[nhist]=new TH1F("histname","histname",100,0,180);
446 fTrackDirLambda2[nhist]->SetName(aux.Data());
447 fTrackDirLambda2[nhist]->SetTitle(aux.Data());
448
449 aux=histnameDirAlpha;
450 aux+=volIDs->At(nhist);
451 fTrackDirAlpha[nhist]=new TH1F("histname","histname",100,-180,180);
452 fTrackDirAlpha[nhist]->SetName(aux.Data());
453 fTrackDirAlpha[nhist]->SetTitle(aux.Data());
454
455 aux=histnameDir;
456 aux+=volIDs->At(nhist);
457 fTrackDir[nhist]=new TH2F("histname","histname",100,-90.,90.,100,-180.,180.);
458 fTrackDir[nhist]->SetName(aux.Data());
459 fTrackDir[nhist]->SetTitle(aux.Data());
460
461 aux=histnameGlob;
462 aux+=volIDs->At(nhist);
463 fResHistGlob[nhist]=new TH1F("histname","histname",400,-0.08,0.08);
464 fResHistGlob[nhist]->SetName(aux.Data());
465 fResHistGlob[nhist]->SetTitle(aux.Data());
466
467 aux=histnameCorrVol;
468 aux+=volIDs->At(nhist);
469 fhistCorrVol[nhist]=new TH2F("histname","histname",50,-3.2,3.2,100,-80,80);
470
471
472 fhistCorrVol[nhist]->SetName(aux.Data());
473 fhistCorrVol[nhist]->SetTitle(aux.Data());
474 fhistCorrVol[nhist]->SetXTitle("Volume #varphi");
475 fhistCorrVol[nhist]->SetYTitle("Volume z ");
476 fVolVolids[nhist]=new TArrayI(100);
477 fVolUsed[nhist]=new TArrayI(1000);
478 fLastVolVolid[nhist]=0;
479
480 }
481 fWriteHist=kFALSE;
482
483 return;
484}
485
486//____________________________________________________________________________
487void AliITSResidualsAnalysis::ListVolUsed(TTree *pointsTree,TArrayI ***arrayIndex,Int_t **lastIndex)
488{
489 //
490 // This is copied from AliAlignmentClass::LoadPoints() method
491 //
bad3394b 492 if(!gGeoManager) AliGeomManager::LoadGeometry(GetFileNameGeometry());
146f76de 493 Int_t volIDalignable,volIDpoint,iModule;
494 AliTrackPoint p;
495 AliTrackPointArray* array = 0;
496 pointsTree->SetBranchAddress("SP", &array);
497
498
499 for(Int_t ivol=0;ivol<fnHist;ivol++){
500 Int_t lastused=0;
501 volIDalignable=fpTrackVolIDs->At(ivol);
502 AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer((UShort_t)volIDalignable,iModule);
503
504 Int_t nArraysId = lastIndex[iLayer-AliGeomManager::kFirstLayer][iModule];
146f76de 505 TArrayI *index = arrayIndex[iLayer-AliGeomManager::kFirstLayer][iModule];
506 for (Int_t iArrayId = 0;iArrayId < nArraysId; iArrayId++) {
507
508 // Get tree entry
509 Int_t entry = (*index)[iArrayId];
510
511 pointsTree->GetEvent(entry);
512 if (!array) {
513 AliWarning("Wrong space point array index!");
514 continue;
515 }
516
517 // Get the space-point array
518 Int_t modnum,nPoints = array->GetNPoints();
519
520 for (Int_t iPoint = 0; iPoint < nPoints; iPoint++) {
bad3394b 521
522
146f76de 523 array->GetPoint(p,iPoint);
524
525 AliGeomManager::ELayerID layer = AliGeomManager::VolUIDToLayer(p.GetVolumeID(),modnum);
526 // check if the layer id is valid
527 if ((layer < AliGeomManager::kFirstLayer) ||
528 (layer >= AliGeomManager::kLastLayer)) {
529 AliError(Form("Layer index is invalid: %d (%d -> %d) !",
530 layer,AliGeomManager::kFirstLayer,AliGeomManager::kLastLayer-1));
531 continue;
532 }
bad3394b 533
534
146f76de 535 if ((modnum >= AliGeomManager::LayerSize(layer)) ||
536 (modnum < 0)) {
537 AliError(Form("Module number inside layer %d is invalid: %d (0 -> %d)",
538 layer,modnum,AliGeomManager::LayerSize(layer)));
539 continue;
540 }
541 if (layer > AliGeomManager::kSSD2) continue; // ITS only
542
bad3394b 543
146f76de 544 volIDpoint=(Int_t)p.GetVolumeID();
bad3394b 545 if (volIDpoint==volIDalignable) continue;
146f76de 546 Int_t size = fVolVolids[ivol]->GetSize();
547 // If needed allocate new size
548 if (fLastVolVolid[ivol]>=size){// Warning: fLAST[NHIST] is useless
549 fVolVolids[ivol]->Set(size + 1000);
550 }
bad3394b 551
552
146f76de 553 fVolVolids[ivol]->AddAt(volIDpoint,fLastVolVolid[ivol]);
554 fLastVolVolid[ivol]++;
555 Bool_t usedVol=kFALSE;
556 for(Int_t used=0;used<lastused;used++){
557 if(fVolUsed[ivol]->At(used)==volIDpoint){
558 usedVol=kTRUE;
559 break;
560 }
561 }
bad3394b 562
563
146f76de 564 if (!usedVol){
565 size = fVolUsed[ivol]->GetSize();
566 // If needed allocate new size
567 if (lastused>= size){
568 fVolUsed[ivol]->Set(size + 1000);
569 }
570 fVolUsed[ivol]->AddAt(volIDpoint,lastused);
571 lastused++;
572 }
573
bad3394b 574 FillVolumeCorrelationHists(ivol,volIDalignable,volIDpoint,usedVol);
575
576 }// end loop
146f76de 577 }
578 }
579 fWriteHist=kTRUE;
580 return;
581}
582
583//____________________________________________________________________________
584void AliITSResidualsAnalysis::FillVolumeCorrelationHists(Int_t ivol,Int_t volIDalignable,Int_t volIDpoint,Bool_t usedVol) const
585{
586 //
587 // Fill the histograms with the correlations between volumes
588 //
589
bad3394b 590
591 if(!gGeoManager) AliGeomManager::LoadGeometry(GetFileNameGeometry());
592 Double_t translGlobal[3];
593 // Double_t radius;
594 Double_t phi;
595 // const char *symname,*volpath;
596 /* TGeoPNEntry *pne;
597 TGeoPhysicalNode *pn;
598 TGeoHMatrix *globMatrix;
599
600
601 symname = AliGeomManager::SymName(volIDalignable);
602 pne = gGeoManager->GetAlignableEntry(symname);
603 volpath=pne->GetTitle();
604 pn=gGeoManager->MakePhysicalNode(volpath);
605 globMatrix=pn->GetMatrix();
606 */
607
608 AliGeomManager::GetOrigTranslation(volIDalignable,translGlobal);
609 // radius=TMath::Sqrt(transGlobal[0]*transGlobal[0]+transGlobal[1]*transGlobal[1]);
610 phi=TMath::ATan2(translGlobal[1],translGlobal[0]);
611 fhistVolNptsUsed->Fill(phi,translGlobal[2]);
612 if(!usedVol){
613 fhistVolUsed->Fill(phi,translGlobal[2]);
614 }
615
616 /* symname = AliGeomManager::SymName(volIDpoint);
617 pne = gGeoManager->GetAlignableEntry(symname);
618 volpath=pne->GetTitle();
619 pn=gGeoManager->MakePhysicalNode(volpath);
620 globMatrix=pn->GetMatrix();
621 transGlobal=globMatrix->GetTranslation();
622 */
623 AliGeomManager::GetOrigTranslation(volIDpoint,translGlobal);
624 // radius=TMath::Sqrt(transGlobal[0]*transGlobal[0]+transGlobal[1]*transGlobal[1]);
625 phi=TMath::ATan2(translGlobal[1],translGlobal[0]);
626
627 fhistCorrVol[ivol]->Fill(phi,translGlobal[2]);
146f76de 628
629 return;
630}
631
146f76de 632//____________________________________________________________________________
bad3394b 633void AliITSResidualsAnalysis::FillResidualsH(AliTrackPointArray *points,
634 AliTrackPointArray *pTrack) const
146f76de 635{
636 //
637 // Method that fills the histograms with the residuals
638 //
639
640 Int_t volIDpoint;
641 Float_t xyz[3],xyz2[3];
bad3394b 642 Double_t xyzD[3],xyz2D[3];
643 Double_t loc[3],loc2[3];
644
645 Float_t resRPHI,resGlob,resZ,resX;
646
647 Double_t pullrphi,sign,phi;
146f76de 648 AliTrackPoint p,pTr;
bad3394b 649
146f76de 650 for(Int_t ipoint=0;ipoint<points->GetNPoints();ipoint++){
bad3394b 651
652 //pTrack->GetPoint(pTr,ipoint);
146f76de 653 points->GetPoint(p,ipoint);
654 volIDpoint=(Int_t)p.GetVolumeID();
655 p.GetXYZ(xyz);
bad3394b 656
146f76de 657 pTrack->GetPoint(pTr,ipoint);
146f76de 658 pTr.GetXYZ(xyz2);
bad3394b 659
660 for(Int_t i=0;i<3;i++){
661 xyzD[i]=xyz[i];
662 xyz2D[i]=xyz2[i];
663 }
664
665 phi = TMath::ATan2(xyz[1],xyz[0]);//<-watch out: phi of the pPoints!
666
667 resZ=xyz2[2]-xyz[2];
668 resX=xyz2[0]-xyz[0];
669
146f76de 670 resRPHI=TMath::Sqrt((xyz2[0]-xyz[0])*(xyz2[0]-xyz[0])+(xyz2[1]-xyz[1])*(xyz2[1]-xyz[1]));
bad3394b 671
146f76de 672 sign=TMath::ATan2(xyz2[1],xyz2[0])-TMath::ATan2(xyz[1],xyz[0]);
673 if(sign!=0.){
674 sign=sign/TMath::Abs(sign);
675 resRPHI=resRPHI*sign;
bad3394b 676
146f76de 677 }
678 else{
679 pullrphi=0.;
680 resRPHI=0.;
681 }
682
146f76de 683 resGlob=TMath::Sqrt((xyz2[0]-xyz[0])*(xyz2[0]-xyz[0])+(xyz2[1]-xyz[1])*(xyz2[1]-xyz[1])+(xyz2[2]-xyz[2])*(xyz2[2]-xyz[2]));
bad3394b 684
146f76de 685 for(Int_t ivolIDs=0;ivolIDs<fpTrackVolIDs->GetSize();ivolIDs++){
686 if(volIDpoint==fpTrackVolIDs->At(ivolIDs)){
bad3394b 687
146f76de 688 fVolResHistRPHI[ivolIDs]->Fill(resRPHI);
689 fResHistZ[ivolIDs]->Fill(resZ);
bad3394b 690 fResHistX[ivolIDs]->Fill(resX);
691 fHistCoordGlobY[ivolIDs]->Fill(xyz[1]);
692
693 Int_t modIndex = -1; // SDD Section
694 if(AliGeomManager::VolUIDToLayer(volIDpoint)==3) modIndex=volIDpoint-6144+240;
695 if(AliGeomManager::VolUIDToLayer(volIDpoint)==4) modIndex=volIDpoint-8192+240+84;
696 if(modIndex>0){
697 AliITSgeomTGeo::GlobalToLocal(modIndex,xyzD,loc); // error here!?
698 AliITSgeomTGeo::GlobalToLocal(modIndex,xyz2D,loc2);
699 Float_t rexloc=loc2[0]-loc[0];
700 //cout<<"Residual: "<<volIDpoint<<" "<<loc[0]<<" -> "<<rexloc<<endl;
701 if(loc[0]<0){
702 fResHistXLocsddR[ivolIDs]->Fill(rexloc);
703 }else{
704 fResHistXLocsddL[ivolIDs]->Fill(rexloc);
705 }
706 }
146f76de 707 fResHistGlob[ivolIDs]->Fill(resGlob);
708
146f76de 709 fTrackDirPhiAll->Fill(phi);
bad3394b 710 fTrackDirPhi[ivolIDs]->Fill(phi);
146f76de 711
146f76de 712 if(fsingleLayer){
713 Int_t binz,binphi;
714 Float_t globalPhi,globalZ;
715 if(kTRUE||(fvolidsToBin[ivolIDs][0]!=volIDpoint)){
716 binphi=GetBinPhiZ((Int_t)volIDpoint,&binz);
717 }
bad3394b 718 else{
719 // This in the case of alignment of one entire layer
720 // (fnHIst=layersize) may reduce iterations:
721 // remind of that fsingleLayer->fnHista<layerSize
146f76de 722 binphi=fvolidsToBin[ivolIDs][1];
723 binz=fvolidsToBin[ivolIDs][2];
724 }
725 globalPhi=fCoordToBinTable[binphi][binz][0];
726 globalZ=fCoordToBinTable[binphi][binz][1];
727
728 fVolNTracks->Fill(globalPhi,globalZ);
729 }
730 else fVolNTracks->Fill(TMath::ATan2(xyz[1],xyz[0]),xyz[2]);
731 }
732 }
733 }
734}
735
146f76de 736//____________________________________________________________________________
bad3394b 737Bool_t AliITSResidualsAnalysis::SaveHists(Int_t minNpoints, TString outname) const
146f76de 738{
739 //
740 // Saves the histograms into a tree and saves the tree into a file
741 //
742
146f76de 743
bad3394b 744 // Output file
745 TFile *hFile=new TFile(outname.Data(),"RECREATE","File containing the Residuals Tree");
746
747 // TTree with the residuals
748 TTree *analysisTree=new TTree("analysisTree","Tree with the residuals");
749
750 // Declares Variables to be stored into the TTree
146f76de 751 TF1 *gauss;
bad3394b 752 Int_t volID,entries,nHistAnalyzed=0;
753 Double_t meanResRPHI,meanResZ,meanResX,rmsResRPHI,rmsResZ,rmsResX,coordVol[3],x,y,z;
754 TH1F *histRPHI = new TH1F();
755 TH1F *histZ = new TH1F();
756 TH1F *histX = new TH1F();
757 TH1F *histXLocsddL = new TH1F();
758 TH1F *histXLocsddR = new TH1F();
759 TH1F *histCoordGlobY = new TH1F();
760 // Note: 0 = RPHI, 1 = Z
761
762
763 // Branching the TTree
146f76de 764 analysisTree->Branch("volID",&volID,"volID/I");
bad3394b 765 analysisTree->Branch("x",&x,"x/D");
766 analysisTree->Branch("y",&y,"y/D");
767 analysisTree->Branch("z",&z,"z/D");
768 analysisTree->Branch("meanResRPHI",&meanResRPHI,"meanResRPHI/D");
769 analysisTree->Branch("meanResZ",&meanResZ,"meanResZ/D");
770 analysisTree->Branch("meanResX",&meanResX,"meanResX/D");
771 analysisTree->Branch("rmsResRPHI",&rmsResRPHI,"rmsResRPHI/D");
772 analysisTree->Branch("rmsResZ",&rmsResZ,"rmsResZ/D");
773
146f76de 774 analysisTree->Branch("histRPHI","TH1F",&histRPHI,128000,0);
146f76de 775 analysisTree->Branch("histZ","TH1F",&histZ,128000,0);
bad3394b 776 analysisTree->Branch("histX","TH1F",&histX,128000,0);
777 analysisTree->Branch("histXLocsddL","TH1F",&histXLocsddL,128000,0);
778 analysisTree->Branch("histXLocsddR","TH1F",&histXLocsddR,128000,0);
779 analysisTree->Branch("histCoordGlobY","TH1F",&histCoordGlobY,128000,0);
780
781 Int_t blimps=0;
782
146f76de 783 for(Int_t j=0;j<fnHist;j++){
bad3394b 784
146f76de 785 volID=fpTrackVolIDs->At(j);
bad3394b 786 AliGeomManager::GetTranslation(volID,coordVol);
787 x=coordVol[0];
788 y=coordVol[1];
789 z=coordVol[2];
790
791 entries=(Int_t)(fResHistGlob[j]->GetEntries());
792 blimps+=entries;
793
794 if(entries>=minNpoints){
146f76de 795 nHistAnalyzed++;
bad3394b 796
797 // Entries
798 //entries=(Int_t)fVolResHistRPHI[j]->GetEntries();
799
800 // Filling the RPHI
146f76de 801 histRPHI=fVolResHistRPHI[j];
bad3394b 802 rmsResRPHI=fVolResHistRPHI[j]->GetRMS();
803 // Fit (for average)
804 gauss=new TF1("gauss","gaus",-3*rmsResRPHI,3*rmsResRPHI);
805 fVolResHistRPHI[j]->Fit("gauss","QRN");
806 meanResRPHI=gauss->GetParameter(1);
807
808 // Filling the Z
146f76de 809 histZ=fResHistZ[j];
bad3394b 810 rmsResZ=fResHistZ[j]->GetRMS();
811 // Fit (for average)
812 gauss=new TF1("gauss","gaus",-3*rmsResZ,3*rmsResZ);
813 fResHistZ[j]->Fit("gauss","QRN");
814 meanResZ=gauss->GetParameter(1);
815
816 // Filling the X
817 histX=fResHistX[j];
818 rmsResX=fResHistX[j]->GetRMS();
819 // Fit (for average)
820 gauss=new TF1("gauss","gaus",-3*rmsResX,3*rmsResX);
821 fResHistX[j]->Fit("gauss","QRN");
822 meanResX=gauss->GetParameter(1);
823
824 histXLocsddL=fResHistXLocsddL[j];
825 histXLocsddR=fResHistXLocsddR[j];
826 histCoordGlobY=fHistCoordGlobY[j];
827
828 analysisTree->Fill();
829 }else{
830
831 // Entries
832 //entries=(Int_t)fVolResHistRPHI[j]->GetEntries();
833
834 // Filling the RPHI
835 histRPHI=fVolResHistRPHI[j];
836 rmsResRPHI=-1.0;
837 meanResRPHI=0.0;
838
839 // Filling the Z
840 histZ=fResHistZ[j];
841 rmsResZ=-1.0;
842 meanResZ=0.0;
843
844 // Filling the X
845 histX=fResHistX[j];
846 rmsResX=-1.0;
847 meanResX=0.0;
848 histXLocsddL=fResHistXLocsddL[j];
849 histXLocsddR=fResHistXLocsddR[j];
850 histCoordGlobY=fHistCoordGlobY[j];
851
146f76de 852 analysisTree->Fill();
bad3394b 853
146f76de 854 }
bad3394b 855
146f76de 856 }
bad3394b 857
858 cout<<"-> Modules Analyzed: "<<nHistAnalyzed<<endl;
859 cout<<" With "<<blimps<<" events"<<endl;
860
861 if(blimps>0){
862 hFile->cd();
863 analysisTree->Write();
146f76de 864 fVolNTracks->Write();
146f76de 865 fhEmpty->Write();
866 if(fWriteHist){
bad3394b 867 //TCanvas *color = new TCanvas("color","fhistVolUsed",800,600);
868 //fhistVolUsed->DrawCopy("COLZ");
146f76de 869 fSigmaVolZ->Write();
870 fhistVolUsed->Write();
bad3394b 871 /* fTrackDirPhiAll->Write();
872 fTrackDirLambdaAll->Write();
873 fTrackDirLambda2All->Write();
874 fTrackDirAlphaAll->Write();
875 fTrackDirAll->Write();
876 fTrackDir2All->Write();
877 fTrackDirXZAll->Write();
878 hFile->Close();*/
146f76de 879 fhistVolNptsUsed->Write();
bad3394b 880 hFile->mkdir("CorrVol");
881 hFile->cd("CorrVol");
882 for(Int_t corr=0;corr<fnHist;corr++)fhistCorrVol[corr]->Write();
146f76de 883 }
bad3394b 884 hFile->cd();
885 // fhistVolNptsUsed->Write();
886 hFile->Close();
146f76de 887 return kTRUE;
bad3394b 888 }else {
146f76de 889 delete analysisTree;
890 delete hFile;
891 return kFALSE;}
892}
893
146f76de 894//____________________________________________________________________________
895void AliITSResidualsAnalysis::DrawHists() const
896{
897 //
898 // Draws the histograms of the residuals and of the number of tracks
899 //
900
901 TString cname;
902 for(Int_t canv=0;canv<fnHist;canv++){
903 cname="canv Residuals";
904 cname+=canv;
905 TCanvas *c=new TCanvas(cname.Data(),cname.Data(),700,700);
906 c->Divide(3,1);
907 c->cd(1);
908 fVolResHistRPHI[canv]->Draw();
909 c->cd(2);
910 fResHistZ[canv]->Draw();
911 c->cd(3);
912 fResHistGlob[canv]->Draw();
913 }
914 cname="canv NVolTracks";
915
916 TCanvas *c2=new TCanvas(cname.Data(),cname.Data(),700,700);
917 c2->cd();
918 fVolNTracks->Draw();
919
920 return;
921}
922
923//____________________________________________________________________________
924Float_t** AliITSResidualsAnalysis::CheckSingleLayer(const TArrayI *volids)
925{
926 //
927 // Checks if volumes array is a single (ITS) layer or not
928 //
929
930 Float_t **binningzphi=new Float_t*[2];
931 Int_t iModule;
932 AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer((UShort_t)volids->At(0),iModule);
933 //Check that one single Layer is going to be aligned
934 for(Int_t nvol=0;nvol<volids->GetSize();nvol++){
935 if(iLayer != AliGeomManager::VolUIDToLayer((UShort_t)volids->At(nvol),iModule)){
936 printf("Wrong Layer! \n %d , %d , %d ,%d \n",(Int_t)AliGeomManager::VolUIDToLayer((UShort_t)volids->At(nvol),iModule),nvol,volids->GetSize(),iModule);
937 fsingleLayer=kFALSE;
bad3394b 938 return binningzphi;
146f76de 939 }
940 }
941
942 //Bool_t used=kFALSE;
943 switch (iLayer) {
944 case AliGeomManager::kSPD1:{
3dbc7836 945 fnPhi=kPhiSPD1;//kPhiSPD1;
946 fnZ=kZSPD1;//nZSPD1;
947 binningzphi[0]=new Float_t[kPhiSPD1+1];
948 binningzphi[1]=new Float_t[kZSPD1+1];
949 fCoordToBinTable=new Double_t**[kPhiSPD1];
146f76de 950 for(Int_t j=0;j<fnPhi;j++){
3dbc7836 951 fCoordToBinTable[j]=new Double_t*[kZSPD1];
146f76de 952 }
953 }; break;
954 case AliGeomManager::kSPD2:{
3dbc7836 955 fnPhi=kPhiSPD2;//kPhiSPD1;
956 fnZ=kZSPD2;//nZSPD1;
957 binningzphi[0]=new Float_t[kPhiSPD2+1];
958 binningzphi[1]=new Float_t[kZSPD2+1];
959 fCoordToBinTable=new Double_t**[kPhiSPD2];
146f76de 960 for(Int_t j=0;j<fnPhi;j++){
3dbc7836 961 fCoordToBinTable[j]=new Double_t*[kZSPD2];
146f76de 962 }
963
964 }; break; case AliGeomManager::kSDD1:{
3dbc7836 965 fnPhi=kPhiSDD1;//kPhiSPD1;
966 fnZ=kZSDD1;//nZSPD1;
967 binningzphi[0]=new Float_t[kPhiSDD1+1];
968 binningzphi[1]=new Float_t[kZSDD1+1];
969 fCoordToBinTable=new Double_t**[kPhiSDD1];
146f76de 970 for(Int_t j=0;j<fnPhi;j++){
3dbc7836 971 fCoordToBinTable[j]=new Double_t*[kZSDD1];
146f76de 972 }
973 }; break; case AliGeomManager::kSDD2:{
3dbc7836 974 fnPhi=kPhiSDD2;//kPhiSPD1;
975 fnZ=kZSDD2;//nZSPD1;
976 binningzphi[0]=new Float_t[kPhiSDD2+1];
977 binningzphi[1]=new Float_t[kZSDD2+1];
978 fCoordToBinTable=new Double_t**[kPhiSDD2];
146f76de 979 for(Int_t j=0;j<fnPhi;j++){
3dbc7836 980 fCoordToBinTable[j]=new Double_t*[kZSDD2];
146f76de 981 }
982 }; break; case AliGeomManager::kSSD1:{
3dbc7836 983 fnPhi=kPhiSSD1;//kPhiSPD1;
984 fnZ=kZSSD1;//nZSPD1;
985 binningzphi[0]=new Float_t[kPhiSSD1+1];
986 binningzphi[1]=new Float_t[kZSSD1+1];
987 fCoordToBinTable=new Double_t**[kPhiSSD1];
146f76de 988 for(Int_t j=0;j<fnPhi;j++){
3dbc7836 989 fCoordToBinTable[j]=new Double_t*[kZSSD1];
146f76de 990 }
991 }; break; case AliGeomManager::kSSD2:{
3dbc7836 992 fnPhi=kPhiSSD2;//kPhiSPD1;
993 fnZ=kZSSD2;//nZSPD1;
994 binningzphi[0]=new Float_t[kPhiSSD2+1];
995 binningzphi[1]=new Float_t[kZSSD2+1];
996 fCoordToBinTable=new Double_t**[kPhiSSD2];
146f76de 997 for(Int_t j=0;j<fnPhi;j++){
3dbc7836 998 fCoordToBinTable[j]=new Double_t*[kZSSD2];
146f76de 999 }
1000 }; break;
1001
1002 default:{
1003 printf("Wrong Layer Label! \n");
1004 fsingleLayer=kFALSE;
5fd5f7a6 1005 //////////////////////
1006 fnPhi=1;//kPhiSPD1;
1007 fnZ=1;//nZSPD1;
1008 binningzphi[0]=new Float_t[1];
1009 binningzphi[1]=new Float_t[1];
1010 fCoordToBinTable=new Double_t**[1];
1011 for(Int_t j=0;j<fnPhi;j++){
1012 fCoordToBinTable[j]=new Double_t*[1];
1013 }
1014 return binningzphi;
1015 /////////////////////
1016 // return 0x0;
146f76de 1017 }
1018 }
1019 fsingleLayer=kTRUE;
1020 return binningzphi;
1021}
1022
1023//____________________________________________________________________________
1024Bool_t AliITSResidualsAnalysis::SetBinning(const TArrayI *volids,Float_t *phiBin,Float_t *zBin)
1025{
1026 //
1027 // Sets the correct binning
1028 //
1029
1030 if(!fsingleLayer)return kFALSE;
bad3394b 1031 // const char *volpath,*symname;
146f76de 1032 Int_t iModule;
1033 Int_t *orderArrayPhi,*orderArrayZ;
1034 UShort_t volID;
bad3394b 1035 Double_t *phiArray,*zArray,*phiArrayOrdered,*zArrayOrdered;
1036 Double_t translGlobal[3];
146f76de 1037 Double_t lastPhimin=-10;
1038 Double_t lastZmin=-99;
1039 Int_t ***orderPhiZ;
bad3394b 1040 /* TGeoPNEntry *pne;
1041 TGeoPhysicalNode *pn;
1042 TGeoHMatrix *globMatrix;*/
146f76de 1043
1044 Bool_t used=kFALSE;
1045
1046 orderPhiZ=new Int_t**[fnPhi];
1047 phiArray=new Double_t[fnPhi];//phiBin[nModulesPhi+1];
1048 zArray=new Double_t[fnZ];//zBin[nModulesZ+1];
1049 phiArrayOrdered=new Double_t[fnPhi];
1050 zArrayOrdered=new Double_t[fnZ];
1051 orderArrayPhi=new Int_t[fnPhi];
1052 orderArrayZ=new Int_t[fnZ];
1053 for(Int_t k=0;k<fnZ;k++){
1054 orderArrayZ[k]=0;
1055 zArray[k]=-1000;
1056 }
1057 for(Int_t k=0;k<fnPhi;k++){
1058 orderArrayPhi[k]=0;
1059 phiArray[k]=-10;
1060 orderPhiZ[k]=new Int_t*[fnZ];
1061 }
1062
1063
1064 AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer((UShort_t)volids->At(0),iModule);
1065
1066 Int_t lastPhi=0,lastZ=0;
1067 for(iModule=0;iModule<AliGeomManager::LayerSize(iLayer);iModule++){
1068 fvolidsToBin[iModule]=new Int_t[3];
1069 volID=AliGeomManager::LayerToVolUID(iLayer,iModule);
1070 fvolidsToBin[iModule][0]=volID;
bad3394b 1071 /* symname = AliGeomManager::SymName(volID);
1072 pne = gGeoManager->GetAlignableEntry(symname);
1073 volpath=pne->GetTitle();
1074 pn=gGeoManager->MakePhysicalNode(volpath);
1075 globMatrix=pn->GetMatrix();
1076 translGlobal=globMatrix->GetTranslation();
1077
1078 */
1079 AliGeomManager::GetOrigTranslation(volID,translGlobal);
146f76de 1080
1081 for(Int_t j=0;j<lastPhi;j++){
1082 used=kFALSE;
bad3394b 1083 if(TMath::Abs(phiArray[j]-TMath::ATan2(translGlobal[1],translGlobal[0]))<2*TMath::Pi()/(10*fnPhi)){//10 is a safety factor but....
146f76de 1084 fvolidsToBin[iModule][1]=j;
1085 used=kTRUE;
1086 break;
1087 }
1088 }
1089 if(!used){
bad3394b 1090 phiArray[lastPhi]=TMath::ATan2(translGlobal[1],translGlobal[0]);
146f76de 1091 fvolidsToBin[iModule][1]=lastPhi;
1092 if(phiArray[lastPhi]<lastPhimin)lastPhimin=phiArray[lastPhi];
1093 lastPhi++;
1094 if(lastPhi>fnPhi){
1095 printf("Wrong Phi! \n");
1096 return kFALSE;}
1097 }
1098
1099 for(Int_t j=0;j<lastZ;j++){
1100 used=kFALSE;
bad3394b 1101 if(TMath::Abs(zArray[j]-translGlobal[2])<0.1){
146f76de 1102 fvolidsToBin[iModule][2]=j;
1103 used=kTRUE;
1104 break;
1105 }
1106 }
1107 if(!used){
1108 fvolidsToBin[iModule][2]=lastZ;
bad3394b 1109 zArray[lastZ]=translGlobal[2];
146f76de 1110 if(zArray[lastZ]<lastZmin)lastZmin=zArray[lastZ];
1111 lastZ++;
1112 if(lastZ>fnZ){
1113 printf("Wrong Z! \n");
1114 return kFALSE;}
1115 }
1116 }
1117
1118
1119 //ORDERING THE ARRAY OF PHI AND Z VALUES
1120 for(Int_t order=0;order<fnPhi;order++){
1121 for(Int_t j=0;j<fnPhi;j++){
1122 if((j!=order)&&(phiArray[order]>phiArray[j]))orderArrayPhi[order]++;
1123 }
1124 }
1125
1126 for(Int_t order=0;order<fnPhi;order++){
1127 for(Int_t j=0;j<fnPhi;j++){
1128 if(orderArrayPhi[j]==order){
1129 phiArrayOrdered[order]=phiArray[j];
1130 break;
1131 }
1132 }
1133 }
1134
1135
1136 for(Int_t order=0;order<fnZ;order++){
1137 for(Int_t j=0;j<fnZ;j++){
1138 if((j!=order)&&(zArray[order]>zArray[j]))orderArrayZ[order]++;
1139 }
1140 }
1141
1142
1143 for(Int_t order=0;order<fnZ;order++){
1144 for(Int_t j=0;j<fnZ;j++){
1145 if(orderArrayZ[j]==order){
1146 zArrayOrdered[order]=zArray[j];
1147 break;
1148 }
1149 }
1150 }
1151
1152
1153 //Filling the fCoordToBinTable
1154 for(Int_t j=0;j<fnPhi;j++){
1155 for(Int_t i=0;i<fnZ;i++){
1156 orderPhiZ[j][i]=new Int_t[2];
1157 orderPhiZ[j][i][0]=orderArrayPhi[j];
1158 orderPhiZ[j][i][1]=orderArrayZ[i];
1159 }
1160 }
1161
1162
1163 for(Int_t j=0;j<fnPhi;j++){
1164 for(Int_t i=0;i<fnZ;i++){
1165 fCoordToBinTable[j][i]=new Double_t[2];
1166 fCoordToBinTable[j][i][0]=phiArrayOrdered[j];
1167 fCoordToBinTable[j][i][1]=zArrayOrdered[i];
bad3394b 1168 printf("Now (binphi,binz)= %d, %d e (phi,z)=%f,%f \n",j,i,fCoordToBinTable[j][i][0],fCoordToBinTable[j][i][1]);
146f76de 1169 }
1170 }
1171 Int_t istar,jstar;
1172 for(iModule=0;iModule<fnPhi*fnZ;iModule++){
1173 istar=fvolidsToBin[iModule][1];
1174 jstar=fvolidsToBin[iModule][2];
1175 fvolidsToBin[iModule][1]=orderPhiZ[istar][jstar][0];
1176 fvolidsToBin[iModule][2]=orderPhiZ[istar][jstar][1];
1177 }
1178
1179
1180 //now constructing the binning
1181 for(Int_t iModPhi=0;iModPhi<fnPhi-1;iModPhi++){
1182 phiBin[iModPhi+1]=(Float_t)phiArrayOrdered[iModPhi]+0.5*(phiArrayOrdered[iModPhi+1]-phiArrayOrdered[iModPhi]);
1183 }
1184
1185 phiBin[0]=(Float_t)phiArrayOrdered[0]-(phiArrayOrdered[1]-phiArrayOrdered[0])/2.;
1186
1187 phiBin[fnPhi]=(Float_t)phiArrayOrdered[fnPhi-1]+(phiArrayOrdered[fnPhi-1]-phiArrayOrdered[fnPhi-2])/2.;
1188 for(Int_t iModPhi=0;iModPhi<fnPhi+1;iModPhi++){
1189 printf("Mean Phi mod %d su %d: %f \n",iModPhi+1,fnPhi,phiBin[iModPhi]);
1190 }
1191
1192 for(Int_t iModZ=0;iModZ<fnZ-1;iModZ++){
1193 zBin[iModZ+1]=(Float_t)zArrayOrdered[iModZ]+0.5*(zArrayOrdered[iModZ+1]-zArrayOrdered[iModZ]);
1194 }
1195 zBin[0]=(Float_t)zArrayOrdered[0]-0.5*(zArrayOrdered[1]-zArrayOrdered[0]);
1196 zBin[fnZ]=(Float_t)zArrayOrdered[fnZ-1]+0.5*(zArrayOrdered[1]-zArrayOrdered[0]);
1197
1198
1199 for(Int_t iModPhi=0;iModPhi<fnZ+1;iModPhi++){
1200 printf("Mean Z mod %d su %d: %f \n",iModPhi+1,fnZ,zBin[iModPhi]);
1201 }
1202 return kTRUE;
1203}
1204
1205
1206//____________________________________________________________________________
1207Int_t AliITSResidualsAnalysis::GetBinPhiZ(const Int_t volID,Int_t *binz) const
1208{
1209 //
1210 // Returns the correct Phi-Z bin
1211 //
1212
1213 if (!fsingleLayer){
1214 printf("No Single Layer reAlignment! \n");
1215 return 100;
1216 }
1217
1218 for(Int_t j=0;j<fnPhi*fnZ;j++){
1219 if(j==fnZ*fnPhi){
1220 printf("Wrong volume UID introduced! fnHist: %d volID: %d iter: %d \n",fnHist,volID,j);
1221 return 100;
1222 }
1223 if(fvolidsToBin[j][0]==volID){
1224
1225 *binz=fvolidsToBin[j][2];
1226 return fvolidsToBin[j][1];
1227 }
1228 }
1229
1230 return 100;
1231
1232}
1233
bad3394b 1234//___________________________________________________________________________
1235Int_t AliITSResidualsAnalysis::WhichSector(Int_t module) const
1236{
1237 //
1238 // This method returns the number of the SPD Sector
1239 // to which belongs the module (Sectors 0-9)
1240
1241 //--->cSect = 0 <---
1242 if(module==2048
1243 || module==2049
1244 || module==2050
1245 || module==2051
1246 || module==2052
1247 || module==2053
1248 || module==2054
1249 || module==2055
1250 || module==4096
1251 || module==4097
1252 || module==4098
1253 || module==4099
1254 || module==4100
1255 || module==4101
1256 || module==4102
1257 || module==4103
1258 || module==4104
1259 || module==4105
1260 || module==4106
1261 || module==4107
1262 || module==4108
1263 || module==4109
1264 || module==4110
1265 || module==4111) return 0;
1266
1267 //--->cSect = 1 <---
1268 if(module==2056
1269 || module==2057
1270 || module==2058
1271 || module==2059
1272 || module==2060
1273 || module==2061
1274 || module==2062
1275 || module==2063
1276 || module==4112
1277 || module==4113
1278 || module==4114
1279 || module==4115
1280 || module==4116
1281 || module==4117
1282 || module==4118
1283 || module==4119
1284 || module==4120
1285 || module==4121
1286 || module==4122
1287 || module==4123
1288 || module==4124
1289 || module==4125
1290 || module==4126
1291 || module==4127) return 1;
1292
1293 //--->cSect = 2 <---
1294 if(module==2064
1295 || module==2065
1296 || module==2066
1297 || module==2067
1298 || module==2068
1299 || module==2069
1300 || module==2070
1301 || module==2071
1302 || module==4128
1303 || module==4129
1304 || module==4130
1305 || module==4131
1306 || module==4132
1307 || module==4133
1308 || module==4134
1309 || module==4135
1310 || module==4136
1311 || module==4137
1312 || module==4138
1313 || module==4139
1314 || module==4140
1315 || module==4141
1316 || module==4142
1317 || module==4143) return 2;
1318
1319 //--->cSect = 3 <---
1320 if(module==2072
1321 || module==2073
1322 || module==2074
1323 || module==2075
1324 || module==2076
1325 || module==2077
1326 || module==2078
1327 || module==2079
1328 || module==4144
1329 || module==4145
1330 || module==4146
1331 || module==4147
1332 || module==4148
1333 || module==4149
1334 || module==4150
1335 || module==4151
1336 || module==4152
1337 || module==4153
1338 || module==4154
1339 || module==4155
1340 || module==4156
1341 || module==4157
1342 || module==4158
1343 || module==4159) return 3;
1344
1345 //--->cSect = 4 <---
1346 if(module==2080
1347 || module==2081
1348 || module==2082
1349 || module==2083
1350 || module==2084
1351 || module==2085
1352 || module==2086
1353 || module==2087
1354 || module==4160
1355 || module==4161
1356 || module==4162
1357 || module==4163
1358 || module==4164
1359 || module==4165
1360 || module==4166
1361 || module==4167
1362 || module==4168
1363 || module==4169
1364 || module==4170
1365 || module==4171
1366 || module==4172
1367 || module==4173
1368 || module==4174
1369 || module==4175) return 4;
1370
1371 //--->cSect = 5 <---
1372 if(module==2088
1373 || module==2089
1374 || module==2090
1375 || module==2091
1376 || module==2092
1377 || module==2093
1378 || module==2094
1379 || module==2095
1380 || module==4176
1381 || module==4177
1382 || module==4178
1383 || module==4179
1384 || module==4180
1385 || module==4181
1386 || module==4182
1387 || module==4183
1388 || module==4184
1389 || module==4185
1390 || module==4186
1391 || module==4187
1392 || module==4188
1393 || module==4189
1394 || module==4190
1395 || module==4191) return 5;
1396
1397 //--->cSect = 6 <---
1398 if(module==2096
1399 || module==2097
1400 || module==2098
1401 || module==2099
1402 || module==2100
1403 || module==2101
1404 || module==2102
1405 || module==2103
1406 || module==4192
1407 || module==4193
1408 || module==4194
1409 || module==4195
1410 || module==4196
1411 || module==4197
1412 || module==4198
1413 || module==4199
1414 || module==4200
1415 || module==4201
1416 || module==4202
1417 || module==4203
1418 || module==4204
1419 || module==4205
1420 || module==4206
1421 || module==4207) return 6;
1422
1423 //--->cSect = 7 <---
1424 if(module==2104
1425 || module==2105
1426 || module==2106
1427 || module==2107
1428 || module==2108
1429 || module==2109
1430 || module==2110
1431 || module==2111
1432 || module==4208
1433 || module==4209
1434 || module==4210
1435 || module==4211
1436 || module==4212
1437 || module==4213
1438 || module==4214
1439 || module==4215
1440 || module==4216
1441 || module==4217
1442 || module==4218
1443 || module==4219
1444 || module==4220
1445 || module==4221
1446 || module==4222
1447 || module==4223) return 7;
1448
1449 //--->cSect = 8 <---
1450 if(module==2112
1451 || module==2113
1452 || module==2114
1453 || module==2115
1454 || module==2116
1455 || module==2117
1456 || module==2118
1457 || module==2119
1458 || module==4224
1459 || module==4225
1460 || module==4226
1461 || module==4227
1462 || module==4228
1463 || module==4229
1464 || module==4230
1465 || module==4231
1466 || module==4232
1467 || module==4233
1468 || module==4234
1469 || module==4235
1470 || module==4236
1471 || module==4237
1472 || module==4238
1473 || module==4239) return 8;
1474
1475 //--->cSect = 9 <---
1476 if(module==2120
1477 || module==2121
1478 || module==2122
1479 || module==2123
1480 || module==2124
1481 || module==2125
1482 || module==2126
1483 || module==2127
1484 || module==4240
1485 || module==4241
1486 || module==4242
1487 || module==4243
1488 || module==4244
1489 || module==4245
1490 || module==4246
1491 || module==4247
1492 || module==4248
1493 || module==4249
1494 || module==4250
1495 || module==4251
1496 || module==4252
1497 || module==4253
1498 || module==4254
1499 || module==4255) return 9;
1500
1501 //printf("Module not belonging to SPD, sorry!");
1502 return -1;
1503
1504}
1505
146f76de 1506//____________________________________________________________________________
bad3394b 1507TArrayI* AliITSResidualsAnalysis::GetSPDSectorsVolids(Int_t sectors[10]) const
146f76de 1508{
1509 //
bad3394b 1510 // This method gets the volID Array for the chosen sectors.
1511 // You have to pass an array with a 1 for each selected sector.
1512 // i.e. sectors[10] = {1,1,0,0,0,0,0,0,1,0} -> Sector 0, 1, 9 selected.
146f76de 1513 //
1514
bad3394b 1515 Int_t nSect=0;
1516 Int_t iModule=0;
1517
146f76de 1518 if(!gGeoManager) AliGeomManager::LoadGeometry(GetFileNameGeometry());
1519
bad3394b 1520 for(Int_t co=0;co<10;co++){ //counts the number of sectors chosen
1521 if(sectors[co]==1) nSect++;
146f76de 1522 }
bad3394b 1523
1524 if(nSect<1){ //if no sector chosen -> exit
1525 Printf("Error! No Sector/s Selected!");
1526 return 0x0;
146f76de 1527 }
1528
bad3394b 1529 TArrayI *volIDs = new TArrayI(nSect*24);
1530
1531 if(sectors[0]==1){ //--->cSect = 0 <---
1532 volIDs->AddAt(2048,iModule); iModule++;
1533 volIDs->AddAt(2049,iModule); iModule++;
1534 volIDs->AddAt(2050,iModule); iModule++;
1535 volIDs->AddAt(2051,iModule); iModule++;
1536 volIDs->AddAt(2052,iModule); iModule++;
1537 volIDs->AddAt(2053,iModule); iModule++;
1538 volIDs->AddAt(2054,iModule); iModule++;
1539 volIDs->AddAt(2055,iModule); iModule++;
1540 volIDs->AddAt(4096,iModule); iModule++;
1541 volIDs->AddAt(4097,iModule); iModule++;
1542 volIDs->AddAt(4098,iModule); iModule++;
1543 volIDs->AddAt(4099,iModule); iModule++;
1544 volIDs->AddAt(4100,iModule); iModule++;
1545 volIDs->AddAt(4101,iModule); iModule++;
1546 volIDs->AddAt(4102,iModule); iModule++;
1547 volIDs->AddAt(4103,iModule); iModule++;
1548 volIDs->AddAt(4104,iModule); iModule++;
1549 volIDs->AddAt(4105,iModule); iModule++;
1550 volIDs->AddAt(4106,iModule); iModule++;
1551 volIDs->AddAt(4107,iModule); iModule++;
1552 volIDs->AddAt(4108,iModule); iModule++;
1553 volIDs->AddAt(4109,iModule); iModule++;
1554 volIDs->AddAt(4110,iModule); iModule++;
1555 volIDs->AddAt(4111,iModule); iModule++;
1556 }
1557 if(sectors[1]==1){ //--->cSect = 1 <//---
1558 volIDs->AddAt(2056,iModule); iModule++;
1559 volIDs->AddAt(2057,iModule); iModule++;
1560 volIDs->AddAt(2058,iModule); iModule++;
1561 volIDs->AddAt(2059,iModule); iModule++;
1562 volIDs->AddAt(2060,iModule); iModule++;
1563 volIDs->AddAt(2061,iModule); iModule++;
1564 volIDs->AddAt(2062,iModule); iModule++;
1565 volIDs->AddAt(2063,iModule); iModule++;
1566 volIDs->AddAt(4112,iModule); iModule++;
1567 volIDs->AddAt(4113,iModule); iModule++;
1568 volIDs->AddAt(4114,iModule); iModule++;
1569 volIDs->AddAt(4115,iModule); iModule++;
1570 volIDs->AddAt(4116,iModule); iModule++;
1571 volIDs->AddAt(4117,iModule); iModule++;
1572 volIDs->AddAt(4118,iModule); iModule++;
1573 volIDs->AddAt(4119,iModule); iModule++;
1574 volIDs->AddAt(4120,iModule); iModule++;
1575 volIDs->AddAt(4121,iModule); iModule++;
1576 volIDs->AddAt(4122,iModule); iModule++;
1577 volIDs->AddAt(4123,iModule); iModule++;
1578 volIDs->AddAt(4124,iModule); iModule++;
1579 volIDs->AddAt(4125,iModule); iModule++;
1580 volIDs->AddAt(4126,iModule); iModule++;
1581 volIDs->AddAt(4127,iModule); iModule++;
1582 }
1583 if(sectors[2]==1){//--->cSect = 2 <//---
1584 volIDs->AddAt(2064,iModule); iModule++;
1585 volIDs->AddAt(2065,iModule); iModule++;
1586 volIDs->AddAt(2066,iModule); iModule++;
1587 volIDs->AddAt(2067,iModule); iModule++;
1588 volIDs->AddAt(2068,iModule); iModule++;
1589 volIDs->AddAt(2069,iModule); iModule++;
1590 volIDs->AddAt(2070,iModule); iModule++;
1591 volIDs->AddAt(2071,iModule); iModule++;
1592 volIDs->AddAt(4128,iModule); iModule++;
1593 volIDs->AddAt(4129,iModule); iModule++;
1594 volIDs->AddAt(4130,iModule); iModule++;
1595 volIDs->AddAt(4131,iModule); iModule++;
1596 volIDs->AddAt(4132,iModule); iModule++;
1597 volIDs->AddAt(4133,iModule); iModule++;
1598 volIDs->AddAt(4134,iModule); iModule++;
1599 volIDs->AddAt(4135,iModule); iModule++;
1600 volIDs->AddAt(4136,iModule); iModule++;
1601 volIDs->AddAt(4137,iModule); iModule++;
1602 volIDs->AddAt(4138,iModule); iModule++;
1603 volIDs->AddAt(4139,iModule); iModule++;
1604 volIDs->AddAt(4140,iModule); iModule++;
1605 volIDs->AddAt(4141,iModule); iModule++;
1606 volIDs->AddAt(4142,iModule); iModule++;
1607 volIDs->AddAt(4143,iModule); iModule++;
1608 }
1609 if(sectors[3]==1){//--->cSect = 3 <//---
1610 volIDs->AddAt(2072,iModule); iModule++;
1611 volIDs->AddAt(2073,iModule); iModule++;
1612 volIDs->AddAt(2074,iModule); iModule++;
1613 volIDs->AddAt(2075,iModule); iModule++;
1614 volIDs->AddAt(2076,iModule); iModule++;
1615 volIDs->AddAt(2077,iModule); iModule++;
1616 volIDs->AddAt(2078,iModule); iModule++;
1617 volIDs->AddAt(2079,iModule); iModule++;
1618 volIDs->AddAt(4144,iModule); iModule++;
1619 volIDs->AddAt(4145,iModule); iModule++;
1620 volIDs->AddAt(4146,iModule); iModule++;
1621 volIDs->AddAt(4147,iModule); iModule++;
1622 volIDs->AddAt(4148,iModule); iModule++;
1623 volIDs->AddAt(4149,iModule); iModule++;
1624 volIDs->AddAt(4150,iModule); iModule++;
1625 volIDs->AddAt(4151,iModule); iModule++;
1626 volIDs->AddAt(4152,iModule); iModule++;
1627 volIDs->AddAt(4153,iModule); iModule++;
1628 volIDs->AddAt(4154,iModule); iModule++;
1629 volIDs->AddAt(4155,iModule); iModule++;
1630 volIDs->AddAt(4156,iModule); iModule++;
1631 volIDs->AddAt(4157,iModule); iModule++;
1632 volIDs->AddAt(4158,iModule); iModule++;
1633 volIDs->AddAt(4159,iModule); iModule++;
1634 }
1635 if(sectors[4]==1){//--->cSect = 4 <//---
1636 volIDs->AddAt(2080,iModule); iModule++;
1637 volIDs->AddAt(2081,iModule); iModule++;
1638 volIDs->AddAt(2082,iModule); iModule++;
1639 volIDs->AddAt(2083,iModule); iModule++;
1640 volIDs->AddAt(2084,iModule); iModule++;
1641 volIDs->AddAt(2085,iModule); iModule++;
1642 volIDs->AddAt(2086,iModule); iModule++;
1643 volIDs->AddAt(2087,iModule); iModule++;
1644 volIDs->AddAt(4160,iModule); iModule++;
1645 volIDs->AddAt(4161,iModule); iModule++;
1646 volIDs->AddAt(4162,iModule); iModule++;
1647 volIDs->AddAt(4163,iModule); iModule++;
1648 volIDs->AddAt(4164,iModule); iModule++;
1649 volIDs->AddAt(4165,iModule); iModule++;
1650 volIDs->AddAt(4166,iModule); iModule++;
1651 volIDs->AddAt(4167,iModule); iModule++;
1652 volIDs->AddAt(4168,iModule); iModule++;
1653 volIDs->AddAt(4169,iModule); iModule++;
1654 volIDs->AddAt(4170,iModule); iModule++;
1655 volIDs->AddAt(4171,iModule); iModule++;
1656 volIDs->AddAt(4172,iModule); iModule++;
1657 volIDs->AddAt(4173,iModule); iModule++;
1658 volIDs->AddAt(4174,iModule); iModule++;
1659 volIDs->AddAt(4175,iModule); iModule++;
1660 }
1661 if(sectors[5]==1){//--->cSect = 5 <//---
1662 volIDs->AddAt(2088,iModule); iModule++;
1663 volIDs->AddAt(2089,iModule); iModule++;
1664 volIDs->AddAt(2090,iModule); iModule++;
1665 volIDs->AddAt(2091,iModule); iModule++;
1666 volIDs->AddAt(2092,iModule); iModule++;
1667 volIDs->AddAt(2093,iModule); iModule++;
1668 volIDs->AddAt(2094,iModule); iModule++;
1669 volIDs->AddAt(2095,iModule); iModule++;
1670 volIDs->AddAt(4176,iModule); iModule++;
1671 volIDs->AddAt(4177,iModule); iModule++;
1672 volIDs->AddAt(4178,iModule); iModule++;
1673 volIDs->AddAt(4179,iModule); iModule++;
1674 volIDs->AddAt(4180,iModule); iModule++;
1675 volIDs->AddAt(4181,iModule); iModule++;
1676 volIDs->AddAt(4182,iModule); iModule++;
1677 volIDs->AddAt(4183,iModule); iModule++;
1678 volIDs->AddAt(4184,iModule); iModule++;
1679 volIDs->AddAt(4185,iModule); iModule++;
1680 volIDs->AddAt(4186,iModule); iModule++;
1681 volIDs->AddAt(4187,iModule); iModule++;
1682 volIDs->AddAt(4188,iModule); iModule++;
1683 volIDs->AddAt(4189,iModule); iModule++;
1684 volIDs->AddAt(4190,iModule); iModule++;
1685 volIDs->AddAt(4191,iModule); iModule++;
1686 }
1687 if(sectors[6]==1){//--->cSect = 6 <//---
1688 volIDs->AddAt(2096,iModule); iModule++;
1689 volIDs->AddAt(2097,iModule); iModule++;
1690 volIDs->AddAt(2098,iModule); iModule++;
1691 volIDs->AddAt(2099,iModule); iModule++;
1692 volIDs->AddAt(2100,iModule); iModule++;
1693 volIDs->AddAt(2101,iModule); iModule++;
1694 volIDs->AddAt(2102,iModule); iModule++;
1695 volIDs->AddAt(2103,iModule); iModule++;
1696 volIDs->AddAt(4192,iModule); iModule++;
1697 volIDs->AddAt(4193,iModule); iModule++;
1698 volIDs->AddAt(4194,iModule); iModule++;
1699 volIDs->AddAt(4195,iModule); iModule++;
1700 volIDs->AddAt(4196,iModule); iModule++;
1701 volIDs->AddAt(4197,iModule); iModule++;
1702 volIDs->AddAt(4198,iModule); iModule++;
1703 volIDs->AddAt(4199,iModule); iModule++;
1704 volIDs->AddAt(4200,iModule); iModule++;
1705 volIDs->AddAt(4201,iModule); iModule++;
1706 volIDs->AddAt(4202,iModule); iModule++;
1707 volIDs->AddAt(4203,iModule); iModule++;
1708 volIDs->AddAt(4204,iModule); iModule++;
1709 volIDs->AddAt(4205,iModule); iModule++;
1710 volIDs->AddAt(4206,iModule); iModule++;
1711 volIDs->AddAt(4207,iModule); iModule++;
1712 }
1713 if(sectors[7]==1){ //--->cSect = 7 <//---
1714 volIDs->AddAt(2104,iModule); iModule++;
1715 volIDs->AddAt(2105,iModule); iModule++;
1716 volIDs->AddAt(2106,iModule); iModule++;
1717 volIDs->AddAt(2107,iModule); iModule++;
1718 volIDs->AddAt(2108,iModule); iModule++;
1719 volIDs->AddAt(2109,iModule); iModule++;
1720 volIDs->AddAt(2110,iModule); iModule++;
1721 volIDs->AddAt(2111,iModule); iModule++;
1722 volIDs->AddAt(4208,iModule); iModule++;
1723 volIDs->AddAt(4209,iModule); iModule++;
1724 volIDs->AddAt(4210,iModule); iModule++;
1725 volIDs->AddAt(4211,iModule); iModule++;
1726 volIDs->AddAt(4212,iModule); iModule++;
1727 volIDs->AddAt(4213,iModule); iModule++;
1728 volIDs->AddAt(4214,iModule); iModule++;
1729 volIDs->AddAt(4215,iModule); iModule++;
1730 volIDs->AddAt(4216,iModule); iModule++;
1731 volIDs->AddAt(4217,iModule); iModule++;
1732 volIDs->AddAt(4218,iModule); iModule++;
1733 volIDs->AddAt(4219,iModule); iModule++;
1734 volIDs->AddAt(4220,iModule); iModule++;
1735 volIDs->AddAt(4221,iModule); iModule++;
1736 volIDs->AddAt(4222,iModule); iModule++;
1737 volIDs->AddAt(4223,iModule); iModule++;
1738 }
1739 if(sectors[8]==1){//--->cSect = 8 <//---
1740 volIDs->AddAt(2112,iModule); iModule++;
1741 volIDs->AddAt(2113,iModule); iModule++;
1742 volIDs->AddAt(2114,iModule); iModule++;
1743 volIDs->AddAt(2115,iModule); iModule++;
1744 volIDs->AddAt(2116,iModule); iModule++;
1745 volIDs->AddAt(2117,iModule); iModule++;
1746 volIDs->AddAt(2118,iModule); iModule++;
1747 volIDs->AddAt(2119,iModule); iModule++;
1748 volIDs->AddAt(4224,iModule); iModule++;
1749 volIDs->AddAt(4225,iModule); iModule++;
1750 volIDs->AddAt(4226,iModule); iModule++;
1751 volIDs->AddAt(4227,iModule); iModule++;
1752 volIDs->AddAt(4228,iModule); iModule++;
1753 volIDs->AddAt(4229,iModule); iModule++;
1754 volIDs->AddAt(4230,iModule); iModule++;
1755 volIDs->AddAt(4231,iModule); iModule++;
1756 volIDs->AddAt(4232,iModule); iModule++;
1757 volIDs->AddAt(4233,iModule); iModule++;
1758 volIDs->AddAt(4234,iModule); iModule++;
1759 volIDs->AddAt(4235,iModule); iModule++;
1760 volIDs->AddAt(4236,iModule); iModule++;
1761 volIDs->AddAt(4237,iModule); iModule++;
1762 volIDs->AddAt(4238,iModule); iModule++;
1763 volIDs->AddAt(4239,iModule); iModule++;
1764 }
1765 if(sectors[9]==1){//--->cSect = 9 <//---
1766 volIDs->AddAt(2120,iModule); iModule++;
1767 volIDs->AddAt(2121,iModule); iModule++;
1768 volIDs->AddAt(2122,iModule); iModule++;
1769 volIDs->AddAt(2123,iModule); iModule++;
1770 volIDs->AddAt(2124,iModule); iModule++;
1771 volIDs->AddAt(2125,iModule); iModule++;
1772 volIDs->AddAt(2126,iModule); iModule++;
1773 volIDs->AddAt(2127,iModule); iModule++;
1774 volIDs->AddAt(4240,iModule); iModule++;
1775 volIDs->AddAt(4241,iModule); iModule++;
1776 volIDs->AddAt(4242,iModule); iModule++;
1777 volIDs->AddAt(4243,iModule); iModule++;
1778 volIDs->AddAt(4244,iModule); iModule++;
1779 volIDs->AddAt(4245,iModule); iModule++;
1780 volIDs->AddAt(4246,iModule); iModule++;
1781 volIDs->AddAt(4247,iModule); iModule++;
1782 volIDs->AddAt(4248,iModule); iModule++;
1783 volIDs->AddAt(4249,iModule); iModule++;
1784 volIDs->AddAt(4250,iModule); iModule++;
1785 volIDs->AddAt(4251,iModule); iModule++;
1786 volIDs->AddAt(4252,iModule); iModule++;
1787 volIDs->AddAt(4253,iModule); iModule++;
1788 volIDs->AddAt(4254,iModule); iModule++;
1789 volIDs->AddAt(4255,iModule); iModule++;
1790 }
1791
146f76de 1792 return volIDs;
bad3394b 1793
1794}
1795
1796//____________________________________________________________________________
1797TArrayI* AliITSResidualsAnalysis::GetITSLayersVolids(Int_t layers[6]) const
1798{
1799 //
1800 // This method gets the volID Array for the chosen layers.
1801 // You have to pass an array with a 1 for each selected layer.
1802 // i.e. layers[6] = {1,1,0,0,1,1} -> SPD + SSD
1803 //
1804
1805 if(!gGeoManager) AliGeomManager::LoadGeometry(GetFileNameGeometry());
1806
1807 Int_t size=0,last=0;
1808
1809 // evaluates the size of the array
1810 for(Int_t i=0;i<6;i++) if(layers[i]==1) size+=AliGeomManager::LayerSize(i+1);
1811
1812 if(size==0){
1813 printf("Error: no layer selected");
1814 return 0x0;
1815 }
1816
1817 TArrayI *volids = new TArrayI(size);
1818
1819 // fills the volId array only for the chosen layers
1820 for(Int_t ilayer=1;ilayer<7;ilayer++){
1821
1822 if(layers[ilayer-1]!=1) continue;
1823
1824 for(Int_t imod=0;imod<AliGeomManager::LayerSize(ilayer);imod++){
1825 volids->AddAt(AliGeomManager::LayerToVolUID(ilayer,imod),last);
1826 last++;
1827 }
1828 }
146f76de 1829
bad3394b 1830 return volids;
1831
146f76de 1832}
1833
1834//____________________________________________________________________________
1835void AliITSResidualsAnalysis::GetTrackDirClusterCov(AliTrackPoint *point,Double_t &phi,Double_t &lambda,Double_t &lambda2,Double_t &alpha,Double_t &xovery,Double_t &zovery) const
1836{
1837 //
1838 // ...
1839 //
1840
1841 TMatrixDSym cov(3);
1842 const Float_t *covvector=point->GetCov();
1843 cov(0,0)=covvector[0];
1844 cov(1,0)=cov(0,1)=covvector[1];
1845 cov(2,0)=cov(0,2)=covvector[2];
1846 cov(1,1)=covvector[3];
1847 cov(1,2)=cov(2,1)=covvector[4];
1848 cov(2,2)=covvector[5];
1849
1850 Double_t determinant=cov.Determinant();
1851 if(determinant!=0.){
1852 TMatrixD vect(3,3);
1853 TVectorD eigenvalues(3);
1854 const TMatrixDSymEigen keigen(cov);
1855 eigenvalues=keigen.GetEigenValues();
1856 vect=keigen.GetEigenVectors();
1857 Double_t mainvect[3];
1858 mainvect[0]=vect(0,0);
1859 mainvect[1]=vect(1,0);
1860 mainvect[2]=vect(2,0);
1861 if(mainvect[1]!=0.){
1862 xovery=mainvect[0]/mainvect[1];
1863 zovery=mainvect[2]/mainvect[1];
1864 }
1865 else {
1866 xovery=9999.;
1867 zovery=9999.;
1868 }
1869 if(mainvect[1]<0.){
1870 mainvect[0]=-1.*mainvect[0];
1871 mainvect[1]=-1.*mainvect[1];
1872 mainvect[2]=-1.*mainvect[2];
1873 }
1874 lambda2=TMath::ATan2(TMath::Sqrt(mainvect[0]*mainvect[0]+mainvect[2]*mainvect[2]),mainvect[1])*TMath::RadToDeg();
1875 lambda=TMath::ATan2(mainvect[2],TMath::Sqrt(mainvect[0]*mainvect[0]+mainvect[1]*mainvect[1]))*TMath::RadToDeg();
1876 phi=TMath::ATan2(mainvect[0],mainvect[2])*TMath::RadToDeg();
1877 alpha=TMath::ATan2(mainvect[1],mainvect[0])*TMath::RadToDeg();
1878 }
1879 else printf("determinant =0!, skip this point \n");
1880
1881 return;
1882}
1883
1884//____________________________________________________________________________
1885void AliITSResidualsAnalysis::CalculateResiduals(const TArrayI *volids,
bad3394b 1886 const TArrayI *volidsfit,
1887 AliGeomManager::ELayerID layerRangeMin,
1888 AliGeomManager::ELayerID layerRangeMax,
1889 TString outname)
146f76de 1890{
1891 // CalculateResiduals for a set of detector volumes.
1892 // Tracks are fitted only within
1893 // the range defined by the user
1894 // (by layerRangeMin and layerRangeMax)
1895 // or within the set of volidsfit
1896 // Repeat the procedure 'iterations' times
1897
146f76de 1898 Int_t nVolIds = volids->GetSize();
bad3394b 1899 if (nVolIds == 0) { AliError("Volume IDs array is empty!"); return; }
146f76de 1900
1901 // Load only the tracks with at least one
1902 // space point in the set of volume (volids)
1903
bad3394b 1904
146f76de 1905 //AliAlignmentTracks::SetPointsFilename(GetFileNameTrackPoints());
1906 AliAlignmentTracks::BuildIndex();
1907
bad3394b 1908 ListVolUsed(fPointsTree,fArrayIndex,fLastIndex);
1909 AliTrackPointArray **points;
146f76de 1910
5fd5f7a6 1911 Int_t pointsDim;
1912 LoadPoints(volids, points,pointsDim);
146f76de 1913
bad3394b 1914 Int_t nArrays = fPointsTree->GetEntries();
146f76de 1915
bad3394b 1916 if (nArrays == 0){ AliError("Points array is empty!"); return; }
1917 AliTrackFitter *fitter = CreateFitter();
146f76de 1918
bad3394b 1919 Int_t ecount=0;
1920 Int_t totcount=0;
146f76de 1921
bad3394b 1922 for (Int_t iArray = 0; iArray < nArrays; iArray++){
146f76de 1923
5fd5f7a6 1924 cout<<"Investigating "<<iArray<<"/"<<nArrays<<endl;
146f76de 1925
bad3394b 1926 if (!points[iArray]){
1927 cout<<" Skipping: "<<iArray<<endl;
1928 continue;
146f76de 1929 }
5fd5f7a6 1930
bad3394b 1931 fitter->SetTrackPointArray(points[iArray],kTRUE); // Watch out, problems
1932 // when few sectors
5fd5f7a6 1933
bad3394b 1934 totcount++;
1935
1936 // *** FITTING ***
1937 if(fitter->Fit(volids,volidsfit,layerRangeMin,layerRangeMax) == kFALSE){
1938 ecount++;
1939 cout<<"->BAD: "<<iArray<<endl;
1940 continue;
1941 } //else cout<<"->GOOD: "<<iArray<<endl;
146f76de 1942
bad3394b 1943 AliTrackPointArray *pVolId,*pTrack;
146f76de 1944
bad3394b 1945 fitter->GetTrackResiduals(pVolId,pTrack);
1946 FillResidualsH(pVolId,pTrack);
146f76de 1947
1948 }
bad3394b 1949
1950 cout<<" -> nVolIds: "<<nVolIds<<endl;
1951 cout<<" -> Non-Fitted tracks: "<<ecount<<"/"<<totcount<<endl;
1952
5fd5f7a6 1953 UnloadPoints(pointsDim,points);
146f76de 1954
bad3394b 1955 SaveHists(3,outname);
1956
146f76de 1957 return;
bad3394b 1958
146f76de 1959}
1960
1961
1962//______________________________________________________________________________
bad3394b 1963void AliITSResidualsAnalysis::ProcessVolumes(Int_t fit,
1964 TArrayI *volIDs,
1965 TArrayI *volIDsFit,
1966 TString misalignmentFile,
1967 TString outname,
1968 Int_t minPoints)
146f76de 1969{
1970 //
bad3394b 1971 // This function process the AliTrackPoints and volID (into residuals)
146f76de 1972 //
bad3394b 1973
1974 // setting up geometry and the trackpoints file
1975 if(!gGeoManager) AliGeomManager::LoadGeometry(GetFileNameGeometry());
1976
146f76de 1977 SetPointsFilename(GetFileNameTrackPoints());
146f76de 1978
bad3394b 1979 // creating some tools
1980 AliTrackFitter *fitter;
146f76de 1981 if(fit==1){
1982 fitter = new AliTrackFitterKalman();
1983 }else fitter = new AliTrackFitterRieman();
1984
bad3394b 1985 fitter->SetMinNPoints(minPoints);
146f76de 1986
bad3394b 1987 SetTrackFitter(fitter);
146f76de 1988
1989 if(misalignmentFile=="")printf("NO FAKE MISALIGNMENT\n");
1990 else {
1991 Bool_t misal=Misalign(misalignmentFile,"ITSAlignObjs");
bad3394b 1992 if(!misal){
1993 printf("PROBLEM WITH FAKE MISALIGNMENT!");
1994 return;
146f76de 1995 }
1996 }
146f76de 1997
bad3394b 1998 CalculateResiduals(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSSD2,outname);
146f76de 1999
146f76de 2000 return;
2001
146f76de 2002}