CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Public Attributes
hcalCalib Class Reference

#include <hcalCalib.h>

Inheritance diagram for hcalCalib:

Public Member Functions

virtual void Begin (TTree *tree)
 
void GetCoefFromMtrxInvOfAve ()
 
virtual Int_t GetEntry (Long64_t entry, Int_t getall=0)
 
virtual TList * GetOutputList () const
 
 hcalCalib (TTree *=0)
 
virtual void Init (TTree *tree)
 
void makeTextFile ()
 
virtual Bool_t Notify ()
 
virtual Bool_t Process (Long64_t entry)
 
Bool_t ReadPhiSymCor ()
 
void SetApplyPhiSymCorFlag (Bool_t b)
 
void SetCalibAbsIEtaMax (Int_t i)
 
void SetCalibAbsIEtaMin (Int_t i)
 
void SetCalibMethod (TString s)
 
void SetCalibType (TString s)
 
void SetCaloGeometry (const CaloGeometry *g, const HcalTopology *topo)
 
void SetCombinePhiFlag (Bool_t b)
 
void SetConeMaxDist (Float_t d)
 
void SetHbClusterSize (Int_t i)
 
void SetHeClusterSize (Int_t i)
 
void SetHistoFileName (TString filename)
 
virtual void SetInputList (TList *input)
 
void SetMaxEOverP (Float_t e)
 
void SetMaxEtThirdJet (Float_t et)
 
void SetMaxProbeJetEmFrac (Float_t f)
 
void SetMaxTagJetAbsEta (Float_t e)
 
void SetMaxTagJetEmFrac (Float_t f)
 
void SetMaxTargetE (Float_t e)
 
void SetMaxTrkEmE (Float_t e)
 
void SetMinCellE (Float_t e)
 
void SetMinDPhiDiJets (Float_t dphi)
 
void SetMinEOverP (Float_t e)
 
void SetMinProbeJetAbsEta (Float_t e)
 
void SetMinTagJetEt (Float_t e)
 
void SetMinTargetE (Float_t e)
 
virtual void SetObject (TObject *obj)
 
virtual void SetOption (const char *option)
 
void SetOutputCorCoefFileName (TString filename)
 
void SetPhiSymCorFileName (TString filename)
 
void SetSumDepthsFlag (Bool_t b)
 
void SetSumSmallDepthsFlag (Bool_t b)
 
void SetUseConeClustering (Bool_t b)
 
virtual void Terminate ()
 
virtual Int_t Version () const
 
virtual ~hcalCalib ()
 

Public Attributes

Bool_t APPLY_PHI_SYM_COR_FLAG
 
TBranch * b_cells
 
TBranch * b_emEnergy
 
TBranch * b_etVetoJet
 
TBranch * b_eventNumber
 
TBranch * b_iEtaHit
 
TBranch * b_iPhiHit
 
TBranch * b_probeJetEmFrac
 
TBranch * b_probeJetP4
 
TBranch * b_runNumber
 
TBranch * b_tagJetEmFrac
 
TBranch * b_tagJetP4
 
TBranch * b_targetE
 
TBranch * b_xTrkEcal
 
TBranch * b_xTrkHcal
 
TBranch * b_yTrkEcal
 
TBranch * b_yTrkHcal
 
TBranch * b_zTrkEcal
 
TBranch * b_zTrkHcal
 
Int_t CALIB_ABS_IETA_MAX
 
Int_t CALIB_ABS_IETA_MIN
 
TString CALIB_METHOD
 
TString CALIB_TYPE
 
std::vector< std::vector
< Float_t > > 
cellEnergies
 
std::vector< std::vector
< UInt_t > > 
cellIds
 
TClonesArray * cells
 
Bool_t COMBINE_PHI
 
Float_t emEnergy
 
Float_t etVetoJet
 
UInt_t eventNumber
 pointer to the analyzed TTree or TChain More...
 
TTree * fChain
 
Int_t HB_CLUSTER_SIZE
 
Int_t HE_CLUSTER_SIZE
 
TString HISTO_FILENAME
 
std::map< Int_t, Float_t > iEtaCoefMap
 
Int_t iEtaHit
 
UInt_t iPhiHit
 
Float_t MAX_CONE_DIST
 
Float_t MAX_EOVERP
 
Float_t MAX_ET_THIRD_JET
 
Float_t MAX_PROBEJET_EMFRAC
 
Float_t MAX_TAGJET_ABSETA
 
Float_t MAX_TAGJET_EMFRAC
 
Float_t MAX_TARGET_E
 
Float_t MAX_TRK_EME
 
Float_t MIN_CELL_E
 
Float_t MIN_DPHI_DIJETS
 
Float_t MIN_EOVERP
 
Float_t MIN_PROBEJET_ABSETA
 
Float_t MIN_TAGJET_ET
 
Float_t MIN_TARGET_E
 
TString OUTPUT_COR_COEF_FILENAME
 
TString PHI_SYM_COR_FILENAME
 
std::map< UInt_t, Float_t > phiSymCor
 
Float_t probeJetEmFrac
 
TLorentzVector * probeJetP4
 
std::vector< std::pair< Int_t,
UInt_t > > 
refIEtaIPhi
 
UInt_t runNumber
 
std::map< UInt_t, Float_t > solution
 
Bool_t SUM_DEPTHS
 
Bool_t SUM_SMALL_DEPTHS
 
Float_t tagJetEmFrac
 
TLorentzVector * tagJetP4
 
Float_t targetE
 
std::vector< Float_t > targetEnergies
 
const CaloGeometrytheCaloGeometry
 
const HcalTopologytopo_
 
Bool_t USE_CONE_CLUSTERING
 
Float_t xTrkEcal
 
Float_t xTrkHcal
 
Float_t yTrkEcal
 
Float_t yTrkHcal
 
Float_t zTrkEcal
 
Float_t zTrkHcal
 

Detailed Description

Definition at line 40 of file hcalCalib.h.

Constructor & Destructor Documentation

hcalCalib::hcalCalib ( TTree *  = 0)
inline

Definition at line 90 of file hcalCalib.h.

90 { }
virtual hcalCalib::~hcalCalib ( )
inlinevirtual

Definition at line 91 of file hcalCalib.h.

91 { }

Member Function Documentation

void hcalCalib::Begin ( TTree *  tree)
virtual

Definition at line 70 of file hcalCalib.cc.

References gather_cfg::cout, cmsRelvalreport::exit, h1_allTrkP, h1_corResp, h1_corRespBarrel, h1_corRespEndcap, h1_corRespIEta, h1_numEventsTwrIEta, h1_rawResp, h1_rawRespBarrel, h1_rawRespEndcap, h1_rawSumE, h1_selTrkP_iEta10, h1_trkP, h2_dHitRefBarrel, h2_dHitRefEndcap, histoFile, i, and nEvents.

