#include <iostream>
#include <fstream>
#include <algorithm>
#include <numeric>
#include "TString.h"
#include "Calibration/HcalCalibAlgos/interface/hcalCalibUtils.h"
Go to the source code of this file.
Functions | |
void | combinePhi (vector< TCell > &selectCells) |
void | combinePhi (vector< TCell > &selectCells, vector< TCell > &combinedCells) |
void | filterCells3x3 (vector< TCell > &selectCells, Int_t iEtaMaxE, UInt_t iPhiMaxE) |
void | filterCells5x5 (vector< TCell > &selectCells, Int_t iEtaMaxE, UInt_t iPhiMaxE) |
void | filterCellsInCone (std::vector< TCell > &selectCells, const GlobalPoint hitPositionHcal, Float_t maxConeDist, const CaloGeometry *theCaloGeometry) |
double | getDistInPlaneSimple (const GlobalPoint caloPoint, const GlobalPoint rechitPoint) |
void | getIEtaIPhiForHighestE (vector< TCell > &selectCells, Int_t &iEtaMostE, UInt_t &iPhiMostE) |
void | sumDepths (vector< TCell > &selectCells) |
void | sumSmallDepths (vector< TCell > &selectCells) |
void combinePhi | ( | vector< TCell > & | selectCells | ) |
Definition at line 91 of file hcalCalibUtils.cc.
{ // Map: NxN -> N cluster // Comine the targetE of cells with the same iEta if (selectCells.size()==0) return; // new container for the TCells // dummy cell id created with iEta; iPhi=1; depth // if combinePhi() is run after combining depths, depth=1 vector<TCell> combinedCells; map<UInt_t, vector<Float_t> > etaSliceE; // keyed by id of cell with iEta and **iPhi=1** // map the cells to the eta ring vector<TCell>::iterator i_it = selectCells.begin(); for (; i_it != selectCells.end(); ++i_it) { DetId id = HcalDetId(i_it->id()); UInt_t thisKey = HcalDetId(HcalDetId(id).subdet(), HcalDetId(id).ieta(), 1, HcalDetId(id).depth() ); etaSliceE[thisKey].push_back(i_it->e()); } map<UInt_t, vector<Float_t> >::iterator m_it = etaSliceE.begin(); for (; m_it != etaSliceE.end(); ++m_it) { combinedCells.push_back(TCell(m_it->first, accumulate(m_it->second.begin(), m_it->second.end(), 0.0) ) ); } // replace the original TCell vector with the new one selectCells = combinedCells; }
Definition at line 124 of file hcalCalibUtils.cc.
{ // Map: NxN -> N cluster // Comine the targetE of cells with the same iEta if (selectCells.size()==0) return; map<UInt_t, vector<Float_t> > etaSliceE; // keyed by id of cell with iEta and **iPhi=1** // map the cells to the eta ring vector<TCell>::iterator i_it = selectCells.begin(); for (; i_it != selectCells.end(); ++i_it) { DetId id = HcalDetId(i_it->id()); UInt_t thisKey = HcalDetId(HcalDetId(id).subdet(), HcalDetId(id).ieta(), 1, HcalDetId(id).depth() ); etaSliceE[thisKey].push_back(i_it->e()); } map<UInt_t, vector<Float_t> >::iterator m_it = etaSliceE.begin(); for (; m_it != etaSliceE.end(); ++m_it) { combinedCells.push_back(TCell(m_it->first, accumulate(m_it->second.begin(), m_it->second.end(), 0.0) ) ); } }
void filterCells3x3 | ( | vector< TCell > & | selectCells, |
Int_t | iEtaMaxE, | ||
UInt_t | iPhiMaxE | ||
) |
Definition at line 182 of file hcalCalibUtils.cc.
{ vector<TCell> filteredCells; Int_t dEta, dPhi; for (vector<TCell>::iterator it=selectCells.begin(); it!=selectCells.end(); ++it) { Bool_t passDEta = false; Bool_t passDPhi = false; dEta = HcalDetId(it->id()).ieta() - iEtaMaxE; dPhi = HcalDetId(it->id()).iphi() - iPhiMaxE; if (dPhi > 36) dPhi -= 72; if (dPhi < -36) dPhi += 72; if (abs(dEta)<=1 || (iEtaMaxE * HcalDetId(it->id()).ieta() == -1)) passDEta = true; if (abs(iEtaMaxE)<=20) { if (abs(HcalDetId(it->id()).ieta())<=20) { if (abs(dPhi)<=1) passDPhi = true; } else { // iPhi is labelled by odd numbers if (iPhiMaxE%2==0){ if (abs(dPhi)<=1) passDPhi = true; } else { if (dPhi== -2 || dPhi==0) passDPhi = true; } } } // if hottest cell with iEta<=20 else { if (abs(HcalDetId(it->id()).ieta())<=20) { if (abs(dPhi)<=1 || dPhi==2) passDPhi = true; } else { if (abs(dPhi)<=2) passDPhi = true; } } // if hottest cell with iEta>20 if (passDEta && passDPhi) filteredCells.push_back(*it); } selectCells = filteredCells; return; }
void filterCells5x5 | ( | vector< TCell > & | selectCells, |
Int_t | iEtaMaxE, | ||
UInt_t | iPhiMaxE | ||
) |
Definition at line 243 of file hcalCalibUtils.cc.
{ vector<TCell> filteredCells; Int_t dEta, dPhi; for (vector<TCell>::iterator it=selectCells.begin(); it!=selectCells.end(); ++it) { dEta = HcalDetId(it->id()).ieta() - iEtaMaxE; dPhi = HcalDetId(it->id()).iphi() - iPhiMaxE; if (dPhi > 36) dPhi -= 72; if (dPhi < -36) dPhi += 72; bool passDPhi = (abs(dPhi)<3); bool passDEta = (abs(dEta)<3 || (iEtaMaxE * HcalDetId(it->id()).ieta() == -2) ); // includes +/- eta boundary if (passDPhi && passDEta) filteredCells.push_back(*it); } selectCells = filteredCells; return; }
void filterCellsInCone | ( | std::vector< TCell > & | selectCells, |
const GlobalPoint | hitPositionHcal, | ||
Float_t | maxConeDist, | ||
const CaloGeometry * | theCaloGeometry | ||
) |
Definition at line 341 of file hcalCalibUtils.cc.
References getDistInPlaneSimple(), CaloGeometry::getPosition(), and diJetCalib::maxConeDist.
Referenced by hcalCalib::Process().
{ vector<TCell> filteredCells; for (vector<TCell>::iterator it=selectCells.begin(); it!=selectCells.end(); ++it) { const GlobalPoint recHitPoint = theCaloGeometry->getPosition(it->id()); if (getDistInPlaneSimple(hitPositionHcal, recHitPoint)<= maxConeDist) filteredCells.push_back(*it); } selectCells = filteredCells; return; }
double getDistInPlaneSimple | ( | const GlobalPoint | caloPoint, |
const GlobalPoint | rechitPoint | ||
) | [inline] |
Definition at line 361 of file hcalCalibUtils.cc.
References Vector3DBase< T, FrameTag >::dot(), PV3DBase< T, PVType, FrameType >::mag(), Vector3DBase< T, FrameTag >::unit(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().
{ // Simplified version of getDistInPlane // Assume track direction is origin -> point of hcal intersection const GlobalVector caloIntersectVector(caloPoint.x(), caloPoint.y(), caloPoint.z()); const GlobalVector caloIntersectUnitVector = caloIntersectVector.unit(); const GlobalVector rechitVector(rechitPoint.x(), rechitPoint.y(), rechitPoint.z()); const GlobalVector rechitUnitVector = rechitVector.unit(); double dotprod = caloIntersectUnitVector.dot(rechitUnitVector); double rechitdist = caloIntersectVector.mag()/dotprod; const GlobalVector effectiveRechitVector = rechitdist*rechitUnitVector; const GlobalPoint effectiveRechitPoint(effectiveRechitVector.x(), effectiveRechitVector.y(), effectiveRechitVector.z()); GlobalVector distance_vector = effectiveRechitPoint-caloPoint; if (dotprod > 0.) { return distance_vector.mag(); } else { return 999999.; } }
void getIEtaIPhiForHighestE | ( | vector< TCell > & | selectCells, |
Int_t & | iEtaMostE, | ||
UInt_t & | iPhiMostE | ||
) |
Definition at line 150 of file hcalCalibUtils.cc.
References diJetCalib::sumDepths.
{ vector<TCell> summedDepthsCells = selectCells; sumDepths(summedDepthsCells); vector<TCell>::iterator highCell = summedDepthsCells.begin(); // sum depths locally to get highest energy tower Float_t highE = -999; for (vector<TCell>::iterator it=summedDepthsCells.begin(); it!=summedDepthsCells.end(); ++it) { if (highE < it->e()) { highCell = it; highE = it->e(); } } iEtaMostE = HcalDetId(highCell->id()).ieta(); iPhiMostE = HcalDetId(highCell->id()).iphi(); return; }
void sumDepths | ( | vector< TCell > & | selectCells | ) |
Definition at line 11 of file hcalCalibUtils.cc.
References abs, gather_cfg::cout, and HcalBarrel.
{ // Assignes teh sum of the energy in cells with the same iEta, iPhi to the cell with depth=1. // All cells with depth>1 are removed form the container. If // the cell at depth=1 is not present: create it and follow the procedure. if (selectCells.size()==0) return; vector<TCell> selectCellsDepth1; vector<TCell> selectCellsHighDepth; // // NB: Here we add depth 3 for iEta==16 in *HE* to the value in the barrel // this approach is reflected in several of the following loops: make sure // to check when making changes. // // In some documents it is described as having depth 1, the mapping in CMSSW uses depth 3. for (vector<TCell>::iterator i_it = selectCells.begin(); i_it != selectCells.end(); ++i_it) { if (HcalDetId(i_it->id()).depth()==1) { selectCellsDepth1.push_back(*i_it); } else { selectCellsHighDepth.push_back(*i_it); } } // case where depth 1 has zero energy, but higher depths with same (iEta, iPhi) have energy. // For iEta<15 there is one depth -> selectCellsHighDepth is empty and we do not get in the loop. for (vector<TCell>::iterator i_it2 = selectCellsHighDepth.begin(); i_it2 != selectCellsHighDepth.end(); ++i_it2) { // protect against corrupt data if (HcalDetId(i_it2->id()).ietaAbs()<15 && HcalDetId(i_it2->id()).depth()>1) { cout << "ERROR!!! there are no HB cells with depth>1 for iEta<15!\n" << "Check the input data..." << endl; cout << "HCalDetId: " << HcalDetId(i_it2->id()) << endl; return; } bool foundDepthOne = false; for (vector<TCell>::iterator i_it = selectCellsDepth1.begin(); i_it != selectCellsDepth1.end(); ++i_it) { if (HcalDetId(i_it->id()).ieta()==HcalDetId(i_it2->id()).ieta() && HcalDetId(i_it->id()).iphi()==HcalDetId(i_it2->id()).iphi()) foundDepthOne = true; continue; } if (!foundDepthOne) { // create entry for depth 1 with 0 energy UInt_t newId; if (abs(HcalDetId(i_it2->id()).ieta())==16) newId = HcalDetId(HcalBarrel, HcalDetId(i_it2->id()).ieta(), HcalDetId(i_it2->id()).iphi(), 1); else newId = HcalDetId(HcalDetId(i_it2->id()).subdet(), HcalDetId(i_it2->id()).ieta(), HcalDetId(i_it2->id()).iphi(), 1); selectCellsDepth1.push_back(TCell(newId, 0.0)); } } for (vector<TCell>::iterator i_it = selectCellsDepth1.begin(); i_it != selectCellsDepth1.end(); ++i_it) { for (vector<TCell>::iterator i_it2 = selectCellsHighDepth.begin(); i_it2 != selectCellsHighDepth.end(); ++i_it2) { if (HcalDetId(i_it->id()).ieta()==HcalDetId(i_it2->id()).ieta() && HcalDetId(i_it->id()).iphi()==HcalDetId(i_it2->id()).iphi()) { i_it->SetE(i_it->e()+i_it2->e()); i_it2->SetE(0.0); // paranoid, aren't we... } } } // replace the original vectors with the new ones selectCells = selectCellsDepth1; return; }
void sumSmallDepths | ( | vector< TCell > & | selectCells | ) |
Definition at line 277 of file hcalCalibUtils.cc.
References spr::find().
{ if (selectCells.size()==0) return; vector<TCell> newCells; // holds unaffected cells to which the modified ones are added vector<TCell> manipulatedCells; // the ones that are combined for (vector<TCell>::iterator i_it = selectCells.begin(); i_it != selectCells.end(); ++i_it) { if ( (HcalDetId(i_it->id()).ietaAbs()==15 && HcalDetId(i_it->id()).depth()<=2) || (HcalDetId(i_it->id()).ietaAbs()==16 && HcalDetId(i_it->id()).depth()<=2)) { manipulatedCells.push_back(*i_it); } else { newCells.push_back(*i_it); } } // if the list is empty there is nothing to manipulate // leave the original vector unchanged if (manipulatedCells.size()<1) { newCells.clear(); return; } // See what cells are needed to hold the combined information: // Make holders for depth=1 for each (iEta,iPhi) // if a cell with these values is present in "manupulatedCells" vector<UInt_t> dummyIds; // to keep track of kreated cells vector<TCell> createdCells; // cells that need to be added or they exists; for (vector<TCell>::iterator i_it = manipulatedCells.begin(); i_it!=manipulatedCells.end(); ++i_it) { UInt_t dummyId = HcalDetId(HcalDetId(i_it->id()).subdet(), HcalDetId(i_it->id()).ieta(), HcalDetId(i_it->id()).iphi(), 1); if (find(dummyIds.begin(), dummyIds.end(), dummyId)==dummyIds.end()) { dummyIds.push_back(dummyId); createdCells.push_back(TCell(dummyId, 0.0)); } } for (vector<TCell>::iterator i_it = createdCells.begin(); i_it!=createdCells.end(); ++i_it) { for (vector<TCell>::iterator i_it2 = manipulatedCells.begin(); i_it2!=manipulatedCells.end(); ++i_it2) { if (HcalDetId(i_it->id()).ieta()==HcalDetId(i_it2->id()).ieta() && HcalDetId(i_it->id()).iphi()==HcalDetId(i_it2->id()).iphi() && HcalDetId(i_it2->id()).depth()<=2) { i_it->SetE(i_it->e()+i_it2->e()); } } } for (vector<TCell>::iterator i_it = createdCells.begin(); i_it!=createdCells.end(); ++i_it) { newCells.push_back(*i_it); } // replace the original vectors with the new ones selectCells = newCells; return; }