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