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