]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/AliITSResidualsAnalysis.cxx
In Open() and GotoEvent() try the ESD operations first, fallback to run-loader.
[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);
bad3394b 384 fVolResHistRPHI[nhist]=new TH1F("histname","histname",4000,-1.0,1.0);
146f76de 385 fVolResHistRPHI[nhist]->SetName(aux.Data());
386 fVolResHistRPHI[nhist]->SetTitle(aux.Data());
387
388 aux=histnameZ;
389 aux+=volIDs->At(nhist);
bad3394b 390 fResHistZ[nhist]=new TH1F("histname","histname",4000,-1.0,1.0);
146f76de 391 fResHistZ[nhist]->SetName(aux.Data());
392 fResHistZ[nhist]->SetTitle(aux.Data());
393
bad3394b 394 aux=histnameX;
395 aux+=volIDs->At(nhist);
396 fResHistX[nhist]=new TH1F("histname","histname",4000,-1.0,1.0);
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");
403 fResHistXLocsddL[nhist]=new TH1F("histname","histname",4000,-1.0,1.0);
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");
410 fResHistXLocsddR[nhist]=new TH1F("histname","histname",4000,-1.0,1.0);
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;
1005 return 0x0;
1006 }
1007 }
1008 fsingleLayer=kTRUE;
1009 return binningzphi;
1010}
1011
1012//____________________________________________________________________________
1013Bool_t AliITSResidualsAnalysis::SetBinning(const TArrayI *volids,Float_t *phiBin,Float_t *zBin)
1014{
1015 //
1016 // Sets the correct binning
1017 //
1018
1019 if(!fsingleLayer)return kFALSE;
bad3394b 1020 // const char *volpath,*symname;
146f76de 1021 Int_t iModule;
1022 Int_t *orderArrayPhi,*orderArrayZ;
1023 UShort_t volID;
bad3394b 1024 Double_t *phiArray,*zArray,*phiArrayOrdered,*zArrayOrdered;
1025 Double_t translGlobal[3];
146f76de 1026 Double_t lastPhimin=-10;
1027 Double_t lastZmin=-99;
1028 Int_t ***orderPhiZ;
bad3394b 1029 /* TGeoPNEntry *pne;
1030 TGeoPhysicalNode *pn;
1031 TGeoHMatrix *globMatrix;*/
146f76de 1032
1033 Bool_t used=kFALSE;
1034
1035 orderPhiZ=new Int_t**[fnPhi];
1036 phiArray=new Double_t[fnPhi];//phiBin[nModulesPhi+1];
1037 zArray=new Double_t[fnZ];//zBin[nModulesZ+1];
1038 phiArrayOrdered=new Double_t[fnPhi];
1039 zArrayOrdered=new Double_t[fnZ];
1040 orderArrayPhi=new Int_t[fnPhi];
1041 orderArrayZ=new Int_t[fnZ];
1042 for(Int_t k=0;k<fnZ;k++){
1043 orderArrayZ[k]=0;
1044 zArray[k]=-1000;
1045 }
1046 for(Int_t k=0;k<fnPhi;k++){
1047 orderArrayPhi[k]=0;
1048 phiArray[k]=-10;
1049 orderPhiZ[k]=new Int_t*[fnZ];
1050 }
1051
1052
1053 AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer((UShort_t)volids->At(0),iModule);
1054
1055 Int_t lastPhi=0,lastZ=0;
1056 for(iModule=0;iModule<AliGeomManager::LayerSize(iLayer);iModule++){
1057 fvolidsToBin[iModule]=new Int_t[3];
1058 volID=AliGeomManager::LayerToVolUID(iLayer,iModule);
1059 fvolidsToBin[iModule][0]=volID;
bad3394b 1060 /* symname = AliGeomManager::SymName(volID);
1061 pne = gGeoManager->GetAlignableEntry(symname);
1062 volpath=pne->GetTitle();
1063 pn=gGeoManager->MakePhysicalNode(volpath);
1064 globMatrix=pn->GetMatrix();
1065 translGlobal=globMatrix->GetTranslation();
1066
1067 */
1068 AliGeomManager::GetOrigTranslation(volID,translGlobal);
146f76de 1069
1070 for(Int_t j=0;j<lastPhi;j++){
1071 used=kFALSE;
bad3394b 1072 if(TMath::Abs(phiArray[j]-TMath::ATan2(translGlobal[1],translGlobal[0]))<2*TMath::Pi()/(10*fnPhi)){//10 is a safety factor but....
146f76de 1073 fvolidsToBin[iModule][1]=j;
1074 used=kTRUE;
1075 break;
1076 }
1077 }
1078 if(!used){
bad3394b 1079 phiArray[lastPhi]=TMath::ATan2(translGlobal[1],translGlobal[0]);
146f76de 1080 fvolidsToBin[iModule][1]=lastPhi;
1081 if(phiArray[lastPhi]<lastPhimin)lastPhimin=phiArray[lastPhi];
1082 lastPhi++;
1083 if(lastPhi>fnPhi){
1084 printf("Wrong Phi! \n");
1085 return kFALSE;}
1086 }
1087
1088 for(Int_t j=0;j<lastZ;j++){
1089 used=kFALSE;
bad3394b 1090 if(TMath::Abs(zArray[j]-translGlobal[2])<0.1){
146f76de 1091 fvolidsToBin[iModule][2]=j;
1092 used=kTRUE;
1093 break;
1094 }
1095 }
1096 if(!used){
1097 fvolidsToBin[iModule][2]=lastZ;
bad3394b 1098 zArray[lastZ]=translGlobal[2];
146f76de 1099 if(zArray[lastZ]<lastZmin)lastZmin=zArray[lastZ];
1100 lastZ++;
1101 if(lastZ>fnZ){
1102 printf("Wrong Z! \n");
1103 return kFALSE;}
1104 }
1105 }
1106
1107
1108 //ORDERING THE ARRAY OF PHI AND Z VALUES
1109 for(Int_t order=0;order<fnPhi;order++){
1110 for(Int_t j=0;j<fnPhi;j++){
1111 if((j!=order)&&(phiArray[order]>phiArray[j]))orderArrayPhi[order]++;
1112 }
1113 }
1114
1115 for(Int_t order=0;order<fnPhi;order++){
1116 for(Int_t j=0;j<fnPhi;j++){
1117 if(orderArrayPhi[j]==order){
1118 phiArrayOrdered[order]=phiArray[j];
1119 break;
1120 }
1121 }
1122 }
1123
1124
1125 for(Int_t order=0;order<fnZ;order++){
1126 for(Int_t j=0;j<fnZ;j++){
1127 if((j!=order)&&(zArray[order]>zArray[j]))orderArrayZ[order]++;
1128 }
1129 }
1130
1131
1132 for(Int_t order=0;order<fnZ;order++){
1133 for(Int_t j=0;j<fnZ;j++){
1134 if(orderArrayZ[j]==order){
1135 zArrayOrdered[order]=zArray[j];
1136 break;
1137 }
1138 }
1139 }
1140
1141
1142 //Filling the fCoordToBinTable
1143 for(Int_t j=0;j<fnPhi;j++){
1144 for(Int_t i=0;i<fnZ;i++){
1145 orderPhiZ[j][i]=new Int_t[2];
1146 orderPhiZ[j][i][0]=orderArrayPhi[j];
1147 orderPhiZ[j][i][1]=orderArrayZ[i];
1148 }
1149 }
1150
1151
1152 for(Int_t j=0;j<fnPhi;j++){
1153 for(Int_t i=0;i<fnZ;i++){
1154 fCoordToBinTable[j][i]=new Double_t[2];
1155 fCoordToBinTable[j][i][0]=phiArrayOrdered[j];
1156 fCoordToBinTable[j][i][1]=zArrayOrdered[i];
bad3394b 1157 printf("Now (binphi,binz)= %d, %d e (phi,z)=%f,%f \n",j,i,fCoordToBinTable[j][i][0],fCoordToBinTable[j][i][1]);
146f76de 1158 }
1159 }
1160 Int_t istar,jstar;
1161 for(iModule=0;iModule<fnPhi*fnZ;iModule++){
1162 istar=fvolidsToBin[iModule][1];
1163 jstar=fvolidsToBin[iModule][2];
1164 fvolidsToBin[iModule][1]=orderPhiZ[istar][jstar][0];
1165 fvolidsToBin[iModule][2]=orderPhiZ[istar][jstar][1];
1166 }
1167
1168
1169 //now constructing the binning
1170 for(Int_t iModPhi=0;iModPhi<fnPhi-1;iModPhi++){
1171 phiBin[iModPhi+1]=(Float_t)phiArrayOrdered[iModPhi]+0.5*(phiArrayOrdered[iModPhi+1]-phiArrayOrdered[iModPhi]);
1172 }
1173
1174 phiBin[0]=(Float_t)phiArrayOrdered[0]-(phiArrayOrdered[1]-phiArrayOrdered[0])/2.;
1175
1176 phiBin[fnPhi]=(Float_t)phiArrayOrdered[fnPhi-1]+(phiArrayOrdered[fnPhi-1]-phiArrayOrdered[fnPhi-2])/2.;
1177 for(Int_t iModPhi=0;iModPhi<fnPhi+1;iModPhi++){
1178 printf("Mean Phi mod %d su %d: %f \n",iModPhi+1,fnPhi,phiBin[iModPhi]);
1179 }
1180
1181 for(Int_t iModZ=0;iModZ<fnZ-1;iModZ++){
1182 zBin[iModZ+1]=(Float_t)zArrayOrdered[iModZ]+0.5*(zArrayOrdered[iModZ+1]-zArrayOrdered[iModZ]);
1183 }
1184 zBin[0]=(Float_t)zArrayOrdered[0]-0.5*(zArrayOrdered[1]-zArrayOrdered[0]);
1185 zBin[fnZ]=(Float_t)zArrayOrdered[fnZ-1]+0.5*(zArrayOrdered[1]-zArrayOrdered[0]);
1186
1187
1188 for(Int_t iModPhi=0;iModPhi<fnZ+1;iModPhi++){
1189 printf("Mean Z mod %d su %d: %f \n",iModPhi+1,fnZ,zBin[iModPhi]);
1190 }
1191 return kTRUE;
1192}
1193
1194
1195//____________________________________________________________________________
1196Int_t AliITSResidualsAnalysis::GetBinPhiZ(const Int_t volID,Int_t *binz) const
1197{
1198 //
1199 // Returns the correct Phi-Z bin
1200 //
1201
1202 if (!fsingleLayer){
1203 printf("No Single Layer reAlignment! \n");
1204 return 100;
1205 }
1206
1207 for(Int_t j=0;j<fnPhi*fnZ;j++){
1208 if(j==fnZ*fnPhi){
1209 printf("Wrong volume UID introduced! fnHist: %d volID: %d iter: %d \n",fnHist,volID,j);
1210 return 100;
1211 }
1212 if(fvolidsToBin[j][0]==volID){
1213
1214 *binz=fvolidsToBin[j][2];
1215 return fvolidsToBin[j][1];
1216 }
1217 }
1218
1219 return 100;
1220
1221}
1222
bad3394b 1223//___________________________________________________________________________
1224Int_t AliITSResidualsAnalysis::WhichSector(Int_t module) const
1225{
1226 //
1227 // This method returns the number of the SPD Sector
1228 // to which belongs the module (Sectors 0-9)
1229
1230 //--->cSect = 0 <---
1231 if(module==2048
1232 || module==2049
1233 || module==2050
1234 || module==2051
1235 || module==2052
1236 || module==2053
1237 || module==2054
1238 || module==2055
1239 || module==4096
1240 || module==4097
1241 || module==4098
1242 || module==4099
1243 || module==4100
1244 || module==4101
1245 || module==4102
1246 || module==4103
1247 || module==4104
1248 || module==4105
1249 || module==4106
1250 || module==4107
1251 || module==4108
1252 || module==4109
1253 || module==4110
1254 || module==4111) return 0;
1255
1256 //--->cSect = 1 <---
1257 if(module==2056
1258 || module==2057
1259 || module==2058
1260 || module==2059
1261 || module==2060
1262 || module==2061
1263 || module==2062
1264 || module==2063
1265 || module==4112
1266 || module==4113
1267 || module==4114
1268 || module==4115
1269 || module==4116
1270 || module==4117
1271 || module==4118
1272 || module==4119
1273 || module==4120
1274 || module==4121
1275 || module==4122
1276 || module==4123
1277 || module==4124
1278 || module==4125
1279 || module==4126
1280 || module==4127) return 1;
1281
1282 //--->cSect = 2 <---
1283 if(module==2064
1284 || module==2065
1285 || module==2066
1286 || module==2067
1287 || module==2068
1288 || module==2069
1289 || module==2070
1290 || module==2071
1291 || module==4128
1292 || module==4129
1293 || module==4130
1294 || module==4131
1295 || module==4132
1296 || module==4133
1297 || module==4134
1298 || module==4135
1299 || module==4136
1300 || module==4137
1301 || module==4138
1302 || module==4139
1303 || module==4140
1304 || module==4141
1305 || module==4142
1306 || module==4143) return 2;
1307
1308 //--->cSect = 3 <---
1309 if(module==2072
1310 || module==2073
1311 || module==2074
1312 || module==2075
1313 || module==2076
1314 || module==2077
1315 || module==2078
1316 || module==2079
1317 || module==4144
1318 || module==4145
1319 || module==4146
1320 || module==4147
1321 || module==4148
1322 || module==4149
1323 || module==4150
1324 || module==4151
1325 || module==4152
1326 || module==4153
1327 || module==4154
1328 || module==4155
1329 || module==4156
1330 || module==4157
1331 || module==4158
1332 || module==4159) return 3;
1333
1334 //--->cSect = 4 <---
1335 if(module==2080
1336 || module==2081
1337 || module==2082
1338 || module==2083
1339 || module==2084
1340 || module==2085
1341 || module==2086
1342 || module==2087
1343 || module==4160
1344 || module==4161
1345 || module==4162
1346 || module==4163
1347 || module==4164
1348 || module==4165
1349 || module==4166
1350 || module==4167
1351 || module==4168
1352 || module==4169
1353 || module==4170
1354 || module==4171
1355 || module==4172
1356 || module==4173
1357 || module==4174
1358 || module==4175) return 4;
1359
1360 //--->cSect = 5 <---
1361 if(module==2088
1362 || module==2089
1363 || module==2090
1364 || module==2091
1365 || module==2092
1366 || module==2093
1367 || module==2094
1368 || module==2095
1369 || module==4176
1370 || module==4177
1371 || module==4178
1372 || module==4179
1373 || module==4180
1374 || module==4181
1375 || module==4182
1376 || module==4183
1377 || module==4184
1378 || module==4185
1379 || module==4186
1380 || module==4187
1381 || module==4188
1382 || module==4189
1383 || module==4190
1384 || module==4191) return 5;
1385
1386 //--->cSect = 6 <---
1387 if(module==2096
1388 || module==2097
1389 || module==2098
1390 || module==2099
1391 || module==2100
1392 || module==2101
1393 || module==2102
1394 || module==2103
1395 || module==4192
1396 || module==4193
1397 || module==4194
1398 || module==4195
1399 || module==4196
1400 || module==4197
1401 || module==4198
1402 || module==4199
1403 || module==4200
1404 || module==4201
1405 || module==4202
1406 || module==4203
1407 || module==4204
1408 || module==4205
1409 || module==4206
1410 || module==4207) return 6;
1411
1412 //--->cSect = 7 <---
1413 if(module==2104
1414 || module==2105
1415 || module==2106
1416 || module==2107
1417 || module==2108
1418 || module==2109
1419 || module==2110
1420 || module==2111
1421 || module==4208
1422 || module==4209
1423 || module==4210
1424 || module==4211
1425 || module==4212
1426 || module==4213
1427 || module==4214
1428 || module==4215
1429 || module==4216
1430 || module==4217
1431 || module==4218
1432 || module==4219
1433 || module==4220
1434 || module==4221
1435 || module==4222
1436 || module==4223) return 7;
1437
1438 //--->cSect = 8 <---
1439 if(module==2112
1440 || module==2113
1441 || module==2114
1442 || module==2115
1443 || module==2116
1444 || module==2117
1445 || module==2118
1446 || module==2119
1447 || module==4224
1448 || module==4225
1449 || module==4226
1450 || module==4227
1451 || module==4228
1452 || module==4229
1453 || module==4230
1454 || module==4231
1455 || module==4232
1456 || module==4233
1457 || module==4234
1458 || module==4235
1459 || module==4236
1460 || module==4237
1461 || module==4238
1462 || module==4239) return 8;
1463
1464 //--->cSect = 9 <---
1465 if(module==2120
1466 || module==2121
1467 || module==2122
1468 || module==2123
1469 || module==2124
1470 || module==2125
1471 || module==2126
1472 || module==2127
1473 || module==4240
1474 || module==4241
1475 || module==4242
1476 || module==4243
1477 || module==4244
1478 || module==4245
1479 || module==4246
1480 || module==4247
1481 || module==4248
1482 || module==4249
1483 || module==4250
1484 || module==4251
1485 || module==4252
1486 || module==4253
1487 || module==4254
1488 || module==4255) return 9;
1489
1490 //printf("Module not belonging to SPD, sorry!");
1491 return -1;
1492
1493}
1494
146f76de 1495//____________________________________________________________________________
bad3394b 1496TArrayI* AliITSResidualsAnalysis::GetSPDSectorsVolids(Int_t sectors[10]) const
146f76de 1497{
1498 //
bad3394b 1499 // This method gets the volID Array for the chosen sectors.
1500 // You have to pass an array with a 1 for each selected sector.
1501 // i.e. sectors[10] = {1,1,0,0,0,0,0,0,1,0} -> Sector 0, 1, 9 selected.
146f76de 1502 //
1503
bad3394b 1504 Int_t nSect=0;
1505 Int_t iModule=0;
1506
146f76de 1507 if(!gGeoManager) AliGeomManager::LoadGeometry(GetFileNameGeometry());
1508
bad3394b 1509 for(Int_t co=0;co<10;co++){ //counts the number of sectors chosen
1510 if(sectors[co]==1) nSect++;
146f76de 1511 }
bad3394b 1512
1513 if(nSect<1){ //if no sector chosen -> exit
1514 Printf("Error! No Sector/s Selected!");
1515 return 0x0;
146f76de 1516 }
1517
bad3394b 1518 TArrayI *volIDs = new TArrayI(nSect*24);
1519
1520 if(sectors[0]==1){ //--->cSect = 0 <---
1521 volIDs->AddAt(2048,iModule); iModule++;
1522 volIDs->AddAt(2049,iModule); iModule++;
1523 volIDs->AddAt(2050,iModule); iModule++;
1524 volIDs->AddAt(2051,iModule); iModule++;
1525 volIDs->AddAt(2052,iModule); iModule++;
1526 volIDs->AddAt(2053,iModule); iModule++;
1527 volIDs->AddAt(2054,iModule); iModule++;
1528 volIDs->AddAt(2055,iModule); iModule++;
1529 volIDs->AddAt(4096,iModule); iModule++;
1530 volIDs->AddAt(4097,iModule); iModule++;
1531 volIDs->AddAt(4098,iModule); iModule++;
1532 volIDs->AddAt(4099,iModule); iModule++;
1533 volIDs->AddAt(4100,iModule); iModule++;
1534 volIDs->AddAt(4101,iModule); iModule++;
1535 volIDs->AddAt(4102,iModule); iModule++;
1536 volIDs->AddAt(4103,iModule); iModule++;
1537 volIDs->AddAt(4104,iModule); iModule++;
1538 volIDs->AddAt(4105,iModule); iModule++;
1539 volIDs->AddAt(4106,iModule); iModule++;
1540 volIDs->AddAt(4107,iModule); iModule++;
1541 volIDs->AddAt(4108,iModule); iModule++;
1542 volIDs->AddAt(4109,iModule); iModule++;
1543 volIDs->AddAt(4110,iModule); iModule++;
1544 volIDs->AddAt(4111,iModule); iModule++;
1545 }
1546 if(sectors[1]==1){ //--->cSect = 1 <//---
1547 volIDs->AddAt(2056,iModule); iModule++;
1548 volIDs->AddAt(2057,iModule); iModule++;
1549 volIDs->AddAt(2058,iModule); iModule++;
1550 volIDs->AddAt(2059,iModule); iModule++;
1551 volIDs->AddAt(2060,iModule); iModule++;
1552 volIDs->AddAt(2061,iModule); iModule++;
1553 volIDs->AddAt(2062,iModule); iModule++;
1554 volIDs->AddAt(2063,iModule); iModule++;
1555 volIDs->AddAt(4112,iModule); iModule++;
1556 volIDs->AddAt(4113,iModule); iModule++;
1557 volIDs->AddAt(4114,iModule); iModule++;
1558 volIDs->AddAt(4115,iModule); iModule++;
1559 volIDs->AddAt(4116,iModule); iModule++;
1560 volIDs->AddAt(4117,iModule); iModule++;
1561 volIDs->AddAt(4118,iModule); iModule++;
1562 volIDs->AddAt(4119,iModule); iModule++;
1563 volIDs->AddAt(4120,iModule); iModule++;
1564 volIDs->AddAt(4121,iModule); iModule++;
1565 volIDs->AddAt(4122,iModule); iModule++;
1566 volIDs->AddAt(4123,iModule); iModule++;
1567 volIDs->AddAt(4124,iModule); iModule++;
1568 volIDs->AddAt(4125,iModule); iModule++;
1569 volIDs->AddAt(4126,iModule); iModule++;
1570 volIDs->AddAt(4127,iModule); iModule++;
1571 }
1572 if(sectors[2]==1){//--->cSect = 2 <//---
1573 volIDs->AddAt(2064,iModule); iModule++;
1574 volIDs->AddAt(2065,iModule); iModule++;
1575 volIDs->AddAt(2066,iModule); iModule++;
1576 volIDs->AddAt(2067,iModule); iModule++;
1577 volIDs->AddAt(2068,iModule); iModule++;
1578 volIDs->AddAt(2069,iModule); iModule++;
1579 volIDs->AddAt(2070,iModule); iModule++;
1580 volIDs->AddAt(2071,iModule); iModule++;
1581 volIDs->AddAt(4128,iModule); iModule++;
1582 volIDs->AddAt(4129,iModule); iModule++;
1583 volIDs->AddAt(4130,iModule); iModule++;
1584 volIDs->AddAt(4131,iModule); iModule++;
1585 volIDs->AddAt(4132,iModule); iModule++;
1586 volIDs->AddAt(4133,iModule); iModule++;
1587 volIDs->AddAt(4134,iModule); iModule++;
1588 volIDs->AddAt(4135,iModule); iModule++;
1589 volIDs->AddAt(4136,iModule); iModule++;
1590 volIDs->AddAt(4137,iModule); iModule++;
1591 volIDs->AddAt(4138,iModule); iModule++;
1592 volIDs->AddAt(4139,iModule); iModule++;
1593 volIDs->AddAt(4140,iModule); iModule++;
1594 volIDs->AddAt(4141,iModule); iModule++;
1595 volIDs->AddAt(4142,iModule); iModule++;
1596 volIDs->AddAt(4143,iModule); iModule++;
1597 }
1598 if(sectors[3]==1){//--->cSect = 3 <//---
1599 volIDs->AddAt(2072,iModule); iModule++;
1600 volIDs->AddAt(2073,iModule); iModule++;
1601 volIDs->AddAt(2074,iModule); iModule++;
1602 volIDs->AddAt(2075,iModule); iModule++;
1603 volIDs->AddAt(2076,iModule); iModule++;
1604 volIDs->AddAt(2077,iModule); iModule++;
1605 volIDs->AddAt(2078,iModule); iModule++;
1606 volIDs->AddAt(2079,iModule); iModule++;
1607 volIDs->AddAt(4144,iModule); iModule++;
1608 volIDs->AddAt(4145,iModule); iModule++;
1609 volIDs->AddAt(4146,iModule); iModule++;
1610 volIDs->AddAt(4147,iModule); iModule++;
1611 volIDs->AddAt(4148,iModule); iModule++;
1612 volIDs->AddAt(4149,iModule); iModule++;
1613 volIDs->AddAt(4150,iModule); iModule++;
1614 volIDs->AddAt(4151,iModule); iModule++;
1615 volIDs->AddAt(4152,iModule); iModule++;
1616 volIDs->AddAt(4153,iModule); iModule++;
1617 volIDs->AddAt(4154,iModule); iModule++;
1618 volIDs->AddAt(4155,iModule); iModule++;
1619 volIDs->AddAt(4156,iModule); iModule++;
1620 volIDs->AddAt(4157,iModule); iModule++;
1621 volIDs->AddAt(4158,iModule); iModule++;
1622 volIDs->AddAt(4159,iModule); iModule++;
1623 }
1624 if(sectors[4]==1){//--->cSect = 4 <//---
1625 volIDs->AddAt(2080,iModule); iModule++;
1626 volIDs->AddAt(2081,iModule); iModule++;
1627 volIDs->AddAt(2082,iModule); iModule++;
1628 volIDs->AddAt(2083,iModule); iModule++;
1629 volIDs->AddAt(2084,iModule); iModule++;
1630 volIDs->AddAt(2085,iModule); iModule++;
1631 volIDs->AddAt(2086,iModule); iModule++;
1632 volIDs->AddAt(2087,iModule); iModule++;
1633 volIDs->AddAt(4160,iModule); iModule++;
1634 volIDs->AddAt(4161,iModule); iModule++;
1635 volIDs->AddAt(4162,iModule); iModule++;
1636 volIDs->AddAt(4163,iModule); iModule++;
1637 volIDs->AddAt(4164,iModule); iModule++;
1638 volIDs->AddAt(4165,iModule); iModule++;
1639 volIDs->AddAt(4166,iModule); iModule++;
1640 volIDs->AddAt(4167,iModule); iModule++;
1641 volIDs->AddAt(4168,iModule); iModule++;
1642 volIDs->AddAt(4169,iModule); iModule++;
1643 volIDs->AddAt(4170,iModule); iModule++;
1644 volIDs->AddAt(4171,iModule); iModule++;
1645 volIDs->AddAt(4172,iModule); iModule++;
1646 volIDs->AddAt(4173,iModule); iModule++;
1647 volIDs->AddAt(4174,iModule); iModule++;
1648 volIDs->AddAt(4175,iModule); iModule++;
1649 }
1650 if(sectors[5]==1){//--->cSect = 5 <//---
1651 volIDs->AddAt(2088,iModule); iModule++;
1652 volIDs->AddAt(2089,iModule); iModule++;
1653 volIDs->AddAt(2090,iModule); iModule++;
1654 volIDs->AddAt(2091,iModule); iModule++;
1655 volIDs->AddAt(2092,iModule); iModule++;
1656 volIDs->AddAt(2093,iModule); iModule++;
1657 volIDs->AddAt(2094,iModule); iModule++;
1658 volIDs->AddAt(2095,iModule); iModule++;
1659 volIDs->AddAt(4176,iModule); iModule++;
1660 volIDs->AddAt(4177,iModule); iModule++;
1661 volIDs->AddAt(4178,iModule); iModule++;
1662 volIDs->AddAt(4179,iModule); iModule++;
1663 volIDs->AddAt(4180,iModule); iModule++;
1664 volIDs->AddAt(4181,iModule); iModule++;
1665 volIDs->AddAt(4182,iModule); iModule++;
1666 volIDs->AddAt(4183,iModule); iModule++;
1667 volIDs->AddAt(4184,iModule); iModule++;
1668 volIDs->AddAt(4185,iModule); iModule++;
1669 volIDs->AddAt(4186,iModule); iModule++;
1670 volIDs->AddAt(4187,iModule); iModule++;
1671 volIDs->AddAt(4188,iModule); iModule++;
1672 volIDs->AddAt(4189,iModule); iModule++;
1673 volIDs->AddAt(4190,iModule); iModule++;
1674 volIDs->AddAt(4191,iModule); iModule++;
1675 }
1676 if(sectors[6]==1){//--->cSect = 6 <//---
1677 volIDs->AddAt(2096,iModule); iModule++;
1678 volIDs->AddAt(2097,iModule); iModule++;
1679 volIDs->AddAt(2098,iModule); iModule++;
1680 volIDs->AddAt(2099,iModule); iModule++;
1681 volIDs->AddAt(2100,iModule); iModule++;
1682 volIDs->AddAt(2101,iModule); iModule++;
1683 volIDs->AddAt(2102,iModule); iModule++;
1684 volIDs->AddAt(2103,iModule); iModule++;
1685 volIDs->AddAt(4192,iModule); iModule++;
1686 volIDs->AddAt(4193,iModule); iModule++;
1687 volIDs->AddAt(4194,iModule); iModule++;
1688 volIDs->AddAt(4195,iModule); iModule++;
1689 volIDs->AddAt(4196,iModule); iModule++;
1690 volIDs->AddAt(4197,iModule); iModule++;
1691 volIDs->AddAt(4198,iModule); iModule++;
1692 volIDs->AddAt(4199,iModule); iModule++;
1693 volIDs->AddAt(4200,iModule); iModule++;
1694 volIDs->AddAt(4201,iModule); iModule++;
1695 volIDs->AddAt(4202,iModule); iModule++;
1696 volIDs->AddAt(4203,iModule); iModule++;
1697 volIDs->AddAt(4204,iModule); iModule++;
1698 volIDs->AddAt(4205,iModule); iModule++;
1699 volIDs->AddAt(4206,iModule); iModule++;
1700 volIDs->AddAt(4207,iModule); iModule++;
1701 }
1702 if(sectors[7]==1){ //--->cSect = 7 <//---
1703 volIDs->AddAt(2104,iModule); iModule++;
1704 volIDs->AddAt(2105,iModule); iModule++;
1705 volIDs->AddAt(2106,iModule); iModule++;
1706 volIDs->AddAt(2107,iModule); iModule++;
1707 volIDs->AddAt(2108,iModule); iModule++;
1708 volIDs->AddAt(2109,iModule); iModule++;
1709 volIDs->AddAt(2110,iModule); iModule++;
1710 volIDs->AddAt(2111,iModule); iModule++;
1711 volIDs->AddAt(4208,iModule); iModule++;
1712 volIDs->AddAt(4209,iModule); iModule++;
1713 volIDs->AddAt(4210,iModule); iModule++;
1714 volIDs->AddAt(4211,iModule); iModule++;
1715 volIDs->AddAt(4212,iModule); iModule++;
1716 volIDs->AddAt(4213,iModule); iModule++;
1717 volIDs->AddAt(4214,iModule); iModule++;
1718 volIDs->AddAt(4215,iModule); iModule++;
1719 volIDs->AddAt(4216,iModule); iModule++;
1720 volIDs->AddAt(4217,iModule); iModule++;
1721 volIDs->AddAt(4218,iModule); iModule++;
1722 volIDs->AddAt(4219,iModule); iModule++;
1723 volIDs->AddAt(4220,iModule); iModule++;
1724 volIDs->AddAt(4221,iModule); iModule++;
1725 volIDs->AddAt(4222,iModule); iModule++;
1726 volIDs->AddAt(4223,iModule); iModule++;
1727 }
1728 if(sectors[8]==1){//--->cSect = 8 <//---
1729 volIDs->AddAt(2112,iModule); iModule++;
1730 volIDs->AddAt(2113,iModule); iModule++;
1731 volIDs->AddAt(2114,iModule); iModule++;
1732 volIDs->AddAt(2115,iModule); iModule++;
1733 volIDs->AddAt(2116,iModule); iModule++;
1734 volIDs->AddAt(2117,iModule); iModule++;
1735 volIDs->AddAt(2118,iModule); iModule++;
1736 volIDs->AddAt(2119,iModule); iModule++;
1737 volIDs->AddAt(4224,iModule); iModule++;
1738 volIDs->AddAt(4225,iModule); iModule++;
1739 volIDs->AddAt(4226,iModule); iModule++;
1740 volIDs->AddAt(4227,iModule); iModule++;
1741 volIDs->AddAt(4228,iModule); iModule++;
1742 volIDs->AddAt(4229,iModule); iModule++;
1743 volIDs->AddAt(4230,iModule); iModule++;
1744 volIDs->AddAt(4231,iModule); iModule++;
1745 volIDs->AddAt(4232,iModule); iModule++;
1746 volIDs->AddAt(4233,iModule); iModule++;
1747 volIDs->AddAt(4234,iModule); iModule++;
1748 volIDs->AddAt(4235,iModule); iModule++;
1749 volIDs->AddAt(4236,iModule); iModule++;
1750 volIDs->AddAt(4237,iModule); iModule++;
1751 volIDs->AddAt(4238,iModule); iModule++;
1752 volIDs->AddAt(4239,iModule); iModule++;
1753 }
1754 if(sectors[9]==1){//--->cSect = 9 <//---
1755 volIDs->AddAt(2120,iModule); iModule++;
1756 volIDs->AddAt(2121,iModule); iModule++;
1757 volIDs->AddAt(2122,iModule); iModule++;
1758 volIDs->AddAt(2123,iModule); iModule++;
1759 volIDs->AddAt(2124,iModule); iModule++;
1760 volIDs->AddAt(2125,iModule); iModule++;
1761 volIDs->AddAt(2126,iModule); iModule++;
1762 volIDs->AddAt(2127,iModule); iModule++;
1763 volIDs->AddAt(4240,iModule); iModule++;
1764 volIDs->AddAt(4241,iModule); iModule++;
1765 volIDs->AddAt(4242,iModule); iModule++;
1766 volIDs->AddAt(4243,iModule); iModule++;
1767 volIDs->AddAt(4244,iModule); iModule++;
1768 volIDs->AddAt(4245,iModule); iModule++;
1769 volIDs->AddAt(4246,iModule); iModule++;
1770 volIDs->AddAt(4247,iModule); iModule++;
1771 volIDs->AddAt(4248,iModule); iModule++;
1772 volIDs->AddAt(4249,iModule); iModule++;
1773 volIDs->AddAt(4250,iModule); iModule++;
1774 volIDs->AddAt(4251,iModule); iModule++;
1775 volIDs->AddAt(4252,iModule); iModule++;
1776 volIDs->AddAt(4253,iModule); iModule++;
1777 volIDs->AddAt(4254,iModule); iModule++;
1778 volIDs->AddAt(4255,iModule); iModule++;
1779 }
1780
146f76de 1781 return volIDs;
bad3394b 1782
1783}
1784
1785//____________________________________________________________________________
1786TArrayI* AliITSResidualsAnalysis::GetITSLayersVolids(Int_t layers[6]) const
1787{
1788 //
1789 // This method gets the volID Array for the chosen layers.
1790 // You have to pass an array with a 1 for each selected layer.
1791 // i.e. layers[6] = {1,1,0,0,1,1} -> SPD + SSD
1792 //
1793
1794 if(!gGeoManager) AliGeomManager::LoadGeometry(GetFileNameGeometry());
1795
1796 Int_t size=0,last=0;
1797
1798 // evaluates the size of the array
1799 for(Int_t i=0;i<6;i++) if(layers[i]==1) size+=AliGeomManager::LayerSize(i+1);
1800
1801 if(size==0){
1802 printf("Error: no layer selected");
1803 return 0x0;
1804 }
1805
1806 TArrayI *volids = new TArrayI(size);
1807
1808 // fills the volId array only for the chosen layers
1809 for(Int_t ilayer=1;ilayer<7;ilayer++){
1810
1811 if(layers[ilayer-1]!=1) continue;
1812
1813 for(Int_t imod=0;imod<AliGeomManager::LayerSize(ilayer);imod++){
1814 volids->AddAt(AliGeomManager::LayerToVolUID(ilayer,imod),last);
1815 last++;
1816 }
1817 }
146f76de 1818
bad3394b 1819 return volids;
1820
146f76de 1821}
1822
1823//____________________________________________________________________________
1824void AliITSResidualsAnalysis::GetTrackDirClusterCov(AliTrackPoint *point,Double_t &phi,Double_t &lambda,Double_t &lambda2,Double_t &alpha,Double_t &xovery,Double_t &zovery) const
1825{
1826 //
1827 // ...
1828 //
1829
1830 TMatrixDSym cov(3);
1831 const Float_t *covvector=point->GetCov();
1832 cov(0,0)=covvector[0];
1833 cov(1,0)=cov(0,1)=covvector[1];
1834 cov(2,0)=cov(0,2)=covvector[2];
1835 cov(1,1)=covvector[3];
1836 cov(1,2)=cov(2,1)=covvector[4];
1837 cov(2,2)=covvector[5];
1838
1839 Double_t determinant=cov.Determinant();
1840 if(determinant!=0.){
1841 TMatrixD vect(3,3);
1842 TVectorD eigenvalues(3);
1843 const TMatrixDSymEigen keigen(cov);
1844 eigenvalues=keigen.GetEigenValues();
1845 vect=keigen.GetEigenVectors();
1846 Double_t mainvect[3];
1847 mainvect[0]=vect(0,0);
1848 mainvect[1]=vect(1,0);
1849 mainvect[2]=vect(2,0);
1850 if(mainvect[1]!=0.){
1851 xovery=mainvect[0]/mainvect[1];
1852 zovery=mainvect[2]/mainvect[1];
1853 }
1854 else {
1855 xovery=9999.;
1856 zovery=9999.;
1857 }
1858 if(mainvect[1]<0.){
1859 mainvect[0]=-1.*mainvect[0];
1860 mainvect[1]=-1.*mainvect[1];
1861 mainvect[2]=-1.*mainvect[2];
1862 }
1863 lambda2=TMath::ATan2(TMath::Sqrt(mainvect[0]*mainvect[0]+mainvect[2]*mainvect[2]),mainvect[1])*TMath::RadToDeg();
1864 lambda=TMath::ATan2(mainvect[2],TMath::Sqrt(mainvect[0]*mainvect[0]+mainvect[1]*mainvect[1]))*TMath::RadToDeg();
1865 phi=TMath::ATan2(mainvect[0],mainvect[2])*TMath::RadToDeg();
1866 alpha=TMath::ATan2(mainvect[1],mainvect[0])*TMath::RadToDeg();
1867 }
1868 else printf("determinant =0!, skip this point \n");
1869
1870 return;
1871}
1872
1873//____________________________________________________________________________
1874void AliITSResidualsAnalysis::CalculateResiduals(const TArrayI *volids,
bad3394b 1875 const TArrayI *volidsfit,
1876 AliGeomManager::ELayerID layerRangeMin,
1877 AliGeomManager::ELayerID layerRangeMax,
1878 TString outname)
146f76de 1879{
1880 // CalculateResiduals for a set of detector volumes.
1881 // Tracks are fitted only within
1882 // the range defined by the user
1883 // (by layerRangeMin and layerRangeMax)
1884 // or within the set of volidsfit
1885 // Repeat the procedure 'iterations' times
1886
146f76de 1887 Int_t nVolIds = volids->GetSize();
bad3394b 1888 if (nVolIds == 0) { AliError("Volume IDs array is empty!"); return; }
146f76de 1889
1890 // Load only the tracks with at least one
1891 // space point in the set of volume (volids)
1892
bad3394b 1893
146f76de 1894 //AliAlignmentTracks::SetPointsFilename(GetFileNameTrackPoints());
1895 AliAlignmentTracks::BuildIndex();
1896
bad3394b 1897 ListVolUsed(fPointsTree,fArrayIndex,fLastIndex);
1898 AliTrackPointArray **points;
146f76de 1899
bad3394b 1900 LoadPoints(volids, points);
146f76de 1901
bad3394b 1902 Int_t nArrays = fPointsTree->GetEntries();
146f76de 1903
bad3394b 1904 if (nArrays == 0){ AliError("Points array is empty!"); return; }
1905 AliTrackFitter *fitter = CreateFitter();
146f76de 1906
bad3394b 1907 Int_t ecount=0;
1908 Int_t totcount=0;
1909
1910 // nArrays=806; // WAAAAAAAAAARNING!
146f76de 1911
bad3394b 1912 Int_t last=0;
146f76de 1913
bad3394b 1914 for (Int_t iArray = 0; iArray < nArrays; iArray++){
146f76de 1915
bad3394b 1916 //cout<<"Investigating "<<iArray<<"/"<<nArrays<<endl;
146f76de 1917
bad3394b 1918 if (!points[iArray]){
1919 cout<<" Skipping: "<<iArray<<endl;
1920 continue;
146f76de 1921 }
146f76de 1922
bad3394b 1923 last++;
1924
1925 fitter->SetTrackPointArray(points[iArray],kTRUE); // Watch out, problems
1926 // when few sectors
1927
1928 totcount++;
1929
1930 // *** FITTING ***
1931 if(fitter->Fit(volids,volidsfit,layerRangeMin,layerRangeMax) == kFALSE){
1932 ecount++;
1933 cout<<"->BAD: "<<iArray<<endl;
1934 continue;
1935 } //else cout<<"->GOOD: "<<iArray<<endl;
146f76de 1936
bad3394b 1937 AliTrackPointArray *pVolId,*pTrack;
146f76de 1938
bad3394b 1939 fitter->GetTrackResiduals(pVolId,pTrack);
1940 FillResidualsH(pVolId,pTrack);
146f76de 1941
1942 }
bad3394b 1943
1944 cout<<" -> nVolIds: "<<nVolIds<<endl;
1945 cout<<" -> Non-Fitted tracks: "<<ecount<<"/"<<totcount<<endl;
1946
1947 UnloadPoints(last, points);
146f76de 1948
bad3394b 1949 SaveHists(3,outname);
1950
146f76de 1951 return;
bad3394b 1952
146f76de 1953}
1954
1955
1956//______________________________________________________________________________
bad3394b 1957void AliITSResidualsAnalysis::ProcessVolumes(Int_t fit,
1958 TArrayI *volIDs,
1959 TArrayI *volIDsFit,
1960 TString misalignmentFile,
1961 TString outname,
1962 Int_t minPoints)
146f76de 1963{
1964 //
bad3394b 1965 // This function process the AliTrackPoints and volID (into residuals)
146f76de 1966 //
bad3394b 1967
1968 // setting up geometry and the trackpoints file
1969 if(!gGeoManager) AliGeomManager::LoadGeometry(GetFileNameGeometry());
1970
146f76de 1971 SetPointsFilename(GetFileNameTrackPoints());
146f76de 1972
bad3394b 1973 // creating some tools
1974 AliTrackFitter *fitter;
146f76de 1975 if(fit==1){
1976 fitter = new AliTrackFitterKalman();
1977 }else fitter = new AliTrackFitterRieman();
1978
bad3394b 1979 fitter->SetMinNPoints(minPoints);
146f76de 1980
bad3394b 1981 SetTrackFitter(fitter);
146f76de 1982
1983 if(misalignmentFile=="")printf("NO FAKE MISALIGNMENT\n");
1984 else {
1985 Bool_t misal=Misalign(misalignmentFile,"ITSAlignObjs");
bad3394b 1986 if(!misal){
1987 printf("PROBLEM WITH FAKE MISALIGNMENT!");
1988 return;
146f76de 1989 }
1990 }
146f76de 1991
bad3394b 1992 CalculateResiduals(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSSD2,outname);
146f76de 1993
146f76de 1994 return;
1995
146f76de 1996}