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