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