70  {
71  TString option = GetOption();
72 
73  nEvents = 0;
74 
76  cout << "\nERROR: Failed to read the phi symmetry corrections." << endl;
77  cout << "Check if the filename is correct. If the corrections are not needed, set the corresponding flag to \"false\"\n" << endl;
78 
79  cout << "\nThe program will be terminated\n" << endl;
80 
81  exit(1);
82 
83  }
84 
85 
86  // cellEnergies.reserve(1000000);
87  // cellIds.reserve(1000000);
88  // targetEnergies.reserve(1000000);
89 
90  histoFile = new TFile(HISTO_FILENAME.Data(), "RECREATE");
91 
92 
93  h1_trkP = new TH1F("h1_trkP", "Track momenta; p_{trk} (GeV); Number of tracks", 100, 0, 200);
94  h1_allTrkP = new TH1F("h1_allTrkP", "Track momenta - all tracks; p_{trk} (GeV); Number of tracks", 100, 0, 200);
95 
96  h1_selTrkP_iEta10 = new TH1F("h1_selTrkP_iEta10", "Track momenta - tracks with |iEta|<10; p_{trk} (GeV); Number of tracks", 100, 0, 200);
97 
98 
99 
100  if (CALIB_TYPE=="ISO_TRACK")
101  h1_rawSumE = new TH1F("h1_rawSumE", "Cluster Energy; E_{cl} (GeV); Number of tracks", 100, 0, 200);
102  else
103  h1_rawSumE = new TH1F("h1_rawSumE", "Cluster Energy; E_{cl} (GeV); Number of tracks", 1000, 0, 2000);
104 
105 
106 
107  h1_rawResp = new TH1F("h1_rawResp", "Uncorrected response: |iEta|<24; E_{had}/p; Number of tracks", 300, 0, 3);
108  h1_corResp = new TH1F("h1_corResp", "Corrected response: |iEta|<24; E_{had}/p; Number of tracks", 300, 0, 3);
109 
110  h1_rawRespBarrel = new TH1F("h1_rawRespBarrel", "Uncorrected response: |iEta|<15; E_{had}/p; Number of tracks", 300, 0, 3);
111  h1_corRespBarrel = new TH1F("h1_corRespBarrel", "Corrected response: |iEta|<15; E_{had}/p; Number of tracks", 300, 0, 3);
112 
113  h1_rawRespEndcap = new TH1F("h1_rawRespEndcap", "Uncorrected response: 17<|iEta|<24; E_{had}/p; Number of tracks", 300, 0, 3);
114  h1_corRespEndcap = new TH1F("h1_corRespEndcap", "Corrected response: 17<|iEta|<24; E_{had}/p; Number of tracks", 300, 0, 3);
115 
116  h1_numEventsTwrIEta = new TH1F("h1_numEventsTwrIEta", "h1_numEventsTwrIEta", 80, -40, 40);
117 
118  h2_dHitRefBarrel = new TH2F("h2_dHitRefBarrel","{#Delta}i{#phi} vs {#Delta}i{#eta} of hit and most energetic tower(|i{#eta}|<16);{#Delta}i{#eta}; {#Delta}i{#phi}", 10, -5, 5, 10, -5, 5);
119  h2_dHitRefEndcap = new TH2F("h2_dHitRefEndcap","{#Delta}i{#phi} vs {#Delta}i{#eta} of hit and most energetic tower (16<|i{#eta}|<25) ;{#Delta}i{#eta}; {#Delta}i{#phi}", 10, -5, 5, 10, -5, 5);
120 
121 
122 
123  TString histoName = "isoTrack_";
124 
125  for (Int_t i=0; i<48; ++i) {
126  Long_t iEta;
127  if (i<24) iEta = i-24;
128  else iEta = i-23;
129  TString hn = histoName + iEta;
130  h1_corRespIEta[i] = new TH1F(hn, hn, 300, 0, 3.0);
131  }
132 
133 
134 } // end of Begin()
TH1F * h1_corRespBarrel
Definition: hcalCalib.cc:57
TFile * histoFile
Definition: hcalCalib.cc:45
TH1F * h1_corRespIEta[48]
Definition: hcalCalib.cc:67
int i
Definition: DBlmapReader.cc:9
TH1F * h1_rawRespBarrel
Definition: hcalCalib.cc:56
TH2F * h2_dHitRefEndcap
Definition: hcalCalib.cc:63
TH2F * h2_dHitRefBarrel
Definition: hcalCalib.cc:62
TString CALIB_TYPE
Definition: hcalCalib.h:138
TH1F * h1_allTrkP
Definition: hcalCalib.cc:49
TH1F * h1_numEventsTwrIEta
Definition: hcalCalib.cc:60
TH1F * h1_rawSumE
Definition: hcalCalib.cc:53
TH1F * h1_selTrkP_iEta10
Definition: hcalCalib.cc:51
Bool_t APPLY_PHI_SYM_COR_FLAG
Definition: hcalCalib.h:142
TH1F * h1_corRespEndcap
Definition: hcalCalib.cc:59
TH1F * h1_rawResp
Definition: hcalCalib.cc:54
TH1F * h1_rawRespEndcap
Definition: hcalCalib.cc:58
TH1F * h1_trkP
Definition: hcalCalib.cc:48
TH1F * h1_corResp
Definition: hcalCalib.cc:55
Bool_t ReadPhiSymCor()
Definition: hcalCalib.cc:615
tuple cout
Definition: gather_cfg.py:121
TString HISTO_FILENAME
Definition: hcalCalib.h:145
UInt_t nEvents
Definition: hcalCalib.cc:43
void hcalCalib::GetCoefFromMtrxInvOfAve ( )

Definition at line 522 of file hcalCalib.cc.

References funct::A, abs, b, begin, end, spr::find(), i, HcalDetId::ieta(), j, findQualityFiles::size, and x.

522  {
523 
524  // these maps are keyed by iEta
525  map<Int_t, Float_t> aveTargetE;
526  map<Int_t, Int_t> nEntries; // count hits
527 
528  // iEtaRef iEtaCell, energy
529  map<Int_t, map<Int_t, Float_t> > aveHitE; // add energies in the loop, normalize after that
530 
531  for (unsigned int i=0; i<cellEnergies.size(); ++i) {
532  Int_t iEtaRef = refIEtaIPhi[i].first;
533  aveTargetE[iEtaRef] += targetEnergies[i];
534  nEntries[iEtaRef]++;
535 
536  // for hybrid method: matrix inv of averages preceeded by L3
537  if (CALIB_METHOD=="L3_AND_MTRX_INV") {
538  for (unsigned int j=0; j<(cellEnergies[i]).size(); ++j) {
539  aveHitE[iEtaRef][HcalDetId(cellIds[i][j]).ieta()] += (solution[cellIds[i][j]] * cellEnergies[i][j]);
540  }
541  }
542  else if (CALIB_METHOD=="MATRIX_INV_OF_ETA_AVE") {
543  for (unsigned int j=0; j<(cellEnergies[i]).size(); ++j) {
544  aveHitE[iEtaRef][HcalDetId(cellIds[i][j]).ieta()] += cellEnergies[i][j];
545  }
546  }
547 
548  } // end of loop of entries
549 
550 
551  // scale by number of entries to get the averages
552  Float_t norm = 1.0;
553  for (map<Int_t, Float_t>::iterator m_it = aveTargetE.begin(); m_it!=aveTargetE.end(); ++m_it){
554  Int_t iEta = m_it->first;
555  norm = (nEntries[iEta] > 0)? 1.0/(nEntries[iEta]) : 1.0;
556  aveTargetE[iEta] *= norm;
557 
558  map<Int_t, Float_t>::iterator n_it = (aveHitE[iEta]).begin();
559 
560  Float_t sumRawE = 0;
561  for (; n_it!=(aveHitE[iEta]).end(); ++n_it) {
562  (n_it->second) *= norm;
563  sumRawE += (n_it->second);
564  }
565 
566  } // end of scaling by number of entries
567 
568  Int_t ONE_SIDE_IETA_RANGE = CALIB_ABS_IETA_MAX - CALIB_ABS_IETA_MIN +1;
569 
570  // conversion from iEta to index for the linear system
571  // contains elements only in the valid range for *matrix inversion*
572  vector<Int_t> iEtaList;
573 
574  for (Int_t i=-CALIB_ABS_IETA_MAX; i<=CALIB_ABS_IETA_MAX; ++i) {
575  if (abs(i)<CALIB_ABS_IETA_MIN) continue;
576  iEtaList.push_back(i);
577  }
578 
579  TMatrixD A(2*ONE_SIDE_IETA_RANGE, 2*ONE_SIDE_IETA_RANGE);
580  TMatrixD b(2*ONE_SIDE_IETA_RANGE, 1);
581  TMatrixD x(2*ONE_SIDE_IETA_RANGE, 1);
582 
583  for (Int_t i=0; i<2*ONE_SIDE_IETA_RANGE; ++i) {
584  for (Int_t j=0; j<2*ONE_SIDE_IETA_RANGE; ++j) {
585  A(i, j) = 0.0;
586  }
587  }
588 
589  for (UInt_t i=0; i<iEtaList.size(); ++i) {
590  b(i, 0) = aveTargetE[iEtaList[i]];
591 
592  map<Int_t, Float_t>::iterator n_it = aveHitE[iEtaList[i]].begin();
593  for (; n_it!=aveHitE[iEtaList[i]].end(); ++n_it){
594 
595  if (fabs(n_it->first)>CALIB_ABS_IETA_MAX || fabs(n_it->first)<CALIB_ABS_IETA_MIN) continue;
596  Int_t j = Int_t(find(iEtaList.begin(),iEtaList.end(), n_it->first) - iEtaList.begin());
597  A(i, j) = n_it->second;
598 
599  }
600 
601  }
602 
603  TMatrixD coef = b;
604  TDecompQRH qrh(A);
605  Bool_t hasSolution = qrh.MultiSolve(coef);
606 
607  if (hasSolution && (CALIB_METHOD=="L3_AND_MTRX_INV" || CALIB_METHOD=="MATRIX_INV_OF_ETA_AVE")) {
608  for (UInt_t i=0; i<iEtaList.size(); ++i) {
609  iEtaCoefMap[iEtaList[i]] = coef(i, 0);
610  }
611  }
612 }
int i
Definition: DBlmapReader.cc:9
std::map< Int_t, Float_t > iEtaCoefMap
Definition: hcalCalib.h:202
#define abs(x)
Definition: mlp_lapack.h:159
std::vector< std::pair< Int_t, UInt_t > > refIEtaIPhi
Definition: hcalCalib.h:196
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:7
Int_t CALIB_ABS_IETA_MIN
Definition: hcalCalib.h:129
std::vector< Float_t > targetEnergies
Definition: hcalCalib.h:197
int ieta() const
get the cell ieta
Definition: HcalDetId.h:38
int j
Definition: DBlmapReader.cc:9
#define end
Definition: vmac.h:38
std::map< UInt_t, Float_t > solution
Definition: hcalCalib.h:201
unsigned int UInt_t
Definition: FUTypes.h:12
TString CALIB_METHOD
Definition: hcalCalib.h:139
double b
Definition: hdecay.h:120
Int_t CALIB_ABS_IETA_MAX
Definition: hcalCalib.h:128
#define begin
Definition: vmac.h:31
Definition: DDAxes.h:10
std::vector< std::vector< Float_t > > cellEnergies
Definition: hcalCalib.h:194
std::vector< std::vector< UInt_t > > cellIds
Definition: hcalCalib.h:195
tuple size
Write out results.
virtual Int_t hcalCalib::GetEntry ( Long64_t  entry,
Int_t  getall = 0 
)
inlinevirtual

