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