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