Definition at line 98 of file hcalCalib.h.

References fChain.

98 { return fChain ? fChain->GetTree()->GetEntry(entry, getall) : 0; }
TTree * fChain
Definition: hcalCalib.h:42
std::pair< std::string, MonitorElement * > entry
Definition: ME_MAP.h:8
virtual TList* hcalCalib::GetOutputList ( ) const
inlinevirtual

Definition at line 102 of file hcalCalib.h.

102 { return fOutput; }
virtual void hcalCalib::Init ( TTree *  tree)
virtual
void hcalCalib::makeTextFile ( )

Definition at line 674 of file hcalCalib.cc.

References abs, alignCSCRings::e, eta(), HcalEndcap, HcalOuter, phi, and sd.

674  {
675 
676  //{ HcalEmpty=0, HcalBarrel=1, HcalEndcap=2, HcalOuter=3, HcalForward=4, HcalTriggerTower=5, HcalOther=7 };
677 
678  TString sdName[8] = {"EMPTY", "HB", "HE", "HO", "HF", "TRITWR", "UNDEF", "OTHER"};
679 
680  FILE *constFile = fopen(OUTPUT_COR_COEF_FILENAME.Data(), "w+");
681 
682  // header of the constants file
683  fprintf(constFile, "%1s%16s%16s%16s%16s%9s%11s\n",
684  "#", "eta", "phi", "depth", "det", "value", "DetId");
685 
686  // Order loops to produce sequence of constants as for phi symmetry
687 
688  for(Int_t sd=1; sd<=4; sd++) {
689 
690  for(Int_t e=1; e<=50; e++) {
691  Int_t eta;
692 
693  for(Int_t side=0; side<2; side++) {
694  eta = (side==0)? -e : e; //ta *= (-1);
695 
696  for(Int_t phi=1; phi<=72; phi++) {
697 
698  for(Int_t d=1; d<5; d++) {
699 
700  if (!topo_->valid(HcalDetId(HcalSubdetector(sd), eta, phi, d))) continue;
701  HcalDetId id(HcalSubdetector(sd), eta, phi, d);
702  Float_t corrFactor = 1.0;
703 
704 
706  // if (abs(eta)>=CALIB_ABS_IETA_MIN && abs(eta)<=22 && HcalSubdetector(sd)!=HcalOuter) {
707 
708 
709  // need some care when depths were summed for iEta=16 =>
710  // the coeficients are saved in depth 1 of HB: affects
711  Int_t subdetInd = sd;
712 
713  if (abs(eta)==16 && HcalSubdetector(sd)==HcalEndcap && SUM_DEPTHS) {
714  subdetInd = 1;
715  }
716 
717  if (CALIB_METHOD=="L3" || CALIB_METHOD=="L3_AND_MTRX_INV") {
718 
719 
720  if (SUM_DEPTHS && COMBINE_PHI) corrFactor = solution[HcalDetId(HcalSubdetector(subdetInd), eta, 1, 1)];
721  else if (SUM_DEPTHS) corrFactor = solution[HcalDetId(HcalSubdetector(subdetInd), eta, phi, 1)];
722  else if (COMBINE_PHI) corrFactor = solution[HcalDetId(HcalSubdetector(sd), eta, 1, d)];
723  else corrFactor = solution[HcalDetId(HcalSubdetector(sd), eta, phi, d)];
724 
725 
726  // Remark: a new case was added (sumSmallDepths) where the first two depths in towers 15,16
727  // are summed and stored in depth 1.
728  // For now we create the correction coef for depth 2 (set equal to depth 1)
729  // after the call to the L3 minimizer so that this case is also handled without modifying the
730  // logic above. Probably it is better to move it here?
731 
732 
733  }// L3
734 
735  else if (CALIB_METHOD=="MATRIX_INV_OF_ETA_AVE") {
736  corrFactor = iEtaCoefMap[eta];
737  }
738 
739  else if (CALIB_METHOD=="ISO_TRK_PHI_SYM") {
740  corrFactor = solution[HcalDetId(HcalSubdetector(sd), eta, phi, d)];
741  }
742 
743 
744  } // if within the calibration range
745 
746  fprintf(constFile, "%17i%16i%16i%16s%9.5f%11X\n",
747  eta, phi, d, sdName[sd].Data(), corrFactor, id.rawId());
748 
749  }
750  }
751  }
752  }
753  }
754 
755  return;
756 }
std::map< Int_t, Float_t > iEtaCoefMap
Definition: hcalCalib.h:202
#define abs(x)
Definition: mlp_lapack.h:159
Bool_t SUM_DEPTHS
Definition: hcalCalib.h:118
T eta() const
Int_t CALIB_ABS_IETA_MIN
Definition: hcalCalib.h:129
HcalSubdetector
Definition: HcalAssistant.h:32
const HcalTopology * topo_
Definition: hcalCalib.h:149
std::map< UInt_t, Float_t > solution
Definition: hcalCalib.h:201
Bool_t COMBINE_PHI
Definition: hcalCalib.h:120
TString CALIB_METHOD
Definition: hcalCalib.h:139
TString OUTPUT_COR_COEF_FILENAME
Definition: hcalCalib.h:144
double sd
Int_t CALIB_ABS_IETA_MAX
Definition: hcalCalib.h:128
virtual bool valid(const DetId &id) const
Definition: HcalTopology.cc:58
Definition: DDAxes.h:10
virtual Bool_t hcalCalib::Notify ( )
virtual
Bool_t hcalCalib::Process ( Long64_t  entry)
virtual

Definition at line 143 of file hcalCalib.cc.

References abs, diJetCalib::combinePhi, gather_cfg::cout, TCell::e(), filterCells3x3(), filterCells5x5(), filterCellsInCone(), getIEtaIPhiForHighestE(), h1_allTrkP, h1_rawSumE, h1_selTrkP_iEta10, h2_dHitRefBarrel, h2_dHitRefEndcap, HcalBarrel, HcalEndcap, HcalForward, HcalOuter, i, TCell::id(), M_PI, nEvents, edm::second(), TCell::SetE(), diJetCalib::sumDepths, and diJetCalib::sumSmallDepths.

