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