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