143  {
144 
145  // fChain->GetTree()->GetEntry(entry);
146  GetEntry(entry);
147 
148 
149  set<UInt_t> uniqueIds; // for testing: check if there are duplicate cells (AA)
150 
151 
152  Bool_t acceptEvent = kTRUE;
153 
154  ++nEvents;
155 
156  if (!(nEvents%100000)) cout << "event: " << nEvents << endl;
157 
158 
159  h1_allTrkP->Fill(targetE);
160 
161  if (targetE < MIN_TARGET_E || targetE > MAX_TARGET_E) return kFALSE;;
162 
163 
164  // make local copy as the cells may be modified due to phi/depth sum, phi corrections etc
165  vector<TCell> selectCells;
166 
167 
168  if (cells->GetSize()==0) return kFALSE;
169 
170  if (CALIB_TYPE=="DI_JET" && probeJetEmFrac > 0.999) return kTRUE;
171 
172 
173  for (Int_t i=0; i<cells->GetSize(); ++i) {
174  TCell* thisCell = (TCell*) cells->At(i);
175 
176  if (HcalDetId(thisCell->id()).subdet()== HcalOuter) continue; // reject HO, make a switch!
177 
178  if (HcalDetId(thisCell->id()).subdet() != HcalBarrel &&
179  HcalDetId(thisCell->id()).subdet() != HcalEndcap &&
180  HcalDetId(thisCell->id()).subdet() != HcalForward) {
181 
182  cout << "Unknown or wrong hcal subdetector: " << HcalDetId(thisCell->id()).subdet() << endl;
183 
184  }
185 
186 
187  // Apply phi symmetry corrections if the flag is set
188  if (APPLY_PHI_SYM_COR_FLAG) thisCell->SetE(phiSymCor[thisCell->id()] * thisCell->e());
189 
190  if (thisCell->e() > MIN_CELL_E) selectCells.push_back(*thisCell);
191  }
192 
193 
194  if (selectCells.size()==0) {
195  cout << "NO CELLS ABOVE THRESHOLD FOUND FOR TARGET!!!" << endl;
196  }
197 
198 
199  if (SUM_DEPTHS) sumDepths(selectCells);
200  else if (SUM_SMALL_DEPTHS) sumSmallDepths(selectCells); // depth 1,2 in twrs 15,16
201 
202 
203 
204  // most energetic tower (IsoTracks) or centroid of probe jet (DiJets)
205  pair<Int_t, UInt_t> refPos;
206 
207  Int_t dEtaHitRef = 999;
208  Int_t dPhiHitRef = 999;
209 
210  if (CALIB_TYPE=="ISO_TRACK") {
211 
212  Int_t iEtaMaxE; // filled by reference in getIEtaIPhiForHighestE
213  UInt_t iPhiMaxE; //
214 
215  getIEtaIPhiForHighestE(selectCells, iEtaMaxE, iPhiMaxE);
216 
217  dEtaHitRef = iEtaMaxE - iEtaHit;
218  dPhiHitRef = iPhiMaxE - iPhiHit;
219 
220  if (dPhiHitRef < -36) dPhiHitRef += 72;
221  if (dPhiHitRef > 36) dPhiHitRef -= 72;
222 
223  if (iEtaHit*iEtaMaxE < 0) {
224  if (dEtaHitRef<0) dEtaHitRef += 1;
225  if (dEtaHitRef>0) dEtaHitRef -= 1;
226  }
227 
228 
229  if (abs(iEtaHit)<16) h2_dHitRefBarrel->Fill(dEtaHitRef, dPhiHitRef);
230  if (abs(iEtaHit)>16 && abs(iEtaHit)<25) h2_dHitRefEndcap->Fill(dEtaHitRef, dPhiHitRef);
231 
232  // --------------------------------------------------
233  // Choice of cluster definition
234  //
235  // fixed size NxN clusters as specified in to config file
236  if (!USE_CONE_CLUSTERING) {
237  if (abs(iEtaMaxE)<16 && HB_CLUSTER_SIZE==3) filterCells3x3(selectCells, iEtaMaxE, iPhiMaxE);
238  if (abs(iEtaMaxE)>15 && HE_CLUSTER_SIZE==3) filterCells3x3(selectCells, iEtaMaxE, iPhiMaxE);
239 
240  if (abs(iEtaMaxE)<16 && HB_CLUSTER_SIZE==5) filterCells5x5(selectCells, iEtaMaxE, iPhiMaxE);
241  if (abs(iEtaMaxE)>15 && HE_CLUSTER_SIZE==5) filterCells5x5(selectCells, iEtaMaxE, iPhiMaxE);
242  }
243  else {
244  // calculate distance at hcal surface
245  const GlobalPoint hitPositionHcal(xTrkHcal, yTrkHcal, zTrkHcal);
246  filterCellsInCone(selectCells, hitPositionHcal, MAX_CONE_DIST, theCaloGeometry);
247  }
248 
249 
250 
251  refPos.first = iEtaMaxE;
252  refPos.second = iPhiMaxE;
253 
254  }
255  else if (CALIB_TYPE=="DI_JET") { // Apply selection cuts on DiJet events here
256 
257  if (etVetoJet>MAX_ET_THIRD_JET) acceptEvent = kFALSE;
258 
259  Float_t jetsDPhi = probeJetP4->DeltaPhi(*tagJetP4);
260  if (fabs(jetsDPhi * 180.0/M_PI) < MIN_DPHI_DIJETS) acceptEvent = kFALSE;
261 
262  if (probeJetEmFrac>MAX_PROBEJET_EMFRAC) acceptEvent = kFALSE;
263  if (fabs(probeJetP4->Eta())<MIN_PROBEJET_ABSETA) acceptEvent = kFALSE;
264  if (fabs(tagJetP4->Eta())>MAX_TAGJET_ABSETA) acceptEvent = kFALSE;
265  if (fabs(tagJetP4->Et())< MIN_TAGJET_ET) acceptEvent = kFALSE;
266 
267  if (acceptEvent) {
268 
269  Int_t iEtaMaxE; // filled by reference in getIEtaIPhiForHighestE
270  UInt_t iPhiMaxE; //
271 
272  getIEtaIPhiForHighestE(selectCells, iEtaMaxE, iPhiMaxE);
273 
274  // The ref position for the jet is not used in the minimization at this time.
275  // It will be needed if we attempt to do matrix inversion: then the question is
276  // which value is better suited: the centroid of the jet or the hottest tower...
277 
278  // refPos.first = iEtaHit;
279  // refPos.second = iPhiHit;
280 
281  refPos.first = iEtaMaxE;
282  refPos.second = iPhiMaxE;
283 
284  if (abs(iEtaMaxE)>40) acceptEvent = kFALSE; // for testing : set as parameter (AA)
285 
286  }
287 
288  }
289 
290 
291  if (COMBINE_PHI) combinePhi(selectCells);
292 
293  // fill the containers for the minimization prcedures
294  vector<Float_t> energies;
295  vector<UInt_t> ids;
296 
297  for (vector<TCell>::iterator i_it = selectCells.begin(); i_it!=selectCells.end(); ++i_it) {
298 
299  // for testing : fill only unique id's
300 
301  if (uniqueIds.insert(i_it->id()).second) {
302  energies.push_back(i_it->e());
303  ids.push_back(i_it->id());
304  }
305 
306  }
307 
308  if (CALIB_TYPE=="ISO_TRACK") {
309  if (accumulate(energies.begin(), energies.end(), 0.0) / targetE < MIN_EOVERP) acceptEvent=kFALSE;
310  if (accumulate(energies.begin(), energies.end(), 0.0) / targetE > MAX_EOVERP) acceptEvent=kFALSE;
311 
312  if (emEnergy > MAX_TRK_EME) acceptEvent=kFALSE;
313 
314  if (abs(dEtaHitRef)>1 || abs(dPhiHitRef)>1) acceptEvent=kFALSE;
315 
316  // Have to check if for |iEta|>20 (and neighboring region) the dPhiHitRef
317  // should be relaxed to 2. The neighboring towers have dPhi=2...
318 
319 
320  }
321 
322 
323  h1_rawSumE->Fill(accumulate(energies.begin(), energies.end(), 0.0));
324 
325 
326  // here we fill the information for the minimization procedure
327  if (acceptEvent) {
328  cellEnergies.push_back(energies);
329  cellIds.push_back(ids);
330  targetEnergies.push_back(targetE);
331  refIEtaIPhi.push_back(refPos);
332 
333  if (abs(refPos.first)<=10)
334  h1_selTrkP_iEta10->Fill(targetE);
335 
336 
337  }
338 
339 
340 
341 
342  // Clean up
343  energies.clear();
344  ids.clear();
345  selectCells.clear();
346 
347  return kTRUE;
348 }
void filterCellsInCone(std::vector< TCell > &selectCells, const GlobalPoint hitPositionHcal, Float_t maxConeDist, const CaloGeometry *theCaloGeometry)
int i
Definition: DBlmapReader.cc:9
Float_t MAX_TARGET_E
Definition: hcalCalib.h:108
tuple sumSmallDepths
Definition: diJetCalib.py:58
TH2F * h2_dHitRefEndcap
Definition: hcalCalib.cc:63
TH2F * h2_dHitRefBarrel
Definition: hcalCalib.cc:62
Float_t xTrkHcal
Definition: hcalCalib.h:53
TString CALIB_TYPE
Definition: hcalCalib.h:138
TLorentzVector * tagJetP4
Definition: hcalCalib.h:60
#define abs(x)
Definition: mlp_lapack.h:159
TH1F * h1_allTrkP
Definition: hcalCalib.cc:49
Float_t zTrkHcal
Definition: hcalCalib.h:55
std::map< UInt_t, Float_t > phiSymCor
Definition: hcalCalib.h:199
Float_t e()
Definition: TCell.h:26
Float_t MAX_EOVERP
Definition: hcalCalib.h:112
Float_t targetE
Definition: hcalCalib.h:50
Bool_t SUM_DEPTHS
Definition: hcalCalib.h:118
Definition: TCell.h:15
Float_t MIN_CELL_E
Definition: hcalCalib.h:110
TLorentzVector * probeJetP4
Definition: hcalCalib.h:61
std::vector< std::pair< Int_t, UInt_t > > refIEtaIPhi
Definition: hcalCalib.h:196
TH1F * h1_rawSumE
Definition: hcalCalib.cc:53
std::vector< Float_t > targetEnergies
Definition: hcalCalib.h:197
U second(std::pair< T, U > const &p)
TH1F * h1_selTrkP_iEta10
Definition: hcalCalib.cc:51
Int_t HE_CLUSTER_SIZE
Definition: hcalCalib.h:123
UInt_t id()
Definition: TCell.h:27
Int_t iEtaHit
Definition: hcalCalib.h:46
Float_t MAX_TRK_EME
Definition: hcalCalib.h:113
Float_t MIN_DPHI_DIJETS
Definition: hcalCalib.h:116
Int_t HB_CLUSTER_SIZE
Definition: hcalCalib.h:122
std::pair< std::string, MonitorElement * > entry
Definition: ME_MAP.h:8
Float_t MIN_TAGJET_ET
Definition: hcalCalib.h:134
Bool_t APPLY_PHI_SYM_COR_FLAG
Definition: hcalCalib.h:142
tuple combinePhi
Definition: diJetCalib.py:54
Float_t emEnergy
Definition: hcalCalib.h:49
Float_t MIN_EOVERP
Definition: hcalCalib.h:111
void SetE(Float_t e)
Definition: TCell.h:29
void filterCells3x3(std::vector< TCell > &selectCells, Int_t iEta, UInt_t iPhi)
Float_t MIN_PROBEJET_ABSETA
Definition: hcalCalib.h:136
Bool_t COMBINE_PHI
Definition: hcalCalib.h:120
tuple sumDepths
Definition: diJetCalib.py:51
Float_t MAX_ET_THIRD_JET
Definition: hcalCalib.h:115
UInt_t iPhiHit
Definition: hcalCalib.h:47
unsigned int UInt_t
Definition: FUTypes.h:12
Bool_t USE_CONE_CLUSTERING
Definition: hcalCalib.h:125
#define M_PI
Definition: BFit3D.cc:3
Float_t etVetoJet
Definition: hcalCalib.h:51
TClonesArray * cells
Definition: hcalCalib.h:48
void filterCells5x5(std::vector< TCell > &selectCells, Int_t iEta, UInt_t iPhi)
Bool_t SUM_SMALL_DEPTHS
Definition: hcalCalib.h:119
tuple cout
Definition: gather_cfg.py:121
Float_t MAX_TAGJET_ABSETA
Definition: hcalCalib.h:133
Float_t MAX_PROBEJET_EMFRAC
Definition: hcalCalib.h:131
const CaloGeometry * theCaloGeometry
Definition: hcalCalib.h:148
std::vector< std::vector< Float_t > > cellEnergies
Definition: hcalCalib.h:194
UInt_t nEvents
Definition: hcalCalib.cc:43
std::vector< std::vector< UInt_t > > cellIds
Definition: hcalCalib.h:195
void getIEtaIPhiForHighestE(std::vector< TCell > &selectCells, Int_t &iEta, UInt_t &iPhi)
Float_t MAX_CONE_DIST
Definition: hcalCalib.h:126
Float_t yTrkHcal
Definition: hcalCalib.h:54
Float_t probeJetEmFrac
Definition: hcalCalib.h:64
virtual Int_t GetEntry(Long64_t entry, Int_t getall=0)
Definition: hcalCalib.h:98
Bool_t hcalCalib::ReadPhiSymCor ( )

