62 _nsamples(iConfig.getUntrackedParameter<unsigned
int>(
"nSamples", 10)),
63 _presample(iConfig.getUntrackedParameter<unsigned
int>(
"nPresamples", 2)),
64 _firstsample(iConfig.getUntrackedParameter<unsigned
int>(
"firstSample", 1)),
65 _lastsample(iConfig.getUntrackedParameter<unsigned
int>(
"lastSample", 2)),
66 _nsamplesPN(iConfig.getUntrackedParameter<unsigned
int>(
"nSamplesPN", 50)),
67 _presamplePN(iConfig.getUntrackedParameter<unsigned
int>(
"nPresamplesPN", 6)),
68 _firstsamplePN(iConfig.getUntrackedParameter<unsigned
int>(
"firstSamplePN", 7)),
69 _lastsamplePN(iConfig.getUntrackedParameter<unsigned
int>(
"lastSamplePN", 8)),
70 _timingcutlow(iConfig.getUntrackedParameter<unsigned
int>(
"timingCutLow", 2)),
71 _timingcuthigh(iConfig.getUntrackedParameter<unsigned
int>(
"timingCutHigh", 9)),
72 _timingquallow(iConfig.getUntrackedParameter<unsigned
int>(
"timingQualLow", 3)),
73 _timingqualhigh(iConfig.getUntrackedParameter<unsigned
int>(
"timingQualHigh", 8)),
74 _ratiomincutlow(iConfig.getUntrackedParameter<double>(
"ratioMinCutLow", 0.4)),
75 _ratiomincuthigh(iConfig.getUntrackedParameter<double>(
"ratioMinCutHigh", 0.95)),
76 _ratiomaxcutlow(iConfig.getUntrackedParameter<double>(
"ratioMaxCutLow", 0.8)),
77 _presamplecut(iConfig.getUntrackedParameter<double>(
"presampleCut", 5.0)),
78 _niter(iConfig.getUntrackedParameter<unsigned
int>(
"nIter", 3)),
79 _fitab(iConfig.getUntrackedParameter<
bool>(
"fitAB",
false)),
80 _alpha(iConfig.getUntrackedParameter<double>(
"alpha", 1.5076494)),
81 _beta(iConfig.getUntrackedParameter<double>(
"beta", 1.5136036)),
82 _nevtmax(iConfig.getUntrackedParameter<unsigned
int>(
"nEvtMax", 200)),
83 _noise(iConfig.getUntrackedParameter<double>(
"noise", 2.0)),
84 _chi2cut(iConfig.getUntrackedParameter<double>(
"chi2cut", 10.0)),
85 _ecalPart(iConfig.getUntrackedParameter<
std::
string>(
"ecalPart",
"EB")),
86 _docorpn(iConfig.getUntrackedParameter<
bool>(
"doCorPN",
false)),
87 _fedid(iConfig.getUntrackedParameter<
int>(
"fedID", -999)),
88 _saveallevents(iConfig.getUntrackedParameter<
bool>(
"saveAllEvents",
false)),
89 _qualpercent(iConfig.getUntrackedParameter<double>(
"qualPercent", 0.2)),
90 _debug(iConfig.getUntrackedParameter<
int>(
"debug", 0)),
146 for (
unsigned int j = 0;
j <
nCrys;
j++) {
158 for (
unsigned int j = 0;
j <
nMod;
j++) {
214 ADCfile +=
"/APDSamplesLaser.root";
217 APDfile +=
"/APDPNLaserAllEvents.root";
221 for (
unsigned int i = 0;
i <
nCrys;
i++) {
222 std::stringstream
name;
223 name <<
"ADCTree" <<
i + 1;
224 ADCtrees[
i] =
new TTree(name.str().c_str(), name.str().c_str());
260 std::stringstream nameabinitfile;
261 nameabinitfile <<
resdir_ <<
"/ABInit.root";
264 std::stringstream nameabfile;
265 nameabfile << resdir_ <<
"/AB.root";
273 if (test ==
nullptr) {
279 TFile* fAB =
nullptr;
280 TTree* ABInit =
nullptr;
282 fAB =
new TFile(nameabinitfile.str().c_str());
285 ABInit = (TTree*)fAB->Get(
"ABCol0");
305 std::stringstream nameapdfile;
306 nameapdfile << resdir_ <<
"/APDPN_LASER.root";
326 DCCHeader = pDCCHeader.
product();
354 std::cout <<
" Wrong ecalPart in cfg file " << std::endl;
375 TheMapping = ecalmapping.
product();
377 std::cerr <<
"Error! can't get the product EcalMappingRcd" << std::endl;
388 int fed = headerItr->fedId();
392 runType = headerItr->getRunType();
393 runNum = headerItr->getRunNumber();
394 event = headerItr->getLV1();
396 dccID = headerItr->getDccInTCCCommand();
397 fedID = headerItr->fedId();
419 if (iter ==
colors.end()) {
442 unsigned int samplemax = 0;
445 std::map<int, std::vector<double> > allPNAmpl;
446 std::map<int, std::vector<double> > allPNGain;
454 std::cout <<
"-- debug -- Inside PNDigi - pnID=" << pnDetId.
iPnId() <<
", dccID=" << pnDetId.
iDCCId()
465 for (
int samId = 0; samId < (*pnItr).size(); samId++) {
466 pn[samId] = (*pnItr).sample(samId).adc();
467 pnG[samId] = (*pnItr).sample(samId).gainId();
475 std::cout <<
"PN gain different from 1" << std::endl;
483 if (chi2pn == 101 || chi2pn == 102 || chi2pn == 103)
500 std::cout <<
"-- debug -- Inside PNDigi - PNampl=" <<
pnAmpl <<
", PNgain=" << pnGain << std::endl;
517 EBDetId id_crystal(digiItr->id());
521 int etaG = id_crystal.ieta();
522 int phiG = id_crystal.iphi();
526 int etaL = LocalCoord.first;
527 int phiL = LocalCoord.second;
530 int xtal = elecid_crystal.
xtalId();
538 unsigned int MyPn0 = pnpair.first;
539 unsigned int MyPn1 = pnpair.second;
543 assert(channel <
nCrys);
545 setGeomEB(etaG, phiG, module, tower, strip, xtal, apdRefTT, channel, lmr);
549 <<
" module:" << module <<
" modules:" <<
modules.size() << std::endl;
556 for (
unsigned int i = 0;
i < (*digiItr).size(); ++
i) {
558 adc[
i] = samp_crystal.adc();
559 adcG[
i] = samp_crystal.gainId();
581 int mem0 =
Mem->
Mem(lmr, 0);
582 int mem1 =
Mem->
Mem(lmr, 1);
584 if (allPNAmpl[mem0].
size() > MyPn0)
585 pn0 = allPNAmpl[mem0][MyPn0];
588 if (allPNAmpl[mem1].
size() > MyPn1)
589 pn1 = allPNAmpl[mem1][MyPn1];
620 EEDetId id_crystal(digiItr->id());
624 int etaG = id_crystal.iy();
625 int phiG = id_crystal.ix();
627 int iX = (phiG - 1) / 5 + 1;
628 int iY = (etaG - 1) / 5 + 1;
634 if (module >= 18 &&
side == 1)
641 unsigned int MyPn0 = pnpair.first;
642 unsigned int MyPn1 = pnpair.second;
650 assert(channel <
nCrys);
652 setGeomEE(etaG, phiG, iX, iY,
iZ, module, tower, ch, apdRefTT, channel, lmr);
656 <<
" module:" << module <<
" modules:" <<
modules.size() << std::endl;
661 if ((*digiItr).size() > 10)
662 std::cout <<
"SAMPLES SIZE > 10!" << (*digiItr).size() << std::endl;
666 for (
unsigned int i = 0;
i < (*digiItr).size(); ++
i) {
668 adc[
i] = samp_crystal.adc();
669 adcG[
i] = samp_crystal.gainId();
692 int mem0 =
Mem->
Mem(lmr, 0);
693 int mem1 =
Mem->
Mem(lmr, 1);
695 if (allPNAmpl[mem0].
size() > MyPn0)
696 pn0 = allPNAmpl[mem0][MyPn0];
699 if (allPNAmpl[mem1].
size() > MyPn1)
700 pn1 = allPNAmpl[mem1][MyPn1];
756 std::cout <<
"\n\t+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+" << std::endl;
757 std::cout <<
"\t+=+ Analyzing data: getting (alpha, beta) +=+" << std::endl;
758 TFile* fAB =
nullptr;
759 TTree* ABInit =
nullptr;
764 ABInit = (TTree*)fAB->Get(
"ABCol0");
767 std::cout <<
"\t+=+ .................................... done +=+" << std::endl;
768 std::cout <<
"\t+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+" << std::endl;
776 std::stringstream del;
778 system(del.str().c_str());
779 std::cout <<
" No Laser Events " << std::endl;
786 double BadGainEvtPercentage = 0.0;
787 double BadTimingEvtPercentage = 0.0;
789 int nChanBadGain = 0;
790 int nChanBadTiming = 0;
792 for (
unsigned int i = 0;
i <
nCrys;
i++) {
807 double BadGainChanPercentage = double(nChanBadGain) / double(nCrys);
808 double BadTimingChanPercentage = double(nChanBadTiming) / double(nCrys);
818 std::cout <<
"\n\t+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+" << std::endl;
819 std::cout <<
"\t+=+ Analyzing laser data: getting APD, PN, APD/PN, PN/PN +=+" << std::endl;
822 std::cout <<
"\t+=+ ............................ WARNING! APD GAIN WAS NOT 1 +=+" << std::endl;
824 std::cout <<
"\t+=+ ............................ WARNING! TIMING WAS BAD +=+" << std::endl;
830 for (
unsigned int i = 0;
i <
nCrys;
i++) {
831 std::stringstream
name;
832 name <<
"APDTree" <<
i + 1;
834 APDtrees[
i] =
new TTree(name.str().c_str(), name.str().c_str());
840 APDtrees[
i]->Branch(
"iphi", &iphi,
"iphi/I");
841 APDtrees[
i]->Branch(
"ieta", &ieta,
"ieta/I");
844 APDtrees[
i]->Branch(
"towerID", &towerID,
"towerID/I");
845 APDtrees[
i]->Branch(
"channelID", &channelID,
"channelID/I");
850 APDtrees[
i]->Branch(
"flag", &flag,
"flag/I");
857 APDtrees[
i]->SetBranchAddress(
"iphi", &iphi);
858 APDtrees[
i]->SetBranchAddress(
"ieta", &ieta);
861 APDtrees[
i]->SetBranchAddress(
"towerID", &towerID);
862 APDtrees[
i]->SetBranchAddress(
"channelID", &channelID);
867 APDtrees[
i]->SetBranchAddress(
"flag", &flag);
873 for (
unsigned int iref = 0; iref <
nRefChan; iref++) {
874 for (
unsigned int imod = 0; imod <
nMod; imod++) {
877 std::stringstream nameref;
878 nameref <<
"refAPDTree" << imod <<
"_" << iref;
880 RefAPDtrees[iref][jmod] =
new TTree(nameref.str().c_str(), nameref.str().c_str());
899 unsigned int nCol =
colors.size();
904 for (
unsigned int iM = 0; iM <
nMod; iM++) {
905 unsigned int iMod =
modules[iM] - 1;
907 for (
unsigned int ich = 0; ich <
nPNPerMod; ich++) {
908 for (
unsigned int icol = 0; icol < nCol; icol++) {
921 for (
unsigned int iCry = 0; iCry <
nCrys; iCry++) {
922 for (
unsigned int icol = 0; icol < nCol; icol++) {
928 std::stringstream
cut;
929 cut <<
"color==" <<
colors[icol];
930 if (
ADCtrees[iCry]->GetEntries(cut.str().c_str()) < 10)
934 unsigned int iMod =
iModule[iCry] - 1;
940 Long64_t nbytes = 0, nb = 0;
941 for (Long64_t jentry = 0; jentry <
ADCtrees[iCry]->GetEntriesFast(); jentry++) {
942 nb =
ADCtrees[iCry]->GetEntry(jentry);
947 unsigned int iCol = 0;
948 for (
unsigned int i = 0;
i < nCol;
i++) {
980 if (chi2 < 0. || chi2 == 102 || chi2 == 101) {
998 if (pn0 < 10 && pn1 > 10) {
1000 }
else if (pn1 < 10 && pn0 > 10) {
1003 pnmean = 0.5 * (
pn0 +
pn1);
1006 std::cout <<
"-- debug test -- pn0=" <<
pn0 <<
", pn1=" <<
pn1 << std::endl;
1012 for (
unsigned int ichan = 0; ichan <
nPNPerMod; ichan++) {
1033 for (
unsigned int ir = 0; ir <
nRefChan; ir++) {
1053 std::stringstream del;
1055 system(del.str().c_str());
1062 for (
unsigned int iColor = 0; iColor < nCol; iColor++) {
1063 std::stringstream nametree;
1064 nametree <<
"APDCol" <<
colors[iColor];
1065 std::stringstream nametree2;
1066 nametree2 <<
"PNCol" << colors[iColor];
1068 restrees[iColor] =
new TTree(nametree.str().c_str(), nametree.str().c_str());
1069 respntrees[iColor] =
new TTree(nametree2.str().c_str(), nametree2.str().c_str());
1071 restrees[iColor]->Branch(
"iphi", &iphi,
"iphi/I");
1072 restrees[iColor]->Branch(
"ieta", &ieta,
"ieta/I");
1076 restrees[iColor]->Branch(
"towerID", &towerID,
"towerID/I");
1077 restrees[iColor]->Branch(
"channelID", &channelID,
"channelID/I");
1078 restrees[iColor]->Branch(
"APD", &
APD,
"APD[6]/D");
1085 restrees[iColor]->Branch(
"flag", &flag,
"flag/I");
1087 respntrees[iColor]->Branch(
"side", &
side,
"side/I");
1088 respntrees[iColor]->Branch(
"moduleID", &
moduleID,
"moduleID/I");
1089 respntrees[iColor]->Branch(
"pnID", &
pnID,
"pnID/I");
1090 respntrees[iColor]->Branch(
"PN", &
PN,
"PN[6]/D");
1091 respntrees[iColor]->Branch(
"PNoPN", &
PNoPN,
"PNoPN[6]/D");
1092 respntrees[iColor]->Branch(
"PNoPNA", &
PNoPNA,
"PNoPNA[6]/D");
1093 respntrees[iColor]->Branch(
"PNoPNB", &
PNoPNB,
"PNoPNB[6]/D");
1095 restrees[iColor]->SetBranchAddress(
"iphi", &iphi);
1096 restrees[iColor]->SetBranchAddress(
"ieta", &ieta);
1100 restrees[iColor]->SetBranchAddress(
"towerID", &towerID);
1101 restrees[iColor]->SetBranchAddress(
"channelID", &channelID);
1109 restrees[iColor]->SetBranchAddress(
"flag", &flag);
1111 respntrees[iColor]->SetBranchAddress(
"side", &
side);
1112 respntrees[iColor]->SetBranchAddress(
"moduleID", &
moduleID);
1113 respntrees[iColor]->SetBranchAddress(
"pnID", &
pnID);
1114 respntrees[iColor]->SetBranchAddress(
"PN",
PN);
1115 respntrees[iColor]->SetBranchAddress(
"PNoPN",
PNoPN);
1116 respntrees[iColor]->SetBranchAddress(
"PNoPNA",
PNoPNA);
1117 respntrees[iColor]->SetBranchAddress(
"PNoPNB",
PNoPNB);
1123 for (
unsigned int iM = 0; iM <
nMod; iM++) {
1124 unsigned int iMod =
modules[iM] - 1;
1126 for (
unsigned int ich = 0; ich <
nPNPerMod; ich++) {
1127 for (
unsigned int icol = 0; icol < nCol; icol++) {
1136 for (
unsigned int imod = 0; imod <
nMod; imod++) {
1147 for (
unsigned int iCry = 0; iCry <
nCrys; iCry++) {
1148 unsigned int iMod =
iModule[iCry] - 1;
1153 for (
unsigned int iCol = 0; iCol < nCol; iCol++) {
1154 std::vector<double> lowcut;
1155 std::vector<double> highcut;
1164 lowcut.push_back(cutMin);
1165 highcut.push_back(cutMax);
1169 lowcut.push_back(cutMin);
1170 highcut.push_back(cutMax);
1190 Long64_t nbytes = 0, nb = 0;
1191 for (Long64_t jentry = 0; jentry <
APDtrees[iCry]->GetEntriesFast(); jentry++) {
1192 nb =
APDtrees[iCry]->GetEntry(jentry);
1196 if (pn0 < 10 && pn1 > 10) {
1198 }
else if (pn1 < 10 && pn0 > 10) {
1201 pnmean = 0.5 * (
pn0 +
pn1);
1206 unsigned int iCol = 0;
1207 for (
unsigned int i = 0;
i < nCol;
i++) {
1218 for (
unsigned int ichan = 0; ichan <
nPNPerMod; ichan++) {
1227 std::cout <<
"-- debug test -- Last Loop event:" <<
event <<
" apdAmpl:" <<
apdAmpl << std::endl;
1231 for (
unsigned int iRef = 0; iRef <
nRefChan; iRef++) {
1237 <<
", event:" <<
event <<
", eventref:" <<
eventref << std::endl;
1253 for (
unsigned int iColor = 0; iColor < nCol; iColor++) {
1254 std::vector<double> apdvec =
APDAnal[iCry][iColor]->
getAPD();
1262 for (
unsigned int i = 0;
i < apdvec.size();
i++) {
1263 APD[
i] = apdvec.at(
i);
1294 for (
unsigned int iM = 0; iM <
nMod; iM++) {
1295 unsigned int iMod =
modules[iM] - 1;
1299 for (
unsigned int ch = 0; ch <
nPNPerMod; ch++) {
1306 for (
unsigned int iColor = 0; iColor < nCol; iColor++) {
1307 std::vector<double> pnvec =
PNAnal[iMod][ch][iColor]->
getPN();
1308 std::vector<double> pnopnvec =
PNAnal[iMod][ch][iColor]->
getPNoPN();
1309 std::vector<double> pnopn0vec =
PNAnal[iMod][ch][iColor]->
getPNoPN0();
1310 std::vector<double> pnopn1vec =
PNAnal[iMod][ch][iColor]->
getPNoPN1();
1312 for (
unsigned int i = 0;
i < pnvec.size();
i++) {
1313 PN[
i] = pnvec.at(
i);
1331 std::stringstream del2;
1333 system(del2.str().c_str());
1346 for (
unsigned int i = 0;
i < nCol;
i++) {
1353 std::cout <<
"\t+=+ .................................................. done +=+" << std::endl;
1354 std::cout <<
"\t+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+" << std::endl;
1358 int etaG,
int phiG,
int module,
int tower,
int strip,
int xtal,
int apdRefTT,
int channel,
int lmr) {
1370 for (
unsigned int iref = 0; iref <
nRefChan; iref++) {
1391 int etaG,
int phiG,
int iX,
int iY,
int iZ,
int module,
int tower,
int ch,
int apdRefTT,
int channel,
int lmr) {
1403 for (
unsigned int iref = 0; iref <
nRefChan; iref++) {
unsigned int _firstsamplePN
void addEntry(double val)
static XYCoord localCoord(int icr)
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
static int lmmod(SuperCrysCoord iX, SuperCrysCoord iY)
TAPD * APDAnal[1700][nColor]
void setAPDCut(double, double)
TPN * PNFirstAnal[22][2][nColor]
int xtalId() const
get the channel id
std::vector< double > getPN()
EcalLaserAnalyzer(const edm::ParameterSet &iConfig)
static int apdRefTower(int ilmr, int ilmmod)
void addEntry(double, double, double)
std::map< int, unsigned int > apdRefMap[2]
int stripId() const
get the tower id
std::string digiPNCollection_
double getDelta(int, int)
std::vector< double > getPNoPN1()
Ecal readout channel identification [32:20] Unused (so far) [19:13] DCC id [12:6] tower [5:3] strip [...
~EcalLaserAnalyzer() override
static int side(SuperCrysCoord iX, SuperCrysCoord iY, int iz)
void setPNCut(double, double)
void addEntry(double, double, double, double, double, double, double)
double getPNCorrectionFactor(double val0, int gain)
std::vector< T >::const_iterator const_iterator
unsigned int _timingcuthigh
void computeShape(std::string namefile, TTree *)
int towerId() const
get the tower id
double * getAdcWithoutPedestal()
const_iterator begin() const
void putAllVals(int, double *, int, int)
static int lmr(EBGlobalCoord ieta, EBGlobalCoord iphi)
unsigned int _timingcutlow
TPN * PNAnal[22][2][nColor]
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
unsigned int _timingquallow
std::vector< double > getAPDoPN0()
void setTimeCut(double, double)
TTree * RefAPDtrees[2][22]
static std::pair< int, int > pn(int ilmmod)
void set2DAPDoAPD1Cut(const std::vector< double > &, const std::vector< double > &)
void setAPDoPN0Cut(double, double)
std::vector< double > getPNoPN0()
std::vector< double > getAPDoPN1()
std::string digiProducer_
int iPnId() const
get the PnId
int hashedIndex(int ieta, int iphi)
EcalElectronicsId getElectronicsId(const DetId &id) const
Get the electronics id for this det id.
static std::pair< int, int > pn(int dee, int ilmod)
std::vector< double > getAPDoAPD0()
#define DEFINE_FWK_MODULE(type)
void set2DAPDoAPD0Cut(const std::vector< double > &, const std::vector< double > &)
TAPD * APDFirstAnal[1700][nColor]
static int electronic_channel(EBLocalCoord ix, EBLocalCoord iy)
unsigned int firstChanMod[22]
std::string digiCollection_
unsigned int isFirstChanModFilled[22]
unsigned int _lastsamplePN
std::string eventHeaderProducer_
double * getAdcWithoutPedestal()
static int apdRefTower(int ilmmod)
unsigned int _firstsample
std::string eventHeaderCollection_
int iDCCId() const
get the DCCId
static int lmmod(EBGlobalCoord ieta, EBGlobalCoord iphi)
std::vector< double > getAPDoAPD1()
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
virtual void init(int, int, int, int, double, double)
std::vector< double > getAPD()
std::map< unsigned int, unsigned int > channelMapEE
unsigned int iModule[1700]
const_iterator end() const
std::string alphainitfile
static int side(EBGlobalCoord ieta, EBGlobalCoord iphi)
T const * product() const
void setGeomEB(int etaG, int phiG, int module, int tower, int strip, int xtal, int apdRefTT, int channel, int lmr)
virtual double doFit(double *)
void analyze(const edm::Event &e, const edm::EventSetup &c) override
void setGeomEE(int etaG, int phiG, int iX, int iY, int iZ, int module, int tower, int ch, int apdRefTT, int channel, int lmr)
int IsThereDataADC[1700][nColor]
boost::transform_iterator< IterHelp, boost::counting_iterator< int > > const_iterator
std::vector< double > getVals(int)
void set_const(int, int, int, int, int, double, double)
const_iterator end() const
unsigned int _timingqualhigh
std::vector< double > getPNoPN()
static int dee(SuperCrysCoord iX, SuperCrysCoord iY, int iz)
unsigned int nevtAB[1700]
double doFit(int, double *)
alpha
zGenParticlesMatch = cms.InputTag(""),
std::vector< double > getAPDoPN()
unsigned int _presamplePN
TTree * respntrees[nColor]
std::vector< double > getTime()
T const * product() const
static std::vector< int > apdRefChannels(ME::LMMid ilmmod, ME::LMRid ilmr)
static std::vector< ME::LMMid > lmmodFromDcc(ME::DCCid idcc)
static int lmr(SuperCrysCoord iX, SuperCrysCoord iY, int iz)
int channelId() const
so far for EndCap only :
const_iterator begin() const
void setAPDoPNCut(double, double)
void setAPDoPN1Cut(double, double)