Removing the tasks from the digitization (Ruben)
[u/mrichter/AliRoot.git] / MUON / AliMUONDigitizerV3.cxx
CommitLineData
92aeef15 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**************************************************************************/
15
16// $Id$
17
9edefa04 18
92aeef15 19#include "AliMUONDigitizerV3.h"
20
4ce497c4 21#include "AliMUON.h"
92aeef15 22#include "AliMUONCalibrationData.h"
92aeef15 23#include "AliMUONConstants.h"
92aeef15 24#include "AliMUONDigit.h"
ad26d4f6 25#include "AliMUONLogger.h"
4f8a3d8b 26#include "AliMUONTriggerElectronics.h"
40e382ae 27#include "AliMUONTriggerStoreV1.h"
05314992 28#include "AliMUONVCalibParam.h"
40e382ae 29#include "AliMUONVDigitStore.h"
331a617a 30#include "AliMUONGeometryTransformer.h" //ADDED for trigger noise
05315e71 31#include "AliMUONRecoParam.h"
d5315275 32#include "AliMUONTriggerChamberEfficiency.h"
33#include "AliMUONTriggerUtilities.h"
331a617a 34
88544f7e 35#include "AliMpCDB.h"
331a617a 36#include "AliMpSegmentation.h"
40e382ae 37#include "AliMpCathodType.h"
38#include "AliMpConstants.h"
4ce497c4 39#include "AliMpDEIterator.h"
92aeef15 40#include "AliMpDEManager.h"
4ce497c4 41#include "AliMpPad.h"
40e382ae 42#include "AliMpStationType.h"
4ce497c4 43#include "AliMpVSegmentation.h"
d5315275 44#include "AliMpDDLStore.h"
331a617a 45
46#include "AliCDBManager.h"
47#include "AliCodeTimer.h"
48#include "AliLog.h"
92aeef15 49#include "AliRun.h"
f21fc003 50#include "AliDigitizationInput.h"
c93255fe 51#include "AliLoader.h"
92aeef15 52#include "AliRunLoader.h"
331a617a 53
866c3232 54#include <Riostream.h>
55#include <TF1.h>
40e382ae 56#include <TFile.h>
866c3232 57#include <TMath.h>
58#include <TRandom.h>
59#include <TString.h>
40e382ae 60#include <TSystem.h>
c93255fe 61#include <TTree.h>
40e382ae 62
3d1463c8 63//-----------------------------------------------------------------------------
9265505b 64/// \class AliMUONDigitizerV3
05315e71 65///
4ce497c4 66/// The digitizer is performing the transformation to go from SDigits (digits
67/// w/o any electronic noise) to Digits (w/ electronic noise, and decalibration)
68///
69/// The decalibration is performed by doing the reverse operation of the
70/// calibration, that is we do (Signal+pedestal)/gain -> ADC
71///
72/// Note also that the digitizer takes care of merging sdigits that belongs
73/// to the same pad, either because we're merging several input sdigit files
74/// or with a single file because the sdigitizer does not merge sdigits itself
75/// (for performance reason mainly, and because anyway we know we have to do it
76/// here, at the digitization level).
77///
05315e71 78/// August 2011. In order to remove the need for specific MC OCDB storages,
79/// we're introducing a dependence of simulation on AliMUONRecoParam (stored
80/// in MUON/Calib/RecoParam in OCDB), which is normally (or conceptually, if
81/// you will) only a reconstruction object. That's not a pretty solution, but,
82/// well, we have to do it...
83/// This dependence comes from the fact that we must know how to decalibrate
84/// the digits, so that the decalibration (done here) - calibration (done during
85/// reco) process is (as much as possible) neutral.
86///
87///
c4ee792d 88/// \author Laurent Aphecetche
05315e71 89///
3d1463c8 90//-----------------------------------------------------------------------------
4ce497c4 91
92namespace
93{
94 AliMUON* muon()
95 {
96 return static_cast<AliMUON*>(gAlice->GetModule("MUON"));
97 }
8c0b5e70 98
99 //ADDED for trigger noise
100 const AliMUONGeometryTransformer* GetTransformer()
101 {
102 return muon()->GetGeometryTransformer();
103 }
4ce497c4 104}
105
5c083cba 106Double_t AliMUONDigitizerV3::fgNSigmas = 4.0;
4ce497c4 107
9265505b 108/// \cond CLASSIMP
92aeef15 109ClassImp(AliMUONDigitizerV3)
9265505b 110/// \endcond
92aeef15 111
112//_____________________________________________________________________________
f21fc003 113AliMUONDigitizerV3::AliMUONDigitizerV3(AliDigitizationInput* digInput,
8c0b5e70 114 Int_t generateNoisyDigits)
f21fc003 115: AliDigitizer(digInput),
92aeef15 116fIsInitialized(kFALSE),
92aeef15 117fCalibrationData(0x0),
118fTriggerProcessor(0x0),
8c0b5e70 119fNoiseFunctionTrig(0x0),
ec9acc85 120fGenerateNoisyDigits(generateNoisyDigits),
121fLogger(new AliMUONLogger(1000)),
40e382ae 122fTriggerStore(new AliMUONTriggerStoreV1),
123fDigitStore(0x0),
9ec6a948 124fOutputDigitStore(0x0),
d5315275 125fInputDigitStores(0x0),
05315e71 126fRecoParam(0x0),
d5315275 127fTriggerEfficiency(0x0),
128fTriggerUtilities(0x0),
129fEfficiencyResponse(2*AliMUONConstants::NTriggerCh()*AliMUONConstants::NTriggerCircuit())
92aeef15 130{
9265505b 131 /// Ctor.
132
f21fc003 133 AliDebug(1,Form("AliDigitizationInput=%p",fDigInput));
ec9acc85 134
92aeef15 135}
136
137//_____________________________________________________________________________
138AliMUONDigitizerV3::~AliMUONDigitizerV3()
139{
9265505b 140 /// Dtor. Note we're the owner of some pointers.
141
92aeef15 142 AliDebug(1,"dtor");
4ce497c4 143
92c23b09 144 // delete fCalibrationData;
92aeef15 145 delete fTriggerProcessor;
8c0b5e70 146 delete fNoiseFunctionTrig;
40e382ae 147 delete fTriggerStore;
148 delete fDigitStore;
149 delete fOutputDigitStore;
ca8c8223 150 delete fInputDigitStores;
d5315275 151 delete fTriggerUtilities;
9ec6a948 152
ad26d4f6 153 AliInfo("Summary of messages");
154 fLogger->Print();
4ce497c4 155
ad26d4f6 156 delete fLogger;
92aeef15 157}
158
159//_____________________________________________________________________________
160void
40e382ae 161AliMUONDigitizerV3::ApplyResponseToTrackerDigit(AliMUONVDigit& digit, Bool_t addNoise)
92aeef15 162{
9265505b 163 /// For tracking digits, starting from an ideal digit's charge, we :
164 ///
cf27231a 165 /// - "divide" by a gain (thus decalibrating the digit)
9265505b 166 /// - add a pedestal (thus decalibrating the digit)
b74ca120 167 /// - add some electronics noise (thus leading to a realistic adc), if requested to do so
9265505b 168 /// - sets the signal to zero if below 3*sigma of the noise
59b91539 169
05315e71 170 Float_t charge = digit.Charge();
171
172 if (!digit.IsChargeInFC())
173 {
174 charge *= AliMUONConstants::DefaultADC2MV()*AliMUONConstants::DefaultA0()*AliMUONConstants::DefaultCapa();
175 fLogger->Log("CHECK ME ! WAS NOT SUPPOSED TO BE HERE !!! ARE YOU RECONSTRUCTING OLD SIMULATIONS ? ");
176 AliError("CHECK ME ! WAS NOT SUPPOSED TO BE HERE !!! ARE YOU RECONSTRUCTING OLD SIMULATIONS ? ");
177 }
cf27231a 178
179 // We set the charge to 0, as the only relevant piece of information
180 // after Digitization is the ADC value.
181 digit.SetCharge(0);
b74ca120 182
4ce497c4 183 Int_t detElemId = digit.DetElemId();
92aeef15 184 Int_t manuId = digit.ManuId();
92aeef15 185
4ce497c4 186 AliMUONVCalibParam* pedestal = fCalibrationData->Pedestals(detElemId,manuId);
92aeef15 187 if (!pedestal)
ad26d4f6 188 {
189 fLogger->Log(Form("%s:%d:Could not get pedestal for DE=%4d manuId=%4d. Disabling.",
190 __FILE__,__LINE__,
191 detElemId,manuId));
ad26d4f6 192 digit.SetADC(0);
193 return;
194 }
59b91539 195
627ba19f 196 Int_t manuChannel = digit.ManuChannel();
197
198 if ( pedestal->ValueAsFloat(manuChannel,0) == AliMUONVCalibParam::InvalidFloatValue() ||
199 pedestal->ValueAsFloat(manuChannel,1) == AliMUONVCalibParam::InvalidFloatValue() )
200 {
201 // protection against invalid pedestal value
202 digit.SetADC(0);
203 return;
204 }
205
4ce497c4 206 AliMUONVCalibParam* gain = fCalibrationData->Gains(detElemId,manuId);
92aeef15 207 if (!gain)
ad26d4f6 208 {
209 fLogger->Log(Form("%s:%d:Could not get gain for DE=%4d manuId=%4d. Disabling.",
210 __FILE__,__LINE__,
211 detElemId,manuId));
ad26d4f6 212 digit.SetADC(0);
213 return;
214 }
cf27231a 215
ec9acc85 216 Int_t adc = DecalibrateTrackerDigit(*pedestal,*gain,manuChannel,charge,addNoise,
05315e71 217 digit.IsNoiseOnly(),
218 fRecoParam->GetCalibrationMode());
59b91539 219
92aeef15 220 digit.SetADC(adc);
221}
222
d5315275 223
224//_____________________________________________________________________________
225void
226AliMUONDigitizerV3::ApplyResponseToTriggerDigit(AliMUONVDigit& digit)
227{
228 /// For trigger digits, starting from an ideal digit, we :
229 ///
230 /// - apply efficiency (on demand)
231 /// - apply trigger masks
232
233 Int_t detElemId = digit.DetElemId();
234 Int_t localCircuit = digit.ManuId();
235 Int_t strip = digit.ManuChannel();
236 Int_t cathode = digit.Cathode();
237 Int_t trigCh = detElemId/100 - 11;
238
239 Int_t arrayIndex = GetArrayIndex(cathode, trigCh, localCircuit);
240
241 // Trigger chamber efficiency
242 if ( fTriggerEfficiency ) {
243 if ( fEfficiencyResponse[arrayIndex] < 0 ) {
244 Bool_t isTrig[2] = {kTRUE, kTRUE};
245 fTriggerEfficiency->IsTriggered(detElemId, localCircuit, isTrig[0], isTrig[1]);
246 Int_t arrayIndexBend = GetArrayIndex(0, trigCh, localCircuit);
247 Int_t arrayIndexNonBend = GetArrayIndex(1, trigCh, localCircuit);
248 fEfficiencyResponse[arrayIndexBend] = isTrig[0];
249 fEfficiencyResponse[arrayIndexNonBend] = isTrig[1];
250 }
251 AliDebug(1,Form("ch %i cath %i board %i strip %i efficiency %i\n", trigCh, cathode, localCircuit, strip, fEfficiencyResponse[arrayIndex]));
252 if ( fEfficiencyResponse[arrayIndex] == 0 ) {
253 digit.SetCharge(0);
254 digit.SetADC(0);
255 //AliDebug(1,Form("ch %i cath %i board %i strip %i NOT efficient\n", trigCh, cathode, localCircuit, strip));
256 return;
257 }
258 }
259
260 // Masked channels
261 Bool_t isMasked = fTriggerUtilities->IsMasked(digit);
262 AliDebug(1,Form("ch %i cath %i board %i strip %i mask %i\n", trigCh, cathode, localCircuit, strip, !isMasked));
263 if ( isMasked ) {
264 digit.SetCharge(0);
265 digit.SetADC(0);
266 //AliDebug(1,Form("ch %i cath %i board %i strip %i masked\n", trigCh, cathode, localCircuit, strip));
267 return;
268 }
269}
270
271
272
92aeef15 273//_____________________________________________________________________________
274void
40e382ae 275AliMUONDigitizerV3::ApplyResponse(const AliMUONVDigitStore& store,
276 AliMUONVDigitStore& filteredStore)
92aeef15 277{
9265505b 278 /// Loop over all chamber digits, and apply the response to them
279 /// Note that this method may remove digits.
280
40e382ae 281 filteredStore.Clear();
282
4ce497c4 283 const Bool_t kAddNoise = kTRUE;
284
40e382ae 285 TIter next(store.CreateIterator());
286 AliMUONVDigit* digit;
287
d5315275 288 if ( fTriggerEfficiency ) fEfficiencyResponse.Reset(-1);
289
40e382ae 290 while ( ( digit = static_cast<AliMUONVDigit*>(next()) ) )
4ce497c4 291 {
40e382ae 292 AliMp::StationType stationType = AliMpDEManager::GetStationType(digit->DetElemId());
293
294 if ( stationType != AliMp::kStationTrigger )
92aeef15 295 {
ac6411cf 296 Bool_t addNoise = kAddNoise;
297 if (digit->IsConverted()) addNoise = kFALSE; // No need to add extra noise to a converted real digit
298 ApplyResponseToTrackerDigit(*digit,addNoise);
92aeef15 299 }
d5315275 300 else {
301 ApplyResponseToTriggerDigit(*digit);
302 }
7d7d22a6 303
cf27231a 304 if ( digit->ADC() > 0 || digit->Charge() > 0 )
40e382ae 305 {
306 filteredStore.Add(*digit,AliMUONVDigitStore::kIgnore);
92aeef15 307 }
92aeef15 308 }
40e382ae 309}
92aeef15 310
311//_____________________________________________________________________________
b74ca120 312Int_t
313AliMUONDigitizerV3::DecalibrateTrackerDigit(const AliMUONVCalibParam& pedestals,
314 const AliMUONVCalibParam& gains,
315 Int_t channel,
316 Float_t charge,
ec9acc85 317 Bool_t addNoise,
05315e71 318 Bool_t noiseOnly,
319 const TString& calibrationMode)
b74ca120 320{
321 /// Decalibrate (i.e. go from charge to adc) a tracker digit, given its
322 /// pedestal and gain parameters.
323 /// Must insure before calling that channel is valid (i.e. between 0 and
324 /// pedestals or gains->GetSize()-1, but also corresponding to a valid channel
325 /// otherwise results are not predictible...)
05315e71 326 ///
327 /// This method is completely tied to what happens in its sister method :
328 /// AliMUONDigitCalibrator::CalibrateDigit, which is doing the reverse work...
329 ///
330
b74ca120 331 static const Int_t kMaxADC = (1<<12)-1; // We code the charge on a 12 bits ADC.
332
05315e71 333 Bool_t nogain = calibrationMode.Contains("NOGAIN");
334
335 Float_t a1(0.0);
336 Int_t thres(4095);
337 Int_t qual(0xF);
338 Float_t capa(AliMUONConstants::DefaultCapa()); // capa = 0.2 and a0 = 1.25
339 Float_t a0(AliMUONConstants::DefaultA0()); // is equivalent to gain = 4 mV/fC
340 Float_t adc2mv(AliMUONConstants::DefaultADC2MV()); // 1 ADC channel = 0.61 mV
341
342 if ( ! nogain )
343 {
344 a0 = gains.ValueAsFloat(channel,0);
345 a1 = gains.ValueAsFloat(channel,1);
346 thres = gains.ValueAsInt(channel,2);
347 qual = gains.ValueAsInt(channel,3);
348 }
349
b74ca120 350 Float_t pedestalMean = pedestals.ValueAsFloat(channel,0);
351 Float_t pedestalSigma = pedestals.ValueAsFloat(channel,1);
352
ec9acc85 353 AliDebugClass(1,Form("DE %04d MANU %04d CH %02d PEDMEAN %7.2f PEDSIGMA %7.2f",
5d950f54 354 pedestals.ID0(),pedestals.ID1(),channel,pedestalMean,pedestalSigma));
ec9acc85 355
b74ca120 356 if ( qual <= 0 ) return 0;
357
358 Float_t chargeThres = a0*thres;
359
360 Float_t padc(0); // (adc - ped) value
361
05315e71 362 if ( nogain || charge <= chargeThres || TMath::Abs(a1) < 1E-12 )
b74ca120 363 {
364 // linear part only
365
366 if ( TMath::Abs(a0) > 1E-12 )
367 {
368 padc = charge/a0;
369 }
370 }
371 else
372 {
05315e71 373 // FIXME: when we'll use capacitances and real gains, revise this part
374 // to take capa properly into account...
375
376 AliWarningClass("YOU PROBABLY NEED TO REVISE THIS PART OF CODE !!!");
377
b74ca120 378 // linear + parabolic part
379 Double_t qt = chargeThres - charge;
380 Double_t delta = a0*a0-4*a1*qt;
381 if ( delta < 0 )
382 {
383 AliErrorClass(Form("delta=%e DE %d Manu %d Channel %d "
384 " charge %e a0 %e a1 %e thres %d ped %e pedsig %e",
385 delta,pedestals.ID0(),pedestals.ID1(),
386 channel, charge, a0, a1, thres, pedestalMean,
387 pedestalSigma));
388 }
389 else
390 {
391 delta = TMath::Sqrt(delta);
392
393 padc = ( ( -a0 + delta ) > 0 ? ( -a0 + delta ) : ( -a0 - delta ) );
394
395 padc /= 2*a1;
396
397 if ( padc < 0 )
398 {
399 if ( TMath::Abs(padc) > 1E-3)
400 {
401 // this is more than a precision problem : let's signal it !
402 AliErrorClass(Form("padc=%e DE %d Manu %d Channel %d "
403 " charge %e a0 %e a1 %e thres %d ped %e pedsig %e delta %e",
404 padc,pedestals.ID0(),pedestals.ID1(),
405 channel, charge, a0, a1, thres, pedestalMean,
406 pedestalSigma,delta));
407 }
408
409 // ok. consider we're just at thres, let it be zero.
410 padc = 0;
411 }
412
413 padc += thres;
414
415 }
416 }
417
05315e71 418 padc /= capa*adc2mv;
419
b74ca120 420 Int_t adc(0);
421
ec9acc85 422 Float_t adcNoise = 0.0;
b74ca120 423
ec9acc85 424 if ( addNoise )
425 {
426 if ( noiseOnly )
427 {
428 adcNoise = NoiseFunction()->GetRandom()*pedestalSigma;
429 }
430 else
431 {
432 adcNoise = gRandom->Gaus(0.0,pedestalSigma);
433 }
434 }
b74ca120 435
ec9acc85 436 adc = TMath::Nint(padc + pedestalMean + adcNoise + 0.5);
437
5c083cba 438 if ( adc < TMath::Nint(pedestalMean + fgNSigmas*pedestalSigma + 0.5) )
ec9acc85 439 {
d537c8e4 440 // this is an error only in specific cases
441 if ( !addNoise || (addNoise && noiseOnly) )
442 {
1b54215c 443 AliDebugClass(1,Form(" DE %04d Manu %04d Channel %02d "
444 " a0 %7.2f a1 %7.2f thres %04d ped %7.2f pedsig %7.2f adcNoise %7.2f "
445 " charge=%7.2f padc=%7.2f adc=%04d ZS=%04d fgNSigmas=%e addNoise %d noiseOnly %d ",
446 pedestals.ID0(),pedestals.ID1(),channel,
447 a0, a1, thres, pedestalMean, pedestalSigma, adcNoise,
448 charge, padc, adc,
449 TMath::Nint(pedestalMean + fgNSigmas*pedestalSigma + 0.5),
450 fgNSigmas,addNoise,noiseOnly));
d537c8e4 451 }
452
ec9acc85 453 adc = 0;
b74ca120 454 }
455
456 // be sure we stick to 12 bits.
457 if ( adc > kMaxADC )
458 {
459 adc = kMaxADC;
460 }
461
462 return adc;
463}
464
465//_____________________________________________________________________________
92aeef15 466void
ca8c8223 467AliMUONDigitizerV3::CreateInputDigitStores()
9ec6a948 468{
ca8c8223 469 /// Create input digit stores
470 ///
9ec6a948 471
ca8c8223 472 if (fInputDigitStores)
473 {
474 AliFatal("Should be called only once !");
475 }
476
477 fInputDigitStores = new TObjArray;
478
479 fInputDigitStores->SetOwner(kTRUE);
9ec6a948 480
f21fc003 481 for ( Int_t iFile = 0; iFile < fDigInput->GetNinputs(); ++iFile )
9ec6a948 482 {
f21fc003 483 AliLoader* inputLoader = GetLoader(fDigInput->GetInputFolderName(iFile));
9ec6a948 484
485 inputLoader->LoadSDigits("READ");
486
487 TTree* iTreeS = inputLoader->TreeS();
488 if (!iTreeS)
489 {
490 AliFatal(Form("Could not get access to input file #%d",iFile));
491 }
492
ca8c8223 493 fInputDigitStores->AddAt(AliMUONVDigitStore::Create(*iTreeS),iFile);
9ec6a948 494 }
495}
496
497//_____________________________________________________________________________
498void
f21fc003 499AliMUONDigitizerV3::Digitize(Option_t*)
92aeef15 500{
9265505b 501 /// Main method.
502 /// We first loop over input files, and merge the sdigits we found there.
503 /// Second, we digitize all the resulting sdigits
504 /// Then we generate noise-only digits (for tracker only)
505 /// And we finally generate the trigger outputs.
4ce497c4 506
99c136e1 507 AliCodeTimerAuto("",0)
b74ca120 508
f21fc003 509 if ( fDigInput->GetNinputs() == 0 )
92aeef15 510 {
511 AliWarning("No input set. Nothing to do.");
512 return;
513 }
514
515 if ( !fIsInitialized )
516 {
517 AliError("Not initialized. Cannot perform the work. Sorry");
518 return;
519 }
520
f21fc003 521 Int_t nInputFiles = fDigInput->GetNinputs();
92aeef15 522
f21fc003 523 AliLoader* outputLoader = GetLoader(fDigInput->GetOutputFolderName());
40e382ae 524
525 outputLoader->MakeDigitsContainer();
526
527 TTree* oTreeD = outputLoader->TreeD();
528
529 if (!oTreeD)
92aeef15 530 {
40e382ae 531 AliFatal("Cannot create output TreeD");
92aeef15 532 }
40e382ae 533
92aeef15 534 // Loop over all the input files, and merge the sdigits found in those
535 // files.
40e382ae 536
92aeef15 537 for ( Int_t iFile = 0; iFile < nInputFiles; ++iFile )
5d950f54 538 {
f21fc003 539 AliLoader* inputLoader = GetLoader(fDigInput->GetInputFolderName(iFile));
40e382ae 540
541 inputLoader->LoadSDigits("READ");
542
543 TTree* iTreeS = inputLoader->TreeS();
544 if (!iTreeS)
92aeef15 545 {
546 AliFatal(Form("Could not get access to input file #%d",iFile));
547 }
9ec6a948 548
ca8c8223 549 if (!fInputDigitStores)
9ec6a948 550 {
ca8c8223 551 CreateInputDigitStores();
9ec6a948 552 }
40e382ae 553
ca8c8223 554 AliMUONVDigitStore* dstore = static_cast<AliMUONVDigitStore*>(fInputDigitStores->At(iFile));
555
556 dstore->Connect(*iTreeS);
40e382ae 557
ca8c8223 558 iTreeS->GetEvent(0);
559
f21fc003 560 MergeWithSDigits(fDigitStore,*dstore,fDigInput->GetMask(iFile));
05314992 561
40e382ae 562 inputLoader->UnloadSDigits();
05314992 563
ca8c8223 564 dstore->Clear();
92aeef15 565 }
5d950f54 566
92aeef15 567
568 // At this point, we do have digit arrays (one per chamber) which contains
569 // the merging of all the sdigits of the input file(s).
570 // We now massage them to apply the detector response, i.e. this
571 // is here that we do the "digitization" work.
572
40e382ae 573 if (!fOutputDigitStore)
574 {
575 fOutputDigitStore = fDigitStore->Create();
576 }
577
8c0b5e70 578 if ( fGenerateNoisyDigits>=2 )
579 {
580 // Generate noise-only digits for trigger.
581 GenerateNoisyDigitsForTrigger(*fDigitStore);
582 }
40e382ae 583 ApplyResponse(*fDigitStore,*fOutputDigitStore);
5d950f54 584
4ce497c4 585 if ( fGenerateNoisyDigits )
586 {
587 // Generate noise-only digits for tracker.
40e382ae 588 GenerateNoisyDigits(*fOutputDigitStore);
4ce497c4 589 }
590
92aeef15 591 // We generate the global and local trigger decisions.
40e382ae 592 fTriggerProcessor->Digits2Trigger(*fOutputDigitStore,*fTriggerStore);
593
594 // Prepare output tree
595 Bool_t okD = fOutputDigitStore->Connect(*oTreeD,kFALSE);
596 Bool_t okT = fTriggerStore->Connect(*oTreeD,kFALSE);
597 if (!okD || !okT)
598 {
599 AliError(Form("Could not make branch : Digit %d Trigger %d",okD,okT));
600 return;
601 }
92aeef15 602
603 // Fill the output treeD
40e382ae 604 oTreeD->Fill();
92aeef15 605
606 // Write to the output tree(D).
607 // Please note that as GlobalTrigger, LocalTrigger and Digits are in the same
608 // tree (=TreeD) in different branches, this WriteDigits in fact writes all of
609 // the 3 branches.
40e382ae 610 outputLoader->WriteDigits("OVERWRITE");
92aeef15 611
40e382ae 612 outputLoader->UnloadDigits();
4ce497c4 613
40e382ae 614 // Finally, we clean up after ourselves.
615 fTriggerStore->Clear();
616 fDigitStore->Clear();
617 fOutputDigitStore->Clear();
92aeef15 618}
619
4ce497c4 620
621//_____________________________________________________________________________
4ce497c4 622void
40e382ae 623AliMUONDigitizerV3::GenerateNoisyDigits(AliMUONVDigitStore& digitStore)
4ce497c4 624{
9265505b 625 /// According to a given probability, generate digits that
626 /// have a signal above the noise cut (ped+n*sigma_ped), i.e. digits
627 /// that are "only noise".
4ce497c4 628
99c136e1 629 AliCodeTimerAuto("",0)
b74ca120 630
4ce497c4 631 for ( Int_t i = 0; i < AliMUONConstants::NTrackingCh(); ++i )
632 {
633 AliMpDEIterator it;
634
635 it.First(i);
636
637 while ( !it.IsDone() )
638 {
639 for ( Int_t cathode = 0; cathode < 2; ++cathode )
640 {
40e382ae 641 GenerateNoisyDigitsForOneCathode(digitStore,it.CurrentDEId(),cathode);
4ce497c4 642 }
643 it.Next();
644 }
645 }
4ce497c4 646}
647
648//_____________________________________________________________________________
649void
40e382ae 650AliMUONDigitizerV3::GenerateNoisyDigitsForOneCathode(AliMUONVDigitStore& digitStore,
651 Int_t detElemId, Int_t cathode)
4ce497c4 652{
9265505b 653 /// Generate noise-only digits for one cathode of one detection element.
654 /// Called by GenerateNoisyDigits()
4ce497c4 655
9265505b 656 const AliMpVSegmentation* seg
866c3232 657 = AliMpSegmentation::Instance()->GetMpSegmentation(detElemId,AliMp::GetCathodType(cathode));
4ce497c4 658 Int_t nofPads = seg->NofPads();
659
660 Int_t maxIx = seg->MaxPadIndexX();
661 Int_t maxIy = seg->MaxPadIndexY();
662
5c083cba 663 static const Double_t kProbToBeOutsideNsigmas = TMath::Erfc(fgNSigmas/TMath::Sqrt(2.0)) / 2. ;
4ce497c4 664
84aac932 665 Int_t nofNoisyPads = TMath::Nint(kProbToBeOutsideNsigmas*nofPads);
4ce497c4 666 if ( !nofNoisyPads ) return;
667
668 nofNoisyPads =
669 TMath::Nint(gRandom->Gaus(nofNoisyPads,
670 nofNoisyPads/TMath::Sqrt(nofNoisyPads)));
671
672 AliDebug(3,Form("DE %d cath %d nofNoisyPads %d",detElemId,cathode,nofNoisyPads));
673
40e382ae 674 for ( Int_t i = 0; i < nofNoisyPads; ++i )
4ce497c4 675 {
676 Int_t ix(-1);
677 Int_t iy(-1);
40e382ae 678 AliMpPad pad;
679
4ce497c4 680 do {
681 ix = gRandom->Integer(maxIx+1);
682 iy = gRandom->Integer(maxIy+1);
168e9c4d 683 pad = seg->PadByIndices(ix,iy,kFALSE);
40e382ae 684 } while ( !pad.IsValid() );
685
168e9c4d 686 Int_t manuId = pad.GetManuId();
687 Int_t manuChannel = pad.GetManuChannel();
40e382ae 688
4ce497c4 689 AliMUONVCalibParam* pedestals = fCalibrationData->Pedestals(detElemId,manuId);
690
ad26d4f6 691 if (!pedestals)
692 {
693 // no pedestal available for this channel, simply give up
40e382ae 694 continue;
ad26d4f6 695 }
696
40e382ae 697 AliMUONVDigit* d = digitStore.CreateDigit(detElemId,manuId,manuChannel,cathode);
698
699 d->SetPadXY(ix,iy);
700
ec9acc85 701 d->SetCharge(0.0); // charge is zero, the ApplyResponseToTrackerDigit will add the noise
40e382ae 702 d->NoiseOnly(kTRUE);
ec9acc85 703 ApplyResponseToTrackerDigit(*d,kTRUE);
cf27231a 704 if ( d->ADC() > 0 )
4ce497c4 705 {
40e382ae 706 Bool_t ok = digitStore.Add(*d,AliMUONVDigitStore::kDeny);
707 // this can happen (that we randomly chose a digit that is
708 // already there). We simply ignore this, but log the occurence
709 // to cross-check that it's not too frequent.
710 if (!ok)
711 {
712 fLogger->Log("Collision while adding noiseOnly digit");
713 }
714 else
715 {
716 fLogger->Log("Added noiseOnly digit");
717 }
4ce497c4 718 }
40e382ae 719 delete d;
4ce497c4 720 }
721}
722
8c0b5e70 723
724//_____________________________________________________________________________
725void
726AliMUONDigitizerV3::GenerateNoisyDigitsForTrigger(AliMUONVDigitStore& digitStore)
727{
728 /// Generate noise-only digits for one cathode of one detection element.
729 /// Called by GenerateNoisyDigits()
730
731 if ( !fNoiseFunctionTrig )
732 {
733 fNoiseFunctionTrig = new TF1("AliMUONDigitizerV3::fNoiseFunctionTrig","landau",
734 50.,270.);
735
736 fNoiseFunctionTrig->SetParameters(3.91070e+02, 9.85026, 9.35881e-02);
737 }
738
739 AliMpPad pad[2];
740 AliMUONVDigit *d[2]={0x0};
741
742 for ( Int_t chamberId = AliMUONConstants::NTrackingCh(); chamberId < AliMUONConstants::NCh(); ++chamberId )
743 {
744
745 Int_t nofNoisyPads = 50;
746
747 Float_t r=-1, fi = 0., gx, gy, x, y, z, xg01, yg01, zg, xg02, yg02;
748 AliMpDEIterator it;
749
750 AliDebug(3,Form("Chamber %d nofNoisyPads %d",chamberId,nofNoisyPads));
751
752 for ( Int_t i = 0; i < nofNoisyPads; ++i )
753 {
754 //printf("Generating noise %i\n",i);
755 Int_t ix(-1);
756 Int_t iy(-1);
757 Bool_t isOk = kFALSE;
758 Int_t detElemId = -1;
759 do {
760 //r = gRandom->Landau(9.85026, 9.35881e-02);
761 r = fNoiseFunctionTrig->GetRandom();
762 fi = 2. * TMath::Pi() * gRandom->Rndm();
763 //printf("r = %f\tfi = %f\n", r, fi);
764 gx = r * TMath::Cos(fi);
765 gy = r * TMath::Sin(fi);
766
767 for ( it.First(chamberId); ! it.IsDone(); it.Next() ){
768 Int_t currDetElemId = it.CurrentDEId();
769 const AliMpVSegmentation* seg
770 = AliMpSegmentation::Instance()->GetMpSegmentation(currDetElemId,AliMp::GetCathodType(0));
771 if (!seg) continue;
6e97fbb8 772 Float_t deltax = seg->GetDimensionX();
773 Float_t deltay = seg->GetDimensionY();
8c0b5e70 774 GetTransformer()->Local2Global(currDetElemId, -deltax, -deltay, 0, xg01, yg01, zg);
775 GetTransformer()->Local2Global(currDetElemId, deltax, deltay, 0, xg02, yg02, zg);
776 Float_t xg1 = xg01, xg2 = xg02, yg1 = yg01, yg2 = yg02;
777 if(xg01>xg02){
778 xg1 = xg02;
779 xg2 = xg01;
780 }
781 if(yg01>yg02){
782 yg1 = yg02;
783 yg2 = yg01;
784 }
785 if(gx>=xg1 && gx<=xg2 && gy>=yg1 && gy<=yg2){
786 detElemId = currDetElemId;
787 GetTransformer()->Global2Local(detElemId, gx, gy, 0, x, y, z);
6e97fbb8 788 pad[0] = seg->PadByPosition(x,y,kFALSE);
8c0b5e70 789 if(!pad[0].IsValid()) continue;
790 isOk = kTRUE;
791 break;
792 }
793 } // loop on slats
794 } while ( !isOk );
795
796 const AliMpVSegmentation* seg1
797 = AliMpSegmentation::Instance()->GetMpSegmentation(detElemId,AliMp::GetCathodType(1));
6e97fbb8 798 pad[1] = seg1->PadByPosition(x,y,kFALSE);
8c0b5e70 799
800 for ( Int_t cathode = 0; cathode < 2; ++cathode ){
168e9c4d 801 Int_t manuId = pad[cathode].GetLocalBoardId(0);
802 Int_t manuChannel = pad[cathode].GetLocalBoardChannel(0);
8c0b5e70 803 d[cathode] = digitStore.CreateDigit(detElemId,manuId,manuChannel,cathode);
168e9c4d 804 ix = pad[cathode].GetIx();
805 iy = pad[cathode].GetIy();
8c0b5e70 806 d[cathode]->SetPadXY(ix,iy);
807 //d[cathode].SetSignal(1);
808 //d[cathode].SetPhysicsSignal(0);
809 d[cathode]->SetCharge(1);
810 d[cathode]->NoiseOnly(kTRUE);
811 AliDebug(3,Form("Adding a pure noise digit :"));
812
813 Bool_t ok = digitStore.Add(*d[cathode],AliMUONVDigitStore::kDeny);
814 if (!ok)
815 {
816 fLogger->Log("Collision while adding TriggerNoise digit");
817 }
818 else
819 {
820 fLogger->Log("Added triggerNoise digit");
821 }
822 } //loop on cathodes
823 } // loop on noisy pads
824 } // loop on chambers
825}
826
827
4ce497c4 828//_____________________________________________________________________________
40e382ae 829AliLoader*
830AliMUONDigitizerV3::GetLoader(const TString& folderName)
92aeef15 831{
40e382ae 832 /// Get a MUON loader
9265505b 833
ad26d4f6 834 AliDebug(2,Form("Getting access to folder %s",folderName.Data()));
40e382ae 835 AliLoader* loader = AliRunLoader::GetDetectorLoader("MUON",folderName.Data());
92aeef15 836 if (!loader)
837 {
838 AliError(Form("Could not get MuonLoader from folder %s",folderName.Data()));
839 return 0x0;
840 }
40e382ae 841 return loader;
92aeef15 842}
843
844//_____________________________________________________________________________
845Bool_t
846AliMUONDigitizerV3::Init()
847{
f21fc003 848 /// Initialization of the digitization :
40e382ae 849 /// a) create the calibrationData, according to run number
850 /// b) create the trigger processing task
9265505b 851
ad26d4f6 852 AliDebug(2,"");
92aeef15 853
854 if ( fIsInitialized )
855 {
856 AliError("Object already initialized.");
857 return kFALSE;
858 }
859
f21fc003 860 if (!fDigInput)
92aeef15 861 {
f21fc003 862 AliError("fDigInput is null !");
92aeef15 863 return kFALSE;
864 }
865
331a617a 866 // Load mapping
867 if ( ! AliMpCDB::LoadDDLStore() ) {
88544f7e 868 AliFatal("Could not access mapping from OCDB !");
869 }
870
92c23b09 871 if (!fCalibrationData)
872 AliFatal("Calibration data object not defined");
873
ad26d4f6 874 if ( !fCalibrationData->Pedestals() )
875 {
876 AliFatal("Could not access pedestals from OCDB !");
877 }
878 if ( !fCalibrationData->Gains() )
879 {
880 AliFatal("Could not access gains from OCDB !");
881 }
92c23b09 882
883
884 AliInfo("Using trigger configuration from CDB");
885
40e382ae 886 fTriggerProcessor = new AliMUONTriggerElectronics(fCalibrationData);
4ce497c4 887
ad26d4f6 888 AliDebug(1, Form("Will %s generate noise-only digits for tracker",
889 (fGenerateNoisyDigits ? "":"NOT")));
d5315275 890
891 fTriggerUtilities = new AliMUONTriggerUtilities(fCalibrationData);
892
893 if ( muon()->GetTriggerEffCells() ) {
894 // Apply trigger efficiency
895 AliDebug(1, "Will apply trigger efficiency");
896 fTriggerEfficiency = new AliMUONTriggerChamberEfficiency(fCalibrationData->TriggerEfficiency());
897 }
ad26d4f6 898
92aeef15 899 fIsInitialized = kTRUE;
900 return kTRUE;
901}
902
d5315275 903
904//_____________________________________________________________________________
905Int_t AliMUONDigitizerV3::GetArrayIndex(Int_t cathode, Int_t trigCh, Int_t localCircuit)
906{
907 /// Get index of array with trigger status map or efficiency
908 return
909 AliMUONConstants::NTriggerCircuit() * AliMUONConstants::NTriggerCh() * cathode +
910 AliMUONConstants::NTriggerCircuit() * trigCh + localCircuit-1;
911}
912
913
92aeef15 914//_____________________________________________________________________________
92aeef15 915void
40e382ae 916AliMUONDigitizerV3::MergeWithSDigits(AliMUONVDigitStore*& outputStore,
917 const AliMUONVDigitStore& input,
918 Int_t mask)
92aeef15 919{
9265505b 920 /// Merge the sdigits in inputData with the digits already present in outputData
92aeef15 921
40e382ae 922 if ( !outputStore ) outputStore = input.Create();
923
924 TIter next(input.CreateIterator());
925 AliMUONVDigit* sdigit;
926
927 while ( ( sdigit = static_cast<AliMUONVDigit*>(next()) ) )
928 {
929 // Update the track references using the mask.
930 // FIXME: this is dirty, for backward compatibility only.
931 // Should re-design all this way of keeping track of MC information...
932 if ( mask ) sdigit->PatchTracks(mask);
933 // Then add or update the digit to the output.
934 AliMUONVDigit* added = outputStore->Add(*sdigit,AliMUONVDigitStore::kMerge);
935 if (!added)
92aeef15 936 {
40e382ae 937 AliError("Could not add digit in merge mode");
92aeef15 938 }
92aeef15 939 }
940}
ec9acc85 941
942//_____________________________________________________________________________
943TF1*
944AliMUONDigitizerV3::NoiseFunction()
945{
946 /// Return noise function
947 static TF1* f = 0x0;
948 if (!f)
949 {
5c083cba 950 f = new TF1("AliMUONDigitizerV3::NoiseFunction","gaus",fgNSigmas,fgNSigmas*10);
ec9acc85 951 f->SetParameters(1,0,1);
952 }
953 return f;
954}
955
05315e71 956//_____________________________________________________________________________
957void AliMUONDigitizerV3::SetCalibrationData(AliMUONCalibrationData* calibrationData,
958 AliMUONRecoParam* recoParam)
959{
960 fCalibrationData = calibrationData;
961 fRecoParam = recoParam;
962 if (!fRecoParam)
963 {
964 AliError("Cannot work (e.g. decalibrate) without recoparams !");
965 }
966}
967