Definition at line 615 of file hcalCalib.cc.

References gather_cfg::cout, HcalBarrel, HcalEndcap, HcalForward, HcalOuter, geometryCSVtoXML::line, sd, AlCaHLTBitMon_QueryRunRegistry::string, and relativeConstraints::value.

615  {
616 
617  ifstream phiSymFile(PHI_SYM_COR_FILENAME.Data());
618 
619  if (!phiSymFile) {
620  cout << "\nERROR: Can not find file with phi symmetry constants \"" << PHI_SYM_COR_FILENAME.Data() << "\"" << endl;
621  return kFALSE;
622  }
623 
624  // assumes the format used in CSA08, first line is a comment
625 
626  Int_t iEta;
627  UInt_t iPhi;
628  UInt_t depth;
629  TString sdName;
630  UInt_t detId;
631 
632  Float_t value;
634 
636 
637  while (getline(phiSymFile, line)) {
638 
639  if(!line.size() || line[0]=='#') continue;
640 
641  std::istringstream linestream(line);
642  linestream >> iEta >> iPhi >> depth >> sdName >> value >> hex >> detId;
643 
644  if (sdName=="HB") sd = HcalBarrel;
645  else if (sdName=="HE") sd = HcalEndcap;
646  else if (sdName=="HO") sd = HcalOuter;
647  else if (sdName=="HF") sd = HcalForward;
648  else {
649  cout << "\nInvalid detector name in phi symmetry constants file: " << sdName.Data() << endl;
650  cout << "Check file and rerun!\n" << endl;
651  return kFALSE;
652  }
653 
654  // check if the data is consistent
655 
656  if (HcalDetId(sd, iEta, iPhi, depth) != HcalDetId(detId)) {
657  cout << "\nInconsistent info in phi symmetry file: subdet, iEta, iPhi, depth do not match rawId!\n" << endl;
658  return kFALSE;
659  }
660  HcalDetId hId(detId);
661  if (!topo_->valid(hId)) {
662  cout << "\nInvalid DetId from: iEta=" << iEta << " iPhi=" << iPhi << " depth=" << depth
663  << " subdet=" << sdName.Data() << " detId=" << detId << endl << endl;
664  return kFALSE;
665  }
666 
667  phiSymCor[HcalDetId(sd, iEta, iPhi, depth)] = value;
668  }
669 
670  return kTRUE;
671 }
TString PHI_SYM_COR_FILENAME
Definition: hcalCalib.h:141
std::map< UInt_t, Float_t > phiSymCor
Definition: hcalCalib.h:199
HcalSubdetector
Definition: HcalAssistant.h:32
const HcalTopology * topo_
Definition: hcalCalib.h:149
unsigned int UInt_t
Definition: FUTypes.h:12
double sd
virtual bool valid(const DetId &id) const
Definition: HcalTopology.cc:58
tuple cout
Definition: gather_cfg.py:121
void hcalCalib::SetApplyPhiSymCorFlag ( Bool_t  b)
inline

Definition at line 173 of file hcalCalib.h.

References APPLY_PHI_SYM_COR_FLAG, and b.

Referenced by HcalCalibrator::endJob().

Bool_t APPLY_PHI_SYM_COR_FLAG
Definition: hcalCalib.h:142
double b
Definition: hdecay.h:120
void hcalCalib::SetCalibAbsIEtaMax ( Int_t  i)
inline

Definition at line 169 of file hcalCalib.h.

References CALIB_ABS_IETA_MAX, and i.

Referenced by HcalCalibrator::endJob().

169 { CALIB_ABS_IETA_MAX = i; }
int i
Definition: DBlmapReader.cc:9
Int_t CALIB_ABS_IETA_MAX
Definition: hcalCalib.h:128
void hcalCalib::SetCalibAbsIEtaMin ( Int_t  i)
inline

Definition at line 170 of file hcalCalib.h.

References CALIB_ABS_IETA_MIN, and i.

Referenced by HcalCalibrator::endJob().

170 { CALIB_ABS_IETA_MIN = i; }
int i
Definition: DBlmapReader.cc:9
Int_t CALIB_ABS_IETA_MIN
Definition: hcalCalib.h:129
void hcalCalib::SetCalibMethod ( TString  s)
inline

Definition at line 162 of file hcalCalib.h.

References CALIB_METHOD, and alignCSCRings::s.

Referenced by HcalCalibrator::endJob().

162 { CALIB_METHOD = s; }
TString CALIB_METHOD
Definition: hcalCalib.h:139
void hcalCalib::SetCalibType ( TString  s)
inline

Definition at line 161 of file hcalCalib.h.

References CALIB_TYPE, and alignCSCRings::s.

Referenced by HcalCalibrator::endJob().

161 { CALIB_TYPE = s; }
TString CALIB_TYPE
Definition: hcalCalib.h:138
void hcalCalib::SetCaloGeometry ( const CaloGeometry g,
const HcalTopology topo 
)
inline

Definition at line 184 of file hcalCalib.h.

References g, theCaloGeometry, and topo_.

Referenced by HcalCalibrator::endJob().

184 { theCaloGeometry = g; topo_=topo; }
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e g
Definition: Activities.doc:4
const HcalTopology * topo_
Definition: hcalCalib.h:149
const CaloGeometry * theCaloGeometry
Definition: hcalCalib.h:148
void hcalCalib::SetCombinePhiFlag ( Bool_t  b)
inline

Definition at line 156 of file hcalCalib.h.

References b, and COMBINE_PHI.

Referenced by HcalCalibrator::endJob().

156 { COMBINE_PHI = b; }
Bool_t COMBINE_PHI
Definition: hcalCalib.h:120
double b
Definition: hdecay.h:120
void hcalCalib::SetConeMaxDist ( Float_t  d)
inline

Definition at line 167 of file hcalCalib.h.

References MAX_CONE_DIST.

Referenced by HcalCalibrator::endJob().

167 { MAX_CONE_DIST = d; }
Float_t MAX_CONE_DIST
Definition: hcalCalib.h:126
void hcalCalib::SetHbClusterSize ( Int_t  i)
inline

Definition at line 163 of file hcalCalib.h.

References HB_CLUSTER_SIZE, and i.

Referenced by HcalCalibrator::endJob().

163 { HB_CLUSTER_SIZE = i; }
int i
Definition: DBlmapReader.cc:9
Int_t HB_CLUSTER_SIZE
Definition: hcalCalib.h:122
void hcalCalib::SetHeClusterSize ( Int_t  i)
inline

Definition at line 164 of file hcalCalib.h.

References HE_CLUSTER_SIZE, and i.

Referenced by HcalCalibrator::endJob().

164 { HE_CLUSTER_SIZE = i; }
int i
Definition: DBlmapReader.cc:9
Int_t HE_CLUSTER_SIZE
Definition: hcalCalib.h:123
void hcalCalib::SetHistoFileName ( TString  filename)
inline

Definition at line 181 of file hcalCalib.h.

References lut2db_cfg::filename, and HISTO_FILENAME.

Referenced by HcalCalibrator::endJob().

tuple filename
Definition: lut2db_cfg.py:20
TString HISTO_FILENAME
Definition: hcalCalib.h:145
virtual void hcalCalib::SetInputList ( TList *  input)
inlinevirtual

Definition at line 101 of file hcalCalib.h.

References LaserDQM_cfg::input.

101 { fInput = input; }
void hcalCalib::SetMaxEOverP ( Float_t  e)
inline

Definition at line 159 of file hcalCalib.h.

References alignCSCRings::e, and MAX_EOVERP.

Referenced by HcalCalibrator::endJob().

