1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
|
#include "Geometry/Records/interface/IdealGeometryRecord.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "Geometry/Records/interface/IdealGeometryRecord.h"
#include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
#include "Geometry/CaloGeometry/interface/CaloGeometry.h"
#include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
#include "CondFormats/DataRecord/interface/EcalIntercalibConstantsRcd.h"
#include "Calibration/Tools/interface/calibXMLwriter.h"
#include "CalibCalorimetry/CaloMiscalibTools/interface/CaloMiscalibTools.h"
#include "CalibCalorimetry/CaloMiscalibTools/interface/CaloMiscalibMapEcal.h"
#include "CalibCalorimetry/CaloMiscalibTools/interface/MiscalibReaderFromXMLEcalBarrel.h"
#include "CalibCalorimetry/CaloMiscalibTools/interface/MiscalibReaderFromXMLEcalEndcap.h"
#include "CondFormats/EcalObjects/interface/EcalIntercalibConstants.h"
#include "FWCore/Framework/interface/EventSetup.h"
#include "FWCore/Framework/interface/ESHandle.h"
//#include "Calibration/EcalAlCaRecoProducers/interface/trivialParser.h"
//#include "Calibration/EcalAlCaRecoProducers/bin/trivialParser.h"
#include "TH2.h"
#include "TH1.h"
#include "TFile.h"
#include "TProfile.h"
#define PI_GRECO 3.14159265
inline int etaShifter(const int etaOld) {
if (etaOld < 0)
return etaOld + 85;
else if (etaOld > 0)
return etaOld + 84;
assert(false);
}
// ------------------------------------------------------------------------
// MF questo file prende due set di coefficienti e li confronta
//
//-------------------------------------------------------------------------
int main(int argc, char* argv[]) {
int EEradStart = 15;
int EEradEnd = 50;
int EEphiStart = 20;
int EEphiEnd = 45;
std::string endcapfile = argv[1];
std::string calibEndcapfile = argv[2];
// std::string endcapfile = "/afs/cern.ch/user/g/govoni/scratch1/CMSSW/CALIB/CMSSW_1_6_0/src/CalibCalorimetry/CaloMiscalibTools/data/ecal_endcap_startup.xml" ;
// std::string calibEndcapfile = "/afs/cern.ch/user/g/govoni/scratch1/CMSSW/CALIB/CMSSW_1_6_0/src/CalibCalorimetry/CaloMiscalibTools/data/inv_ecal_endcap_startup.xml" ;
//PG get the miscalibration files for EB and EE
//PG ------------------------------------------
CaloMiscalibMapEcal EEscalibMap;
EEscalibMap.prefillMap();
MiscalibReaderFromXMLEcalEndcap endcapreader(EEscalibMap);
if (!endcapfile.empty())
endcapreader.parseXMLMiscalibFile(endcapfile);
EcalIntercalibConstants* EEconstants = new EcalIntercalibConstants(EEscalibMap.get());
const EcalIntercalibConstantMap& iEEscalibMap = EEconstants->getMap(); //MF prende i vecchi coeff
//PG get the recalibration files for EB and EE
//PG -----------------------------------------
CaloMiscalibMapEcal EEcalibMap;
EEcalibMap.prefillMap();
MiscalibReaderFromXMLEcalEndcap calibEndcapreader(EEcalibMap);
if (!calibEndcapfile.empty())
calibEndcapreader.parseXMLMiscalibFile(calibEndcapfile);
EcalIntercalibConstants* EECconstants = new EcalIntercalibConstants(EEcalibMap.get());
const EcalIntercalibConstantMap& iEEcalibMap = EECconstants->getMap(); //MF prende i vecchi coeff
//PG fill the histograms
//PG -------------------
TH1F EEPCompareCoeffDistr("EEPCompareCoeffDistr", "EEPCompareCoeffDistr", 5000, 0, 2);
TH2F EEPCompareCoeffMap("EEPCompareCoeffMap", "EEPCompareCoeffMap", 101, 0, 101, 101, 0, 101);
TH2F EEPCompareCoeffEtaTrend("EEPCompareCoeffEtaTrend", "EEPCompareCoeffEtaTrend", 51, 0, 50, 500, 0, 2);
TProfile EEPCompareCoeffEtaProfile("EEPCompareCoeffEtaProfile", "EEPCompareCoeffEtaProfile", 51, 0, 50, 0, 2);
// ECAL endcap +
for (int ix = 1; ix <= 100; ++ix)
for (int iy = 1; iy <= 100; ++iy) {
int rad = static_cast<int>(sqrt((ix - 50) * (ix - 50) + (iy - 50) * (iy - 50)));
if (rad < EEradStart || rad > EEradEnd)
continue;
double phiTemp = atan2(iy - 50, ix - 50);
if (phiTemp < 0)
phiTemp += 2 * PI_GRECO;
int phi = static_cast<int>(phiTemp * 180 / PI_GRECO);
if (phi < EEphiStart || phi > EEphiEnd)
continue;
if (!EEDetId::validDetId(ix, iy, 1))
continue;
EEDetId det = EEDetId(ix, iy, 1, EEDetId::XYMODE);
if (*(iEEscalibMap.find(det.rawId())) == 0)
continue;
double factor = *(iEEcalibMap.find(det.rawId())) / *(iEEscalibMap.find(det.rawId()));
EEPCompareCoeffDistr.Fill(factor);
EEPCompareCoeffMap.Fill(ix, iy, factor);
EEPCompareCoeffEtaTrend.Fill(rad, factor);
EEPCompareCoeffEtaProfile.Fill(rad, factor);
} // ECAL endcap +
// ECAL endcap-
TH1F EEMCompareCoeffDistr("EEMCompareCoeffDistr", "EEMCompareCoeffDistr", 200, 0, 2);
TH2F EEMCompareCoeffMap("EEMCompareCoeffMap", "EEMCompareCoeffMap", 100, 0, 100, 100, 0, 100);
TH2F EEMCompareCoeffEtaTrend("EEMCompareCoeffEtaTrend", "EEMCompareCoeffEtaTrend", 51, 0, 50, 500, 0, 2);
TProfile EEMCompareCoeffEtaProfile("EEMCompareCoeffEtaProfile", "EEMCompareCoeffEtaProfile", 51, 0, 50, 0, 2);
// ECAL endcap -
for (int ix = 1; ix <= 100; ++ix)
for (int iy = 1; iy <= 100; ++iy) {
int rad = static_cast<int>(sqrt((ix - 50) * (ix - 50) + (iy - 50) * (iy - 50)));
if (rad < EEradStart || rad > EEradEnd)
continue;
if (!EEDetId::validDetId(ix, iy, -1))
continue;
EEDetId det = EEDetId(ix, iy, -1, EEDetId::XYMODE);
if (*(iEEscalibMap.find(det.rawId())) == 0)
continue;
double factor = *(iEEcalibMap.find(det.rawId())) / *(iEEscalibMap.find(det.rawId()));
EEMCompareCoeffDistr.Fill(factor);
EEMCompareCoeffMap.Fill(ix, iy, factor);
EEMCompareCoeffEtaTrend.Fill(rad, factor);
EEMCompareCoeffEtaProfile.Fill(rad, factor);
} // ECAL endcap -
std::string filename = "coeffcompareEE.root";
TFile out(filename.c_str(), "recreate");
EEMCompareCoeffMap.Write();
EEMCompareCoeffDistr.Write();
EEMCompareCoeffEtaTrend.Write();
EEMCompareCoeffEtaProfile.Write();
EEPCompareCoeffMap.Write();
EEPCompareCoeffDistr.Write();
EEPCompareCoeffEtaTrend.Write();
EEPCompareCoeffEtaProfile.Write();
out.Close();
}
|