CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
SiStripDigitizerAlgorithm.cc
Go to the documentation of this file.
1 // File: SiStripDigitizerAlgorithm.cc
2 // Description: Steering class for digitization.
3 // Modified 15/May/2013 mark.grimes@bristol.ac.uk - Modified so that the digi-sim link has the correct
4 // index for the sim hits stored. It was previously always set to zero (I won't mention that it was
5 // me who originally wrote that).
6 // Modified on Feb 11, 2015: prolay.kumar.mal@cern.ch & Jean-Laurent.Agram@cern.ch
7 // Added/Modified the individual strip noise in zero suppression
8 // mode from the conditions DB; previously, the digitizer used to
9 // consider the noise value for individual strips inside a module from
10 // the central strip noise value.
12 #include <vector>
13 #include <algorithm>
14 #include <iostream>
31 #include "CLHEP/Random/RandFlat.h"
32 
33 #include <cstring>
34 #include <sstream>
35 
36 #include <boost/algorithm/string.hpp>
37 
39  lorentzAngleName(conf.getParameter<std::string>("LorentzAngle")),
40  theThreshold(conf.getParameter<double>("NoiseSigmaThreshold")),
41  cmnRMStib(conf.getParameter<double>("cmnRMStib")),
42  cmnRMStob(conf.getParameter<double>("cmnRMStob")),
43  cmnRMStid(conf.getParameter<double>("cmnRMStid")),
44  cmnRMStec(conf.getParameter<double>("cmnRMStec")),
45  APVSaturationProbScaling_(conf.getParameter<double>("APVSaturationProbScaling")),
46  makeDigiSimLinks_(conf.getUntrackedParameter<bool>("makeDigiSimLinks", false)),
47  peakMode(conf.getParameter<bool>("APVpeakmode")),
48  noise(conf.getParameter<bool>("Noise")),
49  RealPedestals(conf.getParameter<bool>("RealPedestals")),
50  SingleStripNoise(conf.getParameter<bool>("SingleStripNoise")),
51  CommonModeNoise(conf.getParameter<bool>("CommonModeNoise")),
52  BaselineShift(conf.getParameter<bool>("BaselineShift")),
53  APVSaturationFromHIP(conf.getParameter<bool>("APVSaturationFromHIP")),
54  theFedAlgo(conf.getParameter<int>("FedAlgorithm")),
55  zeroSuppression(conf.getParameter<bool>("ZeroSuppression")),
56  theElectronPerADC(conf.getParameter<double>( peakMode ? "electronPerAdcPeak" : "electronPerAdcDec" )),
57  theTOFCutForPeak(conf.getParameter<double>("TOFCutForPeak")),
58  theTOFCutForDeconvolution(conf.getParameter<double>("TOFCutForDeconvolution")),
59  tofCut(peakMode ? theTOFCutForPeak : theTOFCutForDeconvolution),
60  cosmicShift(conf.getUntrackedParameter<double>("CosmicDelayShift")),
61  inefficiency(conf.getParameter<double>("Inefficiency")),
62  pedOffset((unsigned int)conf.getParameter<double>("PedestalsOffset")),
63  PreMixing_(conf.getParameter<bool>("PreMixingMode")),
64  theSiHitDigitizer(new SiHitDigitizer(conf)),
65  theSiPileUpSignals(new SiPileUpSignals()),
66  theSiNoiseAdder(new SiGaussianTailNoiseAdder(theThreshold)),
67  theSiDigitalConverter(new SiTrivialDigitalConverter(theElectronPerADC, PreMixing_)),
68  theSiZeroSuppress(new SiStripFedZeroSuppression(theFedAlgo)),
69  APVProbabilityFile(conf.getParameter<edm::FileInPath>("APVProbabilityFile")) {
70 
71  if (peakMode) {
72  LogDebug("StripDigiInfo")<<"APVs running in peak mode (poor time resolution)";
73  } else {
74  LogDebug("StripDigiInfo")<<"APVs running in deconvolution mode (good time resolution)";
75  };
76  if(SingleStripNoise) LogDebug("SiStripDigitizerAlgorithm")<<" SingleStripNoise: ON";
77  else LogDebug("SiStripDigitizerAlgorithm")<<" SingleStripNoise: OFF";
78  if(CommonModeNoise) LogDebug("SiStripDigitizerAlgorithm")<<" CommonModeNoise: ON";
79  else LogDebug("SiStripDigitizerAlgorithm")<<" CommonModeNoise: OFF";
80  if(PreMixing_ && APVSaturationFromHIP) throw cms::Exception("PreMixing does not work with HIP loss simulation yet");
83  APVProbaFile.open((APVProbabilityFile.fullPath()).c_str());
84  if (APVProbaFile.is_open()){
85  while ( getline (APVProbaFile,line) ){
86  std::vector<std::string> strs;
87  boost::split(strs,line,boost::is_any_of(" "));
88  if(strs.size()==2){
89  mapOfAPVprobabilities[std::stoi(strs.at(0))]=std::stof(strs.at(1));
90  }
91  }
92  APVProbaFile.close();
93  }else throw cms::Exception("MissingInput")
94  << "It seems that the APV probability list is missing\n";
95  }
96 }
97 
99 }
100 
101 void
103  edm::ESHandle<SiStripBadStrip> deadChannelHandle;
104  iSetup.get<SiStripBadChannelRcd>().get(deadChannelHandle);
105 
106  unsigned int detId = det->geographicalId().rawId();
107  int numStrips = (det->specificTopology()).nstrips();
108 
109  SiStripBadStrip::Range detBadStripRange = deadChannelHandle->getRange(detId);
110  //storing the bad strip of the the module. the module is not removed but just signal put to 0
111  std::vector<bool>& badChannels = allBadChannels[detId];
112  badChannels.clear();
113  badChannels.insert(badChannels.begin(), numStrips, false);
114  for(SiStripBadStrip::ContainerIterator it = detBadStripRange.first; it != detBadStripRange.second; ++it) {
115  SiStripBadStrip::data fs = deadChannelHandle->decode(*it);
116  for(int strip = fs.firstStrip; strip < fs.firstStrip + fs.range; ++strip) {
117  badChannels[strip] = true;
118  }
119  }
120  firstChannelsWithSignal[detId] = numStrips;
121  lastChannelsWithSignal[detId]= 0;
122 
123  // if(APVSaturationFromHIP){
124  // std::bitset<6> &bs=SiStripTrackerAffectedAPVMap[detId];
125  // if(bs.any())theAffectedAPVvector.push_back(std::make_pair(detId,bs));
126  //}
127 }
128 
129 void
131  theSiPileUpSignals->reset();
132  // This should be clear by after all calls to digitize(), but I might as well make sure
133  associationInfoForDetId_.clear();
134 
135  APVSaturationProb_ = APVSaturationProbScaling_; // reset probability
137  FirstLumiCalc_ = true;
138  FirstDigitize_ = true;
139 
140  //get gain noise pedestal lorentzAngle from ES handle
142  iSetup.getData(pdt);
143  setParticleDataTable(&*pdt);
145 }
146 
147 // Run the algorithm for a given module
148 // ------------------------------------
149 
150 void
151 SiStripDigitizerAlgorithm::accumulateSimHits(std::vector<PSimHit>::const_iterator inputBegin,
152  std::vector<PSimHit>::const_iterator inputEnd,
153  size_t inputBeginGlobalIndex,
154  unsigned int tofBin,
155  const StripGeomDetUnit* det,
156  const GlobalVector& bfield,
157  const TrackerTopology *tTopo,
158  CLHEP::HepRandomEngine* engine) {
159  // produce SignalPoints for all SimHits in detector
160  unsigned int detID = det->geographicalId().rawId();
161  int numStrips = (det->specificTopology()).nstrips();
162 
163  size_t thisFirstChannelWithSignal = numStrips;
164  size_t thisLastChannelWithSignal = 0;
165 
166  float langle = (lorentzAngleHandle.isValid()) ? lorentzAngleHandle->getLorentzAngle(detID) : 0.;
167 
168  std::vector<float> locAmpl(numStrips, 0.);
169 
170  // Loop over hits
171 
172  uint32_t detId = det->geographicalId().rawId();
173  // First: loop on the SimHits
174  if(CLHEP::RandFlat::shoot(engine) > inefficiency) {
175  AssociationInfoForChannel* pDetIDAssociationInfo; // I only need this if makeDigiSimLinks_ is true...
176  if( makeDigiSimLinks_ ) pDetIDAssociationInfo=&(associationInfoForDetId_[detId]); // ...so only search the map if that is the case
177  std::vector<float> previousLocalAmplitude; // Only used if makeDigiSimLinks_ is true. Needed to work out the change in amplitude.
178 
179  size_t simHitGlobalIndex=inputBeginGlobalIndex; // This needs to stored to create the digi-sim link later
180  for (std::vector<PSimHit>::const_iterator simHitIter = inputBegin; simHitIter != inputEnd; ++simHitIter, ++simHitGlobalIndex ) {
181  // skip hits not in this detector.
182  if((*simHitIter).detUnitId() != detId) {
183  continue;
184  }
185  // check TOF
186  if (std::fabs(simHitIter->tof() - cosmicShift - det->surface().toGlobal(simHitIter->localPosition()).mag()/30.) < tofCut && simHitIter->energyLoss()>0) {
187  if( makeDigiSimLinks_ ) previousLocalAmplitude=locAmpl; // Not needed except to make the sim link association.
188  size_t localFirstChannel = numStrips;
189  size_t localLastChannel = 0;
190  // process the hit
191  theSiHitDigitizer->processHit(&*simHitIter, *det, bfield, langle, locAmpl, localFirstChannel, localLastChannel, tTopo, engine);
192 
193  if(thisFirstChannelWithSignal > localFirstChannel) thisFirstChannelWithSignal = localFirstChannel;
194  if(thisLastChannelWithSignal < localLastChannel) thisLastChannelWithSignal = localLastChannel;
195 
196  if( makeDigiSimLinks_ ) { // No need to do any of this if truth association was turned off in the configuration
197  for( size_t stripIndex=0; stripIndex<locAmpl.size(); ++stripIndex ) {
198  // Work out the amplitude from this SimHit from the difference of what it was before and what it is now
199  float signalFromThisSimHit=locAmpl[stripIndex]-previousLocalAmplitude[stripIndex];
200  if( signalFromThisSimHit!=0 ) { // If this SimHit had any contribution I need to record it.
201  auto& associationVector=(*pDetIDAssociationInfo)[stripIndex];
202  bool addNewEntry=true;
203  // Make sure the hit isn't in already. I've seen this a few times, it always seems to happen in pairs so I think
204  // it's something to do with the stereo strips.
205  for( auto& associationInfo : associationVector ) {
206  if( associationInfo.trackID==simHitIter->trackId() && associationInfo.eventID==simHitIter->eventId() ) {
207  // The hit is already in, so add this second contribution and move on
208  associationInfo.contributionToADC+=signalFromThisSimHit;
209  addNewEntry=false;
210  break;
211  }
212  } // end of loop over associationVector
213  // If the hit wasn't already in create a new association info structure.
214  if( addNewEntry ) associationVector.push_back( AssociationInfo{ simHitIter->trackId(), simHitIter->eventId(), signalFromThisSimHit, simHitGlobalIndex, tofBin } );
215  } // end of "if( signalFromThisSimHit!=0 )"
216  } // end of loop over locAmpl strips
217  } // end of "if( makeDigiSimLinks_ )"
218  } // end of TOF check
219  } // end for
220  }
221  theSiPileUpSignals->add(detID, locAmpl, thisFirstChannelWithSignal, thisLastChannelWithSignal);
222 
223  if(firstChannelsWithSignal[detID] > thisFirstChannelWithSignal) firstChannelsWithSignal[detID] = thisFirstChannelWithSignal;
224  if(lastChannelsWithSignal[detID] < thisLastChannelWithSignal) lastChannelsWithSignal[detID] = thisLastChannelWithSignal;
225 }
226 
227 //============================================================================
229  //Instlumi scalefactor calculating for dynamic inefficiency
230 
231  if (puInfo && FirstLumiCalc_) {
232 
233  const std::vector<int>& bunchCrossing = puInfo->getMix_bunchCrossing();
234  const std::vector<float>& TrueInteractionList = puInfo->getMix_TrueInteractions();
235  const int bunchSpacing = puInfo->getMix_bunchSpacing();
236 
237  double RevFreq = 11245.;
238  double minBXsec = 70.0E-27; // use 70mb as an approximation
239  double Bunch = 2100.; // 2016 value
240  if (bunchSpacing == 50) Bunch = Bunch/2.;
241 
242  int pui = 0, p = 0;
243  std::vector<int>::const_iterator pu;
244  std::vector<int>::const_iterator pu0 = bunchCrossing.end();
245 
246  for (pu=bunchCrossing.begin(); pu!=bunchCrossing.end(); ++pu) {
247  if (*pu==0) {
248  pu0 = pu;
249  p = pui;
250  }
251  pui++;
252  }
253  if (pu0!=bunchCrossing.end()) { // found the in-time interaction
254  double Tintr = TrueInteractionList.at(p);
255  double instLumi = Bunch*Tintr*RevFreq/minBXsec;
256  APVSaturationProb_ = instLumi/6.0E33;
257  }
258  FirstLumiCalc_ = false;
259  }
260 }
261 
262 //============================================================================
263 
264 
265 void
267  edm::DetSet<SiStripDigi>& outdigi,
268  edm::DetSet<SiStripRawDigi>& outrawdigi,
270  const StripGeomDetUnit *det,
271  edm::ESHandle<SiStripGain> & gainHandle,
272  edm::ESHandle<SiStripThreshold> & thresholdHandle,
273  edm::ESHandle<SiStripNoises> & noiseHandle,
274  edm::ESHandle<SiStripPedestals> & pedestalHandle,
275  std::vector<std::pair<int,std::bitset<6>>> & theAffectedAPVvector,
276  CLHEP::HepRandomEngine* engine) {
277  unsigned int detID = det->geographicalId().rawId();
278  int numStrips = (det->specificTopology()).nstrips();
279 
280  const SiPileUpSignals::SignalMapType* theSignal(theSiPileUpSignals->getSignal(detID));
281 
282  std::vector<float> detAmpl(numStrips, 0.);
283  if(theSignal) {
284  for(const auto& amp : *theSignal) {
285  detAmpl[amp.first] = amp.second;
286  }
287  }
288 
289  //removing signal from the dead (and HIP effected) strips
290  std::vector<bool>& badChannels = allBadChannels[detID];
291  for(int strip =0; strip < numStrips; ++strip) {
292  if(badChannels[strip]) {detAmpl[strip] = 0.;}
293  }
294 
296  //Implementation of the proper charge scaling function. Need consider resaturation effect:
297  //The probability map gives the probability that at least one HIP happened during the last N bunch crossings (cfr APV recovery time).
298  //The impact on the charge depends on the clostest HIP occurance (in terms of bunch crossing).
299  //The function discribing the APV recovery is therefore the weighted average function which takes into account all possibilities of HIP occurances across the last bx's.
300 
301  // do this step here because we now have access to luminosity information
302  if(FirstDigitize_) {
303 
304  for(std::map<int,float>::iterator iter = mapOfAPVprobabilities.begin(); iter != mapOfAPVprobabilities.end(); ++iter){
305  std::bitset<6> bs;
306  for(int Napv=0;Napv<6;Napv++){
307  float cursor=CLHEP::RandFlat::shoot(engine);
308  bs[Napv]=cursor < iter->second*APVSaturationProb_ ? true:false; //APVSaturationProb has been scaled by PU luminosity
309  }
310  SiStripTrackerAffectedAPVMap[iter->first]=bs;
311  }
312 
314  bool HasAtleastOneAffectedAPV=false;
315  while(!HasAtleastOneAffectedAPV){
316  for(int bx=floor(300.0/25.0);bx>0;bx--){ //Reminder: make these numbers not hard coded!!
317  float temp=CLHEP::RandFlat::shoot(engine)<0.5?1:0;
318  if(temp==1 && bx<NumberOfBxBetweenHIPandEvent){
320  HasAtleastOneAffectedAPV=true;
321  }
322  }
323  }
324 
325  FirstDigitize_ = false;
326  }
327 
328  std::bitset<6> & bs=SiStripTrackerAffectedAPVMap[detID];
329 
330  if(bs.any()){
331  // store this information so it can be saved to the event later
332  theAffectedAPVvector.push_back(std::make_pair(detID,bs));
333 
334  if(!PreMixing_) {
335 
336  // Here below is the scaling function which describes the evolution of the baseline (i.e. how the charge is suppressed).
337  // This must be replaced as soon as we have a proper modeling of the baseline evolution from VR runs
338  float Shift=1-NumberOfBxBetweenHIPandEvent/floor(300.0/25.0); //Reminder: make these numbers not hardcoded!!
339  float randomX=CLHEP::RandFlat::shoot(engine);
340  float scalingValue=(randomX-Shift)*10.0/7.0-3.0/7.0;
341 
342  for(int strip =0; strip < numStrips; ++strip) {
343  if(!badChannels[strip] && bs[strip/128]==1){
344  detAmpl[strip] *=scalingValue>0?scalingValue:0.0;
345  }
346  }
347  }
348  }
349  }
350 
351 
352  SiStripNoises::Range detNoiseRange = noiseHandle->getRange(detID);
353  SiStripApvGain::Range detGainRange = gainHandle->getRange(detID);
354  SiStripPedestals::Range detPedestalRange = pedestalHandle->getRange(detID);
355 
356 // -----------------------------------------------------------
357 
358  auto& firstChannelWithSignal = firstChannelsWithSignal[detID];
359  auto& lastChannelWithSignal = lastChannelsWithSignal[detID];
360  auto iAssociationInfoByChannel=associationInfoForDetId_.find(detID); // Use an iterator so that I can easily remove it once finished
361 
362  if(zeroSuppression){
363 
364  //Adding the strip noise
365  //------------------------------------------------------
366  if(noise){
367 
368  if(SingleStripNoise){
369 // std::cout<<"In SSN, detId="<<detID<<std::endl;
370  std::vector<float> noiseRMSv;
371  noiseRMSv.clear();
372  noiseRMSv.insert(noiseRMSv.begin(),numStrips,0.);
373  for(int strip=0; strip< numStrips; ++strip){
374  if(!badChannels[strip]){
375  float gainValue = gainHandle->getStripGain(strip, detGainRange);
376  noiseRMSv[strip] = (noiseHandle->getNoise(strip,detNoiseRange))* theElectronPerADC/gainValue;
377  //std::cout<<"<SiStripDigitizerAlgorithm::digitize>: gainValue: "<<gainValue<<"\tnoiseRMSv["<<strip<<"]: "<<noiseRMSv[strip]<<std::endl;
378  }
379  }
380  theSiNoiseAdder->addNoiseVR(detAmpl, noiseRMSv, engine);
381  } else {
382  int RefStrip = int(numStrips/2.);
383  while(RefStrip<numStrips&&badChannels[RefStrip]){ //if the refstrip is bad, I move up to when I don't find it
384  RefStrip++;
385  }
386  if(RefStrip<numStrips){
387  float RefgainValue = gainHandle->getStripGain(RefStrip, detGainRange);
388  float RefnoiseRMS = noiseHandle->getNoise(RefStrip,detNoiseRange) *theElectronPerADC/RefgainValue;
389 
390  theSiNoiseAdder->addNoise(detAmpl,firstChannelWithSignal,lastChannelWithSignal,numStrips,RefnoiseRMS, engine);
391  //std::cout<<"<SiStripDigitizerAlgorithm::digitize>: RefgainValue: "<<RefgainValue<<"\tRefnoiseRMS: "<<RefnoiseRMS<<std::endl;
392  }
393  }
394  }//if noise
395 
396  DigitalVecType digis;
397  theSiZeroSuppress->suppress(theSiDigitalConverter->convert(detAmpl, gainHandle, detID), digis, detID,noiseHandle,thresholdHandle);
398  // Now do the association to truth. Note that if truth association was turned off in the configuration this map
399  // will be empty and the iterator will always equal associationInfoForDetId_.end().
400  if( iAssociationInfoByChannel!=associationInfoForDetId_.end() ) { // make sure the readings for this DetID aren't completely from noise
401  for( const auto& iDigi : digis ) {
402  auto& associationInfoByChannel=iAssociationInfoByChannel->second;
403  const std::vector<AssociationInfo>& associationInfo=associationInfoByChannel[iDigi.channel()];
404 
405  // Need to find the total from all sim hits, because this might not be the same as the total
406  // digitised due to noise or whatever.
407  float totalSimADC=0;
408  for( const auto& iAssociationInfo : associationInfo ) totalSimADC+=iAssociationInfo.contributionToADC;
409  // Now I know that I can loop again and create the links
410  for( const auto& iAssociationInfo : associationInfo ) {
411  // Note simHitGlobalIndex used to have +1 because TrackerHitAssociator (the only place I can find this value being used)
412  // expected counting to start at 1, not 0. Now changed.
413  outLink.push_back( StripDigiSimLink( iDigi.channel(), iAssociationInfo.trackID, iAssociationInfo.simHitGlobalIndex, iAssociationInfo.tofBin, iAssociationInfo.eventID, iAssociationInfo.contributionToADC/totalSimADC ) );
414  } // end of loop over associationInfo
415  } // end of loop over the digis
416  } // end of check that iAssociationInfoByChannel is a valid iterator
417  outdigi.data = digis;
418  }//if zeroSuppression
419 
420  if(!zeroSuppression){
421  //if(noise){
422  // the constant pedestal offset is needed because
423  // negative adc counts are not allowed in case
424  // Pedestal and CMN subtraction is performed.
425  // The pedestal value read from the conditions
426  // is pedValue and after the pedestal subtraction
427  // the baseline is zero. The Common Mode Noise
428  // is not subtracted from the negative adc counts
429  // channels. Adding pedOffset the baseline is set
430  // to pedOffset after pedestal subtraction and CMN
431  // is subtracted to all the channels since none of
432  // them has negative adc value. The pedOffset is
433  // treated as a constant component in the CMN
434  // estimation and subtracted as CMN.
435 
436 
437  //calculating the charge deposited on each APV and subtracting the shift
438  //------------------------------------------------------
439  if(BaselineShift){
440  theSiNoiseAdder->addBaselineShift(detAmpl, badChannels);
441  }
442 
443  //Adding the strip noise
444  //------------------------------------------------------
445  if(noise){
446  std::vector<float> noiseRMSv;
447  noiseRMSv.clear();
448  noiseRMSv.insert(noiseRMSv.begin(),numStrips,0.);
449 
450  if(SingleStripNoise){
451  for(int strip=0; strip< numStrips; ++strip){
452  if(!badChannels[strip]) noiseRMSv[strip] = (noiseHandle->getNoise(strip,detNoiseRange))* theElectronPerADC;
453  }
454 
455  } else {
456  int RefStrip = 0; //int(numStrips/2.);
457  while(RefStrip<numStrips&&badChannels[RefStrip]){ //if the refstrip is bad, I move up to when I don't find it
458  RefStrip++;
459  }
460  if(RefStrip<numStrips){
461  float noiseRMS = noiseHandle->getNoise(RefStrip,detNoiseRange) *theElectronPerADC;
462  for(int strip=0; strip< numStrips; ++strip){
463  if(!badChannels[strip]) noiseRMSv[strip] = noiseRMS;
464  }
465  }
466  }
467 
468  theSiNoiseAdder->addNoiseVR(detAmpl, noiseRMSv, engine);
469  }
470 
471  //adding the CMN
472  //------------------------------------------------------
473  if(CommonModeNoise){
474  float cmnRMS = 0.;
475  DetId detId(detID);
476  uint32_t SubDet = detId.subdetId();
477  if(SubDet==3){
478  cmnRMS = cmnRMStib;
479  }else if(SubDet==4){
480  cmnRMS = cmnRMStid;
481  }else if(SubDet==5){
482  cmnRMS = cmnRMStob;
483  }else if(SubDet==6){
484  cmnRMS = cmnRMStec;
485  }
486  cmnRMS *= theElectronPerADC;
487  theSiNoiseAdder->addCMNoise(detAmpl, cmnRMS, badChannels, engine);
488  }
489 
490 
491  //Adding the pedestals
492  //------------------------------------------------------
493 
494  std::vector<float> vPeds;
495  vPeds.clear();
496  vPeds.insert(vPeds.begin(),numStrips,0.);
497 
498  if(RealPedestals){
499  for(int strip=0; strip< numStrips; ++strip){
500  if(!badChannels[strip]) vPeds[strip] = (pedestalHandle->getPed(strip,detPedestalRange)+pedOffset)* theElectronPerADC;
501  }
502  } else {
503  for(int strip=0; strip< numStrips; ++strip){
504  if(!badChannels[strip]) vPeds[strip] = pedOffset* theElectronPerADC;
505  }
506  }
507 
508  theSiNoiseAdder->addPedestals(detAmpl, vPeds);
509 
510 
511  //if(!RealPedestals&&!CommonModeNoise&&!noise&&!BaselineShift&&!APVSaturationFromHIP){
512  // edm::LogWarning("SiStripDigitizer")<<"You are running the digitizer without Noise generation and without applying Zero Suppression. ARE YOU SURE???";
513  //}else{
514 
515  DigitalRawVecType rawdigis = theSiDigitalConverter->convertRaw(detAmpl, gainHandle, detID);
516 
517  // Now do the association to truth. Note that if truth association was turned off in the configuration this map
518  // will be empty and the iterator will always equal associationInfoForDetId_.end().
519  if( iAssociationInfoByChannel!=associationInfoForDetId_.end() ) { // make sure the readings for this DetID aren't completely from noise
520  // N.B. For the raw digis the channel is inferred from the position in the vector.
521  // I'VE NOT TESTED THIS YET!!!!!
522  // ToDo Test this properly.
523  for( size_t channel=0; channel<rawdigis.size(); ++channel ) {
524  auto& associationInfoByChannel=iAssociationInfoByChannel->second;
525  const auto iAssociationInfo=associationInfoByChannel.find(channel);
526  if( iAssociationInfo==associationInfoByChannel.end() ) continue; // Skip if there is no sim information for this channel (i.e. it's all noise)
527  const std::vector<AssociationInfo>& associationInfo=iAssociationInfo->second;
528 
529  // Need to find the total from all sim hits, because this might not be the same as the total
530  // digitised due to noise or whatever.
531  float totalSimADC=0;
532  for( const auto& iAssociationInfo : associationInfo ) totalSimADC+=iAssociationInfo.contributionToADC;
533  // Now I know that I can loop again and create the links
534  for( const auto& iAssociationInfo : associationInfo ) {
535  // Note simHitGlobalIndex used to have +1 because TrackerHitAssociator (the only place I can find this value being used)
536  // expected counting to start at 1, not 0. Now changed.
537  outLink.push_back( StripDigiSimLink( channel, iAssociationInfo.trackID, iAssociationInfo.simHitGlobalIndex, iAssociationInfo.tofBin, iAssociationInfo.eventID, iAssociationInfo.contributionToADC/totalSimADC ) );
538  } // end of loop over associationInfo
539  } // end of loop over the digis
540  } // end of check that iAssociationInfoByChannel is a valid iterator
541 
542  outrawdigi.data = rawdigis;
543 
544  //}
545  }
546 
547  // Now that I've finished with this entry in the map of associations, I can remove it.
548  // Note that there might not be an association if the ADC reading is from noise in which
549  // case associationIsValid will be false.
550  if( iAssociationInfoByChannel!=associationInfoForDetId_.end() ) associationInfoForDetId_.erase(iAssociationInfoByChannel);
551 }
#define LogDebug(id)
unsigned short range
std::map< int, float > mapOfAPVprobabilities
GlobalPoint toGlobal(const Point2DBase< Scalar, LocalTag > lp) const
Definition: Surface.h:106
const std::unique_ptr< SiPileUpSignals > theSiPileUpSignals
void digitize(edm::DetSet< SiStripDigi > &outDigis, edm::DetSet< SiStripRawDigi > &outRawDigis, edm::DetSet< StripDigiSimLink > &outLink, const StripGeomDetUnit *stripdet, edm::ESHandle< SiStripGain > &, edm::ESHandle< SiStripThreshold > &, edm::ESHandle< SiStripNoises > &, edm::ESHandle< SiStripPedestals > &, std::vector< std::pair< int, std::bitset< 6 >>> &theAffectedAPVvector, CLHEP::HepRandomEngine *)
edm::ESHandle< SiStripLorentzAngle > lorentzAngleHandle
void push_back(const T &t)
Definition: DetSet.h:68
AssociationInfoForDetId associationInfoForDetId_
Structure that holds the information on the SimTrack contributions. Only filled if makeDigiSimLinks_ ...
const std::vector< float > & getMix_TrueInteractions() const
std::map< int, std::bitset< 6 > > SiStripTrackerAffectedAPVMap
SiDigitalConverter::DigitalRawVecType DigitalRawVecType
SiDigitalConverter::DigitalVecType DigitalVecType
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
const std::vector< int > & getMix_bunchCrossing() const
std::vector< unsigned int >::const_iterator ContainerIterator
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:47
void initializeEvent(const edm::EventSetup &iSetup)
std::pair< ContainerIterator, ContainerIterator > Range
void initializeDetUnit(StripGeomDetUnit const *det, const edm::EventSetup &iSetup)
virtual const StripTopology & specificTopology() const
Returns a reference to the strip proxy topology.
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:42
const std::unique_ptr< const SiGaussianTailNoiseAdder > theSiNoiseAdder
const std::unique_ptr< SiHitDigitizer > theSiHitDigitizer
void getData(T &iHolder) const
Definition: EventSetup.h:98
std::map< unsigned int, size_t > firstChannelsWithSignal
float getPed(const uint16_t &strip, const Range &range) const
std::map< int, Amplitude > SignalMapType
static float getNoise(uint16_t strip, const Range &range)
Definition: SiStripNoises.h:74
float getLorentzAngle(const uint32_t &) const
const std::unique_ptr< SiTrivialDigitalConverter > theSiDigitalConverter
const std::unique_ptr< SiStripFedZeroSuppression > theSiZeroSuppress
static float getStripGain(const uint16_t &strip, const SiStripApvGain::Range &range)
Definition: SiStripGain.h:72
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:41
std::pair< ContainerIterator, ContainerIterator > Range
DetId geographicalId() const
The label of this GeomDet.
Definition: GeomDet.h:79
void calculateInstlumiScale(PileupMixingContent *puInfo)
std::map< unsigned int, std::vector< bool > > allBadChannels
void setParticleDataTable(const ParticleDataTable *pardt)
Definition: DetId.h:18
unsigned short firstStrip
std::map< int, std::vector< AssociationInfo > > AssociationInfoForChannel
RealPedestals
NOTE : turning Noise ON/OFF will make a big change Parameters valid only if Noise = True and ZeroSupp...
SiStripDigitizerAlgorithm(const edm::ParameterSet &conf)
collection_type data
Definition: DetSet.h:78
const Range getRange(const uint32_t detID) const
HLT enums.
const Range getRange(const uint32_t detID) const
std::pair< ContainerIterator, ContainerIterator > Range
std::map< unsigned int, size_t > lastChannelsWithSignal
T get() const
Definition: EventSetup.h:63
void accumulateSimHits(const std::vector< PSimHit >::const_iterator inputBegin, const std::vector< PSimHit >::const_iterator inputEnd, size_t inputBeginGlobalIndex, unsigned int tofBin, const StripGeomDetUnit *stripdet, const GlobalVector &bfield, const TrackerTopology *tTopo, CLHEP::HepRandomEngine *)
const int & getMix_bunchSpacing() const
std::string fullPath() const
Definition: FileInPath.cc:197
const Range getRange(const uint32_t &detID) const
std::pair< ContainerIterator, ContainerIterator > Range
Definition: SiStripNoises.h:50
bool isValid() const
Definition: ESHandle.h:47
double split
Definition: MVATrainer.cc:139
const ParticleDataTable * pdt
data decode(const unsigned int &value) const
const SiStripApvGain::Range getRange(uint32_t detID) const
Definition: SiStripGain.h:70