159 { MAX_EOVERP = e; }
Float_t MAX_EOVERP
Definition: hcalCalib.h:112
void hcalCalib::SetMaxEtThirdJet ( Float_t  et)
inline

Definition at line 171 of file hcalCalib.h.

References MAX_ET_THIRD_JET.

Referenced by HcalCalibrator::endJob().

171 { MAX_ET_THIRD_JET = et; }
Float_t MAX_ET_THIRD_JET
Definition: hcalCalib.h:115
void hcalCalib::SetMaxProbeJetEmFrac ( Float_t  f)
inline

Definition at line 175 of file hcalCalib.h.

References f, and MAX_PROBEJET_EMFRAC.

Referenced by HcalCalibrator::endJob().

175 { MAX_PROBEJET_EMFRAC = f; }
double f[11][100]
Float_t MAX_PROBEJET_EMFRAC
Definition: hcalCalib.h:131
void hcalCalib::SetMaxTagJetAbsEta ( Float_t  e)
inline

Definition at line 177 of file hcalCalib.h.

References alignCSCRings::e, and MAX_TAGJET_ABSETA.

Referenced by HcalCalibrator::endJob().

177 { MAX_TAGJET_ABSETA = e; }
Float_t MAX_TAGJET_ABSETA
Definition: hcalCalib.h:133
void hcalCalib::SetMaxTagJetEmFrac ( Float_t  f)
inline

Definition at line 176 of file hcalCalib.h.

References f, and MAX_TAGJET_EMFRAC.

Referenced by HcalCalibrator::endJob().

176 { MAX_TAGJET_EMFRAC = f; }
double f[11][100]
Float_t MAX_TAGJET_EMFRAC
Definition: hcalCalib.h:132
void hcalCalib::SetMaxTargetE ( Float_t  e)
inline

Definition at line 153 of file hcalCalib.h.

References alignCSCRings::e, and MAX_TARGET_E.

Referenced by HcalCalibrator::endJob().

153 { MAX_TARGET_E = e; }
Float_t MAX_TARGET_E
Definition: hcalCalib.h:108
void hcalCalib::SetMaxTrkEmE ( Float_t  e)
inline

Definition at line 160 of file hcalCalib.h.

References alignCSCRings::e, and MAX_TRK_EME.

Referenced by HcalCalibrator::endJob().

160 { MAX_TRK_EME = e; }
Float_t MAX_TRK_EME
Definition: hcalCalib.h:113
void hcalCalib::SetMinCellE ( Float_t  e)
inline

Definition at line 157 of file hcalCalib.h.

References alignCSCRings::e, and MIN_CELL_E.

Referenced by HcalCalibrator::endJob().

157 { MIN_CELL_E = e; }
Float_t MIN_CELL_E
Definition: hcalCalib.h:110
void hcalCalib::SetMinDPhiDiJets ( Float_t  dphi)
inline

Definition at line 172 of file hcalCalib.h.

References MIN_DPHI_DIJETS.

Referenced by HcalCalibrator::endJob().

172 { MIN_DPHI_DIJETS = dphi; }
Float_t MIN_DPHI_DIJETS
Definition: hcalCalib.h:116
void hcalCalib::SetMinEOverP ( Float_t  e)
inline

Definition at line 158 of file hcalCalib.h.

References alignCSCRings::e, and MIN_EOVERP.

Referenced by HcalCalibrator::endJob().

158 { MIN_EOVERP = e; }
Float_t MIN_EOVERP
Definition: hcalCalib.h:111
void hcalCalib::SetMinProbeJetAbsEta ( Float_t  e)
inline

Definition at line 179 of file hcalCalib.h.

References alignCSCRings::e, and MIN_PROBEJET_ABSETA.

Referenced by HcalCalibrator::endJob().

179 { MIN_PROBEJET_ABSETA = e; }
Float_t MIN_PROBEJET_ABSETA
Definition: hcalCalib.h:136
void hcalCalib::SetMinTagJetEt ( Float_t  e)
inline

Definition at line 178 of file hcalCalib.h.

References alignCSCRings::e, and MIN_TAGJET_ET.

Referenced by HcalCalibrator::endJob().

178 { MIN_TAGJET_ET = e; }
Float_t MIN_TAGJET_ET
Definition: hcalCalib.h:134
void hcalCalib::SetMinTargetE ( Float_t  e)
inline

Definition at line 152 of file hcalCalib.h.

References alignCSCRings::e, and MIN_TARGET_E.

Referenced by HcalCalibrator::endJob().

152 { MIN_TARGET_E = e; }
Float_t MIN_TARGET_E
Definition: hcalCalib.h:107
virtual void hcalCalib::SetObject ( TObject *  obj)
inlinevirtual

Definition at line 100 of file hcalCalib.h.

References getGTfromDQMFile::obj.

100 { fObject = obj; }
virtual void hcalCalib::SetOption ( const char *  option)
inlinevirtual

Definition at line 99 of file hcalCalib.h.

99 { fOption = option; }
void hcalCalib::SetOutputCorCoefFileName ( TString  filename)
inline

Definition at line 180 of file hcalCalib.h.

References lut2db_cfg::filename, and OUTPUT_COR_COEF_FILENAME.

Referenced by HcalCalibrator::endJob().

TString OUTPUT_COR_COEF_FILENAME
Definition: hcalCalib.h:144
tuple filename
Definition: lut2db_cfg.py:20
void hcalCalib::SetPhiSymCorFileName ( TString  filename)
inline

Definition at line 174 of file hcalCalib.h.

References lut2db_cfg::filename, and PHI_SYM_COR_FILENAME.

Referenced by HcalCalibrator::endJob().

TString PHI_SYM_COR_FILENAME
Definition: hcalCalib.h:141
tuple filename
Definition: lut2db_cfg.py:20
void hcalCalib::SetSumDepthsFlag ( Bool_t  b)
inline

Definition at line 154 of file hcalCalib.h.

References b, and SUM_DEPTHS.

Referenced by HcalCalibrator::endJob().

154 { SUM_DEPTHS = b; }
Bool_t SUM_DEPTHS
Definition: hcalCalib.h:118
double b
Definition: hdecay.h:120
void hcalCalib::SetSumSmallDepthsFlag ( Bool_t  b)
inline

Definition at line 155 of file hcalCalib.h.

References b, and SUM_SMALL_DEPTHS.

Referenced by HcalCalibrator::endJob().

155 { SUM_SMALL_DEPTHS = b; }
double b
Definition: hdecay.h:120
Bool_t SUM_SMALL_DEPTHS
Definition: hcalCalib.h:119
void hcalCalib::SetUseConeClustering ( Bool_t  b)
inline

Definition at line 166 of file hcalCalib.h.

References b, and USE_CONE_CLUSTERING.

Referenced by HcalCalibrator::endJob().

166 { USE_CONE_CLUSTERING = b; }
Bool_t USE_CONE_CLUSTERING
Definition: hcalCalib.h:125
double b
Definition: hdecay.h:120
void hcalCalib::Terminate ( )
virtual

Definition at line 352 of file hcalCalib.cc.

References abs, gather_cfg::cout, HcalDetId::depth(), h1_allTrkP, h1_corResp, h1_corRespBarrel, h1_corRespEndcap, h1_corRespIEta, h1_numEventsTwrIEta, h1_rawResp, h1_rawRespBarrel, h1_rawRespEndcap, h1_rawSumE, h1_selTrkP_iEta10, h1_trkP, h2_dHitRefBarrel, h2_dHitRefEndcap, HcalBarrel, histoFile, i, HcalDetId::ieta(), HcalDetId::ietaAbs(), HcalDetId::iphi(), MinL3AlgoUniv< IDdet >::iterate(), j, findQualityFiles::size, HcalDetId::subdet(), and histoStyle::weight.

