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