![]() |
![]() |
#include <FastSimulation/CaloRecHitsProducer/interface/HcalRecHitsMaker.h>
Definition at line 25 of file HcalRecHitsMaker.h.
HcalRecHitsMaker::HcalRecHitsMaker | ( | edm::ParameterSet const & | p, | |
edm::ParameterSet const & | p, | |||
const RandomEngine * | random | |||
) |
Definition at line 32 of file HcalRecHitsMaker.cc.
References doDigis_, doSaturation_, lat::endl(), gains_, edm::ParameterSet::getParameter(), hbhi_, hcalfileinpath_, hcalHotFractionHB_, hcalHotFractionHE_, hcalHotFractionHF_, hcalHotFractionHO_, hehi_, hfhi_, hohi_, inputCol_, maxIndex_, maxIndexDebug_, myGaussianTailGeneratorHB_, myGaussianTailGeneratorHE_, myGaussianTailGeneratorHF_, myGaussianTailGeneratorHO_, noiseFromDb_, noiseHB_, noiseHE_, noiseHF_, noiseHO_, noisesigma_, peds_, random_, refactor_, refactor_mean_, sat_, funct::sqrt(), theDetIds_, thresholdHB_, thresholdHE_, thresholdHF_, and thresholdHO_.
00034 : 00035 initialized_(false), 00036 random_(myrandom), 00037 myGaussianTailGeneratorHB_(0), myGaussianTailGeneratorHE_(0), myGaussianTailGeneratorHO_(0), myGaussianTailGeneratorHF_(0),myHcalSimParameterMap_(0) 00038 { 00039 edm::ParameterSet RecHitsParameters = p1.getParameter<edm::ParameterSet>("HCAL"); 00040 noiseHB_ = RecHitsParameters.getParameter<double>("NoiseHB"); 00041 noiseHE_ = RecHitsParameters.getParameter<double>("NoiseHE"); 00042 noiseHO_ = RecHitsParameters.getParameter<double>("NoiseHO"); 00043 noiseHF_ = RecHitsParameters.getParameter<double>("NoiseHF"); 00044 thresholdHB_ = RecHitsParameters.getParameter<double>("ThresholdHB"); 00045 thresholdHE_ = RecHitsParameters.getParameter<double>("ThresholdHE"); 00046 thresholdHO_ = RecHitsParameters.getParameter<double>("ThresholdHO"); 00047 thresholdHF_ = RecHitsParameters.getParameter<double>("ThresholdHF"); 00048 doSaturation_ = RecHitsParameters.getParameter<bool>("EnableSaturation"); 00049 noiseFromDb_ = RecHitsParameters.getParameter<bool>("NoiseFromDb"); 00050 00051 refactor_ = RecHitsParameters.getParameter<double> ("Refactor"); 00052 refactor_mean_ = RecHitsParameters.getParameter<double> ("Refactor_mean"); 00053 hcalfileinpath_= RecHitsParameters.getParameter<std::string> ("fileNameHcal"); 00054 inputCol_=RecHitsParameters.getParameter<edm::InputTag>("MixedSimHits"); 00055 maxIndex_=0; 00056 maxIndexDebug_=0; 00057 theDetIds_.resize(10000); 00058 00059 peds_.resize(10000); 00060 gains_.resize(10000); 00061 sat_.resize(10000); 00062 noisesigma_.resize(10000); 00063 hbhi_.reserve(2600); 00064 hehi_.reserve(2600); 00065 hohi_.reserve(2200); 00066 hfhi_.reserve(1800); 00067 doDigis_=false; 00068 00069 // edm::ParameterSet hcalparam = p2.getParameter<edm::ParameterSet>("hcalSimParam"); 00070 // myHcalSimParameterMap_ = new HcalSimParameterMap(hcalparam); 00071 00072 // Computes the fraction of HCAL above the threshold 00073 Genfun::Erf myErf; 00074 00075 if(noiseHB_>0.) { 00076 hcalHotFractionHB_ = 0.5-0.5*myErf(thresholdHB_/noiseHB_/sqrt(2.)); 00077 myGaussianTailGeneratorHB_ = new GaussianTail(random_,noiseHB_,thresholdHB_); 00078 if(hcalHotFractionHB_ < 0.3 && noiseFromDb_) 00079 edm::LogWarning("CaloRecHitsProducer") << " WARNING : HCAL Noise simulation, you chose to compute the noise sigma from the database, but the energy threshold is too high, change the value TresholdHB " << std::endl; 00080 } else { 00081 hcalHotFractionHB_ =0.; 00082 } 00083 00084 if(noiseHO_>0.) { 00085 hcalHotFractionHO_ = 0.5-0.5*myErf(thresholdHO_/noiseHO_/sqrt(2.)); 00086 myGaussianTailGeneratorHO_ = new GaussianTail(random_,noiseHO_,thresholdHO_); 00087 if(hcalHotFractionHO_ < 0.3 && noiseFromDb_) 00088 edm::LogWarning("CaloRecHitsProducer") << " WARNING : HCAL Noise simulation, you chose to compute the noise sigma from the database, but the energy threshold is too high, change the value TresholdHO " << std::endl; 00089 } else { 00090 hcalHotFractionHO_ =0.; 00091 } 00092 00093 if(noiseHE_>0.) { 00094 hcalHotFractionHE_ = 0.5-0.5*myErf(thresholdHE_/noiseHE_/sqrt(2.)); 00095 myGaussianTailGeneratorHE_ = new GaussianTail(random_,noiseHE_,thresholdHE_); 00096 if(hcalHotFractionHE_ < 0.3 && noiseFromDb_) 00097 edm::LogWarning("CaloRecHitsProducer") << " WARNING : HCAL Noise simulation, you chose to compute the noise sigma from the database, but the energy threshold is too high, change the value TresholdHE " << std::endl; 00098 } else { 00099 hcalHotFractionHE_ =0.; 00100 } 00101 00102 if(noiseHF_>0.) { 00103 hcalHotFractionHF_ = 0.5-0.5*myErf(thresholdHF_/noiseHF_/sqrt(2.)); 00104 myGaussianTailGeneratorHF_ = new GaussianTail(random_,noiseHF_,thresholdHF_); 00105 if(hcalHotFractionHF_ < 0.3 && noiseFromDb_) 00106 edm::LogWarning("CaloRecHitsProducer") << " WARNING : HCAL Noise simulation, you chose to compute the noise sigma from the database, but the energy threshold is too high, change the value TresholdHF " << std::endl; 00107 } else { 00108 hcalHotFractionHF_ =0.; 00109 } 00110 00111 }
HcalRecHitsMaker::~HcalRecHitsMaker | ( | ) |
Definition at line 113 of file HcalRecHitsMaker.cc.
References clean(), hbhi_, hehi_, hfhi_, hohi_, myGaussianTailGeneratorHB_, myGaussianTailGeneratorHE_, myGaussianTailGeneratorHF_, myGaussianTailGeneratorHO_, myHcalSimParameterMap_, and theDetIds_.
00114 { 00115 clean(); 00116 if(myGaussianTailGeneratorHB_) delete myGaussianTailGeneratorHB_; 00117 if(myGaussianTailGeneratorHE_) delete myGaussianTailGeneratorHE_; 00118 if(myGaussianTailGeneratorHO_) delete myGaussianTailGeneratorHO_; 00119 if(myGaussianTailGeneratorHF_) delete myGaussianTailGeneratorHF_; 00120 if(myHcalSimParameterMap_) delete myHcalSimParameterMap_; 00121 theDetIds_.clear(); 00122 hbhi_.clear(); 00123 hehi_.clear(); 00124 hohi_.clear(); 00125 hfhi_.clear(); 00126 00127 }
void HcalRecHitsMaker::clean | ( | ) | [private] |
Definition at line 633 of file HcalRecHitsMaker.cc.
References cleanSubDet(), firedCellsHB_, firedCellsHE_, firedCellsHF_, firedCellsHO_, and hcalRecHits_.
Referenced by loadPCaloHits(), and ~HcalRecHitsMaker().
00634 { 00635 cleanSubDet(hcalRecHits_,firedCellsHB_); 00636 cleanSubDet(hcalRecHits_,firedCellsHE_); 00637 cleanSubDet(hcalRecHits_,firedCellsHO_); 00638 cleanSubDet(hcalRecHits_,firedCellsHF_); 00639 }
void HcalRecHitsMaker::cleanSubDet | ( | std::vector< float > & | hits, | |
std::vector< int > & | cells | |||
) | [private] |
Definition at line 641 of file HcalRecHitsMaker.cc.
References size.
Referenced by clean().
00642 { 00643 unsigned size=cells.size(); 00644 // Reset the energies 00645 for(unsigned ic=0;ic<size;++ic) 00646 { 00647 hits[cells[ic]] = 0.; 00648 } 00649 // Clear the list of fired cells 00650 cells.clear(); 00651 }
unsigned HcalRecHitsMaker::createVectorOfSubdetectorCells | ( | const CaloGeometry & | cg, | |
int | subdetn, | |||
std::vector< int > & | cellsvec | |||
) | [private] |
Definition at line 490 of file HcalRecHitsMaker.cc.
References CaloGeometry::getSubdetectorGeometry(), CaloSubdetectorGeometry::getValidDetIds(), HcalDetId::hashed_index(), DetId::Hcal, i, maxIndex_, and theDetIds_.
Referenced by createVectorsOfCells().
00491 { 00492 const CaloSubdetectorGeometry* geom=cg.getSubdetectorGeometry(DetId::Hcal,subdetn); 00493 const std::vector<DetId>& ids=geom->getValidDetIds(DetId::Hcal,subdetn); 00494 for (std::vector<DetId>::const_iterator i=ids.begin(); i!=ids.end(); i++) 00495 { 00496 HcalDetId myDetId(*i); 00497 // std::cout << myDetId << myHcalSimParameterMap_->simParameters(myDetId).simHitToPhotoelectrons() << std::endl;; 00498 unsigned hi=myDetId.hashed_index(); 00499 theDetIds_[hi]=myDetId; 00500 // std::cout << myDetId << " " << hi << std::endl; 00501 cellsvec.push_back(hi); 00502 00503 if(hi>maxIndex_) 00504 maxIndex_=hi; 00505 } 00506 return cellsvec.size(); 00507 }
unsigned HcalRecHitsMaker::createVectorsOfCells | ( | const edm::EventSetup & | es | ) | [private] |
Definition at line 478 of file HcalRecHitsMaker.cc.
References createVectorOfSubdetectorCells(), edm::EventSetup::get(), hbhi_, HcalBarrel, HcalEndcap, HcalForward, HcalOuter, hehi_, hfhi_, hohi_, nhbcells_, nhecells_, nhfcells_, and nhocells_.
Referenced by init().
00479 { 00480 edm::ESHandle<CaloGeometry> pG; 00481 es.get<CaloGeometryRecord>().get(pG); 00482 nhbcells_ = createVectorOfSubdetectorCells(*pG, HcalBarrel, hbhi_); 00483 nhecells_ = createVectorOfSubdetectorCells(*pG, HcalEndcap, hehi_); 00484 nhocells_ = createVectorOfSubdetectorCells(*pG, HcalOuter, hohi_); 00485 nhfcells_ = createVectorOfSubdetectorCells(*pG, HcalForward, hfhi_); 00486 return nhbcells_+nhecells_+nhocells_+nhfcells_; 00487 }
int HcalRecHitsMaker::fCtoAdc | ( | double | fc | ) | const [private] |
Definition at line 654 of file HcalRecHitsMaker.cc.
References fctoadc_.
Referenced by loadHcalRecHits().
00655 { 00656 if(fc<0.) return 0; 00657 if(fc>9985.) return 127; 00658 return fctoadc_[(unsigned)floor(fc)]; 00659 }
void HcalRecHitsMaker::Fill | ( | int | id, | |
float | energy, | |||
std::vector< int > & | myHits, | |||
float | noise | |||
) | [private] |
Definition at line 510 of file HcalRecHitsMaker.cc.
References doMiscalib_, RandomEngine::gaussShoot(), hcalRecHits_, miscalib_, noiseFromDb_, noisesigma_, and random_.
Referenced by loadPCaloHits().
00511 { 00512 if(doMiscalib_) 00513 energy*=miscalib_[id]; 00514 00515 if(noiseFromDb_) 00516 noise=noisesigma_[id]; 00517 00518 // Check if the RecHit exists 00519 if(hcalRecHits_[id]>0.) 00520 hcalRecHits_[id]+=energy; 00521 else 00522 { 00523 // the noise is injected only the first time 00524 hcalRecHits_[id]=energy + random_->gaussShoot(0.,noise); 00525 theHits.push_back(id); 00526 } 00527 }
void HcalRecHitsMaker::init | ( | const edm::EventSetup & | es, | |
bool | dodigis, | |||
bool | domiscalib | |||
) |
Definition at line 129 of file HcalRecHitsMaker.cc.
References GenMuonPlsPt100GeV_cfg::cout, createVectorsOfCells(), doDigis_, doMiscalib_, doSaturation_, lat::endl(), fctoadc_, edm::FileInPath::fullPath(), gains_, edm::EventSetup::get(), CaloMiscalibMapHcal::get(), HcalDbService::getGain(), HcalDbService::getPedestal(), HcalPedestal::getValue(), HcalGain::getValue(), HcalDetId::hashed_index(), hbhi_, hcalfileinpath_, hcalRecHits_, hehi_, hfhi_, hohi_, i, in, index, initialized_, int, it, maxIndex_, miscalib_, nhbcells_, nhecells_, nhfcells_, nhocells_, noiseInfCfromDB(), noisesigma_, MiscalibReaderFromXML::parseXMLMiscalibFile(), peds_, CaloMiscalibMapHcal::prefillMap(), CaloMiscalibMapHcal::print(), edm::ESHandle< T >::product(), refactor_, refactor_mean_, sat_, size, theDetIds_, TPGFactor_, x, and y.
Referenced by CaloRecHitsProducer::beginRun().
00130 { 00131 doDigis_=doDigis; 00132 doMiscalib_=doMiscalib; 00133 if(initialized_) return; 00134 00135 00136 unsigned ncells=createVectorsOfCells(es); 00137 edm::LogInfo("CaloRecHitsProducer") << "Total number of cells in HCAL " << ncells << std::endl; 00138 hcalRecHits_.resize(maxIndex_+1,0.); 00139 edm::LogInfo("CaloRecHitsProducer") << "Largest HCAL hashedindex" << maxIndex_ << std::endl; 00140 00141 if(doMiscalib_) 00142 { 00143 miscalib_.resize(maxIndex_+1,1.); 00144 // Read from file ( a la HcalRecHitsRecalib.cc) 00145 // here read them from xml (particular to HCAL) 00146 CaloMiscalibMapHcal mapHcal; 00147 mapHcal.prefillMap(); 00148 00149 edm::FileInPath hcalfiletmp("CalibCalorimetry/CaloMiscalibTools/data/"+hcalfileinpath_); 00150 std::string hcalfile=hcalfiletmp.fullPath(); 00151 MiscalibReaderFromXMLHcal hcalreader_(mapHcal); 00152 if(!hcalfile.empty()) 00153 { 00154 hcalreader_.parseXMLMiscalibFile(hcalfile); 00155 mapHcal.print(); 00156 std::map<uint32_t,float>::const_iterator it=mapHcal.get().begin(); 00157 std::map<uint32_t,float>::const_iterator itend=mapHcal.get().end(); 00158 for(;it!=itend;++it) 00159 { 00160 HcalDetId myDetId(it->first); 00161 float icalconst=it->second; 00162 miscalib_[myDetId.hashed_index()]=refactor_mean_+(icalconst-refactor_mean_)*refactor_; 00163 } 00164 } 00165 } 00166 00167 00168 // Will be needed for the DB-based miscalibration 00169 edm::ESHandle<HcalDbService> conditions; 00170 es.get<HcalDbRecord>().get(conditions); 00171 const HcalDbService * theDbService=conditions.product(); 00172 // myHcalSimParameterMap_->setDbService(theDbService); 00173 // Open the histogram for the fC to ADC conversion 00174 gROOT->cd(); 00175 edm::FileInPath myDataFile("FastSimulation/CaloRecHitsProducer/data/adcvsfc.root"); 00176 TFile * myFile = new TFile(myDataFile.fullPath().c_str(),"READ"); 00177 TGraph * myGraf = (TGraph*)myFile->Get("adcvsfc"); 00178 unsigned size=myGraf->GetN(); 00179 fctoadc_.resize(10000); 00180 unsigned p_index=0; 00181 fctoadc_[0]=0; 00182 int prev_nadc=0; 00183 int nadc=0; 00184 for(unsigned ibin=0;ibin<size;++ibin) 00185 { 00186 double x,y; 00187 myGraf->GetPoint(ibin,x,y); 00188 int index=(int)floor(x); 00189 if(index<0||index>=10000) continue; 00190 prev_nadc=nadc; 00191 nadc=(int)y; 00192 // Now fills the vector 00193 for(unsigned ivec=p_index;ivec<(unsigned)index;++ivec) 00194 { 00195 fctoadc_[ivec] = prev_nadc; 00196 } 00197 p_index = index; 00198 } 00199 myFile->Close(); 00200 gROOT->cd(); 00201 edm::FileInPath myTPGFilePath("CalibCalorimetry/HcalTPGAlgos/data/RecHit-TPG-calib.dat"); 00202 TPGFactor_.resize(87,1.2); 00203 std::ifstream myTPGFile(myTPGFilePath.fullPath().c_str(),ifstream::in); 00204 if(myTPGFile) 00205 { 00206 float gain; 00207 myTPGFile >> gain; 00208 for(unsigned i=0;i<86;++i) 00209 { 00210 myTPGFile >> TPGFactor_[i] ; 00211 // std::cout << TPGFactor_[i] << std::endl; 00212 } 00213 } 00214 else 00215 { 00216 std::cout << " Unable to open CalibCalorimetry/HcalTPGAlgos/data/RecHit-TPG-calib.dat" << std::endl; 00217 std::cout << " Using a constant 1.2 factor " << std::endl; 00218 } 00219 00220 //HB 00221 for(unsigned ic=0;ic<nhbcells_;++ic) 00222 { 00223 float gain = theDbService->getGain(theDetIds_[hbhi_[ic]])->getValue(0); 00224 float mgain=0.; 00225 for(unsigned ig=0;ig<4;++ig) 00226 mgain+=theDbService->getGain(theDetIds_[hbhi_[ic]])->getValue(ig); 00227 noisesigma_[hbhi_[ic]]=noiseInfCfromDB(theDbService,theDetIds_[hbhi_[ic]])*mgain*0.25; 00228 // std::cout << " NOISEHB " << theDetIds_[hbhi_[ic]].ieta() << " " << noisesigma_[hbhi_[ic]] << " "<< std::endl; 00229 mgain*=2500.;// 10000 (ADC scale) / 4. (to compute the mean) 00230 sat_[hbhi_[ic]]=(doSaturation_)?mgain:99999.; 00231 00232 peds_[hbhi_[ic]]=theDbService->getPedestal(theDetIds_[hbhi_[ic]])->getValue(0); 00233 int ieta=theDetIds_[hbhi_[ic]].ieta(); 00234 gain*=TPGFactor_[(ieta>0)?ieta+43:-ieta]; 00235 gains_[hbhi_[ic]]=gain; 00236 } 00237 //HE 00238 for(unsigned ic=0;ic<nhecells_;++ic) 00239 { 00240 float gain= theDbService->getGain(theDetIds_[hehi_[ic]])->getValue(0); 00241 int ieta=theDetIds_[hehi_[ic]].ieta(); 00242 float mgain=0.; 00243 for(unsigned ig=0;ig<4;++ig) 00244 mgain+=theDbService->getGain(theDetIds_[hehi_[ic]])->getValue(ig); 00245 noisesigma_[hehi_[ic]]=noiseInfCfromDB(theDbService,theDetIds_[hehi_[ic]])*mgain*0.25; 00246 00247 // std::cout << " NOISEHE " << theDetIds_[hehi_[ic]].ieta() << " " << noisesigma_[hehi_[ic]] << " "<< std::endl; 00248 mgain*=2500.; 00249 sat_[hehi_[ic]]=(doSaturation_)?mgain:99999.; 00250 00251 gain*=TPGFactor_[(ieta>0)?ieta+44:-ieta+1]; 00252 // if(abs(ieta)>=21&&abs(ieta)<=26) 00253 // gain*=2.; 00254 // if(abs(ieta)>=27&&abs(ieta)<=29) 00255 // gain*=5.; 00256 peds_[hehi_[ic]]=theDbService->getPedestal(theDetIds_[hehi_[ic]])->getValue(0); 00257 gains_[hehi_[ic]]=gain; 00258 } 00259 //HO 00260 for(unsigned ic=0;ic<nhocells_;++ic) 00261 { 00262 float ped=theDbService->getPedestal(theDetIds_[hohi_[ic]])->getValue(0); 00263 float gain=theDbService->getGain(theDetIds_[hohi_[ic]])->getValue(0); 00264 float mgain=0.; 00265 for(unsigned ig=0;ig<4;++ig) 00266 mgain+=theDbService->getGain(theDetIds_[hohi_[ic]])->getValue(ig); 00267 noisesigma_[hohi_[ic]]=noiseInfCfromDB(theDbService,theDetIds_[hohi_[ic]])*mgain*0.25; 00268 // std::cout << " NOISEHO " << theDetIds_[hohi_[ic]].ieta() << " " << noisesigma_[hohi_[ic]] << " "<< std::endl; 00269 00270 mgain*=2500.; 00271 sat_[hohi_[ic]]=(doSaturation_)?mgain:99999.; 00272 int ieta=HcalDetId(hohi_[ic]).ieta(); 00273 gain*=TPGFactor_[(ieta>0)?ieta+43:-ieta]; 00274 peds_[hohi_[ic]]=ped; 00275 gains_[hohi_[ic]]=gain; 00276 } 00277 //HF 00278 for(unsigned ic=0;ic<nhfcells_;++ic) 00279 { 00280 float ped=theDbService->getPedestal(theDetIds_[hfhi_[ic]])->getValue(0); 00281 float gain=theDbService->getGain(theDetIds_[hfhi_[ic]])->getValue(0); 00282 float mgain=0.; 00283 for(unsigned ig=0;ig<4;++ig) 00284 mgain+=theDbService->getGain(theDetIds_[hfhi_[ic]])->getValue(ig); 00285 // additional 1/2 factor for the HF (Salavat) 00286 noisesigma_[hfhi_[ic]]=noiseInfCfromDB(theDbService,theDetIds_[hfhi_[ic]])*mgain*0.25; 00287 // std::cout << " NOISEHF " << theDetIds_[hfhi_[ic]].ieta() << " " << noisesigma_[hfhi_[ic]] << " "<< std::endl; 00288 00289 mgain*=2500.; 00290 sat_[hfhi_[ic]]=(doSaturation_)?mgain:99999.; 00291 int ieta=theDetIds_[hfhi_[ic]].ieta(); 00292 gain*=TPGFactor_[(ieta>0)?ieta+45:-ieta+2]; 00293 peds_[hfhi_[ic]]=ped; 00294 gains_[hfhi_[ic]]=gain; 00295 } 00296 00297 initialized_=true; 00298 }
void HcalRecHitsMaker::loadHcalRecHits | ( | edm::Event & | iEvent, | |
HBHERecHitCollection & | hbheHits, | |||
HORecHitCollection & | hoHits, | |||
HFRecHitCollection & | hfHits, | |||
HBHEDigiCollection & | hbheDigis, | |||
HODigiCollection & | hoDigis, | |||
HFDigiCollection & | hfDigis | |||
) |
Definition at line 356 of file HcalRecHitsMaker.cc.
References doDigis_, relval_parameters_module::energy, fCtoAdc(), firedCellsHB_, firedCellsHE_, firedCellsHF_, firedCellsHO_, gains_, hcalRecHits_, loadPCaloHits(), noisify(), peds_, edm::SortedCollection< T, SORT >::push_back(), edm::SortedCollection< T, SORT >::reserve(), sat_, HBHEDataFrame::setSample(), HFDataFrame::setSample(), HBHEDataFrame::setSize(), HFDataFrame::setSize(), theDetIds_, thresholdHB_, thresholdHE_, thresholdHF_, and thresholdHO_.
Referenced by CaloRecHitsProducer::produce().
00357 { 00358 loadPCaloHits(iEvent); 00359 noisify(); 00360 hbheHits.reserve(firedCellsHB_.size()+firedCellsHE_.size()); 00361 hoHits.reserve(firedCellsHO_.size()); 00362 hfHits.reserve(firedCellsHF_.size()); 00363 if(doDigis_) 00364 { 00365 hbheDigis.reserve(firedCellsHB_.size()+firedCellsHE_.size()); 00366 hfDigis.reserve(firedCellsHF_.size()); 00367 hoDigis.reserve(firedCellsHO_.size()); 00368 } 00369 static HcalQIESample zeroSample(0,0,0,0); 00370 unsigned nhits=firedCellsHB_.size(); 00371 // HB 00372 for(unsigned ihit=0;ihit<nhits;++ihit) 00373 { 00374 unsigned cellhashedindex=firedCellsHB_[ihit]; 00375 // Check if it is above the threshold 00376 if(hcalRecHits_[cellhashedindex]<thresholdHB_) continue; 00377 float energy=hcalRecHits_[cellhashedindex]; 00378 // poor man saturation 00379 if(energy>sat_[cellhashedindex]) 00380 { 00381 // std::cout << " Saturation " << energy << " " << sat_[cellhashedindex] << " HB " << std::endl; 00382 energy=sat_[cellhashedindex]; 00383 } 00384 00385 const HcalDetId& detid = theDetIds_[cellhashedindex]; 00386 hbheHits.push_back(HBHERecHit(detid,energy,0.)); 00387 if(doDigis_) 00388 { 00389 HBHEDataFrame myDataFrame(detid); 00390 myDataFrame.setSize(2); 00391 double nfc=hcalRecHits_[cellhashedindex]/gains_[cellhashedindex]+peds_[cellhashedindex]; 00392 int nadc=fCtoAdc(nfc); 00393 HcalQIESample qie(nadc, 0, 0, 0) ; 00394 myDataFrame.setSample(0,qie); 00395 myDataFrame.setSample(1,zeroSample); 00396 hbheDigis.push_back(myDataFrame); 00397 } 00398 } 00399 00400 // HE 00401 nhits=firedCellsHE_.size(); 00402 for(unsigned ihit=0;ihit<nhits;++ihit) 00403 { 00404 unsigned cellhashedindex=firedCellsHE_[ihit]; 00405 // Check if it is above the threshold 00406 if(hcalRecHits_[cellhashedindex]<thresholdHE_) continue; 00407 float energy=hcalRecHits_[cellhashedindex]; 00408 // poor man saturation 00409 if(energy>sat_[cellhashedindex]) 00410 { 00411 // std::cout << " Op saturation " << energy << " " << sat_[cellhashedindex] << " HE " << std::endl; 00412 energy=sat_[cellhashedindex]; 00413 } 00414 00415 const HcalDetId & detid= theDetIds_[cellhashedindex]; 00416 hbheHits.push_back(HBHERecHit(detid,energy,0.)); 00417 if(doDigis_) 00418 { 00419 HBHEDataFrame myDataFrame(detid); 00420 myDataFrame.setSize(2); 00421 double nfc=hcalRecHits_[cellhashedindex]/gains_[cellhashedindex]+peds_[cellhashedindex]; 00422 int nadc=fCtoAdc(nfc); 00423 HcalQIESample qie(nadc, 0, 0, 0) ; 00424 myDataFrame.setSample(0,qie); 00425 myDataFrame.setSample(1,zeroSample); 00426 hbheDigis.push_back(myDataFrame); 00427 } 00428 } 00429 00430 00431 // HO 00432 nhits=firedCellsHO_.size(); 00433 for(unsigned ihit=0;ihit<nhits;++ihit) 00434 { 00435 00436 unsigned cellhashedindex=firedCellsHO_[ihit]; 00437 // Check if it is above the threshold 00438 if(hcalRecHits_[cellhashedindex]<thresholdHO_) continue; 00439 00440 float energy=hcalRecHits_[cellhashedindex]; 00441 // poor man saturation 00442 if(energy>sat_[cellhashedindex]) energy=sat_[cellhashedindex]; 00443 00444 const HcalDetId& detid=theDetIds_[cellhashedindex]; 00445 00446 hoHits.push_back(HORecHit(detid,energy,0)); 00447 } 00448 00449 // HF 00450 nhits=firedCellsHF_.size(); 00451 for(unsigned ihit=0;ihit<nhits;++ihit) 00452 { 00453 unsigned cellhashedindex=firedCellsHF_[ihit]; 00454 // Check if it is above the threshold 00455 if(hcalRecHits_[cellhashedindex]<thresholdHF_) continue; 00456 00457 float energy=hcalRecHits_[cellhashedindex]; 00458 // poor man saturation 00459 if(energy>sat_[cellhashedindex]) energy=sat_[cellhashedindex]; 00460 00461 const HcalDetId & detid=theDetIds_[cellhashedindex]; 00462 hfHits.push_back(HFRecHit(detid,energy,0.)); 00463 if(doDigis_) 00464 { 00465 HFDataFrame myDataFrame(detid); 00466 myDataFrame.setSize(1); 00467 double nfc=hcalRecHits_[cellhashedindex]/gains_[cellhashedindex]+peds_[cellhashedindex]; 00468 int nadc=fCtoAdc(nfc/2.6); 00469 HcalQIESample qie(nadc, 0, 0, 0) ; 00470 myDataFrame.setSample(0,qie); 00471 hfDigis.push_back(myDataFrame); 00472 } 00473 } 00474 }
void HcalRecHitsMaker::loadPCaloHits | ( | const edm::Event & | iEvent | ) | [private] |
Definition at line 303 of file HcalRecHitsMaker.cc.
References clean(), counter(), Fill(), firedCellsHB_, firedCellsHE_, firedCellsHF_, firedCellsHO_, edm::Event::getByLabel(), HcalBarrel, HcalEndcap, HcalForward, HcalOuter, inputCol_, it, sistrip::extrainfo::noise_, noiseHB_, noiseHE_, noiseHF_, noiseHO_, and edm::Handle< T >::product().
Referenced by loadHcalRecHits().
00304 { 00305 00306 clean(); 00307 00308 edm::Handle<CrossingFrame<PCaloHit> > cf; 00309 iEvent.getByLabel(inputCol_,cf); 00310 std::auto_ptr<MixCollection<PCaloHit> > colcalo(new MixCollection<PCaloHit>(cf.product(),std::pair<int,int>(0,0) )); 00311 00312 MixCollection<PCaloHit>::iterator it=colcalo->begin();; 00313 MixCollection<PCaloHit>::iterator itend=colcalo->end(); 00314 unsigned counter=0; 00315 for(;it!=itend;++it) 00316 { 00317 HcalDetId detid(it->id()); 00318 int hashedindex=detid.hashed_index(); 00319 double noise_ = 0.; 00320 switch(detid.subdet()) 00321 { 00322 case HcalBarrel: 00323 { 00324 noise_ = noiseHB_; 00325 Fill(hashedindex,it->energy(),firedCellsHB_,noise_); 00326 } 00327 break; 00328 case HcalEndcap: 00329 { 00330 noise_ = noiseHE_; 00331 Fill(hashedindex,it->energy(),firedCellsHE_,noise_); 00332 } 00333 break; 00334 case HcalOuter: 00335 { 00336 // std::cout << detid << " " << it->energy() << std::endl; 00337 noise_ = noiseHO_; 00338 Fill(hashedindex,it->energy(),firedCellsHO_,noise_); 00339 } 00340 break; 00341 case HcalForward: 00342 { 00343 noise_ = noiseHF_; 00344 Fill(hashedindex,it->energy(),firedCellsHF_,noise_); 00345 } 00346 break; 00347 default: 00348 edm::LogWarning("CaloRecHitsProducer") << "RecHit not registered\n"; 00349 ; 00350 } 00351 ++counter; 00352 } 00353 }
double HcalRecHitsMaker::noiseInfCfromDB | ( | const HcalDbService * | conditions, | |
const HcalDetId & | detId | |||
) | [private] |
Definition at line 661 of file HcalRecHitsMaker.cc.
References DeDxTools::beta(), f, HcalPedestalWidth::getSigma(), funct::sqrt(), and HcalDetId::subdet().
Referenced by init().
00662 { 00663 // method from Salavat 00664 // fetch pedestal widths (noise sigmas for all 4 CapID) 00665 const HcalPedestalWidth* pedWidth = 00666 conditions-> getPedestalWidth(detId); 00667 double ssqq_1 = pedWidth->getSigma(0,0); 00668 double ssqq_2 = pedWidth->getSigma(1,1); 00669 double ssqq_3 = pedWidth->getSigma(2,2); 00670 double ssqq_4 = pedWidth->getSigma(3,3); 00671 00672 // correction factors (hb,he,ho,hf) 00673 static float corrfac[4]={1.25,1.20,1.40,0.67}; 00674 00675 int sub = detId.subdet(); 00676 00677 // effective RecHits (for this particular detId) noise calculation : 00678 00679 double sig_sq_mean = 0.25 * ( ssqq_1 + ssqq_2 + ssqq_3 + ssqq_4); 00680 00681 // f - external parameter (currently set to 0.5 in the FullSim) !!! 00682 double f=0.5; 00683 00684 double term = sqrt (1. + sqrt(1. - 2.*f*f)); 00685 double alpha = sqrt (0.5) * term; 00686 double beta = sqrt (0.5) * f / term; 00687 00688 double RMS1 = sqrt(sig_sq_mean) * sqrt (2.*beta*beta + alpha*alpha) ; 00689 double RMS4 = sqrt(sig_sq_mean) * sqrt (2.*beta*beta + 2.*(alpha-beta)*(alpha-beta) + 2.*(alpha-2.*beta)*(alpha-2.*beta)) ; 00690 00691 double noise_rms_fC; 00692 if(sub == 4) noise_rms_fC = RMS1; 00693 else noise_rms_fC = RMS4; 00694 00695 noise_rms_fC *= corrfac[sub-1]; 00696 00697 // to convert from above fC to GeV - multiply by gain (GeV/fC) 00698 // const HcalGain* gain = conditions->getGain(detId); 00699 // double noise_rms_GeV = noise_rms_fC * gain->getValue(0); // Noise RMS (GeV) for detId channel 00700 return noise_rms_fC; 00701 }
void HcalRecHitsMaker::noisify | ( | ) | [private] |
Definition at line 529 of file HcalRecHitsMaker.cc.
References lat::endl(), firedCellsHB_, firedCellsHE_, firedCellsHF_, firedCellsHO_, hbhi_, hcalHotFractionHB_, hcalHotFractionHE_, hcalHotFractionHF_, hcalHotFractionHO_, hcalRecHits_, hehi_, hfhi_, hohi_, myGaussianTailGeneratorHB_, myGaussianTailGeneratorHE_, myGaussianTailGeneratorHF_, myGaussianTailGeneratorHO_, nhbcells_, nhecells_, nhfcells_, nhocells_, noiseHB_, noiseHE_, noiseHF_, noiseHO_, noisifySubdet(), thresholdHB_, thresholdHE_, thresholdHF_, and thresholdHO_.
Referenced by loadHcalRecHits().
00530 { 00531 unsigned total=0; 00532 if(noiseHB_ > 0.) { 00533 if(firedCellsHB_.size()<nhbcells_) 00534 { 00535 // No need to do it anymore. The noise on the signal has been added 00536 // when loading the PCaloHits 00537 // noisifySignal(hbheRecHits_); 00538 total+=noisifySubdet(hcalRecHits_,firedCellsHB_,hbhi_,nhbcells_,hcalHotFractionHB_,myGaussianTailGeneratorHB_,noiseHB_,thresholdHB_); 00539 } 00540 else 00541 edm::LogWarning("CaloRecHitsProducer") << "All HCAL (HB) cells on ! " << std::endl; 00542 } 00543 00544 00545 if(noiseHE_ > 0.) { 00546 if(firedCellsHE_.size()<nhecells_) 00547 { 00548 // No need to do it anymore. The noise on the signal has been added 00549 // when loading the PCaloHits 00550 // noisifySignal(hbheRecHits_); 00551 total+=noisifySubdet(hcalRecHits_,firedCellsHE_,hehi_,nhecells_,hcalHotFractionHE_,myGaussianTailGeneratorHE_,noiseHE_,thresholdHE_); 00552 } 00553 else 00554 edm::LogWarning("CaloRecHitsProducer") << "All HCAL (HE) cells on ! " << std::endl; 00555 } 00556 00557 if(noiseHO_ > 0.) { 00558 if( firedCellsHO_.size()<nhocells_) 00559 { 00560 // noisifySignal(hoRecHits_); 00561 total+=noisifySubdet(hcalRecHits_,firedCellsHO_,hohi_,nhocells_,hcalHotFractionHO_,myGaussianTailGeneratorHO_,noiseHO_,thresholdHO_); 00562 } 00563 else 00564 edm::LogWarning("CaloRecHitsProducer") << "All HCAL(HO) cells on ! " << std::endl; 00565 } 00566 00567 if(noiseHF_ > 0.) { 00568 if(firedCellsHF_.size()<nhfcells_) 00569 { 00570 // noisifySignal(hfRecHits_); 00571 total+=noisifySubdet(hcalRecHits_,firedCellsHF_,hfhi_,nhfcells_,hcalHotFractionHF_,myGaussianTailGeneratorHF_,noiseHF_,thresholdHF_); 00572 } 00573 else 00574 edm::LogWarning("CaloRecHitsProducer") << "All HCAL(HF) cells on ! " << std::endl; 00575 } 00576 00577 edm::LogInfo("CaloRecHitsProducer") << "CaloRecHitsProducer : added noise in "<< total << " HCAL cells " << std::endl; 00578 }
unsigned HcalRecHitsMaker::noisifySubdet | ( | std::vector< float > & | theMap, | |
std::vector< int > & | theHits, | |||
const std::vector< int > & | thecells, | |||
unsigned | ncells, | |||
double | hcalHotFraction_, | |||
const GaussianTail * | myGT, | |||
double | sigma, | |||
double | threshold | |||
) | [private] |
Definition at line 580 of file HcalRecHitsMaker.cc.
References RandomEngine::flatShoot(), RandomEngine::gaussShoot(), hcalRecHits_, mean(), noiseFromDb_, noisesigma_, RandomEngine::poissonShoot(), random_, and GaussianTail::shoot().
Referenced by noisify().
00581 { 00582 // If the fraction of "hot " is small, use an optimized method to inject noise only in noisy cells. The 30% has not been tuned 00583 00584 if(hcalHotFraction<0.3) 00585 { 00586 double mean = (double)(ncells-theHits.size())*hcalHotFraction; 00587 unsigned nhcal = random_->poissonShoot(mean); 00588 00589 unsigned ncell=0; 00590 unsigned cellindex=0; 00591 uint32_t cellhashedindex=0; 00592 00593 while(ncell < nhcal) 00594 { 00595 cellindex = (unsigned)(random_->flatShoot()*ncells); 00596 cellhashedindex = thecells[cellindex]; 00597 if(hcalRecHits_[cellhashedindex]==0.) // new cell 00598 { 00599 hcalRecHits_[cellhashedindex]=myGT->shoot(); 00600 theHits.push_back(cellhashedindex); 00601 ++ncell; 00602 } 00603 } 00604 return ncell; 00605 } 00606 else // otherwise, inject noise everywhere 00607 { 00608 uint32_t cellhashedindex=0; 00609 unsigned nhcal=thecells.size(); 00610 00611 00612 for(unsigned ncell=0;ncell<nhcal;++ncell) 00613 { 00614 cellhashedindex = thecells[ncell]; 00615 if(hcalRecHits_[cellhashedindex]==0.) // new cell 00616 { 00617 // if noise from db activated 00618 if(noiseFromDb_) sigma=noisesigma_[cellhashedindex]; 00619 00620 double noise =random_->gaussShoot(0.,sigma); 00621 if(noise>threshold) 00622 { 00623 hcalRecHits_[cellhashedindex]=noise; 00624 theHits.push_back(cellhashedindex); 00625 } 00626 } 00627 } 00628 return nhcal; 00629 } 00630 return 0; 00631 }
bool HcalRecHitsMaker::doDigis_ [private] |
Definition at line 58 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), init(), and loadHcalRecHits().
bool HcalRecHitsMaker::doMiscalib_ [private] |
bool HcalRecHitsMaker::doSaturation_ [private] |
std::vector<int> HcalRecHitsMaker::fctoadc_ [private] |
std::vector<int> HcalRecHitsMaker::firedCellsHB_ [private] |
Definition at line 68 of file HcalRecHitsMaker.h.
Referenced by clean(), loadHcalRecHits(), loadPCaloHits(), and noisify().
std::vector<int> HcalRecHitsMaker::firedCellsHE_ [private] |
Definition at line 69 of file HcalRecHitsMaker.h.
Referenced by clean(), loadHcalRecHits(), loadPCaloHits(), and noisify().
std::vector<int> HcalRecHitsMaker::firedCellsHF_ [private] |
Definition at line 71 of file HcalRecHitsMaker.h.
Referenced by clean(), loadHcalRecHits(), loadPCaloHits(), and noisify().
std::vector<int> HcalRecHitsMaker::firedCellsHO_ [private] |
Definition at line 70 of file HcalRecHitsMaker.h.
Referenced by clean(), loadHcalRecHits(), loadPCaloHits(), and noisify().
std::vector<float> HcalRecHitsMaker::gains_ [private] |
Definition at line 80 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), init(), and loadHcalRecHits().
std::vector<int> HcalRecHitsMaker::hbhi_ [private] |
Definition at line 88 of file HcalRecHitsMaker.h.
Referenced by createVectorsOfCells(), HcalRecHitsMaker(), init(), noisify(), and ~HcalRecHitsMaker().
std::string HcalRecHitsMaker::hcalfileinpath_ [private] |
double HcalRecHitsMaker::hcalHotFractionHB_ [private] |
double HcalRecHitsMaker::hcalHotFractionHE_ [private] |
double HcalRecHitsMaker::hcalHotFractionHF_ [private] |
double HcalRecHitsMaker::hcalHotFractionHO_ [private] |
std::vector<float> HcalRecHitsMaker::hcalRecHits_ [private] |
Definition at line 66 of file HcalRecHitsMaker.h.
Referenced by clean(), Fill(), init(), loadHcalRecHits(), noisify(), and noisifySubdet().
std::vector<int> HcalRecHitsMaker::hehi_ [private] |
Definition at line 89 of file HcalRecHitsMaker.h.
Referenced by createVectorsOfCells(), HcalRecHitsMaker(), init(), noisify(), and ~HcalRecHitsMaker().
std::vector<int> HcalRecHitsMaker::hfhi_ [private] |
Definition at line 91 of file HcalRecHitsMaker.h.
Referenced by createVectorsOfCells(), HcalRecHitsMaker(), init(), noisify(), and ~HcalRecHitsMaker().
std::vector<int> HcalRecHitsMaker::hohi_ [private] |
Definition at line 90 of file HcalRecHitsMaker.h.
Referenced by createVectorsOfCells(), HcalRecHitsMaker(), init(), noisify(), and ~HcalRecHitsMaker().
bool HcalRecHitsMaker::initialized_ [private] |
edm::InputTag HcalRecHitsMaker::inputCol_ [private] |
Definition at line 56 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), and loadPCaloHits().
unsigned HcalRecHitsMaker::maxIndex_ [private] |
Definition at line 86 of file HcalRecHitsMaker.h.
Referenced by createVectorOfSubdetectorCells(), HcalRecHitsMaker(), and init().
unsigned HcalRecHitsMaker::maxIndexDebug_ [private] |
std::vector<float> HcalRecHitsMaker::miscalib_ [private] |
const HcalTPGCoder* HcalRecHitsMaker::myCoder_ [private] |
Definition at line 103 of file HcalRecHitsMaker.h.
const GaussianTail* HcalRecHitsMaker::myGaussianTailGeneratorHB_ [private] |
Definition at line 98 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), noisify(), and ~HcalRecHitsMaker().
const GaussianTail* HcalRecHitsMaker::myGaussianTailGeneratorHE_ [private] |
Definition at line 99 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), noisify(), and ~HcalRecHitsMaker().
const GaussianTail* HcalRecHitsMaker::myGaussianTailGeneratorHF_ [private] |
Definition at line 101 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), noisify(), and ~HcalRecHitsMaker().
const GaussianTail* HcalRecHitsMaker::myGaussianTailGeneratorHO_ [private] |
Definition at line 100 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), noisify(), and ~HcalRecHitsMaker().
unsigned HcalRecHitsMaker::nhbcells_ [private] |
Definition at line 92 of file HcalRecHitsMaker.h.
Referenced by createVectorsOfCells(), init(), and noisify().
unsigned HcalRecHitsMaker::nhecells_ [private] |
Definition at line 93 of file HcalRecHitsMaker.h.
Referenced by createVectorsOfCells(), init(), and noisify().
unsigned HcalRecHitsMaker::nhfcells_ [private] |
Definition at line 95 of file HcalRecHitsMaker.h.
Referenced by createVectorsOfCells(), init(), and noisify().
unsigned HcalRecHitsMaker::nhocells_ [private] |
Definition at line 94 of file HcalRecHitsMaker.h.
Referenced by createVectorsOfCells(), init(), and noisify().
bool HcalRecHitsMaker::noiseFromDb_ [private] |
Definition at line 61 of file HcalRecHitsMaker.h.
Referenced by Fill(), HcalRecHitsMaker(), and noisifySubdet().
float HcalRecHitsMaker::noiseHB_ [private] |
Definition at line 52 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), loadPCaloHits(), and noisify().
float HcalRecHitsMaker::noiseHE_ [private] |
Definition at line 52 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), loadPCaloHits(), and noisify().
float HcalRecHitsMaker::noiseHF_ [private] |
Definition at line 52 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), loadPCaloHits(), and noisify().
float HcalRecHitsMaker::noiseHO_ [private] |
Definition at line 52 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), loadPCaloHits(), and noisify().
std::vector<float> HcalRecHitsMaker::noisesigma_ [private] |
Definition at line 82 of file HcalRecHitsMaker.h.
Referenced by Fill(), HcalRecHitsMaker(), init(), and noisifySubdet().
std::vector<float> HcalRecHitsMaker::peds_ [private] |
Definition at line 79 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), init(), and loadHcalRecHits().
const RandomEngine* HcalRecHitsMaker::random_ [private] |
Definition at line 97 of file HcalRecHitsMaker.h.
Referenced by Fill(), HcalRecHitsMaker(), and noisifySubdet().
double HcalRecHitsMaker::refactor_ [private] |
double HcalRecHitsMaker::refactor_mean_ [private] |
std::vector<float> HcalRecHitsMaker::sat_ [private] |
Definition at line 81 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), init(), and loadHcalRecHits().
std::vector<HcalDetId> HcalRecHitsMaker::theDetIds_ [private] |
Definition at line 73 of file HcalRecHitsMaker.h.
Referenced by createVectorOfSubdetectorCells(), HcalRecHitsMaker(), init(), loadHcalRecHits(), and ~HcalRecHitsMaker().
float HcalRecHitsMaker::thresholdHB_ [private] |
Definition at line 51 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), loadHcalRecHits(), and noisify().
float HcalRecHitsMaker::thresholdHE_ [private] |
Definition at line 51 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), loadHcalRecHits(), and noisify().
float HcalRecHitsMaker::thresholdHF_ [private] |
Definition at line 51 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), loadHcalRecHits(), and noisify().
float HcalRecHitsMaker::thresholdHO_ [private] |
Definition at line 51 of file HcalRecHitsMaker.h.
Referenced by HcalRecHitsMaker(), loadHcalRecHits(), and noisify().
std::vector<float> HcalRecHitsMaker::TPGFactor_ [private] |