352  {
353 
354  cout << "\n\nFinished reading the events.\n";
355  cout << "Number of input objects: " << cellIds.size() << endl;
356  cout << "Performing minimization: depending on selected method can take some time...\n\n";
357 
358  for (vector<pair<Int_t, UInt_t> >::iterator it_rp=refIEtaIPhi.begin(); it_rp!=refIEtaIPhi.end(); ++it_rp ) {
359  Float_t weight = (abs(it_rp->first)<21)? 1.0/72.0 : 1.0/36.0;
360  h1_numEventsTwrIEta->Fill(it_rp->first, weight);
361  }
362 
363  if (CALIB_METHOD=="MATRIX_INV_OF_ETA_AVE") {
365  }
366  else if (CALIB_METHOD=="L3" || CALIB_METHOD=="L3_AND_MTRX_INV") {
367  int eventWeight = 2; // 2 is the default (try at some point 0,1,2,3)
368  MinL3AlgoUniv<UInt_t>* thisL3Algo = new MinL3AlgoUniv<UInt_t>(eventWeight);
369  int numIterations = 10; // default 10
370 
371  solution = thisL3Algo->iterate(cellEnergies, cellIds, targetEnergies, numIterations);
372 
373  // in order to handle the case where sumDepths="false", but the flag to sum depths 1,2 in HB towers 15, 16
374  // is set (sumSmallDepths) we create entries in "solution" to create equal correction here
375  // for each encountered coef in depth one.
376 
377  if (!SUM_DEPTHS && SUM_SMALL_DEPTHS) {
378  vector<UInt_t> idForSummedCells;
379 
380  for (map<UInt_t, Float_t>::iterator m_it = solution.begin(); m_it != solution.end(); ++m_it) {
381  if (HcalDetId(m_it->first).ietaAbs()!=15 && HcalDetId(m_it->first).ietaAbs()!=16) continue;
382  if (HcalDetId(m_it->first).subdet()!=HcalBarrel) continue;
383  if (HcalDetId(m_it->first).depth()==1)
384  idForSummedCells.push_back(HcalDetId(m_it->first));
385  }
386 
387  for (vector<UInt_t>::iterator v_it=idForSummedCells.begin(); v_it!=idForSummedCells.end(); ++v_it) {
388  UInt_t addCoefId = HcalDetId(HcalBarrel, HcalDetId(*v_it).ieta(), HcalDetId(*v_it).iphi(), 2);
389  solution[addCoefId] = solution[*v_it];
390  }
391 
392  } // end of special treatment for "sumSmallDepths" mode
393 
394 
395  if (CALIB_METHOD=="L3_AND_MTRX_INV") {
397 
398  // loop over the solution from L3 and multiply by the additional factor from
399  // the matrix inversion. Set coef outside of the valid calibration region =1.
400  for (map<UInt_t, Float_t>::iterator it_s=solution.begin(); it_s!=solution.end(); ++it_s) {
401  Int_t iEtaSol = HcalDetId(it_s->first).ieta();
402  if (abs(iEtaSol)<CALIB_ABS_IETA_MIN || abs(iEtaSol)> CALIB_ABS_IETA_MAX) it_s->second = 1.0;
403  else it_s->second *= iEtaCoefMap[iEtaSol];
404  }
405 
406  } // if CALIB_METHOD=="L3_AND_MTRX_INV"
407 
408  } // end of L3 or L3 + mtrx inversion minimization
409 
410  // done getting the constants -> write the formatted file
411  makeTextFile();
412 
413  // fill some histograms
414  Float_t rawResp = 0;
415  Float_t corResp = 0;
416  Int_t maxIEta = 0;
417  Int_t minIEta = 999;
418 
419 
420  for (unsigned int i=0; i<cellEnergies.size(); ++i) {
421  Int_t iEta;
422  for (unsigned int j=0; j<(cellEnergies[i]).size(); ++j) {
423  iEta = HcalDetId(cellIds[i][j]).ieta();
424  rawResp += (cellEnergies[i])[j];
425 
426  if (CALIB_METHOD=="L3_AND_MTRX_INV") {
427  corResp += (solution[cellIds[i][j]] * cellEnergies[i][j]);
428  }
429  else if (CALIB_METHOD=="L3") {
430  corResp += (solution[cellIds[i][j]] * (cellEnergies[i][j]));
431  }
432  else if (CALIB_METHOD=="MATRIX_INV_OF_ETA_AVE") {
433  corResp += (iEtaCoefMap[iEta]* cellEnergies[i][j]);
434  }
435 
436  if (maxIEta < abs(iEta) ) maxIEta = abs(iEta);
437  if (minIEta > abs(iEta) ) minIEta = abs(iEta);
438  }
439 
440  rawResp /= targetEnergies[i];
441  corResp /= targetEnergies[i];
442 
443 
444  // fill histograms based on iEta on ref point of the cluster (for now: hottest tower)
445  // expect range |iEta|<=24 (to do: add flexibility for arbitrary range)
446 
447  if (CALIB_TYPE=="ISO_TRACK") {
448  Int_t ind = refIEtaIPhi[i].first;
449  (ind<0) ? (ind +=24) : (ind +=23);
450  if (ind>=0 && ind<48) {
451  h1_corRespIEta[ind]->Fill(corResp);
452  }
453 
454  // fill histograms for cases where all towers are in the barrel or endcap
455  if (maxIEta<25) {
456  h1_rawResp->Fill(rawResp);
457  h1_corResp->Fill(corResp);
458  }
459  if (maxIEta<15) {
460  h1_rawRespBarrel->Fill(rawResp);
461  h1_corRespBarrel->Fill(corResp);
462  }
463  else if (maxIEta<25 && minIEta>16) {
464  h1_rawRespEndcap->Fill(rawResp);
465  h1_corRespEndcap->Fill(corResp);
466  }
467  } // histograms for isotrack calibration
468 
469 
470  else {
471 
472  // put jet plots here
473 
474  }
475 
476 
477 
478  rawResp = 0;
479  corResp = 0;
480  maxIEta = 0;
481  minIEta = 999;
482  }
483 
484 
485  // save the histograms
486  h1_trkP->Write();
487  h1_allTrkP->Write();
488 
489 h1_selTrkP_iEta10->Write();
490 
491  h1_rawSumE->Write();
492  h1_rawResp->Write();
493  h1_corResp->Write();
494  h1_rawRespBarrel->Write();
495  h1_corRespBarrel->Write();
496  h1_rawRespEndcap->Write();
497  h1_corRespEndcap->Write();
498  h1_numEventsTwrIEta->Write();
499  h2_dHitRefBarrel->Write();
500  h2_dHitRefEndcap->Write();
501  for (Int_t i=0; i<48; ++i) {
502  h1_corRespIEta[i]->Write();
503  }
504 
505 
506 
507 
508 
509  histoFile->Write();
510  histoFile->Close();
511 
512  cout << "\n Finished calibration.\n " << endl;
513 
514 
515 } // end of Terminate()
TH1F * h1_corRespBarrel
Definition: hcalCalib.cc:57
TFile * histoFile
Definition: hcalCalib.cc:45
TH1F * h1_corRespIEta[48]
Definition: hcalCalib.cc:67
int i
Definition: DBlmapReader.cc:9
std::map< Int_t, Float_t > iEtaCoefMap
Definition: hcalCalib.h:202
TH1F * h1_rawRespBarrel
Definition: hcalCalib.cc:56
HcalSubdetector subdet() const
get the subdetector
Definition: HcalDetId.h:32
void GetCoefFromMtrxInvOfAve()
Definition: hcalCalib.cc:522
TH2F * h2_dHitRefEndcap
Definition: hcalCalib.cc:63
TH2F * h2_dHitRefBarrel
Definition: hcalCalib.cc:62
TString CALIB_TYPE
Definition: hcalCalib.h:138
#define abs(x)
Definition: mlp_lapack.h:159
TH1F * h1_allTrkP
Definition: hcalCalib.cc:49
Bool_t SUM_DEPTHS
Definition: hcalCalib.h:118
TH1F * h1_numEventsTwrIEta
Definition: hcalCalib.cc:60
std::vector< std::pair< Int_t, UInt_t > > refIEtaIPhi
Definition: hcalCalib.h:196
TH1F * h1_rawSumE
Definition: hcalCalib.cc:53
Int_t CALIB_ABS_IETA_MIN
Definition: hcalCalib.h:129
std::vector< Float_t > targetEnergies
Definition: hcalCalib.h:197
TH1F * h1_selTrkP_iEta10
Definition: hcalCalib.cc:51
int depth() const
get the tower depth
Definition: HcalDetId.h:42
int ieta() const
get the cell ieta
Definition: HcalDetId.h:38
int j
Definition: DBlmapReader.cc:9
std::map< UInt_t, Float_t > solution
Definition: hcalCalib.h:201
int ietaAbs() const
get the absolute value of the cell ieta
Definition: HcalDetId.h:36
TH1F * h1_corRespEndcap
Definition: hcalCalib.cc:59
TH1F * h1_rawResp
Definition: hcalCalib.cc:54
int iphi() const
get the cell iphi
Definition: HcalDetId.h:40
unsigned int UInt_t
Definition: FUTypes.h:12
TString CALIB_METHOD
Definition: hcalCalib.h:139
TH1F * h1_rawRespEndcap
Definition: hcalCalib.cc:58
TH1F * h1_trkP
Definition: hcalCalib.cc:48
TH1F * h1_corResp
Definition: hcalCalib.cc:55
Int_t CALIB_ABS_IETA_MAX
Definition: hcalCalib.h:128
void makeTextFile()
Definition: hcalCalib.cc:674
Bool_t SUM_SMALL_DEPTHS
Definition: hcalCalib.h:119
tuple cout
Definition: gather_cfg.py:121
IDmap iterate(const std::vector< std::vector< float > > &eventMatrix, const std::vector< std::vector< IDdet > > &idMatrix, const std::vector< float > &energyVector, const int &nIter, const bool &normalizeFlag=false)
Definition: MinL3AlgoUniv.h:83
std::vector< std::vector< Float_t > > cellEnergies
Definition: hcalCalib.h:194
int weight
Definition: histoStyle.py:50
std::vector< std::vector< UInt_t > > cellIds
Definition: hcalCalib.h:195
tuple size
Write out results.
virtual Int_t hcalCalib::Version ( ) const
inlinevirtual

Definition at line 92 of file hcalCalib.h.

92 { return 2; }

Member Data Documentation

Bool_t hcalCalib::APPLY_PHI_SYM_COR_FLAG

Definition at line 142 of file hcalCalib.h.

Referenced by SetApplyPhiSymCorFlag().

