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