]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TRD/AliTRDCalibraFillHisto.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
- include header suppression in raw data simulation (as used in real
[u/mrichter/AliRoot.git] / TRD / AliTRDCalibraFillHisto.cxx
CommitLineData
64942b85 1//**************************************************************************
2// * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3// * *
4// * Author: The ALICE Off-line Project. *
5// * Contributors are mentioned in the code where appropriate. *
6// * *
7// * Permission to use, copy, modify and distribute this software and its *
8// * documentation strictly for non-commercial purposes is hereby granted *
9// * without fee, provided that the above copyright notice appears in all *
10// * copies and that both the copyright notice and this permission notice *
11//* appear in the supporting documentation. The authors make no claims *
12//* about the suitability of this software for any purpose. It is *
13//* provided "as is" without express or implied warranty. *
14//**************************************************************************/
55a288e5 15
16/* $Id$ */
17
18/////////////////////////////////////////////////////////////////////////////////
170c35f1 19//
55a288e5 20// AliTRDCalibraFillHisto
21//
22// This class is for the TRD calibration of the relative gain factor, the drift velocity,
23// the time 0 and the pad response function. It fills histos or vectors.
24// It can be used for the calibration per chamber but also per group of pads and eventually per pad.
170c35f1 25// The user has to choose with the functions SetNz and SetNrphi the precision of the calibration (see AliTRDCalibraMode).
55a288e5 26// 2D Histograms (Histo2d) or vectors (Vector2d), then converted in Trees, will be filled
27// from RAW DATA in a run or from reconstructed TRD tracks during the offline tracking
28// in the function "FollowBackProlongation" (AliTRDtracker)
29// Per default the functions to fill are off.
30//
1785640c 31// Authors:
32// R. Bailhache (R.Bailhache@gsi.de, rbailhache@ikf.uni-frankfurt.de)
33// J. Book (jbook@ikf.uni-frankfurt.de)
55a288e5 34//
35//////////////////////////////////////////////////////////////////////////////////////
36
55a288e5 37#include <TProfile2D.h>
38#include <TProfile.h>
39#include <TFile.h>
55a288e5 40#include <TStyle.h>
41#include <TCanvas.h>
55a288e5 42#include <TObjArray.h>
d0569428 43#include <TObject.h>
55a288e5 44#include <TH1F.h>
45#include <TH2I.h>
46#include <TH2.h>
47#include <TStopwatch.h>
48#include <TMath.h>
49#include <TDirectory.h>
170c35f1 50#include <TTreeStream.h>
51#include <TVectorD.h>
6bbdc11a 52#include <TLinearFitter.h>
55a288e5 53
54#include "AliLog.h"
55a288e5 55
56#include "AliTRDCalibraFillHisto.h"
57#include "AliTRDCalibraMode.h"
58#include "AliTRDCalibraVector.h"
3a0f6479 59#include "AliTRDCalibraVdriftLinearFit.h"
55a288e5 60#include "AliTRDcalibDB.h"
61#include "AliTRDCommonParam.h"
55a288e5 62#include "AliTRDpadPlane.h"
63#include "AliTRDcluster.h"
f3d3af1b 64#include "AliTRDtrack.h"
bcb6fb78 65#include "AliTRDtrackV1.h"
b70c68da 66#include "AliTRDrawStreamBase.h"
3a0f6479 67#include "AliRawReader.h"
68#include "AliRawReaderDate.h"
f162af62 69#include "AliTRDgeometry.h"
d0569428 70#include "./Cal/AliTRDCalROC.h"
71#include "./Cal/AliTRDCalDet.h"
170c35f1 72
1785640c 73#include "AliTRDrawFastStream.h"
74#include "AliTRDdigitsManager.h"
75#include "AliTRDdigitsParam.h"
76#include "AliTRDSignalIndex.h"
77#include "AliTRDarrayADC.h"
78
170c35f1 79#ifdef ALI_DATE
80#include "event.h"
81#endif
55a288e5 82
3a0f6479 83
55a288e5 84ClassImp(AliTRDCalibraFillHisto)
85
86AliTRDCalibraFillHisto* AliTRDCalibraFillHisto::fgInstance = 0;
87Bool_t AliTRDCalibraFillHisto::fgTerminated = kFALSE;
88
89//_____________singleton implementation_________________________________________________
90AliTRDCalibraFillHisto *AliTRDCalibraFillHisto::Instance()
91{
92 //
93 // Singleton implementation
94 //
95
96 if (fgTerminated != kFALSE) {
97 return 0;
98 }
99
100 if (fgInstance == 0) {
101 fgInstance = new AliTRDCalibraFillHisto();
102 }
103
104 return fgInstance;
105
106}
107
108//______________________________________________________________________________________
109void AliTRDCalibraFillHisto::Terminate()
110{
111 //
112 // Singleton implementation
113 // Deletes the instance of this class
114 //
115
116 fgTerminated = kTRUE;
117
118 if (fgInstance != 0) {
119 delete fgInstance;
120 fgInstance = 0;
121 }
122
123}
3a0f6479 124
55a288e5 125//______________________________________________________________________________________
126AliTRDCalibraFillHisto::AliTRDCalibraFillHisto()
127 :TObject()
f162af62 128 ,fGeo(0)
01239968 129 ,fCalibDB(0)
64942b85 130 ,fIsHLT(kFALSE)
55a288e5 131 ,fCH2dOn(kFALSE)
132 ,fPH2dOn(kFALSE)
133 ,fPRF2dOn(kFALSE)
134 ,fHisto2d(kFALSE)
135 ,fVector2d(kFALSE)
170c35f1 136 ,fLinearFitterOn(kFALSE)
137 ,fLinearFitterDebugOn(kFALSE)
55a288e5 138 ,fRelativeScale(0)
170c35f1 139 ,fThresholdClusterPRF2(15.0)
b70c68da 140 ,fLimitChargeIntegration(kFALSE)
f9ed5395 141 ,fFillWithZero(kFALSE)
64942b85 142 ,fNormalizeNbOfCluster(kFALSE)
143 ,fMaxCluster(0)
144 ,fNbMaxCluster(0)
170c35f1 145 ,fCalibraMode(new AliTRDCalibraMode())
146 ,fDebugStreamer(0)
147 ,fDebugLevel(0)
55a288e5 148 ,fDetectorPreviousTrack(-1)
e4db522f 149 ,fMCMPrevious(-1)
150 ,fROBPrevious(-1)
64942b85 151 ,fNumberClusters(1)
bcb6fb78 152 ,fNumberClustersf(30)
6aafa7ea 153 ,fNumberClustersProcent(0.5)
154 ,fThresholdClustersDAQ(120.0)
155 ,fNumberRowDAQ(2)
156 ,fNumberColDAQ(4)
170c35f1 157 ,fProcent(6.0)
158 ,fDifference(17)
159 ,fNumberTrack(0)
160 ,fTimeMax(0)
161 ,fSf(10.0)
e526983e 162 ,fNumberBinCharge(50)
163 ,fNumberBinPRF(10)
164 ,fNgroupprf(3)
55a288e5 165 ,fAmpTotal(0x0)
166 ,fPHPlace(0x0)
167 ,fPHValue(0x0)
170c35f1 168 ,fGoodTracklet(kTRUE)
1ca79a00 169 ,fLinearFitterTracklet(0x0)
170c35f1 170 ,fEntriesCH(0x0)
171 ,fEntriesLinearFitter(0x0)
172 ,fCalibraVector(0x0)
55a288e5 173 ,fPH2d(0x0)
174 ,fPRF2d(0x0)
175 ,fCH2d(0x0)
d0569428 176 ,fLinearFitterArray(540)
3a0f6479 177 ,fLinearVdriftFit(0x0)
d0569428 178 ,fCalDetGain(0x0)
179 ,fCalROCGain(0x0)
55a288e5 180{
181 //
182 // Default constructor
183 //
184
170c35f1 185 //
186 // Init some default values
187 //
55a288e5 188
170c35f1 189 fNumberUsedCh[0] = 0;
190 fNumberUsedCh[1] = 0;
191 fNumberUsedPh[0] = 0;
192 fNumberUsedPh[1] = 0;
3a0f6479 193
f162af62 194 fGeo = new AliTRDgeometry();
01239968 195 fCalibDB = AliTRDcalibDB::Instance();
55a288e5 196}
f162af62 197
55a288e5 198//______________________________________________________________________________________
199AliTRDCalibraFillHisto::AliTRDCalibraFillHisto(const AliTRDCalibraFillHisto &c)
f162af62 200 :TObject(c)
201 ,fGeo(0)
01239968 202 ,fCalibDB(0)
64942b85 203 ,fIsHLT(c.fIsHLT)
f162af62 204 ,fCH2dOn(c.fCH2dOn)
205 ,fPH2dOn(c.fPH2dOn)
206 ,fPRF2dOn(c.fPRF2dOn)
207 ,fHisto2d(c.fHisto2d)
208 ,fVector2d(c.fVector2d)
209 ,fLinearFitterOn(c.fLinearFitterOn)
210 ,fLinearFitterDebugOn(c.fLinearFitterDebugOn)
211 ,fRelativeScale(c.fRelativeScale)
212 ,fThresholdClusterPRF2(c.fThresholdClusterPRF2)
b70c68da 213 ,fLimitChargeIntegration(c.fLimitChargeIntegration)
37b0cf5e 214 ,fFillWithZero(c.fFillWithZero)
64942b85 215 ,fNormalizeNbOfCluster(c.fNormalizeNbOfCluster)
216 ,fMaxCluster(c.fMaxCluster)
217 ,fNbMaxCluster(c.fNbMaxCluster)
f162af62 218 ,fCalibraMode(0x0)
219 ,fDebugStreamer(0)
220 ,fDebugLevel(c.fDebugLevel)
f162af62 221 ,fDetectorPreviousTrack(c.fDetectorPreviousTrack)
e4db522f 222 ,fMCMPrevious(c.fMCMPrevious)
223 ,fROBPrevious(c.fROBPrevious)
f162af62 224 ,fNumberClusters(c.fNumberClusters)
bcb6fb78 225 ,fNumberClustersf(c.fNumberClustersf)
6aafa7ea 226 ,fNumberClustersProcent(c.fNumberClustersProcent)
227 ,fThresholdClustersDAQ(c.fThresholdClustersDAQ)
228 ,fNumberRowDAQ(c.fNumberRowDAQ)
1785640c 229 ,fNumberColDAQ(c.fNumberColDAQ)
f162af62 230 ,fProcent(c.fProcent)
231 ,fDifference(c.fDifference)
232 ,fNumberTrack(c.fNumberTrack)
233 ,fTimeMax(c.fTimeMax)
234 ,fSf(c.fSf)
235 ,fNumberBinCharge(c.fNumberBinCharge)
236 ,fNumberBinPRF(c.fNumberBinPRF)
237 ,fNgroupprf(c.fNgroupprf)
3a0f6479 238 ,fAmpTotal(0x0)
239 ,fPHPlace(0x0)
240 ,fPHValue(0x0)
f162af62 241 ,fGoodTracklet(c.fGoodTracklet)
1ca79a00 242 ,fLinearFitterTracklet(0x0)
3a0f6479 243 ,fEntriesCH(0x0)
244 ,fEntriesLinearFitter(0x0)
f162af62 245 ,fCalibraVector(0x0)
246 ,fPH2d(0x0)
247 ,fPRF2d(0x0)
248 ,fCH2d(0x0)
d0569428 249 ,fLinearFitterArray(540)
3a0f6479 250 ,fLinearVdriftFit(0x0)
d0569428 251 ,fCalDetGain(0x0)
252 ,fCalROCGain(0x0)
55a288e5 253{
254 //
255 // Copy constructor
256 //
170c35f1 257 if(c.fCalibraMode) fCalibraMode = new AliTRDCalibraMode(*c.fCalibraMode);
258 if(c.fCalibraVector) fCalibraVector = new AliTRDCalibraVector(*c.fCalibraVector);
259 if(c.fPH2d) {
260 fPH2d = new TProfile2D(*c.fPH2d);
261 fPH2d->SetDirectory(0);
262 }
263 if(c.fPRF2d) {
264 fPRF2d = new TProfile2D(*c.fPRF2d);
265 fPRF2d->SetDirectory(0);
266 }
267 if(c.fCH2d) {
268 fCH2d = new TH2I(*c.fCH2d);
269 fCH2d->SetDirectory(0);
270 }
3a0f6479 271 if(c.fLinearVdriftFit){
272 fLinearVdriftFit = new AliTRDCalibraVdriftLinearFit(*c.fLinearVdriftFit);
170c35f1 273 }
3a0f6479 274
d0569428 275 if(c.fCalDetGain) fCalDetGain = new AliTRDCalDet(*c.fCalDetGain);
d0569428 276 if(c.fCalROCGain) fCalROCGain = new AliTRDCalROC(*c.fCalROCGain);
d0569428 277
f162af62 278 if (fGeo) {
279 delete fGeo;
280 }
281 fGeo = new AliTRDgeometry();
01239968 282 fCalibDB = AliTRDcalibDB::Instance();
55a288e5 283}
f162af62 284
55a288e5 285//____________________________________________________________________________________
286AliTRDCalibraFillHisto::~AliTRDCalibraFillHisto()
287{
288 //
289 // AliTRDCalibraFillHisto destructor
290 //
291
292 ClearHistos();
170c35f1 293 if ( fDebugStreamer ) delete fDebugStreamer;
f162af62 294
d0569428 295 if ( fCalDetGain ) delete fCalDetGain;
d0569428 296 if ( fCalROCGain ) delete fCalROCGain;
d0569428 297
1ca79a00 298 if( fLinearFitterTracklet ) { delete fLinearFitterTracklet; }
299
300 delete [] fPHPlace;
301 delete [] fPHValue;
302 delete [] fEntriesCH;
303 delete [] fEntriesLinearFitter;
304 delete [] fAmpTotal;
305
306 for(Int_t idet=0; idet<AliTRDgeometry::kNdet; idet++){
307 TLinearFitter *f = (TLinearFitter*)fLinearFitterArray.At(idet);
308 if(f) { delete f;}
309 }
46c4871e 310 if(fLinearVdriftFit) delete fLinearVdriftFit;
f162af62 311 if (fGeo) {
312 delete fGeo;
313 }
55a288e5 314
315}
55a288e5 316//_____________________________________________________________________________
317void AliTRDCalibraFillHisto::Destroy()
318{
319 //
320 // Delete instance
321 //
322
323 if (fgInstance) {
324 delete fgInstance;
325 fgInstance = 0x0;
326 }
55a288e5 327}
bcb6fb78 328//_____________________________________________________________________________
329void AliTRDCalibraFillHisto::DestroyDebugStreamer()
330{
331 //
332 // Delete DebugStreamer
333 //
55a288e5 334
bcb6fb78 335 if ( fDebugStreamer ) delete fDebugStreamer;
64942b85 336 fDebugStreamer = 0x0;
bcb6fb78 337
338}
55a288e5 339//_____________________________________________________________________________
340void AliTRDCalibraFillHisto::ClearHistos()
341{
342 //
343 // Delete the histos
344 //
345
346 if (fPH2d) {
347 delete fPH2d;
348 fPH2d = 0x0;
349 }
350 if (fCH2d) {
351 delete fCH2d;
352 fCH2d = 0x0;
353 }
354 if (fPRF2d) {
355 delete fPRF2d;
356 fPRF2d = 0x0;
357 }
3a0f6479 358
55a288e5 359}
bcb6fb78 360//////////////////////////////////////////////////////////////////////////////////
361// calibration with AliTRDtrackV1: Init, Update
362//////////////////////////////////////////////////////////////////////////////////
55a288e5 363//____________Functions for initialising the AliTRDCalibraFillHisto in the code_________
9ff7be6d 364Bool_t AliTRDCalibraFillHisto::Init2Dhistos(Int_t nboftimebin)
55a288e5 365{
366 //
bcb6fb78 367 // Init the histograms and stuff to be filled
d0569428 368 //
369
55a288e5 370 // DB Setting
371 // Get cal
372 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
373 if (!cal) {
374 AliInfo("Could not get calibDB");
375 return kFALSE;
376 }
377 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
378 if (!parCom) {
379 AliInfo("Could not get CommonParam");
380 return kFALSE;
381 }
382
383 // Some parameters
9ff7be6d 384 if(nboftimebin > 0) fTimeMax = nboftimebin;
c698a3d9 385 else fTimeMax = cal->GetNumberOfTimeBinsDCS();
5496295d 386 if(fTimeMax <= 0) fTimeMax = 30;
09ea3770 387 printf("////////////////////////////////////////////\n");
388 printf("Number of time bins in calibration component %d\n",fTimeMax);
389 printf("////////////////////////////////////////////\n");
d95484e4 390 fSf = parCom->GetSamplingFrequency();
64942b85 391 if(!fNormalizeNbOfCluster) fRelativeScale = 20.0;
392 else fRelativeScale = 1.18;
d95484e4 393 fNumberClustersf = fTimeMax;
6aafa7ea 394 fNumberClusters = (Int_t)(fNumberClustersProcent*fTimeMax);
d0569428 395
1ca79a00 396 // Init linear fitter
397 if(!fLinearFitterTracklet) {
398 fLinearFitterTracklet = new TLinearFitter(2,"pol1");
64942b85 399 fLinearFitterTracklet->StoreData(kTRUE);
1ca79a00 400 }
170c35f1 401
1ca79a00 402 //calib object from database used for reconstruction
403 if( fCalDetGain ){
404 fCalDetGain->~AliTRDCalDet();
405 new(fCalDetGain) AliTRDCalDet(*(cal->GetGainFactorDet()));
406 }else fCalDetGain = new AliTRDCalDet(*(cal->GetGainFactorDet()));
407
3a0f6479 408 // Calcul Xbins Chambd0, Chamb2
0c349049 409 Int_t ntotal0 = CalculateTotalNumberOfBins(0);
410 Int_t ntotal1 = CalculateTotalNumberOfBins(1);
411 Int_t ntotal2 = CalculateTotalNumberOfBins(2);
3a0f6479 412
413 // If vector method On initialised all the stuff
414 if(fVector2d){
55a288e5 415 fCalibraVector = new AliTRDCalibraVector();
170c35f1 416 fCalibraVector->SetNumberBinCharge(fNumberBinCharge);
417 fCalibraVector->SetTimeMax(fTimeMax);
418 if(fNgroupprf != 0) {
419 fCalibraVector->SetNumberBinPRF(2*fNgroupprf*fNumberBinPRF);
420 fCalibraVector->SetPRFRange((Float_t)(3.0*fNgroupprf));
421 }
422 else {
423 fCalibraVector->SetNumberBinPRF(fNumberBinPRF);
424 fCalibraVector->SetPRFRange(1.5);
425 }
3a0f6479 426 for(Int_t k = 0; k < 3; k++){
427 fCalibraVector->SetDetCha0(k,fCalibraMode->GetDetChamb0(k));
428 fCalibraVector->SetDetCha2(k,fCalibraMode->GetDetChamb2(k));
429 }
e526983e 430 fCalibraVector->SetNzNrphi(0,fCalibraMode->GetNz(0),fCalibraMode->GetNrphi(0));
431 fCalibraVector->SetNzNrphi(1,fCalibraMode->GetNz(1),fCalibraMode->GetNrphi(1));
432 fCalibraVector->SetNzNrphi(2,fCalibraMode->GetNz(2),fCalibraMode->GetNrphi(2));
433 fCalibraVector->SetNbGroupPRF(fNgroupprf);
55a288e5 434 }
170c35f1 435
55a288e5 436 // Create the 2D histos corresponding to the pad groupCalibration mode
437 if (fCH2dOn) {
438
439 AliInfo(Form("The pad calibration mode for the relative gain calibration: Nz %d, and Nrphi %d"
440 ,fCalibraMode->GetNz(0)
441 ,fCalibraMode->GetNrphi(0)));
442
55a288e5 443 // Create the 2D histo
444 if (fHisto2d) {
0c349049 445 CreateCH2d(ntotal0);
55a288e5 446 }
55a288e5 447 // Variable
448 fAmpTotal = new Float_t[TMath::Max(fCalibraMode->GetDetChamb2(0),fCalibraMode->GetDetChamb0(0))];
449 for (Int_t k = 0; k < TMath::Max(fCalibraMode->GetDetChamb2(0),fCalibraMode->GetDetChamb0(0)); k++) {
450 fAmpTotal[k] = 0.0;
451 }
170c35f1 452 //Statistics
0c349049 453 fEntriesCH = new Int_t[ntotal0];
454 for(Int_t k = 0; k < ntotal0; k++){
170c35f1 455 fEntriesCH[k] = 0;
456 }
3a0f6479 457
55a288e5 458 }
55a288e5 459 if (fPH2dOn) {
460
461 AliInfo(Form("The pad calibration mode for the drift velocity calibration: Nz %d, and Nrphi %d"
462 ,fCalibraMode->GetNz(1)
463 ,fCalibraMode->GetNrphi(1)));
464
55a288e5 465 // Create the 2D histo
466 if (fHisto2d) {
0c349049 467 CreatePH2d(ntotal1);
55a288e5 468 }
55a288e5 469 // Variable
470 fPHPlace = new Short_t[fTimeMax];
471 for (Int_t k = 0; k < fTimeMax; k++) {
472 fPHPlace[k] = -1;
473 }
474 fPHValue = new Float_t[fTimeMax];
475 for (Int_t k = 0; k < fTimeMax; k++) {
476 fPHValue[k] = 0.0;
477 }
170c35f1 478 }
479 if (fLinearFitterOn) {
46c4871e 480 if(fLinearFitterDebugOn) {
481 fLinearFitterArray.SetName("ArrayLinearFitters");
482 fEntriesLinearFitter = new Int_t[540];
483 for(Int_t k = 0; k < 540; k++){
484 fEntriesLinearFitter[k] = 0;
485 }
170c35f1 486 }
3a0f6479 487 fLinearVdriftFit = new AliTRDCalibraVdriftLinearFit();
55a288e5 488 }
489
490 if (fPRF2dOn) {
491
492 AliInfo(Form("The pad calibration mode for the PRF calibration: Nz %d, and Nrphi %d"
493 ,fCalibraMode->GetNz(2)
494 ,fCalibraMode->GetNrphi(2)));
55a288e5 495 // Create the 2D histo
496 if (fHisto2d) {
0c349049 497 CreatePRF2d(ntotal2);
55a288e5 498 }
55a288e5 499 }
500
501 return kTRUE;
502
503}
d0569428 504//____________Offline tracking in the AliTRDtracker____________________________
6bbdc11a 505Bool_t AliTRDCalibraFillHisto::UpdateHistograms(const AliTRDtrack *t)
d0569428 506{
507 //
bcb6fb78 508 // Use AliTRDtrack for the calibration
d0569428 509 //
510
d0569428 511
bcb6fb78 512 AliTRDcluster *cl = 0x0; // pointeur to cluster
513 Int_t index0 = 0; // index of the first cluster in the current chamber
514 Int_t index1 = 0; // index of the current cluster in the current chamber
d0569428 515
bcb6fb78 516 //////////////////////////////
d0569428 517 // loop over the clusters
bcb6fb78 518 ///////////////////////////////
d0569428 519 while((cl = t->GetCluster(index1))){
520
bcb6fb78 521 /////////////////////////////////////////////////////////////////////////
522 // Fill the infos for the previous clusters if not the same detector
523 ////////////////////////////////////////////////////////////////////////
d0569428 524 Int_t detector = cl->GetDetector();
d0569428 525 if ((detector != fDetectorPreviousTrack) &&
526 (index0 != index1)) {
527
528 fNumberTrack++;
529
d0569428 530 //If the same track, then look if the previous detector is in
531 //the same plane, if yes: not a good track
053767a4 532 if ((GetLayer(detector) == GetLayer(fDetectorPreviousTrack))) {
bcb6fb78 533 return kFALSE;
d0569428 534 }
535
536 // Fill only if the track doesn't touch a masked pad or doesn't
bcb6fb78 537 if (fGoodTracklet) {
d0569428 538
539 // drift velocity unables to cut bad tracklets
bcb6fb78 540 Bool_t pass = FindP1TrackPHtracklet(t,index0,index1);
d0569428 541
542 // Gain calibration
543 if (fCH2dOn) {
bcb6fb78 544 FillTheInfoOfTheTrackCH(index1-index0);
d0569428 545 }
546
547 // PH calibration
548 if (fPH2dOn) {
549 FillTheInfoOfTheTrackPH();
550 }
551
bcb6fb78 552 if(pass && fPRF2dOn) HandlePRFtracklet(t,index0,index1);
d0569428 553
554
bcb6fb78 555 } // if a good tracklet
d0569428 556
557 // reset stuff
bcb6fb78 558 ResetfVariablestracklet();
d0569428 559 index0 = index1;
560
561 } // Fill at the end the charge
bcb6fb78 562
563
564 /////////////////////////////////////////////////////////////
565 // Localise and take the calib gain object for the detector
566 ////////////////////////////////////////////////////////////
d0569428 567 if (detector != fDetectorPreviousTrack) {
568
569 //Localise the detector
570 LocalisationDetectorXbins(detector);
571
572 // Get cal
573 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
574 if (!cal) {
575 AliInfo("Could not get calibDB");
576 return kFALSE;
577 }
578
579 // Get calib objects
1ca79a00 580 if( fCalROCGain ){
64942b85 581 if(!fIsHLT){
582 fCalROCGain->~AliTRDCalROC();
583 new(fCalROCGain) AliTRDCalROC(*(cal->GetGainFactorROC(detector)));
584 }
1ca79a00 585 }else fCalROCGain = new AliTRDCalROC(*(cal->GetGainFactorROC(detector)));
586
d0569428 587 }
588
589 // Reset the detectbjobsor
590 fDetectorPreviousTrack = detector;
591
bcb6fb78 592 ///////////////////////////////////////
593 // Store the info of the cluster
594 ///////////////////////////////////////
595 Int_t row = cl->GetPadRow();
596 Int_t col = cl->GetPadCol();
d95484e4 597 CheckGoodTrackletV1(cl);
d0569428 598 Int_t group[2] = {0,0};
bcb6fb78 599 if(fCH2dOn) group[0] = CalculateCalibrationGroup(0,row,col);
600 if(fPH2dOn) group[1] = CalculateCalibrationGroup(1,row,col);
601 StoreInfoCHPHtrack(cl,t->GetClusterdQdl(index1),group,row,col);
d0569428 602
603 index1++;
604
605 } // while on clusters
606
bcb6fb78 607 ///////////////////////////
d0569428 608 // Fill the last plane
bcb6fb78 609 //////////////////////////
d0569428 610 if( index0 != index1 ){
d0569428 611
612 fNumberTrack++;
613
bcb6fb78 614 if (fGoodTracklet) {
d0569428 615
616 // drift velocity unables to cut bad tracklets
bcb6fb78 617 Bool_t pass = FindP1TrackPHtracklet(t,index0,index1);
d0569428 618
619 // Gain calibration
620 if (fCH2dOn) {
bcb6fb78 621 FillTheInfoOfTheTrackCH(index1-index0);
d0569428 622 }
623
624 // PH calibration
625 if (fPH2dOn) {
626 FillTheInfoOfTheTrackPH();
627 }
628
bcb6fb78 629 if(pass && fPRF2dOn) HandlePRFtracklet(t,index0,index1);
d0569428 630
bcb6fb78 631 } // if a good tracklet
d0569428 632
633 }
170c35f1 634
d0569428 635 // reset stuff
bcb6fb78 636 ResetfVariablestracklet();
d0569428 637
638 return kTRUE;
639
170c35f1 640}
55a288e5 641//____________Offline tracking in the AliTRDtracker____________________________
6bbdc11a 642Bool_t AliTRDCalibraFillHisto::UpdateHistogramsV1(const AliTRDtrackV1 *t)
55a288e5 643{
644 //
bcb6fb78 645 // Use AliTRDtrackV1 for the calibration
55a288e5 646 //
647
bcb6fb78 648
649 const AliTRDseedV1 *tracklet = 0x0; // tracklet per plane
4aad967c 650 AliTRDcluster *cl = 0x0; // cluster attached now to the tracklet
1785640c 651 AliTRDcluster *cls = 0x0; // shared cluster attached now to the tracklet
4aad967c 652 Bool_t newtr = kTRUE; // new track
bcb6fb78 653
1ca79a00 654 // Get cal
01239968 655 // AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
656 /*
1ca79a00 657 if (!cal) {
658 AliInfo("Could not get calibDB");
659 return kFALSE;
660 }
01239968 661*/
662 if (!fCalibDB) {
663 AliInfo("Could not get calibDB");
664 return kFALSE;
665 }
666
bcb6fb78 667
668 ///////////////////////////
669 // loop over the tracklet
670 ///////////////////////////
671 for(Int_t itr = 0; itr < 6; itr++){
55a288e5 672
bcb6fb78 673 if(!(tracklet = t->GetTracklet(itr))) continue;
674 if(!tracklet->IsOK()) continue;
675 fNumberTrack++;
676 ResetfVariablestracklet();
677
678 //////////////////////////////////////////
679 // localisation of the tracklet and dqdl
680 //////////////////////////////////////////
053767a4 681 Int_t layer = tracklet->GetPlane();
bcb6fb78 682 Int_t ic = 0;
683 while(!(cl = tracklet->GetClusters(ic++))) continue;
684 Int_t detector = cl->GetDetector();
685 if (detector != fDetectorPreviousTrack) {
4aad967c 686 // if not a new track
687 if(!newtr){
688 // don't use the rest of this track if in the same plane
689 if (layer == GetLayer(fDetectorPreviousTrack)) {
690 //printf("bad tracklet, same layer for detector %d\n",detector);
691 break;
692 }
bcb6fb78 693 }
694 //Localise the detector bin
695 LocalisationDetectorXbins(detector);
bcb6fb78 696 // Get calib objects
1ca79a00 697 if( fCalROCGain ){
64942b85 698 if(!fIsHLT){
699 fCalROCGain->~AliTRDCalROC();
01239968 700 new(fCalROCGain) AliTRDCalROC(*(fCalibDB->GetGainFactorROC(detector)));
64942b85 701 }
01239968 702 }else fCalROCGain = new AliTRDCalROC(*(fCalibDB->GetGainFactorROC(detector)));
1ca79a00 703
bcb6fb78 704 // reset
705 fDetectorPreviousTrack = detector;
55a288e5 706 }
4aad967c 707 newtr = kFALSE;
55a288e5 708
bcb6fb78 709 ////////////////////////////
710 // loop over the clusters
711 ////////////////////////////
712 Int_t nbclusters = 0;
e3cf3d02 713 for(int jc=0; jc<AliTRDseedV1::kNtb; jc++){
7bce990c 714 if(!(cl = tracklet->GetClusters(jc))) continue;
bcb6fb78 715 nbclusters++;
716
717 // Store the info bis of the tracklet
718 Int_t row = cl->GetPadRow();
719 Int_t col = cl->GetPadCol();
d95484e4 720 CheckGoodTrackletV1(cl);
bcb6fb78 721 Int_t group[2] = {0,0};
722 if(fCH2dOn) group[0] = CalculateCalibrationGroup(0,row,col);
723 if(fPH2dOn) group[1] = CalculateCalibrationGroup(1,row,col);
1785640c 724 // Add the charge if shared cluster
725 cls = tracklet->GetClusters(jc+AliTRDseedV1::kNtb);
726 //
727 StoreInfoCHPHtrack(cl, tracklet->GetdQdl(jc),group,row,col,cls);
bcb6fb78 728 }
729
730 ////////////////////////////////////////
731 // Fill the stuffs if a good tracklet
732 ////////////////////////////////////////
733 if (fGoodTracklet) {
e526983e 734
bcb6fb78 735 // drift velocity unables to cut bad tracklets
736 Bool_t pass = FindP1TrackPHtrackletV1(tracklet, nbclusters);
737
e526983e 738 //printf("pass %d and nbclusters %d\n",pass,nbclusters);
739
55a288e5 740 // Gain calibration
741 if (fCH2dOn) {
bcb6fb78 742 FillTheInfoOfTheTrackCH(nbclusters);
55a288e5 743 }
bcb6fb78 744
55a288e5 745 // PH calibration
746 if (fPH2dOn) {
747 FillTheInfoOfTheTrackPH();
748 }
bcb6fb78 749
750 if(pass && fPRF2dOn) HandlePRFtrackletV1(tracklet,nbclusters);
751
752 } // if a good tracklet
bcb6fb78 753 }
4aad967c 754
55a288e5 755 return kTRUE;
756
757}
bcb6fb78 758///////////////////////////////////////////////////////////////////////////////////
759// Routine inside the update with AliTRDtrack
760///////////////////////////////////////////////////////////////////////////////////
761//____________Offine tracking in the AliTRDtracker_____________________________
6bbdc11a 762Bool_t AliTRDCalibraFillHisto::FindP1TrackPHtracklet(const AliTRDtrack *t, Int_t index0, Int_t index1)
55a288e5 763{
764 //
bcb6fb78 765 // Drift velocity calibration:
766 // Fit the clusters with a straight line
767 // From the slope find the drift velocity
55a288e5 768 //
55a288e5 769
bcb6fb78 770 //Number of points: if less than 3 return kFALSE
771 Int_t npoints = index1-index0;
772 if(npoints <= 2) return kFALSE;
55a288e5 773
bcb6fb78 774 /////////////////
775 //Variables
776 ////////////////
777 Int_t detector = ((AliTRDcluster *) t->GetCluster(index0))->GetDetector(); //detector
778 // parameters of the track
053767a4 779 Double_t snp = t->GetSnpPlane(GetLayer(detector)); // sin angle in the plan yx track
780 Double_t tgl = t->GetTglPlane(GetLayer(detector)); // dz/dl and not dz/dx!
bcb6fb78 781 Double_t tnp = 0.0; // tan angle in the plan xy track
782 if( TMath::Abs(snp) < 1.){
60e55aee 783 tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp));
bcb6fb78 784 }
785 Float_t dzdx = tgl*TMath::Sqrt(1+tnp*tnp); // dz/dx now from dz/dl
786 // tilting pad and cross row
787 Int_t crossrow = 0; // if it crosses a pad row
788 Int_t rowp = -1; // if it crosses a pad row
053767a4 789 AliTRDpadPlane *padplane = fGeo->GetPadPlane(GetLayer(detector),GetStack(detector));
bcb6fb78 790 Double_t tiltingangle = padplane->GetTiltingAngle(); // tiltingangle of the pad
791 Float_t tnt = TMath::Tan(tiltingangle/180.*TMath::Pi()); // tan tiltingangle
792 // linear fit
1ca79a00 793 fLinearFitterTracklet->ClearPoints();
bcb6fb78 794 Double_t dydt = 0.0; // dydt tracklet after straight line fit
795 Double_t errorpar = 0.0; // error after straight line fit on dy/dt
796 Double_t pointError = 0.0; // error after straight line fit
797 Int_t nbli = 0; // number linear fitter points
bcb6fb78 798
799 //////////////////////////////
800 // loop over clusters
801 ////////////////////////////
802 for(Int_t k = 0; k < npoints; k++){
803
804 AliTRDcluster *cl = (AliTRDcluster *) t->GetCluster(k+index0);
b70c68da 805 if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue;
bcb6fb78 806 Double_t ycluster = cl->GetY();
807 Int_t time = cl->GetPadTime();
808 Double_t timeis = time/fSf;
809 //Double_t q = cl->GetQ();
810 //See if cross two pad rows
811 Int_t row = cl->GetPadRow();
812 if(k==0) rowp = row;
813 if(row != rowp) crossrow = 1;
55a288e5 814
1ca79a00 815 fLinearFitterTracklet->AddPoint(&timeis,ycluster,1);
bcb6fb78 816 nbli++;
55a288e5 817
bcb6fb78 818 }
55a288e5 819
bcb6fb78 820 //////////////////////////////
821 // linear fit
822 //////////////////////////////
1ca79a00 823 if(nbli <= 2){
824 fLinearFitterTracklet->ClearPoints();
825 return kFALSE;
826 }
bcb6fb78 827 TVectorD pars;
1ca79a00 828 fLinearFitterTracklet->Eval();
829 fLinearFitterTracklet->GetParameters(pars);
830 pointError = TMath::Sqrt(fLinearFitterTracklet->GetChisquare()/(nbli-2));
831 errorpar = fLinearFitterTracklet->GetParError(1)*pointError;
bcb6fb78 832 dydt = pars[1];
1ca79a00 833 fLinearFitterTracklet->ClearPoints();
834
bcb6fb78 835 /////////////////////////////
836 // debug
837 ////////////////////////////
838 if(fDebugLevel > 0){
839 if ( !fDebugStreamer ) {
840 //debug stream
841 TDirectory *backup = gDirectory;
842 fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root");
843 if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer
844 }
55a288e5 845
55a288e5 846
bcb6fb78 847 (* fDebugStreamer) << "FindP1TrackPHtracklet0"<<
848 //"snpright="<<snpright<<
849 "npoints="<<npoints<<
850 "nbli="<<nbli<<
851 "detector="<<detector<<
852 "snp="<<snp<<
853 "tnp="<<tnp<<
854 "tgl="<<tgl<<
855 "tnt="<<tnt<<
856 "dydt="<<dydt<<
857 "dzdx="<<dzdx<<
858 "crossrow="<<crossrow<<
859 "errorpar="<<errorpar<<
860 "pointError="<<pointError<<
861 "\n";
170c35f1 862
170c35f1 863
bcb6fb78 864 Int_t nbclusters = index1-index0;
053767a4 865 Int_t layer = GetLayer(fDetectorPreviousTrack);
170c35f1 866
bcb6fb78 867 (* fDebugStreamer) << "FindP1TrackPHtracklet1"<<
868 //"snpright="<<snpright<<
869 "nbclusters="<<nbclusters<<
870 "detector="<<fDetectorPreviousTrack<<
053767a4 871 "layer="<<layer<<
bcb6fb78 872 "\n";
170c35f1 873
bcb6fb78 874 }
875
876 ///////////////////////////
877 // quality cuts
878 ///////////////////////////
879 if(npoints < fNumberClusters) return kFALSE;
880 if(npoints > fNumberClustersf) return kFALSE;
881 if(pointError >= 0.1) return kFALSE;
882 if(crossrow == 1) return kFALSE;
170c35f1 883
bcb6fb78 884 ////////////////////////////
885 // fill
886 ////////////////////////////
887 if(fLinearFitterOn){
888 //Add to the linear fitter of the detector
889 if( TMath::Abs(snp) < 1.){
890 Double_t x = tnp-dzdx*tnt;
bcb6fb78 891 if(fLinearFitterDebugOn) {
46c4871e 892 (GetLinearFitter(detector,kTRUE))->AddPoint(&x,dydt);
893 fEntriesLinearFitter[detector]++;
bcb6fb78 894 }
46c4871e 895 fLinearVdriftFit->Update(detector,x,pars[1]);
bcb6fb78 896 }
897 }
170c35f1 898
bcb6fb78 899 return kTRUE;
170c35f1 900}
bcb6fb78 901//____________Offine tracking in the AliTRDtracker_____________________________
902Bool_t AliTRDCalibraFillHisto::FindP1TrackPHtrackletV1(const AliTRDseedV1 *tracklet, Int_t nbclusters)
170c35f1 903{
904 //
bcb6fb78 905 // Drift velocity calibration:
906 // Fit the clusters with a straight line
907 // From the slope find the drift velocity
170c35f1 908 //
909
bcb6fb78 910 ////////////////////////////////////////////////
911 //Number of points: if less than 3 return kFALSE
912 /////////////////////////////////////////////////
913 if(nbclusters <= 2) return kFALSE;
170c35f1 914
bcb6fb78 915 ////////////
916 //Variables
917 ////////////
bcb6fb78 918 // results of the linear fit
919 Double_t dydt = 0.0; // dydt tracklet after straight line fit
920 Double_t errorpar = 0.0; // error after straight line fit on dy/dt
921 Double_t pointError = 0.0; // error after straight line fit
922 // pad row problemes: avoid tracklet that cross pad rows, tilting angle in the constant
923 Int_t crossrow = 0; // if it crosses a pad row
924 Int_t rowp = -1; // if it crosses a pad row
925 Float_t tnt = tracklet->GetTilt(); // tan tiltingangle
1ca79a00 926 fLinearFitterTracklet->ClearPoints();
927
170c35f1 928
bcb6fb78 929 ///////////////////////////////////////////
930 // Take the parameters of the track
931 //////////////////////////////////////////
932 // take now the snp, tnp and tgl from the track
933 Double_t snp = tracklet->GetSnp(); // sin dy/dx at the end of the chamber
934 Double_t tnp = 0.0; // dy/dx at the end of the chamber
935 if( TMath::Abs(snp) < 1.){
60e55aee 936 tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp));
bcb6fb78 937 }
938 Double_t tgl = tracklet->GetTgl(); // dz/dl
939 Double_t dzdx = tgl*TMath::Sqrt(1+tnp*tnp); // dz/dx calculated from dz/dl
940 // at the entrance
941 //Double_t tnp = tracklet->GetYref(1); // dy/dx at the entrance of the chamber
942 //Double_t tgl = tracklet->GetZref(1); // dz/dl at the entrance of the chamber
943 //Double_t dzdx = tgl; //*TMath::Sqrt(1+tnp*tnp); // dz/dx from dz/dl
944 // at the end with correction due to linear fit
945 //Double_t tnp = tracklet->GetYfit(1); // dy/dx at the end of the chamber after fit correction
946 //Double_t tgl = tracklet->GetZfit(1); // dz/dl at the end of the chamber after fit correction
947
948
949 ////////////////////////////
950 // loop over the clusters
951 ////////////////////////////
952 Int_t nbli = 0;
953 AliTRDcluster *cl = 0x0;
6aafa7ea 954 //////////////////////////////
955 // Check no shared clusters
956 //////////////////////////////
957 for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
958 if((cl = tracklet->GetClusters(icc))) crossrow = 1;
959 }
960 //////////////////////////////////
961 // Loop clusters
962 //////////////////////////////////
e3cf3d02 963 for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
bcb6fb78 964 if(!(cl = tracklet->GetClusters(ic))) continue;
b70c68da 965 if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue;
bcb6fb78 966
967 Double_t ycluster = cl->GetY();
968 Int_t time = cl->GetPadTime();
969 Double_t timeis = time/fSf;
970 //See if cross two pad rows
971 Int_t row = cl->GetPadRow();
972 if(rowp==-1) rowp = row;
973 if(row != rowp) crossrow = 1;
974
1ca79a00 975 fLinearFitterTracklet->AddPoint(&timeis,ycluster,1);
bcb6fb78 976 nbli++;
170c35f1 977
1785640c 978
55a288e5 979 }
980
bcb6fb78 981 ////////////////////////////////////
982 // Do the straight line fit now
983 ///////////////////////////////////
1ca79a00 984 if(nbli <= 2){
985 fLinearFitterTracklet->ClearPoints();
986 return kFALSE;
987 }
bcb6fb78 988 TVectorD pars;
1ca79a00 989 fLinearFitterTracklet->Eval();
990 fLinearFitterTracklet->GetParameters(pars);
991 pointError = TMath::Sqrt(fLinearFitterTracklet->GetChisquare()/(nbli-2));
992 errorpar = fLinearFitterTracklet->GetParError(1)*pointError;
bcb6fb78 993 dydt = pars[1];
64942b85 994 //printf("chis %f, nbli %d, pointError %f, parError %f, errorpar %f\n",fLinearFitterTracklet->GetChisquare(),nbli,pointError,fLinearFitterTracklet->GetParError(1),errorpar);
1ca79a00 995 fLinearFitterTracklet->ClearPoints();
170c35f1 996
bcb6fb78 997 ////////////////////////////////
998 // Debug stuff
999 ///////////////////////////////
1000
1001
1002 if(fDebugLevel > 0){
1003 if ( !fDebugStreamer ) {
1004 //debug stream
1005 TDirectory *backup = gDirectory;
1006 fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root");
1007 if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer
1008 }
1009
1010
053767a4 1011 Int_t layer = GetLayer(fDetectorPreviousTrack);
64942b85 1012
bcb6fb78 1013 (* fDebugStreamer) << "FindP1TrackPHtrackletV1"<<
1014 //"snpright="<<snpright<<
64942b85 1015 "nbli="<<nbli<<
bcb6fb78 1016 "nbclusters="<<nbclusters<<
1017 "detector="<<fDetectorPreviousTrack<<
053767a4 1018 "layer="<<layer<<
bcb6fb78 1019 "snp="<<snp<<
1020 "tnp="<<tnp<<
1021 "tgl="<<tgl<<
1022 "tnt="<<tnt<<
1023 "dydt="<<dydt<<
1024 "dzdx="<<dzdx<<
1025 "crossrow="<<crossrow<<
1026 "errorpar="<<errorpar<<
1027 "pointError="<<pointError<<
1028 "\n";
1029
170c35f1 1030 }
1031
bcb6fb78 1032 /////////////////////////
1033 // Cuts quality
1034 ////////////////////////
170c35f1 1035
bcb6fb78 1036 if(nbclusters < fNumberClusters) return kFALSE;
1037 if(nbclusters > fNumberClustersf) return kFALSE;
64942b85 1038 if(pointError >= 0.3) return kFALSE;
bcb6fb78 1039 if(crossrow == 1) return kFALSE;
170c35f1 1040
bcb6fb78 1041 ///////////////////////
1042 // Fill
1043 //////////////////////
170c35f1 1044
bcb6fb78 1045 if(fLinearFitterOn){
1046 //Add to the linear fitter of the detector
1047 if( TMath::Abs(snp) < 1.){
1048 Double_t x = tnp-dzdx*tnt;
bcb6fb78 1049 if(fLinearFitterDebugOn) {
46c4871e 1050 (GetLinearFitter(fDetectorPreviousTrack,kTRUE))->AddPoint(&x,dydt);
1051 fEntriesLinearFitter[fDetectorPreviousTrack]++;
bcb6fb78 1052 }
46c4871e 1053 fLinearVdriftFit->Update(fDetectorPreviousTrack,x,pars[1]);
bcb6fb78 1054 }
170c35f1 1055 }
1056
1057 return kTRUE;
170c35f1 1058}
bcb6fb78 1059//____________Offine tracking in the AliTRDtracker_____________________________
6bbdc11a 1060Bool_t AliTRDCalibraFillHisto::HandlePRFtracklet(const AliTRDtrack *t, Int_t index0, Int_t index1)
170c35f1 1061{
1062 //
bcb6fb78 1063 // PRF width calibration
1064 // Assume a Gaussian shape: determinate the position of the three pad clusters
1065 // Fit with a straight line
1066 // Take the fitted values for all the clusters (3 or 2 pad clusters)
1067 // Fill the PRF as function of angle of the track
1068 //
170c35f1 1069 //
d0569428 1070
1071
bcb6fb78 1072 //////////////////////////
1073 // variables
1074 /////////////////////////
1075 Int_t npoints = index1-index0; // number of total points
1076 Int_t nb3pc = 0; // number of three pads clusters used for fit
1077 Int_t detector = ((AliTRDcluster *) t->GetCluster(index0))->GetDetector(); // detector
1078 // To see the difference due to the fit
1079 Double_t *padPositions;
1080 padPositions = new Double_t[npoints];
1081 for(Int_t k = 0; k < npoints; k++){
1082 padPositions[k] = 0.0;
1083 }
1084 // Take the tgl and snp with the track t now
053767a4 1085 Double_t tgl = t->GetTglPlane(GetLayer(detector)); //dz/dl and not dz/dx
1086 Double_t snp = t->GetSnpPlane(GetLayer(detector)); // sin angle in xy plan
bcb6fb78 1087 Float_t dzdx = 0.0; // dzdx
1088 Float_t tnp = 0.0;
1089 if(TMath::Abs(snp) < 1.0){
60e55aee 1090 tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp));
bcb6fb78 1091 dzdx = tgl*TMath::Sqrt(1+tnp*tnp);
170c35f1 1092 }
bcb6fb78 1093 // linear fitter
1ca79a00 1094 fLinearFitterTracklet->ClearPoints();
170c35f1 1095
bcb6fb78 1096 ///////////////////////////
1097 // calcul the tnp group
1098 ///////////////////////////
1099 Bool_t echec = kFALSE;
1100 Double_t shift = 0.0;
1101 //Calculate the shift in x coresponding to this tnp
1102 if(fNgroupprf != 0.0){
1103 shift = -3.0*(fNgroupprf-1)-1.5;
1104 Double_t limithigh = -0.2*(fNgroupprf-1);
1105 if((tnp < (-0.2*fNgroupprf)) || (tnp > (0.2*fNgroupprf))) echec = kTRUE;
1106 else{
1107 while(tnp > limithigh){
1108 limithigh += 0.2;
1109 shift += 3.0;
1110 }
1111 }
d0569428 1112 }
1ca79a00 1113 if(echec) {
1114 delete [] padPositions;
1115 return kFALSE;
1116 }
d0569428 1117
bcb6fb78 1118 //////////////////////
1119 // loop clusters
1120 /////////////////////
1121 for(Int_t k = 0; k < npoints; k++){
1122 //Take the cluster
1123 AliTRDcluster *cl = (AliTRDcluster *) t->GetCluster(k+index0);
e526983e 1124 if(!cl) continue;
bcb6fb78 1125 Short_t *signals = cl->GetSignals();
e526983e 1126 Double_t time = cl->GetPadTime();
bcb6fb78 1127 //Calculate x if possible
1128 Float_t xcenter = 0.0;
7bce990c 1129 Bool_t echec1 = kTRUE;
bcb6fb78 1130 if((time<=7) || (time>=21)) continue;
1131 // Center 3 balanced: position with the center of the pad
1132 if ((((Float_t) signals[3]) > 0.0) &&
1133 (((Float_t) signals[2]) > 0.0) &&
1134 (((Float_t) signals[4]) > 0.0)) {
7bce990c 1135 echec1 = kFALSE;
bcb6fb78 1136 // Security if the denomiateur is 0
1137 if ((((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3]))) /
1138 ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))) != 1.0) {
1139 xcenter = 0.5 * (TMath::Log((Float_t) (((Float_t) signals[4]) / ((Float_t) signals[2]))))
1140 / (TMath::Log(((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3])))
1141 / ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))));
1142 }
1143 else {
7bce990c 1144 echec1 = kTRUE;
bcb6fb78 1145 }
1146 }
1147 if(TMath::Abs(xcenter) > 0.5) echec = kTRUE;
1148 if(echec) continue;
1149 //if no echec: calculate with the position of the pad
1150 // Position of the cluster
1151 Double_t padPosition = xcenter + cl->GetPadCol();
1152 padPositions[k] = padPosition;
1153 nb3pc++;
1ca79a00 1154 fLinearFitterTracklet->AddPoint(&time, padPosition,1);
bcb6fb78 1155 }//clusters loop
d0569428 1156
d0569428 1157
bcb6fb78 1158 /////////////////////////////
1159 // fit
1160 ////////////////////////////
1ca79a00 1161 if(nb3pc < 3) {
1162 delete [] padPositions;
1163 fLinearFitterTracklet->ClearPoints();
1164 return kFALSE;
1165 }
1166 fLinearFitterTracklet->Eval();
bcb6fb78 1167 TVectorD line(2);
1ca79a00 1168 fLinearFitterTracklet->GetParameters(line);
bcb6fb78 1169 Float_t pointError = -1.0;
1ca79a00 1170 if( fLinearFitterTracklet->GetChisquare()>=0.0) {
1171 pointError = TMath::Sqrt( fLinearFitterTracklet->GetChisquare()/(nb3pc-2));
37b0cf5e 1172 }
1ca79a00 1173 fLinearFitterTracklet->ClearPoints();
170c35f1 1174
bcb6fb78 1175 /////////////////////////////////////////////////////
1176 // Now fill the PRF: second loop over clusters
1177 /////////////////////////////////////////////////////
1178 for(Int_t k = 0; k < npoints; k++){
1179 //Take the cluster
1180 AliTRDcluster *cl = (AliTRDcluster *) t->GetCluster(k+index0);
2e32a5ae 1181 Short_t *signals = cl->GetSignals(); // signal
1182 Double_t time = cl->GetPadTime(); // time bin
bcb6fb78 1183 Float_t padPosTracklet = line[0]+line[1]*time; // reconstruct position from fit
1184 Float_t padPos = cl->GetPadCol(); // middle pad
1185 Double_t dpad = padPosTracklet - padPos; // reconstruct position relative to middle pad from fit
1186 Float_t ycenter = 0.0; // relative center charge
1187 Float_t ymin = 0.0; // relative left charge
1188 Float_t ymax = 0.0; // relative right charge
e4db522f 1189
bcb6fb78 1190 //Requiere simply two pads clusters at least
1191 if(((((Float_t) signals[3]) > 0.0) && (((Float_t) signals[2]) > 0.0)) ||
1192 ((((Float_t) signals[3]) > 0.0) && (((Float_t) signals[4]) > 0.0))){
1193 Float_t sum = ((Float_t) signals[2]) + ((Float_t) signals[3]) + ((Float_t) signals[4]);
1194 if(sum > 0.0) ycenter = ((Float_t) signals[3])/ sum;
1195 if(sum > 0.0) ymin = ((Float_t) signals[2])/ sum;
1196 if(sum > 0.0) ymax = ((Float_t) signals[4])/ sum;
e4db522f 1197 }
bcb6fb78 1198
1199 //calibration group
1200 Int_t rowcl = cl->GetPadRow(); // row of cluster
1201 Int_t colcl = cl->GetPadCol(); // col of cluster
1202 Int_t grouplocal = CalculateCalibrationGroup(2,rowcl,colcl); // calcul the corresponding group
1203 Int_t caligroup = fCalibraMode->GetXbins(2)+ grouplocal; // calcul the corresponding group
1204 Float_t xcl = cl->GetY(); // y cluster
1205 Float_t qcl = cl->GetQ(); // charge cluster
053767a4 1206 Int_t layer = GetLayer(detector); // layer
1207 Int_t stack = GetStack(detector); // stack
bcb6fb78 1208 Double_t xdiff = dpad; // reconstructed position constant
1209 Double_t x = dpad; // reconstructed position moved
1210 Float_t ep = pointError; // error of fit
1211 Float_t signal1 = (Float_t)signals[1]; // signal at the border
1212 Float_t signal3 = (Float_t)signals[3]; // signal
1213 Float_t signal2 = (Float_t)signals[2]; // signal
1214 Float_t signal4 = (Float_t)signals[4]; // signal
1215 Float_t signal5 = (Float_t)signals[5]; // signal at the border
1216
1217 //////////////////////////////
1218 // debug
1219 /////////////////////////////
1220 if(fDebugLevel > 0){
1221 if ( !fDebugStreamer ) {
1222 //debug stream
1223 TDirectory *backup = gDirectory;
1224 fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root");
1225 if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer
1226 }
1227
3a0f6479 1228
bcb6fb78 1229 x = xdiff;
1230 Int_t type=0;
1231 Float_t y = ycenter;
1232 (* fDebugStreamer) << "HandlePRFtracklet"<<
1233 "caligroup="<<caligroup<<
1234 "detector="<<detector<<
053767a4 1235 "layer="<<layer<<
1236 "stack="<< stack <<
bcb6fb78 1237 "npoints="<<npoints<<
1238 "Np="<<nb3pc<<
1239 "ep="<<ep<<
1240 "type="<<type<<
1241 "snp="<<snp<<
1242 "tnp="<<tnp<<
1243 "tgl="<<tgl<<
1244 "dzdx="<<dzdx<<
1245 "padPos="<<padPos<<
1246 "padPosition="<<padPositions[k]<<
1247 "padPosTracklet="<<padPosTracklet<<
1248 "x="<<x<<
1249 "y="<<y<<
1250 "xcl="<<xcl<<
1251 "qcl="<<qcl<<
1252 "signal1="<<signal1<<
1253 "signal2="<<signal2<<
1254 "signal3="<<signal3<<
1255 "signal4="<<signal4<<
1256 "signal5="<<signal5<<
1257 "time="<<time<<
1258 "\n";
1259 x=-(xdiff+1);
1260 y = ymin;
1261 type=-1;
1262 (* fDebugStreamer) << "HandlePRFtracklet"<<
1263 "caligroup="<<caligroup<<
1264 "detector="<<detector<<
053767a4 1265 "layer="<<layer<<
1266 "stack="<<stack<<
bcb6fb78 1267 "npoints="<<npoints<<
1268 "Np="<<nb3pc<<
1269 "ep="<<ep<<
1270 "type="<<type<<
1271 "snp="<<snp<<
1272 "tnp="<<tnp<<
1273 "tgl="<<tgl<<
1274 "dzdx="<<dzdx<<
1275 "padPos="<<padPos<<
1276 "padPosition="<<padPositions[k]<<
1277 "padPosTracklet="<<padPosTracklet<<
1278 "x="<<x<<
1279 "y="<<y<<
1280 "xcl="<<xcl<<
1281 "qcl="<<qcl<<
1282 "signal1="<<signal1<<
1283 "signal2="<<signal2<<
1284 "signal3="<<signal3<<
1285 "signal4="<<signal4<<
1286 "signal5="<<signal5<<
1287 "time="<<time<<
1288 "\n";
1289 x=1-xdiff;
1290 y = ymax;
1291 type=1;
1292 (* fDebugStreamer) << "HandlePRFtracklet"<<
1293 "caligroup="<<caligroup<<
1294 "detector="<<detector<<
053767a4 1295 "layer="<<layer<<
1296 "stack="<<stack<<
bcb6fb78 1297 "npoints="<<npoints<<
1298 "Np="<<nb3pc<<
1299 "ep="<<ep<<
1300 "type="<<type<<
1301 "snp="<<snp<<
1302 "tnp="<<tnp<<
1303 "tgl="<<tgl<<
1304 "dzdx="<<dzdx<<
1305 "padPos="<<padPos<<
1306 "padPosition="<<padPositions[k]<<
1307 "padPosTracklet="<<padPosTracklet<<
1308 "x="<<x<<
1309 "y="<<y<<
1310 "xcl="<<xcl<<
1311 "qcl="<<qcl<<
1312 "signal1="<<signal1<<
1313 "signal2="<<signal2<<
1314 "signal3="<<signal3<<
1315 "signal4="<<signal4<<
1316 "signal5="<<signal5<<
1317 "time="<<time<<
1318 "\n";
e4db522f 1319
bcb6fb78 1320 }
e4db522f 1321
bcb6fb78 1322 ////////////////////////////
1323 // quality cuts
1324 ///////////////////////////
1325 if(npoints < fNumberClusters) continue;
1326 if(npoints > fNumberClustersf) continue;
1327 if(nb3pc <= 5) continue;
1328 if((time >= 21) || (time < 7)) continue;
1329 if(TMath::Abs(snp) >= 1.0) continue;
1330 if(TMath::Abs(qcl) < 80) continue;
1331
1332 ////////////////////////////
1333 // Fill
1334 ///////////////////////////
1335 if (fHisto2d) {
1336 if(TMath::Abs(dpad) < 1.5) {
1337 fPRF2d->Fill(shift+dpad,(caligroup+0.5),ycenter);
1338 fPRF2d->Fill(shift-dpad,(caligroup+0.5),ycenter);
3a0f6479 1339 }
bcb6fb78 1340 if((ymin > 0.0) && (TMath::Abs(dpad+1.0) < 1.5)) {
1341 fPRF2d->Fill(shift-(dpad+1.0),(caligroup+0.5),ymin);
1342 fPRF2d->Fill(shift+(dpad+1.0),(caligroup+0.5),ymin);
1343 }
1344 if((ymax > 0.0) && (TMath::Abs(dpad-1.0) < 1.5)) {
1345 fPRF2d->Fill(shift+1.0-dpad,(caligroup+0.5),ymax);
1346 fPRF2d->Fill(shift-1.0+dpad,(caligroup+0.5),ymax);
3a0f6479 1347 }
1348 }
bcb6fb78 1349 if (fVector2d) {
1350 if(TMath::Abs(dpad) < 1.5) {
1351 fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+dpad,ycenter);
1352 fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-dpad,ycenter);
e4db522f 1353 }
bcb6fb78 1354 if((ymin > 0.0) && (TMath::Abs(dpad+1.0) < 1.5)) {
1355 fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-(dpad+1.0),ymin);
1356 fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+(dpad+1.0),ymin);
e4db522f 1357 }
bcb6fb78 1358 if((ymax > 0.0) && (TMath::Abs(dpad-1.0) < 1.5)) {
1359 fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+1.0-dpad,ymax);
1360 fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-1.0+dpad,ymax);
e4db522f 1361 }
bcb6fb78 1362 }
1363 }
1ca79a00 1364 delete [] padPositions;
bcb6fb78 1365 return kTRUE;
1366
1367}
1368//____________Offine tracking in the AliTRDtracker_____________________________
1369Bool_t AliTRDCalibraFillHisto::HandlePRFtrackletV1(const AliTRDseedV1 *tracklet, Int_t nbclusters)
1370{
1371 //
1372 // PRF width calibration
1373 // Assume a Gaussian shape: determinate the position of the three pad clusters
1374 // Fit with a straight line
1375 // Take the fitted values for all the clusters (3 or 2 pad clusters)
1376 // Fill the PRF as function of angle of the track
1377 //
1378 //
1379
4aad967c 1380 //printf("begin\n");
bcb6fb78 1381 ///////////////////////////////////////////
1382 // Take the parameters of the track
1383 //////////////////////////////////////////
1384 // take now the snp, tnp and tgl from the track
1385 Double_t snp = tracklet->GetSnp(); // sin dy/dx at the end of the chamber
1386 Double_t tnp = 0.0; // dy/dx at the end of the chamber
1387 if( TMath::Abs(snp) < 1.){
60e55aee 1388 tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp));
bcb6fb78 1389 }
1390 Double_t tgl = tracklet->GetTgl(); // dz/dl
1391 Double_t dzdx = tgl*TMath::Sqrt(1+tnp*tnp); // dz/dx calculated from dz/dl
1392 // at the entrance
1393 //Double_t tnp = tracklet->GetYref(1); // dy/dx at the entrance of the chamber
1394 //Double_t tgl = tracklet->GetZref(1); // dz/dl at the entrance of the chamber
1395 //Double_t dzdx = tgl; //*TMath::Sqrt(1+tnp*tnp); // dz/dx from dz/dl
1396 // at the end with correction due to linear fit
1397 //Double_t tnp = tracklet->GetYfit(1); // dy/dx at the end of the chamber after fit correction
1398 //Double_t tgl = tracklet->GetZfit(1); // dz/dl at the end of the chamber after fit correction
1399
1400 ///////////////////////////////
1401 // Calculate tnp group shift
1402 ///////////////////////////////
1403 Bool_t echec = kFALSE;
1404 Double_t shift = 0.0;
1405 //Calculate the shift in x coresponding to this tnp
1406 if(fNgroupprf != 0.0){
1407 shift = -3.0*(fNgroupprf-1)-1.5;
1408 Double_t limithigh = -0.2*(fNgroupprf-1);
1409 if((tnp < (-0.2*fNgroupprf)) || (tnp > (0.2*fNgroupprf))) echec = kTRUE;
1410 else{
1411 while(tnp > limithigh){
1412 limithigh += 0.2;
1413 shift += 3.0;
170c35f1 1414 }
1415 }
170c35f1 1416 }
bcb6fb78 1417 // do nothing if out of tnp range
4aad967c 1418 //printf("echec %d\n",(Int_t)echec);
bcb6fb78 1419 if(echec) return kFALSE;
3a0f6479 1420
bcb6fb78 1421 ///////////////////////
1422 // Variables
1423 //////////////////////
3a0f6479 1424
bcb6fb78 1425 Int_t nb3pc = 0; // number of three pads clusters used for fit
e3cf3d02 1426 // to see the difference between the fit and the 3 pad clusters position
1427 Double_t padPositions[AliTRDseedV1::kNtb];
1428 memset(padPositions, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
1ca79a00 1429 fLinearFitterTracklet->ClearPoints();
1430
4aad967c 1431 //printf("loop clusters \n");
bcb6fb78 1432 ////////////////////////////
1433 // loop over the clusters
1434 ////////////////////////////
1435 AliTRDcluster *cl = 0x0;
e3cf3d02 1436 for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
6cb3ebda 1437 // reject shared clusters on pad row
e526983e 1438 if((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNclusters) {
1439 if((cl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) continue;
1440 }
1441 if(!(cl = tracklet->GetClusters(ic))) continue;
1442 Double_t time = cl->GetPadTime();
bcb6fb78 1443 if((time<=7) || (time>=21)) continue;
1444 Short_t *signals = cl->GetSignals();
1445 Float_t xcenter = 0.0;
7bce990c 1446 Bool_t echec1 = kTRUE;
e4db522f 1447
bcb6fb78 1448 /////////////////////////////////////////////////////////////
1449 // Center 3 balanced: position with the center of the pad
1450 /////////////////////////////////////////////////////////////
1451 if ((((Float_t) signals[3]) > 0.0) &&
1452 (((Float_t) signals[2]) > 0.0) &&
1453 (((Float_t) signals[4]) > 0.0)) {
7bce990c 1454 echec1 = kFALSE;
bcb6fb78 1455 // Security if the denomiateur is 0
1456 if ((((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3]))) /
1457 ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))) != 1.0) {
1458 xcenter = 0.5 * (TMath::Log((Float_t) (((Float_t) signals[4]) / ((Float_t) signals[2]))))
1459 / (TMath::Log(((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3])))
1460 / ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))));
3a0f6479 1461 }
bcb6fb78 1462 else {
7bce990c 1463 echec1 = kTRUE;
bcb6fb78 1464 }
1465 }
7bce990c 1466 if(TMath::Abs(xcenter) > 0.5) echec1 = kTRUE;
1467 if(echec1) continue;
e4db522f 1468
bcb6fb78 1469 ////////////////////////////////////////////////////////
7bce990c 1470 //if no echec1: calculate with the position of the pad
bcb6fb78 1471 // Position of the cluster
1472 // fill the linear fitter
1473 ///////////////////////////////////////////////////////
1474 Double_t padPosition = xcenter + cl->GetPadCol();
1475 padPositions[ic] = padPosition;
1476 nb3pc++;
1ca79a00 1477 fLinearFitterTracklet->AddPoint(&time, padPosition,1);
e4db522f 1478
e4db522f 1479
bcb6fb78 1480 }//clusters loop
1481
4aad967c 1482 //printf("Fin loop clusters \n");
bcb6fb78 1483 //////////////////////////////
1484 // fit with a straight line
1485 /////////////////////////////
1ca79a00 1486 if(nb3pc < 3){
1ca79a00 1487 fLinearFitterTracklet->ClearPoints();
1ca79a00 1488 return kFALSE;
1489 }
1490 fLinearFitterTracklet->Eval();
bcb6fb78 1491 TVectorD line(2);
1ca79a00 1492 fLinearFitterTracklet->GetParameters(line);
bcb6fb78 1493 Float_t pointError = -1.0;
1ca79a00 1494 if( fLinearFitterTracklet->GetChisquare()>=0.0) {
1495 pointError = TMath::Sqrt( fLinearFitterTracklet->GetChisquare()/(nb3pc-2));
4aad967c 1496 }
1ca79a00 1497 fLinearFitterTracklet->ClearPoints();
1498
4aad967c 1499 //printf("PRF second loop \n");
bcb6fb78 1500 ////////////////////////////////////////////////
1501 // Fill the PRF: Second loop over clusters
1502 //////////////////////////////////////////////
e3cf3d02 1503 for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
6cb3ebda 1504 // reject shared clusters on pad row
8d2bec9e 1505 if(((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNclusters) && (cl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) continue;
6cb3ebda 1506 //
bcb6fb78 1507 if(!(cl = tracklet->GetClusters(ic))) continue;
162a3e73 1508
bcb6fb78 1509 Short_t *signals = cl->GetSignals(); // signal
e526983e 1510 Double_t time = cl->GetPadTime(); // time bin
bcb6fb78 1511 Float_t padPosTracklet = line[0]+line[1]*time; // reconstruct position from fit
1512 Float_t padPos = cl->GetPadCol(); // middle pad
1513 Double_t dpad = padPosTracklet - padPos; // reconstruct position relative to middle pad from fit
1514 Float_t ycenter = 0.0; // relative center charge
1515 Float_t ymin = 0.0; // relative left charge
1516 Float_t ymax = 0.0; // relative right charge
162a3e73 1517
bcb6fb78 1518 ////////////////////////////////////////////////////////////////
1519 // Calculate ycenter, ymin and ymax even for two pad clusters
1520 ////////////////////////////////////////////////////////////////
1521 if(((((Float_t) signals[3]) > 0.0) && (((Float_t) signals[2]) > 0.0)) ||
1522 ((((Float_t) signals[3]) > 0.0) && (((Float_t) signals[4]) > 0.0))){
1523 Float_t sum = ((Float_t) signals[2]) + ((Float_t) signals[3]) + ((Float_t) signals[4]);
1524 if(sum > 0.0) ycenter = ((Float_t) signals[3])/ sum;
1525 if(sum > 0.0) ymin = ((Float_t) signals[2])/ sum;
1526 if(sum > 0.0) ymax = ((Float_t) signals[4])/ sum;
3a0f6479 1527 }
1528
bcb6fb78 1529 /////////////////////////
1530 // Calibration group
1531 ////////////////////////
1532 Int_t rowcl = cl->GetPadRow(); // row of cluster
1533 Int_t colcl = cl->GetPadCol(); // col of cluster
1534 Int_t grouplocal = CalculateCalibrationGroup(2,rowcl,colcl); // calcul the corresponding group
1535 Int_t caligroup = fCalibraMode->GetXbins(2)+ grouplocal; // calcul the corresponding group
1536 Float_t xcl = cl->GetY(); // y cluster
1537 Float_t qcl = cl->GetQ(); // charge cluster
053767a4 1538 Int_t layer = GetLayer(fDetectorPreviousTrack); // layer
1539 Int_t stack = GetStack(fDetectorPreviousTrack); // stack
bcb6fb78 1540 Double_t xdiff = dpad; // reconstructed position constant
1541 Double_t x = dpad; // reconstructed position moved
1542 Float_t ep = pointError; // error of fit
1543 Float_t signal1 = (Float_t)signals[1]; // signal at the border
1544 Float_t signal3 = (Float_t)signals[3]; // signal
1545 Float_t signal2 = (Float_t)signals[2]; // signal
1546 Float_t signal4 = (Float_t)signals[4]; // signal
1547 Float_t signal5 = (Float_t)signals[5]; // signal at the border
1548
1ca79a00 1549
1550
bcb6fb78 1551 /////////////////////
1552 // Debug stuff
1553 ////////////////////
1554
1555 if(fDebugLevel > 0){
1556 if ( !fDebugStreamer ) {
1557 //debug stream
1558 TDirectory *backup = gDirectory;
1559 fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root");
1560 if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer
1561 }
1562
1563 x = xdiff;
1564 Int_t type=0;
1565 Float_t y = ycenter;
1566 (* fDebugStreamer) << "HandlePRFtrackletV1"<<
1567 "caligroup="<<caligroup<<
1568 "detector="<<fDetectorPreviousTrack<<
053767a4 1569 "layer="<<layer<<
1570 "stack="<<stack<<
bcb6fb78 1571 "npoints="<<nbclusters<<
1572 "Np="<<nb3pc<<
1573 "ep="<<ep<<
1574 "type="<<type<<
1575 "snp="<<snp<<
1576 "tnp="<<tnp<<
1577 "tgl="<<tgl<<
1578 "dzdx="<<dzdx<<
1579 "padPos="<<padPos<<
1580 "padPosition="<<padPositions[ic]<<
1581 "padPosTracklet="<<padPosTracklet<<
1582 "x="<<x<<
1583 "y="<<y<<
1584 "xcl="<<xcl<<
1585 "qcl="<<qcl<<
1586 "signal1="<<signal1<<
1587 "signal2="<<signal2<<
1588 "signal3="<<signal3<<
1589 "signal4="<<signal4<<
1590 "signal5="<<signal5<<
1591 "time="<<time<<
1592 "\n";
1593 x=-(xdiff+1);
1594 y = ymin;
1595 type=-1;
1596 (* fDebugStreamer) << "HandlePRFtrackletV1"<<
1597 "caligroup="<<caligroup<<
1598 "detector="<<fDetectorPreviousTrack<<
053767a4 1599 "layer="<<layer<<
1600 "stack="<<stack<<
bcb6fb78 1601 "npoints="<<nbclusters<<
1602 "Np="<<nb3pc<<
1603 "ep="<<ep<<
1604 "type="<<type<<
1605 "snp="<<snp<<
1606 "tnp="<<tnp<<
1607 "tgl="<<tgl<<
1608 "dzdx="<<dzdx<<
1609 "padPos="<<padPos<<
1610 "padPosition="<<padPositions[ic]<<
1611 "padPosTracklet="<<padPosTracklet<<
1612 "x="<<x<<
1613 "y="<<y<<
1614 "xcl="<<xcl<<
1615 "qcl="<<qcl<<
1616 "signal1="<<signal1<<
1617 "signal2="<<signal2<<
1618 "signal3="<<signal3<<
1619 "signal4="<<signal4<<
1620 "signal5="<<signal5<<
1621 "time="<<time<<
1622 "\n";
1623 x=1-xdiff;
1624 y = ymax;
1625 type=1;
1626 (* fDebugStreamer) << "HandlePRFtrackletV1"<<
1627 "caligroup="<<caligroup<<
1628 "detector="<<fDetectorPreviousTrack<<
053767a4 1629 "layer="<<layer<<
1630 "stack="<<stack<<
bcb6fb78 1631 "npoints="<<nbclusters<<
1632 "Np="<<nb3pc<<
1633 "ep="<<ep<<
1634 "type="<<type<<
1635 "snp="<<snp<<
1636 "tnp="<<tnp<<
1637 "tgl="<<tgl<<
1638 "dzdx="<<dzdx<<
1639 "padPos="<<padPos<<
1640 "padPosition="<<padPositions[ic]<<
1641 "padPosTracklet="<<padPosTracklet<<
1642 "x="<<x<<
1643 "y="<<y<<
1644 "xcl="<<xcl<<
1645 "qcl="<<qcl<<
1646 "signal1="<<signal1<<
1647 "signal2="<<signal2<<
1648 "signal3="<<signal3<<
1649 "signal4="<<signal4<<
1650 "signal5="<<signal5<<
1651 "time="<<time<<
1652 "\n";
e4db522f 1653
bcb6fb78 1654 }
1655
1656 /////////////////////
1657 // Cuts quality
1658 /////////////////////
1659 if(nbclusters < fNumberClusters) continue;
1660 if(nbclusters > fNumberClustersf) continue;
1661 if(nb3pc <= 5) continue;
1662 if((time >= 21) || (time < 7)) continue;
1663 if(TMath::Abs(qcl) < 80) continue;
1664 if( TMath::Abs(snp) > 1.) continue;
1665
1666
1667 ////////////////////////
1668 // Fill the histos
1669 ///////////////////////
1670 if (fHisto2d) {
1671 if(TMath::Abs(dpad) < 1.5) {
1672 fPRF2d->Fill(shift+dpad,(caligroup+0.5),ycenter);
1673 fPRF2d->Fill(shift-dpad,(caligroup+0.5),ycenter);
1674 //printf("place %f, ycenter %f\n",(shift+dpad),ycenter);
e4db522f 1675 }
bcb6fb78 1676 if((ymin > 0.0) && (TMath::Abs(dpad+1.0) < 1.5)) {
1677 fPRF2d->Fill(shift-(dpad+1.0),(caligroup+0.5),ymin);
1678 fPRF2d->Fill(shift+(dpad+1.0),(caligroup+0.5),ymin);
e4db522f 1679 }
bcb6fb78 1680 if((ymax > 0.0) && (TMath::Abs(dpad-1.0) < 1.5)) {
1681 fPRF2d->Fill(shift+1.0-dpad,(caligroup+0.5),ymax);
1682 fPRF2d->Fill(shift-1.0+dpad,(caligroup+0.5),ymax);
e4db522f 1683 }
bcb6fb78 1684 }
1685 // vector method
1686 if (fVector2d) {
1687 if(TMath::Abs(dpad) < 1.5) {
1688 fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+dpad,ycenter);
1689 fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-dpad,ycenter);
1690 }
1691 if((ymin > 0.0) && (TMath::Abs(dpad+1.0) < 1.5)) {
1692 fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-(dpad+1.0),ymin);
1693 fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+(dpad+1.0),ymin);
1694 }
1695 if((ymax > 0.0) && (TMath::Abs(dpad-1.0) < 1.5)) {
1696 fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+1.0-dpad,ymax);
1697 fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-1.0+dpad,ymax);
170c35f1 1698 }
3a0f6479 1699 }
bcb6fb78 1700 } // second loop over clusters
1701
1702
1703 return kTRUE;
170c35f1 1704}
bcb6fb78 1705///////////////////////////////////////////////////////////////////////////////////////
1706// Pad row col stuff: see if masked or not
1707///////////////////////////////////////////////////////////////////////////////////////
1708//_____________________________________________________________________________
6bbdc11a 1709void AliTRDCalibraFillHisto::CheckGoodTrackletV1(const AliTRDcluster *cl)
d95484e4 1710{
1711 //
1712 // See if we are not near a masked pad
1713 //
1714
1715 if(cl->IsMasked()) fGoodTracklet = kFALSE;
1716
1717
1718}
1719//_____________________________________________________________________________
6bbdc11a 1720void AliTRDCalibraFillHisto::CheckGoodTrackletV0(const Int_t detector,const Int_t row,const Int_t col)
170c35f1 1721{
1722 //
bcb6fb78 1723 // See if we are not near a masked pad
170c35f1 1724 //
1725
bcb6fb78 1726 if (!IsPadOn(detector, col, row)) {
1727 fGoodTracklet = kFALSE;
1728 }
170c35f1 1729
bcb6fb78 1730 if (col > 0) {
1731 if (!IsPadOn(detector, col-1, row)) {
1732 fGoodTracklet = kFALSE;
1733 }
1734 }
170c35f1 1735
bcb6fb78 1736 if (col < 143) {
1737 if (!IsPadOn(detector, col+1, row)) {
1738 fGoodTracklet = kFALSE;
1739 }
1740 }
1741
170c35f1 1742}
bcb6fb78 1743//_____________________________________________________________________________
1744Bool_t AliTRDCalibraFillHisto::IsPadOn(Int_t detector, Int_t row, Int_t col) const
170c35f1 1745{
1746 //
bcb6fb78 1747 // Look in the choosen database if the pad is On.
1748 // If no the track will be "not good"
170c35f1 1749 //
170c35f1 1750
bcb6fb78 1751 // Get the parameter object
1752 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
1753 if (!cal) {
1754 AliInfo("Could not get calibDB");
1755 return kFALSE;
1756 }
1757
1758 if (!cal->IsChamberInstalled(detector) ||
1759 cal->IsChamberMasked(detector) ||
1760 cal->IsPadMasked(detector,col,row)) {
1761 return kFALSE;
1762 }
1763 else {
1764 return kTRUE;
1765 }
1766
170c35f1 1767}
bcb6fb78 1768///////////////////////////////////////////////////////////////////////////////////////
1769// Calibration groups: calculate the number of groups, localise...
1770////////////////////////////////////////////////////////////////////////////////////////
1771//_____________________________________________________________________________
1772Int_t AliTRDCalibraFillHisto::CalculateCalibrationGroup(Int_t i, Int_t row, Int_t col) const
170c35f1 1773{
1774 //
bcb6fb78 1775 // Calculate the calibration group number for i
170c35f1 1776 //
170c35f1 1777
bcb6fb78 1778 // Row of the cluster and position in the pad groups
1779 Int_t posr = 0;
1780 if (fCalibraMode->GetNnZ(i) != 0) {
1781 posr = (Int_t) row / fCalibraMode->GetNnZ(i);
1782 }
1783
1784
1785 // Col of the cluster and position in the pad groups
1786 Int_t posc = 0;
1787 if (fCalibraMode->GetNnRphi(i) != 0) {
1788 posc = (Int_t) col / fCalibraMode->GetNnRphi(i);
1789 }
170c35f1 1790
bcb6fb78 1791 return posc*fCalibraMode->GetNfragZ(i)+posr;
170c35f1 1792
1793}
bcb6fb78 1794//____________________________________________________________________________________
1795Int_t AliTRDCalibraFillHisto::CalculateTotalNumberOfBins(Int_t i)
170c35f1 1796{
3a0f6479 1797 //
bcb6fb78 1798 // Calculate the total number of calibration groups
3a0f6479 1799 //
1800
bcb6fb78 1801 Int_t ntotal = 0;
64942b85 1802
1803 // All together
1804 if((fCalibraMode->GetNz(i)==100) && (fCalibraMode->GetNrphi(i)==100)){
1805 ntotal = 1;
1806 AliInfo(Form("Total number of Xbins: %d for i %d",ntotal,i));
1807 return ntotal;
1808 }
1809
1810 // Per Supermodule
1811 if((fCalibraMode->GetNz(i)==10) && (fCalibraMode->GetNrphi(i)==10)){
1812 ntotal = 18;
1813 AliInfo(Form("Total number of Xbins: %d for i %d",ntotal,i));
1814 return ntotal;
1815 }
1816
1817 // More
bcb6fb78 1818 fCalibraMode->ModePadCalibration(2,i);
1819 fCalibraMode->ModePadFragmentation(0,2,0,i);
1820 fCalibraMode->SetDetChamb2(i);
1821 ntotal += 6 * 18 * fCalibraMode->GetDetChamb2(i);
1822 fCalibraMode->ModePadCalibration(0,i);
1823 fCalibraMode->ModePadFragmentation(0,0,0,i);
1824 fCalibraMode->SetDetChamb0(i);
1825 ntotal += 6 * 4 * 18 * fCalibraMode->GetDetChamb0(i);
1826 AliInfo(Form("Total number of Xbins: %d for i %d",ntotal,i));
1827 return ntotal;
1828
1829}
1830//_____________________________________________________________________________
1831void AliTRDCalibraFillHisto::SetNz(Int_t i, Short_t Nz)
1832{
1833 //
1834 // Set the mode of calibration group in the z direction for the parameter i
1835 //
1836
1837 if ((Nz >= 0) &&
1838 (Nz < 5)) {
1839 fCalibraMode->SetNz(i, Nz);
1840 }
1841 else {
1842 AliInfo("You have to choose between 0 and 4");
55a288e5 1843 }
170c35f1 1844
bcb6fb78 1845}
1846//_____________________________________________________________________________
1847void AliTRDCalibraFillHisto::SetNrphi(Int_t i, Short_t Nrphi)
1848{
1849 //
1850 // Set the mode of calibration group in the rphi direction for the parameter i
1851 //
1852
1853 if ((Nrphi >= 0) &&
1854 (Nrphi < 7)) {
1855 fCalibraMode->SetNrphi(i ,Nrphi);
1856 }
1857 else {
1858 AliInfo("You have to choose between 0 and 6");
3a0f6479 1859 }
55a288e5 1860
bcb6fb78 1861}
64942b85 1862
1863//_____________________________________________________________________________
1864void AliTRDCalibraFillHisto::SetAllTogether(Int_t i)
1865{
1866 //
1867 // Set the mode of calibration group all together
1868 //
1869 if(fVector2d == kTRUE) {
1870 AliInfo("Can not work with the vector method");
1871 return;
1872 }
1873 fCalibraMode->SetAllTogether(i);
1874
1875}
1876
1877//_____________________________________________________________________________
1878void AliTRDCalibraFillHisto::SetPerSuperModule(Int_t i)
1879{
1880 //
1881 // Set the mode of calibration group per supermodule
1882 //
1883 if(fVector2d == kTRUE) {
1884 AliInfo("Can not work with the vector method");
1885 return;
1886 }
1887 fCalibraMode->SetPerSuperModule(i);
1888
1889}
1890
bcb6fb78 1891//____________Set the pad calibration variables for the detector_______________
1892Bool_t AliTRDCalibraFillHisto::LocalisationDetectorXbins(Int_t detector)
1893{
1894 //
1895 // For the detector calcul the first Xbins and set the number of row
1896 // and col pads per calibration groups, the number of calibration
1897 // groups in the detector.
1898 //
1899
1900 // first Xbins of the detector
1901 if (fCH2dOn) {
1902 fCalibraMode->CalculXBins(detector,0);
55a288e5 1903 }
bcb6fb78 1904 if (fPH2dOn) {
1905 fCalibraMode->CalculXBins(detector,1);
55a288e5 1906 }
170c35f1 1907 if (fPRF2dOn) {
bcb6fb78 1908 fCalibraMode->CalculXBins(detector,2);
55a288e5 1909 }
bcb6fb78 1910
1911 // fragmentation of idect
1912 for (Int_t i = 0; i < 3; i++) {
053767a4 1913 fCalibraMode->ModePadCalibration((Int_t) GetStack(detector),i);
1914 fCalibraMode->ModePadFragmentation((Int_t) GetLayer(detector)
1915 , (Int_t) GetStack(detector)
bcb6fb78 1916 , (Int_t) GetSector(detector),i);
170c35f1 1917 }
55a288e5 1918
bcb6fb78 1919 return kTRUE;
1920
55a288e5 1921}
bcb6fb78 1922//_____________________________________________________________________________
1923void AliTRDCalibraFillHisto::SetNumberGroupsPRF(Short_t numberGroupsPRF)
55a288e5 1924{
1925 //
bcb6fb78 1926 // Should be between 0 and 6
55a288e5 1927 //
bcb6fb78 1928
1929 if ((numberGroupsPRF < 0) || (numberGroupsPRF > 6)) {
1930 AliInfo("The number of groups must be between 0 and 6!");
1931 }
1932 else {
1933 fNgroupprf = numberGroupsPRF;
1934 }
55a288e5 1935
bcb6fb78 1936}
1937///////////////////////////////////////////////////////////////////////////////////////////
1938// Per tracklet: store or reset the info, fill the histos with the info
1939//////////////////////////////////////////////////////////////////////////////////////////
1940//_____________________________________________________________________________
1785640c 1941void AliTRDCalibraFillHisto::StoreInfoCHPHtrack(const AliTRDcluster *cl,const Double_t dqdl,const Int_t *group,const Int_t row,const Int_t col,const AliTRDcluster *cls)
bcb6fb78 1942{
1943 //
1944 // Store the infos in fAmpTotal, fPHPlace and fPHValue
1945 // Correct from the gain correction before
1785640c 1946 // cls is shared cluster if any
bcb6fb78 1947 //
1948
09ea3770 1949 //printf("StoreInfoCHPHtrack\n");
1950
bcb6fb78 1951 // time bin of the cluster not corrected
1952 Int_t time = cl->GetPadTime();
e526983e 1953 Float_t charge = TMath::Abs(cl->GetQ());
09ea3770 1954 if(cls) {
1955 charge += TMath::Abs(cls->GetQ());
1956 //printf("AliTRDCalibraFillHisto::Add the cluster charge");
1957 }
e526983e 1958
1959 //printf("Store::time %d, amplitude %f\n",time,dqdl);
1960
bcb6fb78 1961 //Correct for the gain coefficient used in the database for reconstruction
64942b85 1962 Float_t correctthegain = 1.0;
1963 if(fIsHLT) correctthegain = fCalDetGain->GetValue(fDetectorPreviousTrack);
1964 else correctthegain = fCalDetGain->GetValue(fDetectorPreviousTrack)*fCalROCGain->GetValue(col,row);
bcb6fb78 1965 Float_t correction = 1.0;
1966 Float_t normalisation = 6.67;
1967 // we divide with gain in AliTRDclusterizer::Transform...
1968 if( correctthegain > 0 ) normalisation /= correctthegain;
55a288e5 1969
55a288e5 1970
bcb6fb78 1971 // take dd/dl corrected from the angle of the track
1972 correction = dqdl / (normalisation);
1973
55a288e5 1974
bcb6fb78 1975 // Fill the fAmpTotal with the charge
b70c68da 1976 if (fCH2dOn) {
e526983e 1977 if((!fLimitChargeIntegration) || (cl->IsInChamber())) {
1978 //printf("Store::group %d, amplitude %f\n",group[0],correction);
1979 fAmpTotal[(Int_t) group[0]] += correction;
1980 }
bcb6fb78 1981 }
55a288e5 1982
bcb6fb78 1983 // Fill the fPHPlace and value
1984 if (fPH2dOn) {
1985 fPHPlace[time] = group[1];
e526983e 1986 fPHValue[time] = charge;
bcb6fb78 1987 }
1988
1989}
1990//____________Offine tracking in the AliTRDtracker_____________________________
1991void AliTRDCalibraFillHisto::ResetfVariablestracklet()
1992{
1993 //
1994 // Reset values per tracklet
1995 //
1996
1997 //Reset good tracklet
1998 fGoodTracklet = kTRUE;
1999
2000 // Reset the fPHValue
2001 if (fPH2dOn) {
2002 //Reset the fPHValue and fPHPlace
2003 for (Int_t k = 0; k < fTimeMax; k++) {
2004 fPHValue[k] = 0.0;
2005 fPHPlace[k] = -1;
55a288e5 2006 }
bcb6fb78 2007 }
55a288e5 2008
bcb6fb78 2009 // Reset the fAmpTotal where we put value
2010 if (fCH2dOn) {
2011 for (Int_t k = 0; k < fCalibraMode->GetNfragZ(0)*fCalibraMode->GetNfragRphi(0); k++) {
2012 fAmpTotal[k] = 0.0;
55a288e5 2013 }
55a288e5 2014 }
55a288e5 2015}
bcb6fb78 2016//____________Offine tracking in the AliTRDtracker_____________________________
2017void AliTRDCalibraFillHisto::FillTheInfoOfTheTrackCH(Int_t nbclusters)
170c35f1 2018{
2019 //
bcb6fb78 2020 // For the offline tracking or mcm tracklets
2021 // This function will be called in the functions UpdateHistogram...
2022 // to fill the info of a track for the relativ gain calibration
170c35f1 2023 //
bcb6fb78 2024
2025 Int_t nb = 0; // Nombre de zones traversees
2026 Int_t fd = -1; // Premiere zone non nulle
2027 Float_t totalcharge = 0.0; // Total charge for the supermodule histo
170c35f1 2028
e526983e 2029 //printf("CH2d nbclusters %d, fNumberClusters %d, fNumberClustersf %d\n",nbclusters,fNumberClusters,fNumberClustersf);
2030
bcb6fb78 2031 if(nbclusters < fNumberClusters) return;
2032 if(nbclusters > fNumberClustersf) return;
64942b85 2033
2034
2035 // Normalize with the number of clusters
2036 Double_t normalizeCst = fRelativeScale;
2037 if(fNormalizeNbOfCluster) normalizeCst = normalizeCst*nbclusters;
e526983e 2038
2039 //printf("Number of groups in one detector %d\n",fCalibraMode->GetNfragZ(0)*fCalibraMode->GetNfragRphi(0));
bcb6fb78 2040
2041 // See if the track goes through different zones
2042 for (Int_t k = 0; k < fCalibraMode->GetNfragZ(0)*fCalibraMode->GetNfragRphi(0); k++) {
e526983e 2043 //printf("fAmpTotal %f for %d\n",fAmpTotal[k],k);
bcb6fb78 2044 if (fAmpTotal[k] > 0.0) {
2045 totalcharge += fAmpTotal[k];
2046 nb++;
2047 if (nb == 1) {
2048 fd = k;
2049 }
170c35f1 2050 }
170c35f1 2051 }
170c35f1 2052
e526983e 2053 //printf("CH2d: nb %d, fd %d, calibration group %d, amplitude %f, detector %d\n",nb,fd,fCalibraMode->GetXbins(0),fAmpTotal[fd]/normalizeCst,fDetectorPreviousTrack);
bcb6fb78 2054
2055 switch (nb)
2056 {
2057 case 1:
2058 fNumberUsedCh[0]++;
2059 fEntriesCH[fCalibraMode->GetXbins(0)+fd]++;
2060 if (fHisto2d) {
64942b85 2061 FillCH2d(fCalibraMode->GetXbins(0)+fd,fAmpTotal[fd]/normalizeCst);
2062 //fCH2d->Fill(fAmpTotal[fd]/normalizeCst,fCalibraMode->GetXbins(0)+fd+0.5);
bcb6fb78 2063 }
2064 if (fVector2d) {
64942b85 2065 fCalibraVector->UpdateVectorCH(fDetectorPreviousTrack,fd,fAmpTotal[fd]/normalizeCst);
bcb6fb78 2066 }
2067 break;
2068 case 2:
2069 if ((fAmpTotal[fd] > 0.0) &&
2070 (fAmpTotal[fd+1] > 0.0)) {
2071 // One of the two very big
2072 if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+1]) {
2073 if (fHisto2d) {
64942b85 2074 FillCH2d(fCalibraMode->GetXbins(0)+fd,fAmpTotal[fd]/normalizeCst);
2075 //fCH2d->Fill(fAmpTotal[fd]/normalizeCst,fCalibraMode->GetXbins(0)+fd+0.5);
bcb6fb78 2076 }
2077 if (fVector2d) {
64942b85 2078 fCalibraVector->UpdateVectorCH(fDetectorPreviousTrack,fd,fAmpTotal[fd]/normalizeCst);
bcb6fb78 2079 }
2080 fNumberUsedCh[1]++;
2081 fEntriesCH[fCalibraMode->GetXbins(0)+fd]++;
2082 }
2083 if (fAmpTotal[fd+1] > fProcent*fAmpTotal[fd]) {
2084 if (fHisto2d) {
64942b85 2085 FillCH2d(fCalibraMode->GetXbins(0)+fd+1,fAmpTotal[fd+1]/normalizeCst);
2086 //fCH2d->Fill(fAmpTotal[fd+1]/normalizeCst,fCalibraMode->GetXbins(0)+fd+1.5);
bcb6fb78 2087 }
2088 if (fVector2d) {
64942b85 2089 fCalibraVector->UpdateVectorCH(fDetectorPreviousTrack,fd+1,fAmpTotal[fd+1]/normalizeCst);
bcb6fb78 2090 }
2091 fNumberUsedCh[1]++;
2092 fEntriesCH[fCalibraMode->GetXbins(0)+fd+1]++;
2093 }
2094 }
2095 if (fCalibraMode->GetNfragZ(0) > 1) {
2096 if (fAmpTotal[fd] > 0.0) {
2097 if ((fd+fCalibraMode->GetNfragZ(0)) < (fCalibraMode->GetNfragZ(0)*fCalibraMode->GetNfragRphi(0))) {
2098 if (fAmpTotal[fd+fCalibraMode->GetNfragZ(0)] > 0.0) {
2099 // One of the two very big
2100 if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+fCalibraMode->GetNfragZ(0)]) {
2101 if (fHisto2d) {
64942b85 2102 FillCH2d(fCalibraMode->GetXbins(0)+fd,fAmpTotal[fd]/normalizeCst);
2103 //fCH2d->Fill(fAmpTotal[fd]/normalizeCst,fCalibraMode->GetXbins(0)+fd+0.5);
bcb6fb78 2104 }
2105 if (fVector2d) {
64942b85 2106 fCalibraVector->UpdateVectorCH(fDetectorPreviousTrack,fd,fAmpTotal[fd]/normalizeCst);
bcb6fb78 2107 }
2108 fNumberUsedCh[1]++;
2109 fEntriesCH[fCalibraMode->GetXbins(0)+fd]++;
2110 }
2111 if (fAmpTotal[fd+fCalibraMode->GetNfragZ(0)] > fProcent*fAmpTotal[fd]) {
2112 if (fHisto2d) {
64942b85 2113 FillCH2d(fCalibraMode->GetXbins(0)+fd+fCalibraMode->GetNfragZ(0),fAmpTotal[fd+fCalibraMode->GetNfragZ(0)]/normalizeCst);
2114 //fCH2d->Fill(fAmpTotal[fd+fCalibraMode->GetNfragZ(0)]/normalizeCst,fCalibraMode->GetXbins(0)+fd+fCalibraMode->GetNfragZ(0)+0.5);
bcb6fb78 2115 }
2116 fNumberUsedCh[1]++;
2117 fEntriesCH[fCalibraMode->GetXbins(0)+fd+fCalibraMode->GetNfragZ(0)]++;
2118 if (fVector2d) {
64942b85 2119 fCalibraVector->UpdateVectorCH(fDetectorPreviousTrack,fd+fCalibraMode->GetNfragZ(0),fAmpTotal[fd+fCalibraMode->GetNfragZ(0)]/normalizeCst);
bcb6fb78 2120 }
2121 }
2122 }
2123 }
2124 }
2125 }
2126 break;
2127 default: break;
2128 }
170c35f1 2129}
bcb6fb78 2130//____________Offine tracking in the AliTRDtracker_____________________________
2131void AliTRDCalibraFillHisto::FillTheInfoOfTheTrackPH()
170c35f1 2132{
2133 //
bcb6fb78 2134 // For the offline tracking or mcm tracklets
2135 // This function will be called in the functions UpdateHistogram...
2136 // to fill the info of a track for the drift velocity calibration
170c35f1 2137 //
bcb6fb78 2138
2139 Int_t nb = 1; // Nombre de zones traversees 1, 2 ou plus de 3
2140 Int_t fd1 = -1; // Premiere zone non nulle
2141 Int_t fd2 = -1; // Deuxieme zone non nulle
2142 Int_t k1 = -1; // Debut de la premiere zone
2143 Int_t k2 = -1; // Debut de la seconde zone
2144 Int_t nbclusters = 0; // number of clusters
170c35f1 2145
170c35f1 2146
170c35f1 2147
bcb6fb78 2148 // See if the track goes through different zones
2149 for (Int_t k = 0; k < fTimeMax; k++) {
2150 if (fPHValue[k] > 0.0) {
2151 nbclusters++;
2152 if (fd1 == -1) {
2153 fd1 = fPHPlace[k];
2154 k1 = k;
2155 }
2156 if (fPHPlace[k] != fd1) {
2157 if (fd2 == -1) {
2158 k2 = k;
2159 fd2 = fPHPlace[k];
2160 nb = 2;
2161 }
2162 if (fPHPlace[k] != fd2) {
2163 nb = 3;
2164 }
2165 }
55a288e5 2166 }
2167 }
55a288e5 2168
64942b85 2169 // See if noise before and after
2170 if(fMaxCluster > 0) {
2171 if(fPHValue[0] > fMaxCluster) return;
2172 if(fTimeMax > fNbMaxCluster) {
2173 for(Int_t k = (fTimeMax-fNbMaxCluster); k < fTimeMax; k++){
2174 if(fPHValue[k] > fMaxCluster) return;
2175 }
2176 }
2177 }
2178
4aad967c 2179 //printf("nbclusters %d, low limit %d, high limit %d\n",nbclusters,fNumberClusters,fNumberClustersf);
2180
bcb6fb78 2181 if(nbclusters < fNumberClusters) return;
2182 if(nbclusters > fNumberClustersf) return;
3a0f6479 2183
bcb6fb78 2184 switch(nb)
2185 {
170c35f1 2186 case 1:
bcb6fb78 2187 fNumberUsedPh[0]++;
2188 for (Int_t i = 0; i < fTimeMax; i++) {
2189 if (fHisto2d) {
37b0cf5e 2190 if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2191 else {
2192 if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2193 }
bcb6fb78 2194 //printf("Fill the time bin %d with %f\n",i,fPHValue[i]);
2195 }
2196 if (fVector2d) {
37b0cf5e 2197 if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2198 else {
2199 if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2200 }
bcb6fb78 2201 }
55a288e5 2202 }
170c35f1 2203 break;
2204 case 2:
bcb6fb78 2205 if ((fd1 == fd2+1) ||
2206 (fd2 == fd1+1)) {
2207 // One of the two fast all the think
2208 if (k2 > (k1+fDifference)) {
2209 //we choose to fill the fd1 with all the values
2210 fNumberUsedPh[1]++;
2211 for (Int_t i = 0; i < fTimeMax; i++) {
2212 if (fHisto2d) {
37b0cf5e 2213 if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2214 else {
2215 if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2216 }
bcb6fb78 2217 }
2218 if (fVector2d) {
37b0cf5e 2219 if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2220 else {
2221 if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2222 }
bcb6fb78 2223 }
55a288e5 2224 }
55a288e5 2225 }
bcb6fb78 2226 if ((k2+fDifference) < fTimeMax) {
2227 //we choose to fill the fd2 with all the values
2228 fNumberUsedPh[1]++;
2229 for (Int_t i = 0; i < fTimeMax; i++) {
2230 if (fHisto2d) {
37b0cf5e 2231 if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) fPHValue[i]);
2232 else {
2233 if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) fPHValue[i]);
2234 }
bcb6fb78 2235 }
55a288e5 2236 if (fVector2d) {
37b0cf5e 2237 if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd2,i,fPHValue[i]);
2238 else {
2239 if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd2,i,fPHValue[i]);
2240 }
bcb6fb78 2241 }
55a288e5 2242 }
55a288e5 2243 }
2244 }
bcb6fb78 2245 // Two zones voisines sinon rien!
2246 if (fCalibraMode->GetNfragZ(1) > 1) {
2247 // Case 2
2248 if ((fd1+fCalibraMode->GetNfragZ(1)) < (fCalibraMode->GetNfragZ(1)*fCalibraMode->GetNfragRphi(1))) {
2249 if (fd2 == (fd1+fCalibraMode->GetNfragZ(1))) {
2250 // One of the two fast all the think
2251 if (k2 > (k1+fDifference)) {
2252 //we choose to fill the fd1 with all the values
2253 fNumberUsedPh[1]++;
2254 for (Int_t i = 0; i < fTimeMax; i++) {
55a288e5 2255 if (fHisto2d) {
37b0cf5e 2256 if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2257 else {
2258 if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2259 }
55a288e5 2260 }
2261 if (fVector2d) {
37b0cf5e 2262 if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2263 else {
2264 if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2265 }
55a288e5 2266 }
55a288e5 2267 }
bcb6fb78 2268 }
2269 if ((k2+fDifference) < fTimeMax) {
2270 //we choose to fill the fd2 with all the values
2271 fNumberUsedPh[1]++;
2272 for (Int_t i = 0; i < fTimeMax; i++) {
55a288e5 2273 if (fHisto2d) {
37b0cf5e 2274 if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) fPHValue[i]);
2275 else {
2276 if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) fPHValue[i]);
2277 }
170c35f1 2278 }
2279 if (fVector2d) {
37b0cf5e 2280 if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd2,i,fPHValue[i]);
2281 else {
2282 if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd2,i,fPHValue[i]);
2283 }
170c35f1 2284 }
55a288e5 2285 }
2286 }
2287 }
2288 }
170c35f1 2289 // Two zones voisines sinon rien!
2290 // Case 3
2291 if ((fd1 - fCalibraMode->GetNfragZ(1)) >= 0) {
2292 if (fd2 == (fd1 - fCalibraMode->GetNfragZ(1))) {
2293 // One of the two fast all the think
2294 if (k2 > (k1 + fDifference)) {
2295 //we choose to fill the fd1 with all the values
2296 fNumberUsedPh[1]++;
2297 for (Int_t i = 0; i < fTimeMax; i++) {
2298 if (fHisto2d) {
37b0cf5e 2299 if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2300 else {
2301 if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2302 }
170c35f1 2303 }
2304 if (fVector2d) {
37b0cf5e 2305 if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2306 else {
2307 if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2308 }
170c35f1 2309 }
55a288e5 2310 }
2311 }
170c35f1 2312 if ((k2+fDifference) < fTimeMax) {
2313 //we choose to fill the fd2 with all the values
2314 fNumberUsedPh[1]++;
2315 for (Int_t i = 0; i < fTimeMax; i++) {
2316 if (fHisto2d) {
37b0cf5e 2317 if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) fPHValue[i]);
2318 else {
2319 if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) fPHValue[i]);
2320 }
170c35f1 2321 }
2322 if (fVector2d) {
37b0cf5e 2323 if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd2,i,fPHValue[i]);
2324 else {
2325 if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd2,i,fPHValue[i]);
2326 }
170c35f1 2327 }
55a288e5 2328 }
2329 }
2330 }
2331 }
2332 }
170c35f1 2333 break;
2334 default: break;
2335 }
55a288e5 2336}
bcb6fb78 2337//////////////////////////////////////////////////////////////////////////////////////////
2338// DAQ process functions
2339/////////////////////////////////////////////////////////////////////////////////////////
2340//_____________________________________________________________________
b70c68da 2341Int_t AliTRDCalibraFillHisto::ProcessEventDAQ(AliTRDrawStreamBase *rawStream, Bool_t nocheck)
55a288e5 2342{
2343 //
b70c68da 2344 // Event Processing loop - AliTRDrawStreamBase
bcb6fb78 2345 // TestBeam 2007 version
2346 // 0 timebin problem
2347 // 1 no input
2348 // 2 input
bcb6fb78 2349 // Same algorithm as TestBeam but different reader
55a288e5 2350 //
0ca1720e 2351
2352 rawStream->SetSharedPadReadout(kFALSE);
170c35f1 2353
bcb6fb78 2354 Int_t withInput = 1;
2355
2356 Double_t phvalue[16][144][36];
2357 for(Int_t k = 0; k < 36; k++){
2358 for(Int_t j = 0; j < 16; j++){
2359 for(Int_t c = 0; c < 144; c++){
2360 phvalue[j][c][k] = 0.0;
2361 }
170c35f1 2362 }
55a288e5 2363 }
170c35f1 2364
bcb6fb78 2365 fDetectorPreviousTrack = -1;
2366 fMCMPrevious = -1;
2367 fROBPrevious = -1;
2368 Int_t nbtimebin = 0;
d95484e4 2369 Int_t baseline = 10;
1785640c 2370 //printf("------------Detector\n");
2371
bcb6fb78 2372 if(!nocheck){
170c35f1 2373
bcb6fb78 2374 fTimeMax = 0;
2375
2376 while (rawStream->Next()) {
2377
2378 Int_t idetector = rawStream->GetDet(); // current detector
2379 Int_t imcm = rawStream->GetMCM(); // current MCM
2380 Int_t irob = rawStream->GetROB(); // current ROB
d95484e4 2381
2382 //printf("Detector %d\n",idetector);
2383
bcb6fb78 2384 if((fDetectorPreviousTrack != idetector) && (fDetectorPreviousTrack != -1)){
2385
2386 // Fill
2387 withInput = TMath::Max(FillDAQ(phvalue),withInput);
55a288e5 2388
1785640c 2389
bcb6fb78 2390 // reset
2391 for(Int_t k = 0; k < 36; k++){
2392 for(Int_t j = 0; j < 16; j++){
2393 for(Int_t c = 0; c < 144; c++){
2394 phvalue[j][c][k] = 0.0;
2395 }
2396 }
2397 }
2398 }
170c35f1 2399
bcb6fb78 2400 fDetectorPreviousTrack = idetector;
2401 fMCMPrevious = imcm;
2402 fROBPrevious = irob;
170c35f1 2403
bcb6fb78 2404 nbtimebin = rawStream->GetNumberOfTimeBins(); // number of time bins read from data
2405 if(nbtimebin == 0) return 0;
2406 if((fTimeMax != 0) && (nbtimebin != fTimeMax)) return 0;
2407 fTimeMax = nbtimebin;
170c35f1 2408
d95484e4 2409 //baseline = rawStream->GetCommonAdditive(); // common additive baseline
2410 fNumberClustersf = fTimeMax;
6aafa7ea 2411 fNumberClusters = (Int_t)(fNumberClustersProcent*fTimeMax);
d95484e4 2412
2413
bcb6fb78 2414 Int_t *signal = rawStream->GetSignals(); // current ADC signal
2415 Int_t col = rawStream->GetCol();
2416 Int_t row = rawStream->GetRow();
2417
2418
1785640c 2419 // printf("detector %d, signal[0] %d, signal[1] %d, signal[2] %d, baseline %d\n",idetector,signal[0],signal[1],signal[2], baseline);
bcb6fb78 2420
2421
d95484e4 2422 for(Int_t itime = 0; itime < nbtimebin; itime++){
2423 phvalue[row][col][itime] = signal[itime]-baseline;
170c35f1 2424 }
2425 }
bcb6fb78 2426
2427 // fill the last one
2428 if(fDetectorPreviousTrack != -1){
170c35f1 2429
bcb6fb78 2430 // Fill
2431 withInput = TMath::Max(FillDAQ(phvalue),withInput);
2432
2433 // reset
2434 for(Int_t k = 0; k < 36; k++){
2435 for(Int_t j = 0; j < 16; j++){
2436 for(Int_t c = 0; c < 144; c++){
2437 phvalue[j][c][k] = 0.0;
2438 }
2439 }
2440 }
2441 }
2442
2443 }
2444 else{
170c35f1 2445
1785640c 2446 while (rawStream->Next()) { //iddetecte
170c35f1 2447
bcb6fb78 2448 Int_t idetector = rawStream->GetDet(); // current detector
2449 Int_t imcm = rawStream->GetMCM(); // current MCM
2450 Int_t irob = rawStream->GetROB(); // current ROB
170c35f1 2451
d95484e4 2452 //printf("Detector %d\n",idetector);
2453
bcb6fb78 2454 if((fDetectorPreviousTrack != idetector) && (fDetectorPreviousTrack != -1)){
55a288e5 2455
bcb6fb78 2456 // Fill
2457 withInput = TMath::Max(FillDAQ(phvalue),withInput);
2458
2459 // reset
2460 for(Int_t k = 0; k < 36; k++){
2461 for(Int_t j = 0; j < 16; j++){
2462 for(Int_t c = 0; c < 144; c++){
2463 phvalue[j][c][k] = 0.0;
2464 }
2465 }
2466 }
2467 }
170c35f1 2468
bcb6fb78 2469 fDetectorPreviousTrack = idetector;
2470 fMCMPrevious = imcm;
2471 fROBPrevious = irob;
2472
d95484e4 2473 //baseline = rawStream->GetCommonAdditive(); // common baseline
170c35f1 2474
bcb6fb78 2475 fTimeMax = rawStream->GetNumberOfTimeBins(); // number of time bins read from data
d95484e4 2476 fNumberClustersf = fTimeMax;
6aafa7ea 2477 fNumberClusters = (Int_t)(fNumberClustersProcent*fTimeMax);
bcb6fb78 2478 Int_t *signal = rawStream->GetSignals(); // current ADC signal
2479 Int_t col = rawStream->GetCol();
2480 Int_t row = rawStream->GetRow();
2481
d95484e4 2482
bcb6fb78 2483 //printf("detector %d, signal[0] %d, signal[1] %d, signal[2] %d, baseline %d\n",idetector,signal[0],signal[1],signal[2], baseline);
170c35f1 2484
d95484e4 2485 for(Int_t itime = 0; itime < fTimeMax; itime++){
2486 phvalue[row][col][itime] = signal[itime]-baseline;
1785640c 2487 /*if(phvalue[row][col][itime] >= 20) {
2488 printf("----------> phvalue[%d][%d][%d] %d baseline %d \n",
2489 row,
2490 col,
2491 itime,
2492 signal[itime],
2493 baseline);
2494 }*/
bcb6fb78 2495 }
170c35f1 2496 }
2497
bcb6fb78 2498 // fill the last one
2499 if(fDetectorPreviousTrack != -1){
2500
2501 // Fill
2502 withInput = TMath::Max(FillDAQ(phvalue),withInput);
2503
2504 // reset
2505 for(Int_t k = 0; k < 36; k++){
2506 for(Int_t j = 0; j < 16; j++){
2507 for(Int_t c = 0; c < 144; c++){
2508 phvalue[j][c][k] = 0.0;
2509 }
170c35f1 2510 }
2511 }
2512 }
bcb6fb78 2513 }
2514
2515 return withInput;
2516
bcb6fb78 2517}
2518//_____________________________________________________________________
2519Int_t AliTRDCalibraFillHisto::ProcessEventDAQ(AliRawReader *rawReader, Bool_t nocheck)
2520{
2521 //
2522 // Event processing loop - AliRawReader
2523 // Testbeam 2007 version
2524 //
2525
b70c68da 2526 AliTRDrawStreamBase rawStream(rawReader);
bcb6fb78 2527
2528 rawReader->Select("TRD");
2529
2530 return ProcessEventDAQ(&rawStream, nocheck);
2531}
2532
2533//_________________________________________________________________________
2534Int_t AliTRDCalibraFillHisto::ProcessEventDAQ(
2535#ifdef ALI_DATE
6bbdc11a 2536 const eventHeaderStruct *event,
bcb6fb78 2537 Bool_t nocheck
2538#else
6bbdc11a 2539 const eventHeaderStruct* /*event*/,
bcb6fb78 2540 Bool_t /*nocheck*/
2541
2542#endif
2543 )
2544{
2545 //
2546 // process date event
2547 // Testbeam 2007 version
2548 //
2549#ifdef ALI_DATE
2550 AliRawReader *rawReader = new AliRawReaderDate((void*)event);
2551 Int_t result=ProcessEventDAQ(rawReader, nocheck);
2552 delete rawReader;
2553 return result;
2554#else
2555 Fatal("AliTRDCalibraFillHisto", "this class was compiled without DATE");
2556 return 0;
2557#endif
2558
bcb6fb78 2559}
1785640c 2560//_____________________________________________________________________
2561Int_t AliTRDCalibraFillHisto::ProcessEventDAQ2(AliRawReader *rawReader)
2562 { //main
2563 //
2564 // Event Processing loop - AliTRDrawFastStream
2565 //
2566 // 0 timebin problem
2567 // 1 no input
2568 // 2 input
2569 // Same algorithm as TestBeam but different reader
2570 //
2571
2572 // AliTRDrawFastStream *rawStream = AliTRDrawFastStream::GetRawStream(rawReader);
2573 AliTRDrawFastStream *rawStream = new AliTRDrawFastStream(rawReader);
09ea3770 2574 rawStream->SetNoErrorWarning();
1785640c 2575 rawStream->SetSharedPadReadout(kFALSE);
2576
2577 AliTRDdigitsManager *digitsManager = new AliTRDdigitsManager(kTRUE);
2578 digitsManager->CreateArrays();
2579
2580 Int_t withInput = 1;
2581
2582 Double_t phvalue[16][144][36];
2583 for(Int_t k = 0; k < 36; k++){
2584 for(Int_t j = 0; j < 16; j++){
2585 for(Int_t c = 0; c < 144; c++){
2586 phvalue[j][c][k] = 0.0;
2587 }
2588 }
2589 }
2590
2591 fDetectorPreviousTrack = -1;
2592 fMCMPrevious = -1;
2593 fROBPrevious = -1;
2594
2595 Int_t nbtimebin = 0;
2596 Int_t baseline = 10;
2597
2598
2599 fTimeMax = 0;
2600
2601 Int_t det = 0;
2602 while ((det = rawStream->NextChamber(digitsManager, NULL, NULL)) >= 0) { //idetector
2603
2604 if (digitsManager->GetIndexes(det)->HasEntry()) {//QA
2605 // printf("there is ADC data on this chamber!\n");
2606
2607 AliTRDarrayADC *digits = (AliTRDarrayADC *) digitsManager->GetDigits(det); //mod
2608 if (digits->HasData()) { //array
2609
2610 AliTRDSignalIndex *indexes = digitsManager->GetIndexes(det);
2611 if (indexes->IsAllocated() == kFALSE) {
2612 AliError("Indexes do not exist!");
2613 break;
2614 }
2615 Int_t iRow = 0;
2616 Int_t iCol = 0;
2617 indexes->ResetCounters();
2618
2619 while (indexes->NextRCIndex(iRow, iCol)) { //column,row
2620 //printf(" det %d \t row %d \t col %d \t digit\n",det,iRow,iCol);
2621 //while (rawStream->Next()) {
2622
2623 Int_t idetector = det; // current detector
2624 //Int_t imcm = rawStream->GetMCM(); // current MCM
2625 //Int_t irob = rawStream->GetROB(); // current ROB
2626
2627
2628 if((fDetectorPreviousTrack != idetector) && (fDetectorPreviousTrack != -1)) {
2629 // Fill
2630 withInput = TMath::Max(FillDAQ(phvalue),withInput);
2631
2632 // reset
2633 for(Int_t k = 0; k < 36; k++){
2634 for(Int_t j = 0; j < 16; j++){
2635 for(Int_t c = 0; c < 144; c++){
2636 phvalue[j][c][k] = 0.0;
2637 }
2638 }
2639 }
2640 }
2641
2642 fDetectorPreviousTrack = idetector;
2643 //fMCMPrevious = imcm;
2644 //fROBPrevious = irob;
2645
2646 // nbtimebin = rawStream->GetNumberOfTimeBins(); // number of time bins read from data
2647 AliTRDdigitsParam *digitParam = (AliTRDdigitsParam *)digitsManager->GetDigitsParam();
a0446ff1 2648 nbtimebin = digitParam->GetNTimeBins(det); // number of time bins read from data
2649 baseline = digitParam->GetADCbaseline(det); // baseline
1785640c 2650
2651 if(nbtimebin == 0) return 0;
2652 if((fTimeMax != 0) && (nbtimebin != fTimeMax)) return 0;
2653 fTimeMax = nbtimebin;
2654
2655 fNumberClustersf = fTimeMax;
2656 fNumberClusters = (Int_t)(fNumberClustersProcent*fTimeMax);
2657
2658
2659 for(Int_t itime = 0; itime < nbtimebin; itime++) {
2660 // phvalue[row][col][itime] = signal[itime]-baseline;
2661 phvalue[iRow][iCol][itime] = (Short_t)(digits->GetData(iRow,iCol,itime) - baseline);
2662 /*if(phvalue[iRow][iCol][itime] >= 20) {
2663 printf("----------> phvalue[%d][%d][%d] %d baseline %d \n",
2664 iRow,
2665 iCol,
2666 itime,
2667 (Short_t)(digits->GetData(iRow,iCol,itime)),
2668 baseline);
2669 }*/
2670 }
2671
2672 }//column,row
2673
2674 // fill the last one
2675 if(fDetectorPreviousTrack != -1){
2676
2677 // Fill
2678 withInput = TMath::Max(FillDAQ(phvalue),withInput);
2679 // printf("\n ---> withinput %d\n\n",withInput);
2680 // reset
2681 for(Int_t k = 0; k < 36; k++){
2682 for(Int_t j = 0; j < 16; j++){
2683 for(Int_t c = 0; c < 144; c++){
2684 phvalue[j][c][k] = 0.0;
2685 }
2686 }
2687 }
2688 }
2689
2690 }//array
2691 }//QA
2692 digitsManager->ClearArrays(det);
2693 }//idetector
2694 delete digitsManager;
2695
2696 delete rawStream;
2697 return withInput;
2698 }//main
2699//_____________________________________________________________________
bcb6fb78 2700//////////////////////////////////////////////////////////////////////////////
2701// Routine inside the DAQ process
2702/////////////////////////////////////////////////////////////////////////////
2703//_______________________________________________________________________
bcb6fb78 2704Int_t AliTRDCalibraFillHisto::FillDAQ(Double_t phvalue[16][144][36]){
2705
2706 //
2707 // Look for the maximum by collapsing over the time
2708 // Sum over four pad col and two pad row
2709 //
2710
2711 Int_t used = 0;
2712
2713
2714 Int_t idect = fDetectorPreviousTrack;
d95484e4 2715 //printf("Enter Detector %d\n",fDetectorPreviousTrack);
bcb6fb78 2716 Double_t sum[36];
2717 for(Int_t tb = 0; tb < 36; tb++){
2718 sum[tb] = 0.0;
2719 }
2720
d95484e4 2721 //fGoodTracklet = kTRUE;
2722 //fDetectorPreviousTrack = 0;
bcb6fb78 2723
2724
2725 ///////////////////////////
2726 // look for maximum
2727 /////////////////////////
2728
2729 Int_t imaxRow = 0;
2730 Int_t imaxCol = 0;
2731 Double_t integralMax = -1;
2732
2733 for (Int_t ir = 1; ir <= 15; ir++)
2734 {
2735 for (Int_t ic = 2; ic <= 142; ic++)
2736 {
2737 Double_t integral = 0;
6aafa7ea 2738 for (Int_t ishiftR = 0; ishiftR < fNumberRowDAQ; ishiftR++)
bcb6fb78 2739 {
6aafa7ea 2740 for (Int_t ishiftC = -fNumberColDAQ; ishiftC < fNumberColDAQ; ishiftC++)
bcb6fb78 2741 {
2742 if (ir + ishiftR >= 1 && ir + ishiftR <= 16 &&
2743 ic + ishiftC >= 1 && ic + ishiftC <= 144)
2744 {
2745
2746 for(Int_t tb = 0; tb< fTimeMax; tb++){
2747 integral += phvalue[ir + ishiftR-1][ic + ishiftC-1][tb];
2748 }// addtb
2749 } //addsignal
2750 } //shiftC
2751 } // shiftR
bcb6fb78 2752 if (integralMax < integral)
2753 {
2754 imaxRow = ir;
2755 imaxCol = ic;
2756 integralMax = integral;
1785640c 2757
bcb6fb78 2758 } // check max integral
2759 } //ic
2760 } // ir
2761
1785640c 2762 // printf("imaxRow %d, imaxCol %d, fTimeMax %d, integralMax %f\n",imaxRow,imaxCol,fTimeMax, integralMax);
6aafa7ea 2763 //if((imaxRow == 0) || (imaxRow >= 15) || (imaxCol <= 3) || (imaxCol >= 140)) {
2764 // used=1;
2765 // return used;
2766 // }
1785640c 2767
6aafa7ea 2768 if(((imaxRow + fNumberRowDAQ) > 16) || (imaxRow == 0) || ((imaxCol - fNumberColDAQ) <= 1) || ((imaxCol + fNumberColDAQ) >= 144)) {
27b46d95 2769 used=1;
2770 return used;
2771 }
d95484e4 2772 //CheckGoodTrackletV0(fDetectorPreviousTrack,imaxRow,imaxCol);
2773 //if(!fGoodTracklet) used = 1;;
bcb6fb78 2774
2775 // /////////////////////////////////////////////////////
2776 // sum ober 2 row and 4 pad cols for each time bins
2777 // ////////////////////////////////////////////////////
2778
2779
1785640c 2780
6aafa7ea 2781 for (Int_t ishiftR = 0; ishiftR < fNumberRowDAQ; ishiftR++)
bcb6fb78 2782 {
6aafa7ea 2783 for (Int_t ishiftC = -fNumberColDAQ; ishiftC < fNumberColDAQ; ishiftC++)
bcb6fb78 2784 {
6aafa7ea 2785 if (imaxRow + ishiftR >= 1 && imaxRow + ishiftR <= 16 &&
2786 imaxCol + ishiftC >= 1 && imaxCol + ishiftC <= 144)
2787 {
2788 for(Int_t it = 0; it < fTimeMax; it++){
2789 sum[it] += phvalue[imaxRow + ishiftR-1][imaxCol + ishiftC-1][it];
2790 }
2791 }
2792 } // col shift
2793 }// row shift
bcb6fb78 2794
2795 Int_t nbcl = 0;
2796 Double_t sumcharge = 0.0;
2797 for(Int_t it = 0; it < fTimeMax; it++){
2798 sumcharge += sum[it];
1785640c 2799 if(sum[it] > fThresholdClustersDAQ) nbcl++;
bcb6fb78 2800 }
2801
2802
2803 /////////////////////////////////////////////////////////
2804 // Debug
2805 ////////////////////////////////////////////////////////
2806 if(fDebugLevel > 0){
2807 if ( !fDebugStreamer ) {
2808 //debug stream
2809 TDirectory *backup = gDirectory;
2810 fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root");
2811 if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer
2812 }
2813
2814 Double_t amph0 = sum[0];
2815 Double_t amphlast = sum[fTimeMax-1];
2816 Double_t rms = TMath::RMS(fTimeMax,sum);
2817 Int_t goodtracklet = (Int_t) fGoodTracklet;
2818 for(Int_t it = 0; it < fTimeMax; it++){
2819 Double_t clustera = sum[it];
2820
2821 (* fDebugStreamer) << "FillDAQa"<<
2822 "ampTotal="<<sumcharge<<
2823 "row="<<imaxRow<<
2824 "col="<<imaxCol<<
2825 "detector="<<idect<<
2826 "amph0="<<amph0<<
2827 "amphlast="<<amphlast<<
2828 "goodtracklet="<<goodtracklet<<
2829 "clustera="<<clustera<<
2830 "it="<<it<<
2831 "rms="<<rms<<
6aafa7ea 2832 "nbcl="<<nbcl<<
bcb6fb78 2833 "\n";
2834 }
2835 }
2836
2837 ////////////////////////////////////////////////////////
2838 // fill
2839 ///////////////////////////////////////////////////////
6aafa7ea 2840 //printf("fNumberClusters %d, fNumberClustersf %d\n",fNumberClusters,fNumberClustersf);
bcb6fb78 2841 if(sum[0] > 100.0) used = 1;
2842 if(nbcl < fNumberClusters) used = 1;
2843 if(nbcl > fNumberClustersf) used = 1;
2844
2845 //if(fDetectorPreviousTrack == 15){
2846 // printf("rms %f and first time bin %f\n",TMath::RMS(fTimeMax,sum),sum[0]);
2847 //}
2848 //if((TMath::RMS(fTimeMax,sum) <= 10.0) && (sum[0] > 200.0)) return 1;
2849 if(used == 0){
2850 for(Int_t it = 0; it < fTimeMax; it++){
37b0cf5e 2851 if(fFillWithZero) UpdateDAQ(fDetectorPreviousTrack,0,0,it,sum[it],fTimeMax);
2852 else{
1785640c 2853 if(sum[it] > 0.0) UpdateDAQ(fDetectorPreviousTrack,0,0,it,sum[it],fTimeMax);
37b0cf5e 2854 }
1785640c 2855 //if(fFillWithZero) UpdateDAQ(0,0,0,it,sum[it],fTimeMax);
6aafa7ea 2856 //else{
1785640c 2857 // if(sum[it] > 0.0) UpdateDAQ(0,0,0,it,sum[it],fTimeMax);
2858 //}
bcb6fb78 2859 }
2860
2861
2862 //((TH2I *)GetCH2d()->Fill(sumcharge/30.0,fDetectorPreviousTrack));
2863 used = 2;
d95484e4 2864 //printf("Pass Detector %d\n",fDetectorPreviousTrack);
2865
bcb6fb78 2866 }
2867
2868 return used;
2869
2870}
2871//____________Online trackling in AliTRDtrigger________________________________
2872Bool_t AliTRDCalibraFillHisto::UpdateDAQ(Int_t det, Int_t /*row*/, Int_t /*col*/, Int_t timebin, Float_t signal, Int_t nbtimebins)
2873{
2874 //
2875 // For the DAQ
2876 // Fill a simple average pulse height
2877 //
2878
2879
2880 ((TProfile2D *)GetPH2d(nbtimebins,fSf))->Fill((Float_t) timebin/fSf,det+0.5,(Float_t) signal);
2881
2882
2883 return kTRUE;
2884
2885}
2886//____________Write_____________________________________________________
2887//_____________________________________________________________________
2888void AliTRDCalibraFillHisto::Write2d(const Char_t *filename, Bool_t append)
2889{
2890 //
2891 // Write infos to file
2892 //
2893
2894 //For debugging
2895 if ( fDebugStreamer ) {
2896 delete fDebugStreamer;
2897 fDebugStreamer = 0x0;
2898 }
2899
2900 AliInfo(Form("Numbertrack: %d Numberusedch[0]: %d, Numberusedch[1]: %d Numberusedph[0]: %d, Numberusedph[1]: %d"
2901 ,fNumberTrack
2902 ,fNumberUsedCh[0]
2903 ,fNumberUsedCh[1]
2904 ,fNumberUsedPh[0]
2905 ,fNumberUsedPh[1]));
2906
2907 TDirectory *backup = gDirectory;
2908 TString option;
2909
2910 if ( append )
2911 option = "update";
2912 else
2913 option = "recreate";
2914
2915 TFile f(filename,option.Data());
2916
2917 TStopwatch stopwatch;
2918 stopwatch.Start();
2919 if(fVector2d) {
2920 f.WriteTObject(fCalibraVector);
2921 }
2922
2923 if (fCH2dOn ) {
2924 if (fHisto2d) {
2925 f.WriteTObject(fCH2d);
2926 }
2927 }
2928 if (fPH2dOn ) {
2929 if (fHisto2d) {
2930 f.WriteTObject(fPH2d);
2931 }
2932 }
2933 if (fPRF2dOn) {
2934 if (fHisto2d) {
2935 f.WriteTObject(fPRF2d);
2936 }
2937 }
2938 if(fLinearFitterOn){
46c4871e 2939 if(fLinearFitterDebugOn) AnalyseLinearFitter();
bcb6fb78 2940 f.WriteTObject(fLinearVdriftFit);
2941 }
2942
2943 f.Close();
2944
2945 if ( backup ) backup->cd();
2946
2947 AliInfo(Form("Execution time Write2d: R:%.2fs C:%.2fs"
2948 ,stopwatch.RealTime(),stopwatch.CpuTime()));
2949}
2950//////////////////////////////////////////////////////////////////////////////////////////////////////////////
2951// Stats stuff
2952//////////////////////////////////////////////////////////////////////////////////////////////////////////////
2953//___________________________________________probe the histos__________________________________________________
2954Double_t *AliTRDCalibraFillHisto::StatH(TH2 *h, Int_t i)
2955{
2956 //
2957 // Check the number of stats in h, 0 is TH2I 1 is TProfile2D
2958 // debug mode with 2 for TH2I and 3 for TProfile2D
2959 // It gives a pointer to a Double_t[7] with the info following...
2960 // [0] : number of calibration groups with entries
2961 // [1] : minimal number of entries found
2962 // [2] : calibration group number of the min
2963 // [3] : maximal number of entries found
2964 // [4] : calibration group number of the max
2965 // [5] : mean number of entries found
2966 // [6] : mean relative error
2967 //
2968
2969 Double_t *info = new Double_t[7];
2970
2971 // Number of Xbins (detectors or groups of pads)
2972 Int_t nbins = h->GetNbinsY(); //number of calibration groups
2973 Int_t nxbins = h->GetNbinsX(); //number of bins per histo
2974
2975 // Initialise
2976 Double_t nbwe = 0; //number of calibration groups with entries
2977 Double_t minentries = 0; //minimal number of entries found
2978 Double_t maxentries = 0; //maximal number of entries found
2979 Double_t placemin = 0; //calibration group number of the min
2980 Double_t placemax = -1; //calibration group number of the max
2981 Double_t meanstats = 0.0; //mean number of entries over the calibration group with at least ome entry
2982 Double_t meanrelativerror = 0.0; //mean relativ error in the TProfile2D
2983
2984 Double_t counter = 0;
2985
2986 //Debug
2987 TH1F *nbEntries = 0x0;//distribution of the number of entries
2988 TH1F *nbEntriesPerGroup = 0x0;//Number of entries per group
2989 TProfile *nbEntriesPerSp = 0x0;//Number of entries for one supermodule
2990
2991 // Beginning of the loop over the calibration groups
2992 for (Int_t idect = 0; idect < nbins; idect++) {
2993
2994 TH1I *projch = (TH1I *) h->ProjectionX("projch",idect+1,idect+1,(Option_t *)"e");
2995 projch->SetDirectory(0);
2996
2997 // Number of entries for this calibration group
2998 Double_t nentries = 0.0;
2999 if((i%2) == 0){
3000 for (Int_t k = 0; k < nxbins; k++) {
3001 nentries += h->GetBinContent(h->GetBin(k+1,idect+1));
3a0f6479 3002 }
bcb6fb78 3003 }
3004 else{
3005 for (Int_t k = 0; k < nxbins; k++) {
3006 nentries += ((TProfile2D *)h)->GetBinEntries(h->GetBin(k+1,idect+1));
3007 if(h->GetBinContent(h->GetBin(k+1,idect+1)) != 0) {
3008 meanrelativerror += (h->GetBinError(h->GetBin(k+1,idect+1))/(TMath::Abs(h->GetBinContent(h->GetBin(k+1,idect+1)))));
3009 counter++;
3010 }
3011 }
3012 }
3013
3014 //Debug
3015 if(i > 1){
3016 if((!((Bool_t)nbEntries)) && (nentries > 0)){
3017 nbEntries = new TH1F("Number of entries","Number of entries"
3018 ,100,(Int_t)nentries/2,nentries*2);
3019 nbEntries->SetDirectory(0);
3020 nbEntriesPerGroup = new TH1F("Number of entries per group","Number of entries per group"
3021 ,nbins,0,nbins);
3022 nbEntriesPerGroup->SetDirectory(0);
3023 nbEntriesPerSp = new TProfile("Number of entries per supermodule","Number of entries per supermodule"
3024 ,(Int_t)(nbins/18),0,(Int_t)(nbins/18));
3025 nbEntriesPerSp->SetDirectory(0);
3026 }
3027 if(nbEntries){
3028 if(nentries > 0) nbEntries->Fill(nentries);
3029 nbEntriesPerGroup->Fill(idect+0.5,nentries);
3030 nbEntriesPerSp->Fill((idect%((Int_t)(nbins/18)))+0.5,nentries);
170c35f1 3031 }
3032 }
bcb6fb78 3033
3034 //min amd max
3035 if(nentries > maxentries){
3036 maxentries = nentries;
3037 placemax = idect;
3038 }
3039 if(idect == 0) {
3040 minentries = nentries;
3041 }
3042 if(nentries < minentries){
3043 minentries = nentries;
3044 placemin = idect;
3045 }
3046 //nbwe
3047 if(nentries > 0) {
3048 nbwe++;
3049 meanstats += nentries;
3050 }
3051 }//calibration groups loop
3052
3053 if(nbwe > 0) meanstats /= nbwe;
3054 if(counter > 0) meanrelativerror /= counter;
3055
3056 AliInfo(Form("There are %f calibration groups with entries",nbwe));
3057 AliInfo(Form("The minimum number of entries is %f for the group %f",minentries,placemin));
3058 AliInfo(Form("The maximum number of entries is %f for the group %f",maxentries,placemax));
3059 AliInfo(Form("The mean number of entries is %f",meanstats));
3060 if((i%2) == 1) AliInfo(Form("The mean relative error is %f",meanrelativerror));
3061
3062 info[0] = nbwe;
3063 info[1] = minentries;
3064 info[2] = placemin;
3065 info[3] = maxentries;
3066 info[4] = placemax;
3067 info[5] = meanstats;
3068 info[6] = meanrelativerror;
3069
3070 if(i > 1){
3071 gStyle->SetPalette(1);
3072 gStyle->SetOptStat(1111);
3073 gStyle->SetPadBorderMode(0);
3074 gStyle->SetCanvasColor(10);
3075 gStyle->SetPadLeftMargin(0.13);
3076 gStyle->SetPadRightMargin(0.01);
3077 TCanvas *stat = new TCanvas("stat","",50,50,600,800);
3078 stat->Divide(2,1);
3079 stat->cd(1);
3080 nbEntries->Draw("");
3081 stat->cd(2);
3082 nbEntriesPerSp->SetStats(0);
3083 nbEntriesPerSp->Draw("");
3084 TCanvas *stat1 = new TCanvas("stat1","",50,50,600,800);
3085 stat1->cd();
3086 nbEntriesPerGroup->SetStats(0);
3087 nbEntriesPerGroup->Draw("");
3088 }
3089
3090 return info;
3091
3092}
3093//____________________________________________________________________________
3094Double_t *AliTRDCalibraFillHisto::GetMeanMedianRMSNumberCH()
3095{
3096 //
3097 // Return a Int_t[4] with:
3098 // 0 Mean number of entries
3099 // 1 median of number of entries
3100 // 2 rms of number of entries
3101 // 3 number of group with entries
3102 //
3103
3104 Double_t *stat = new Double_t[4];
3105 stat[3] = 0.0;
3106
3107 Int_t nbofgroups = CalculateTotalNumberOfBins(0);
3108 Double_t *weight = new Double_t[nbofgroups];
3109 Int_t *nonul = new Int_t[nbofgroups];
3110
3111 for(Int_t k = 0; k < nbofgroups; k++){
3112 if(fEntriesCH[k] > 0) {
3113 weight[k] = 1.0;
3114 nonul[(Int_t)stat[3]] = fEntriesCH[k];
3115 stat[3]++;
3116 }
3117 else weight[k] = 0.0;
170c35f1 3118 }
bcb6fb78 3119 stat[0] = TMath::Mean(nbofgroups,fEntriesCH,weight);
3120 stat[1] = TMath::Median(nbofgroups,fEntriesCH,weight);
3121 stat[2] = TMath::RMS((Int_t)stat[3],nonul);
3122
3123 return stat;
3124
170c35f1 3125}
bcb6fb78 3126//____________________________________________________________________________
3127Double_t *AliTRDCalibraFillHisto::GetMeanMedianRMSNumberLinearFitter() const
55a288e5 3128{
3129 //
bcb6fb78 3130 // Return a Int_t[4] with:
3131 // 0 Mean number of entries
3132 // 1 median of number of entries
3133 // 2 rms of number of entries
3134 // 3 number of group with entries
55a288e5 3135 //
3136
bcb6fb78 3137 Double_t *stat = new Double_t[4];
3138 stat[3] = 0.0;
55a288e5 3139
bcb6fb78 3140 Int_t nbofgroups = 540;
3141 Double_t *weight = new Double_t[nbofgroups];
3142 Int_t *nonul = new Int_t[nbofgroups];
55a288e5 3143
bcb6fb78 3144 for(Int_t k = 0; k < nbofgroups; k++){
3145 if(fEntriesLinearFitter[k] > 0) {
3146 weight[k] = 1.0;
3147 nonul[(Int_t) stat[3]] = fEntriesLinearFitter[k];
3148 stat[3]++;
d0569428 3149 }
bcb6fb78 3150 else weight[k] = 0.0;
d0569428 3151 }
bcb6fb78 3152 stat[0] = TMath::Mean(nbofgroups,fEntriesLinearFitter,weight);
3153 stat[1] = TMath::Median(nbofgroups,fEntriesLinearFitter,weight);
3154 stat[2] = TMath::RMS((Int_t)stat[3],nonul);
55a288e5 3155
bcb6fb78 3156 return stat;
170c35f1 3157
bcb6fb78 3158}
3159//////////////////////////////////////////////////////////////////////////////////////
3160// Create Histos
3161//////////////////////////////////////////////////////////////////////////////////////
3162//_____________________________________________________________________________
3163void AliTRDCalibraFillHisto::CreatePRF2d(Int_t nn)
3164{
3165 //
3166 // Create the 2D histos: here we have 2*fNgroupprf bins in tnp of 0.2 amplitude each
3167 // If fNgroupprf is zero then no binning in tnp
3168 //
d0569428 3169
bcb6fb78 3170 TString name("Nz");
3171 name += fCalibraMode->GetNz(2);
3172 name += "Nrphi";
3173 name += fCalibraMode->GetNrphi(2);
3174 name += "Ngp";
3175 name += fNgroupprf;
d0569428 3176
bcb6fb78 3177 if(fNgroupprf != 0){
3178
3179 fPRF2d = new TProfile2D("PRF2d",(const Char_t *) name
3180 ,2*fNgroupprf*fNumberBinPRF,-3.0*fNgroupprf,3.0*fNgroupprf,nn,0,nn);
3181 fPRF2d->SetYTitle("Det/pad groups");
3182 fPRF2d->SetXTitle("Position x/W [pad width units]");
3183 fPRF2d->SetZTitle("Q_{i}/Q_{total}");
3184 fPRF2d->SetStats(0);
d0569428 3185 }
bcb6fb78 3186 else{
3187 fPRF2d = new TProfile2D("PRF2d",(const Char_t *) name
3188 ,fNumberBinPRF,-1.5,1.5,nn,0,nn);
3189 fPRF2d->SetYTitle("Det/pad groups");
3190 fPRF2d->SetXTitle("Position x/W [pad width units]");
3191 fPRF2d->SetZTitle("Q_{i}/Q_{total}");
3192 fPRF2d->SetStats(0);
d0569428 3193 }
d0569428 3194
3195}
bcb6fb78 3196
3197//_____________________________________________________________________________
3198void AliTRDCalibraFillHisto::CreatePH2d(Int_t nn)
d0569428 3199{
3200 //
bcb6fb78 3201 // Create the 2D histos
d0569428 3202 //
3203
bcb6fb78 3204 TString name("Nz");
3205 name += fCalibraMode->GetNz(1);
3206 name += "Nrphi";
3207 name += fCalibraMode->GetNrphi(1);
d0569428 3208
bcb6fb78 3209 fPH2d = new TProfile2D("PH2d",(const Char_t *) name
3210 ,fTimeMax,-0.5/fSf,(Float_t) (fTimeMax-0.5)/fSf
3211 ,nn,0,nn);
3212 fPH2d->SetYTitle("Det/pad groups");
3213 fPH2d->SetXTitle("time [#mus]");
3214 fPH2d->SetZTitle("<PH> [a.u.]");
3215 fPH2d->SetStats(0);
d0569428 3216
bcb6fb78 3217}
3218//_____________________________________________________________________________
3219void AliTRDCalibraFillHisto::CreateCH2d(Int_t nn)
3220{
3221 //
3222 // Create the 2D histos
3223 //
d0569428 3224
bcb6fb78 3225 TString name("Nz");
3226 name += fCalibraMode->GetNz(0);
3227 name += "Nrphi";
3228 name += fCalibraMode->GetNrphi(0);
d0569428 3229
bcb6fb78 3230 fCH2d = new TH2I("CH2d",(const Char_t *) name
3231 ,fNumberBinCharge,0,300,nn,0,nn);
3232 fCH2d->SetYTitle("Det/pad groups");
3233 fCH2d->SetXTitle("charge deposit [a.u]");
3234 fCH2d->SetZTitle("counts");
3235 fCH2d->SetStats(0);
3236 fCH2d->Sumw2();
d0569428 3237
bcb6fb78 3238}
3239//////////////////////////////////////////////////////////////////////////////////
3240// Set relative scale
3241/////////////////////////////////////////////////////////////////////////////////
3242//_____________________________________________________________________________
3243void AliTRDCalibraFillHisto::SetRelativeScale(Float_t RelativeScale)
3244{
3245 //
3246 // Set the factor that will divide the deposited charge
3247 // to fit in the histo range [0,300]
3248 //
3249
3250 if (RelativeScale > 0.0) {
3251 fRelativeScale = RelativeScale;
3252 }
3253 else {
3254 AliInfo("RelativeScale must be strict positif!");
d0569428 3255 }
bcb6fb78 3256
3257}
3258//////////////////////////////////////////////////////////////////////////////////
3259// Quick way to fill a histo
3260//////////////////////////////////////////////////////////////////////////////////
3261//_____________________________________________________________________
3262void AliTRDCalibraFillHisto::FillCH2d(Int_t x, Float_t y)
3263{
3264 //
3265 // FillCH2d: Marian style
3266 //
3267
3268 //skip simply the value out of range
3269 if((y>=300.0) || (y<0.0)) return;
3270
3271 //Calcul the y place
3272 Int_t yplace = (Int_t) (fNumberBinCharge*y/300.0)+1;
3273 Int_t place = (fNumberBinCharge+2)*(x+1)+yplace;
d0569428 3274
bcb6fb78 3275 //Fill
3276 fCH2d->GetArray()[place]++;
3277
d0569428 3278}
bcb6fb78 3279
3280//////////////////////////////////////////////////////////////////////////////////
3281// Geometrical functions
3282///////////////////////////////////////////////////////////////////////////////////
d0569428 3283//_____________________________________________________________________________
053767a4 3284Int_t AliTRDCalibraFillHisto::GetLayer(Int_t d) const
d0569428 3285{
3286 //
053767a4 3287 // Reconstruct the layer number from the detector number
d0569428 3288 //
3289
3290 return ((Int_t) (d % 6));
3291
3292}
3293
3294//_____________________________________________________________________________
053767a4 3295Int_t AliTRDCalibraFillHisto::GetStack(Int_t d) const
d0569428 3296{
3297 //
053767a4 3298 // Reconstruct the stack number from the detector number
d0569428 3299 //
053767a4 3300 const Int_t kNlayer = 6;
d0569428 3301
053767a4 3302 return ((Int_t) (d % 30) / kNlayer);
d0569428 3303
3304}
3305
3306//_____________________________________________________________________________
3307Int_t AliTRDCalibraFillHisto::GetSector(Int_t d) const
3308{
3309 //
3310 // Reconstruct the sector number from the detector number
3311 //
3312 Int_t fg = 30;
3313
3314 return ((Int_t) (d / fg));
3315
3316}
bcb6fb78 3317///////////////////////////////////////////////////////////////////////////////////
3318// Getter functions for DAQ of the CH2d and the PH2d
3319//////////////////////////////////////////////////////////////////////////////////
d0569428 3320//_____________________________________________________________________
3321TProfile2D* AliTRDCalibraFillHisto::GetPH2d(Int_t nbtimebin, Float_t samplefrequency)
3322{
3323 //
3324 // return pointer to fPH2d TProfile2D
3325 // create a new TProfile2D if it doesn't exist allready
3326 //
3327 if ( fPH2d )
3328 return fPH2d;
3329
3330 // Some parameters
3331 fTimeMax = nbtimebin;
3332 fSf = samplefrequency;
3333
170c35f1 3334 CreatePH2d(540);
3335
3336 return fPH2d;
3337}
3338//_____________________________________________________________________
bcb6fb78 3339TH2I* AliTRDCalibraFillHisto::GetCH2d()
3340{
3341 //
3342 // return pointer to fCH2d TH2I
3343 // create a new TH2I if it doesn't exist allready
3344 //
3345 if ( fCH2d )
3346 return fCH2d;
3347
3348 CreateCH2d(540);
3349
3350 return fCH2d;
3351}
3352////////////////////////////////////////////////////////////////////////////////////////////
3353// Drift velocity calibration
3354///////////////////////////////////////////////////////////////////////////////////////////
3355//_____________________________________________________________________
170c35f1 3356TLinearFitter* AliTRDCalibraFillHisto::GetLinearFitter(Int_t detector, Bool_t force)
3357{
3358 //
3359 // return pointer to TLinearFitter Calibration
3360 // if force is true create a new TLinearFitter if it doesn't exist allready
3361 //
170c35f1 3362
d0569428 3363 if ((!force) || (fLinearFitterArray.UncheckedAt(detector))){
3364 return (TLinearFitter*)fLinearFitterArray.UncheckedAt(detector);
3365 }
3366
3367 // if we are forced and TLinearFitter doesn't yet exist create it
170c35f1 3368
d0569428 3369 // new TLinearFitter
3370 TLinearFitter *linearfitter = new TLinearFitter(2,"pol1");
3371 fLinearFitterArray.AddAt(linearfitter,detector);
3372 return linearfitter;
170c35f1 3373}
170c35f1 3374
3a0f6479 3375//____________________________________________________________________________
3376void AliTRDCalibraFillHisto::AnalyseLinearFitter()
3377{
3378 //
3379 // Analyse array of linear fitter because can not be written
3380 // Store two arrays: one with the param the other one with the error param + number of entries
3381 //
3382
3383 for(Int_t k = 0; k < 540; k++){
3384 TLinearFitter *linearfitter = GetLinearFitter(k);
3385 if((linearfitter!=0) && (fEntriesLinearFitter[k]>10)){
3386 TVectorD *par = new TVectorD(2);
3387 TVectorD pare = TVectorD(2);
3388 TVectorD *parE = new TVectorD(3);
3389 linearfitter->Eval();
3390 linearfitter->GetParameters(*par);
3391 linearfitter->GetErrors(pare);
3392 Float_t ppointError = TMath::Sqrt(TMath::Abs(linearfitter->GetChisquare())/fEntriesLinearFitter[k]);
3393 (*parE)[0] = pare[0]*ppointError;
3394 (*parE)[1] = pare[1]*ppointError;
3395 (*parE)[2] = (Double_t) fEntriesLinearFitter[k];
3396 ((TObjArray *)fLinearVdriftFit->GetPArray())->AddAt(par,k);
3397 ((TObjArray *)fLinearVdriftFit->GetEArray())->AddAt(parE,k);
3398 }
3399 }
1f3c89d8 3400}
1785640c 3401
3402