TBranch* hcalCalib::b_cells

Definition at line 71 of file hcalCalib.h.

TBranch* hcalCalib::b_emEnergy

Definition at line 72 of file hcalCalib.h.

TBranch* hcalCalib::b_etVetoJet

Definition at line 74 of file hcalCalib.h.

TBranch* hcalCalib::b_eventNumber

Definition at line 67 of file hcalCalib.h.

TBranch* hcalCalib::b_iEtaHit

Definition at line 69 of file hcalCalib.h.

TBranch* hcalCalib::b_iPhiHit

Definition at line 70 of file hcalCalib.h.

TBranch* hcalCalib::b_probeJetEmFrac

Definition at line 84 of file hcalCalib.h.

TBranch* hcalCalib::b_probeJetP4

Definition at line 87 of file hcalCalib.h.

TBranch* hcalCalib::b_runNumber

Definition at line 68 of file hcalCalib.h.

TBranch* hcalCalib::b_tagJetEmFrac

Definition at line 83 of file hcalCalib.h.

TBranch* hcalCalib::b_tagJetP4

Definition at line 86 of file hcalCalib.h.

TBranch* hcalCalib::b_targetE

Definition at line 73 of file hcalCalib.h.

TBranch* hcalCalib::b_xTrkEcal

Definition at line 79 of file hcalCalib.h.

TBranch* hcalCalib::b_xTrkHcal

Definition at line 76 of file hcalCalib.h.

TBranch* hcalCalib::b_yTrkEcal

Definition at line 80 of file hcalCalib.h.

TBranch* hcalCalib::b_yTrkHcal

Definition at line 77 of file hcalCalib.h.

TBranch* hcalCalib::b_zTrkEcal

Definition at line 81 of file hcalCalib.h.

TBranch* hcalCalib::b_zTrkHcal

Definition at line 78 of file hcalCalib.h.

Int_t hcalCalib::CALIB_ABS_IETA_MAX

Definition at line 128 of file hcalCalib.h.

Referenced by SetCalibAbsIEtaMax().

Int_t hcalCalib::CALIB_ABS_IETA_MIN

Definition at line 129 of file hcalCalib.h.

Referenced by SetCalibAbsIEtaMin().

TString hcalCalib::CALIB_METHOD

Definition at line 139 of file hcalCalib.h.

Referenced by SetCalibMethod().

TString hcalCalib::CALIB_TYPE

Definition at line 138 of file hcalCalib.h.

Referenced by SetCalibType().

std::vector< std::vector<Float_t> > hcalCalib::cellEnergies

Definition at line 194 of file hcalCalib.h.

std::vector< std::vector<UInt_t> > hcalCalib::cellIds

Definition at line 195 of file hcalCalib.h.

TClonesArray* hcalCalib::cells

Definition at line 48 of file hcalCalib.h.

Bool_t hcalCalib::COMBINE_PHI

Definition at line 120 of file hcalCalib.h.

Referenced by SetCombinePhiFlag().

Float_t hcalCalib::emEnergy

Definition at line 49 of file hcalCalib.h.

Float_t hcalCalib::etVetoJet

Definition at line 51 of file hcalCalib.h.

UInt_t hcalCalib::eventNumber

pointer to the analyzed TTree or TChain

Definition at line 44 of file hcalCalib.h.

Referenced by Vispa.Plugins.EdmBrowser.EdmDataAccessor.EdmDataAccessor::setFilterBranches().

TTree* hcalCalib::fChain

Definition at line 42 of file hcalCalib.h.

Referenced by GetEntry().

Int_t hcalCalib::HB_CLUSTER_SIZE

Definition at line 122 of file hcalCalib.h.

Referenced by SetHbClusterSize().

Int_t hcalCalib::HE_CLUSTER_SIZE

Definition at line 123 of file hcalCalib.h.

Referenced by SetHeClusterSize().

TString hcalCalib::HISTO_FILENAME

Definition at line 145 of file hcalCalib.h.

Referenced by SetHistoFileName().

std::map<Int_t, Float_t> hcalCalib::iEtaCoefMap

Definition at line 202 of file hcalCalib.h.

Int_t hcalCalib::iEtaHit

Definition at line 46 of file hcalCalib.h.

UInt_t hcalCalib::iPhiHit

Definition at line 47 of file hcalCalib.h.

Float_t hcalCalib::MAX_CONE_DIST

Definition at line 126 of file hcalCalib.h.

Referenced by SetConeMaxDist().

Float_t hcalCalib::MAX_EOVERP

Definition at line 112 of file hcalCalib.h.

Referenced by SetMaxEOverP().

Float_t hcalCalib::MAX_ET_THIRD_JET

Definition at line 115 of file hcalCalib.h.

Referenced by SetMaxEtThirdJet().

Float_t hcalCalib::MAX_PROBEJET_EMFRAC

Definition at line 131 of file hcalCalib.h.

Referenced by SetMaxProbeJetEmFrac().

Float_t hcalCalib::MAX_TAGJET_ABSETA

Definition at line 133 of file hcalCalib.h.

Referenced by SetMaxTagJetAbsEta().

Float_t hcalCalib::MAX_TAGJET_EMFRAC

Definition at line 132 of file hcalCalib.h.

Referenced by SetMaxTagJetEmFrac().

Float_t hcalCalib::MAX_TARGET_E

Definition at line 108 of file hcalCalib.h.

Referenced by SetMaxTargetE().

Float_t hcalCalib::MAX_TRK_EME

Definition at line 113 of file hcalCalib.h.

Referenced by SetMaxTrkEmE().

Float_t hcalCalib::MIN_CELL_E

Definition at line 110 of file hcalCalib.h.

Referenced by SetMinCellE().

Float_t hcalCalib::MIN_DPHI_DIJETS

Definition at line 116 of file hcalCalib.h.

Referenced by SetMinDPhiDiJets().

Float_t hcalCalib::MIN_EOVERP

Definition at line 111 of file hcalCalib.h.

Referenced by SetMinEOverP().

Float_t hcalCalib::MIN_PROBEJET_ABSETA

Definition at line 136 of file hcalCalib.h.

Referenced by SetMinProbeJetAbsEta().

Float_t hcalCalib::MIN_TAGJET_ET

Definition at line 134 of file hcalCalib.h.

Referenced by SetMinTagJetEt().

Float_t hcalCalib::MIN_TARGET_E

Definition at line 107 of file hcalCalib.h.

Referenced by SetMinTargetE().

TString hcalCalib::OUTPUT_COR_COEF_FILENAME

Definition at line 144 of file hcalCalib.h.

Referenced by SetOutputCorCoefFileName().

TString hcalCalib::PHI_SYM_COR_FILENAME

Definition at line 141 of file hcalCalib.h.

Referenced by SetPhiSymCorFileName().

std::map<UInt_t, Float_t> hcalCalib::phiSymCor

Definition at line 199 of file hcalCalib.h.

Float_t hcalCalib::probeJetEmFrac

Definition at line 64 of file hcalCalib.h.

TLorentzVector* hcalCalib::probeJetP4

Definition at line 61 of file hcalCalib.h.

std::vector< std::pair<Int_t, UInt_t> > hcalCalib::refIEtaIPhi

Definition at line 196 of file hcalCalib.h.

UInt_t hcalCalib::runNumber

Definition at line 45 of file hcalCalib.h.

std::map<UInt_t, Float_t> hcalCalib::solution

Definition at line 201 of file hcalCalib.h.

Bool_t hcalCalib::SUM_DEPTHS

Definition at line 118 of file hcalCalib.h.

Referenced by SetSumDepthsFlag().

Bool_t hcalCalib::SUM_SMALL_DEPTHS

Definition at line 119 of file hcalCalib.h.

Referenced by SetSumSmallDepthsFlag().

Float_t hcalCalib::tagJetEmFrac

Definition at line 63 of file hcalCalib.h.

TLorentzVector* hcalCalib::tagJetP4

Definition at line 60 of file hcalCalib.h.

Float_t hcalCalib::targetE

Definition at line 50 of file hcalCalib.h.

std::vector< Float_t> hcalCalib::targetEnergies

Definition at line 197 of file hcalCalib.h.

const CaloGeometry* hcalCalib::theCaloGeometry

Definition at line 148 of file hcalCalib.h.

Referenced by SetCaloGeometry().

const HcalTopology* hcalCalib::topo_

Definition at line 149 of file hcalCalib.h.

Referenced by SetCaloGeometry().

Bool_t hcalCalib::USE_CONE_CLUSTERING

Definition at line 125 of file hcalCalib.h.

Referenced by SetUseConeClustering().

Float_t hcalCalib::xTrkEcal

Definition at line 56 of file hcalCalib.h.

Float_t hcalCalib::xTrkHcal

Definition at line 53 of file hcalCalib.h.

Float_t hcalCalib::yTrkEcal

Definition at line 57 of file hcalCalib.h.

Float_t hcalCalib::yTrkHcal

Definition at line 54 of file hcalCalib.h.

Float_t hcalCalib::zTrkEcal

Definition at line 58 of file hcalCalib.h.

Float_t hcalCalib::zTrkHcal

Definition at line 55 of file hcalCalib.h.