From 1b322caa5a782c817e5b618e9eca1e10bace8f67 Mon Sep 17 00:00:00 2001 From: Sunanda Date: Wed, 20 May 2020 17:03:20 +0200 Subject: [PATCH 1/2] Update material budget studies for HCAL --- Validation/Geometry/src/MaterialBudgetHcal.cc | 6 +- .../Geometry/src/MaterialBudgetHcalHistos.cc | 109 ++++++--------- Validation/Geometry/test/MatBudgetHcal.C | 128 +++++++++++------- Validation/Geometry/test/runP_Calo_cfg.py | 5 +- .../Geometry/test/runP_Castor_Debug_cfg.py | 24 +--- Validation/Geometry/test/runP_Castor_cfg.py | 1 + Validation/Geometry/test/runP_ECAL_cfg.py | 1 + Validation/Geometry/test/runP_Forward_cfg.py | 1 + .../Geometry/test/runP_H2TB_Debug_cfg.py | 1 + .../Geometry/test/runP_HCAL_Debug_cfg.py | 26 +--- Validation/Geometry/test/runP_HCAL_cfg.py | 28 ++-- .../Geometry/test/runP_InFrontOfECAL_cfg.py | 1 + 12 files changed, 147 insertions(+), 184 deletions(-) diff --git a/Validation/Geometry/src/MaterialBudgetHcal.cc b/Validation/Geometry/src/MaterialBudgetHcal.cc index e01257896fd4f..e7ca6573a07a0 100644 --- a/Validation/Geometry/src/MaterialBudgetHcal.cc +++ b/Validation/Geometry/src/MaterialBudgetHcal.cc @@ -22,9 +22,7 @@ MaterialBudgetHcal::MaterialBudgetHcal(const edm::ParameterSet& p) : theHistoHca rMax = m_p.getUntrackedParameter("RMax", 4.5) * CLHEP::m; zMax = m_p.getUntrackedParameter("ZMax", 13.0) * CLHEP::m; bool doHcal = m_p.getUntrackedParameter("DoHCAL", true); - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcal initialized with rMax " << rMax << " mm and zMax " << zMax - << " mm" - << " doHcal is set to " << doHcal; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcal initialized with rMax " << rMax << " mm and zMax " << zMax << " mm doHcal is set to " << doHcal; if (doHcal) theHistoHcal = new MaterialBudgetHcalHistos(m_p); else @@ -82,7 +80,7 @@ bool MaterialBudgetHcal::stopAfter(const G4Step* aStep) { double zz = std::abs(hitPoint.z()); if (rr > rMax || zz > zMax) { - LogDebug("MaterialBudget") << " MaterialBudgetHcal::StopAfter R = " << rr << " and Z = " << zz; + edm::LogVerbatim("MaterialBudget") << " MaterialBudgetHcal::StopAfter R = " << rr << " and Z = " << zz; return true; } else { return false; diff --git a/Validation/Geometry/src/MaterialBudgetHcalHistos.cc b/Validation/Geometry/src/MaterialBudgetHcalHistos.cc index 04cd339f1e1d1..0f1b8c3a5246a 100644 --- a/Validation/Geometry/src/MaterialBudgetHcalHistos.cc +++ b/Validation/Geometry/src/MaterialBudgetHcalHistos.cc @@ -14,6 +14,8 @@ #include +//#define EDM_ML_DEBUG + MaterialBudgetHcalHistos::MaterialBudgetHcalHistos(const edm::ParameterSet& p) { binEta = p.getUntrackedParameter("NBinEta", 260); binPhi = p.getUntrackedParameter("NBinPhi", 180); @@ -22,10 +24,7 @@ MaterialBudgetHcalHistos::MaterialBudgetHcalHistos(const edm::ParameterSet& p) { etaHigh = p.getUntrackedParameter("EtaHigh", 5.2); fillHistos = p.getUntrackedParameter("FillHisto", true); printSum = p.getUntrackedParameter("PrintSummary", false); - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos: FillHisto : " << fillHistos << " PrintSummary " - << printSum << " == Eta plot: NX " << binEta << " Range " << -maxEta << ":" << maxEta - << " Phi plot: NX " << binPhi << " Range " << -pi << ":" << pi << " (Eta limit " - << etaLow << ":" << etaHigh << ")"; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: FillHisto : " << fillHistos << " PrintSummary " << printSum << " == Eta plot: NX " << binEta << " Range " << -maxEta << ":" << maxEta << " Phi plot: NX " << binPhi << " Range " << -pi << ":" << pi << " (Eta limit " << etaLow << ":" << etaHigh << ")"; if (fillHistos) book(); } @@ -37,11 +36,9 @@ void MaterialBudgetHcalHistos::fillBeginJob(const DDCompactView& cpv) { DDSpecificsMatchesValueFilter filter1{DDValue(attribute, value, 0)}; DDFilteredView fv1(cpv, filter1); sensitives = getNames(fv1); - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos: Names to be " - << "tested for " << attribute << " = " << value << " has " << sensitives.size() - << " elements"; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Names to be tested for " << attribute << " = " << value << " has " << sensitives.size() << " elements"; for (unsigned int i = 0; i < sensitives.size(); i++) - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos: sensitives[" << i << "] = " << sensitives[i]; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: sensitives[" << i << "] = " << sensitives[i]; attribute = "Volume"; value = "HF"; @@ -51,14 +48,11 @@ void MaterialBudgetHcalHistos::fillBeginJob(const DDCompactView& cpv) { fv2.firstChild(); DDsvalues_type sv(fv2.mergedSpecifics()); std::vector temp = getDDDArray("Levels", sv); - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos: Names to be " - << "tested for " << attribute << " = " << value << " has " << hfNames.size() - << " elements"; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Names to be tested for " << attribute << " = " << value << " has " << hfNames.size() << " elements"; for (unsigned int i = 0; i < hfNames.size(); i++) { int level = static_cast(temp[i]); hfLevels.push_back(level); - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos: HF[" << i << "] = " << hfNames[i] << " at level " - << hfLevels[i]; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: HF[" << i << "] = " << hfNames[i] << " at level " << hfLevels[i]; } std::string ecalRO[2] = {"EcalHitsEB", "EcalHitsEE"}; @@ -68,14 +62,12 @@ void MaterialBudgetHcalHistos::fillBeginJob(const DDCompactView& cpv) { DDSpecificsMatchesValueFilter filter3{DDValue(attribute, value, 0)}; DDFilteredView fv3(cpv, filter3); std::vector senstmp = getNames(fv3); - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos: Names to be" - << " tested for " << attribute << " = " << value << " has " << senstmp.size() - << " elements"; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Names to be tested for " << attribute << " = " << value << " has " << senstmp.size() << " elements"; for (unsigned int i = 0; i < senstmp.size(); i++) sensitiveEC.push_back(senstmp[i]); } for (unsigned int i = 0; i < sensitiveEC.size(); i++) - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos:sensitiveEC[" << i << "] = " << sensitiveEC[i]; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos:sensitiveEC[" << i << "] = " << sensitiveEC[i]; } } @@ -101,9 +93,7 @@ void MaterialBudgetHcalHistos::fillStartTrack(const G4Track* aTrack) { intLength.clear(); } - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos: Track " << aTrack->GetTrackID() << " Code " << theID - << " Energy " << theEnergy / GeV << " GeV; Eta " << eta << " Phi " << phi / deg - << " PT " << dir.perp() / GeV << " GeV *****"; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Track " << aTrack->GetTrackID() << " Code " << theID << " Energy " << theEnergy / CLHEP::GeV << " GeV; Eta " << eta << " Phi " << phi / CLHEP::deg << " PT " << dir.perp() / CLHEP::GeV << " GeV *****"; } void MaterialBudgetHcalHistos::fillPerStep(const G4Step* aStep) { @@ -134,17 +124,9 @@ void MaterialBudgetHcalHistos::fillPerStep(const G4Step* aStep) { radLength.push_back(step / radl); intLength.push_back(step / intl); } - edm::LogInfo("MaterialBudget") << name << " " << step << " " << matName << " " << stepLen << " " << step / radl - << " " << radLen << " " << step / intl << " " << intLen; + edm::LogVerbatim("MaterialBudget") << "Volume " << name << " id " << id << ":" << idOld << " Step " << step << " Material " << matName << " Old Length " << stepLen << " X0 " << step / radl << ":" << radLen << " Lambda " << step / intl << ":" << intLen; } else { - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos: Step at " << name << " Length " << step << " in " - << matName << " of density " << density << " g/cc; Radiation Length " << radl - << " mm;" - << " Interaction Length " << intl << " mm\n" - << " Position " << aStep->GetPreStepPoint()->GetPosition() - << " Cylindrical R " << aStep->GetPreStepPoint()->GetPosition().perp() - << " Length (so far) " << stepLen << " L/X0 " << step / radl << "/" << radLen - << " L/Lambda " << step / intl << "/" << intLen; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Step at " << name << " id " << id << ":" << idOld << " Length " << step << " in " << matName << " of density " << density << " g/cc; Radiation Length " << radl << " mm; Interaction Length " << intl << " mm\n" << " Position " << aStep->GetPreStepPoint()->GetPosition() << " Cylindrical R " << aStep->GetPreStepPoint()->GetPosition().perp() << " Length (so far) " << stepLen << " L/X0 " << step / radl << "/" << radLen << " L/Lambda " << step / intl << "/" << intLen; } int det = 0, lay = 0; @@ -183,8 +165,9 @@ void MaterialBudgetHcalHistos::fillPerStep(const G4Step* aStep) { nlayHB++; } } - LogDebug("MaterialBudget") << "MaterialBudgetHcalHistos: Det " << det << " Layer " << lay << " Eta " << eta - << " Phi " << phi / deg; +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Det " << det << " Layer " << lay << " Eta " << eta << " Phi " << phi / CLHEP::deg; +#endif } else if (layer == 1) { det = -1; lay = 2; @@ -198,38 +181,39 @@ void MaterialBudgetHcalHistos::fillPerStep(const G4Step* aStep) { } if (id > idOld) { - // edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos: Step at " << name; +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Step at " << name << " calls filHisto with " << (id - 1); +#endif fillHisto(id - 1); } } stepLen += step; - radLen += step / radl; - intLen += step / intl; + radLen += (step / radl); + intLen += (step / intl); if (fillHistos) { - if (layer == 21 && det == 5) { + if (id == 21) { if (!isItHF(aStep->GetPostStepPoint()->GetTouchable())) { - LogDebug("MaterialBudget") << "MaterialBudgetHcalHistos: After HF in " - << aStep->GetPostStepPoint()->GetTouchable()->GetVolume(0)->GetName(); - fillHisto(id); - id++; - layer = 0; +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: After HF in " << name << ":" << aStep->GetPostStepPoint()->GetTouchable()->GetVolume(0)->GetName() << " calls fillHisto with " << id; +#endif + fillHisto(idOld); + id++; + layer = 0; } } } } void MaterialBudgetHcalHistos::fillEndTrack() { - edm::LogInfo("MaterialBudget") << "Number of layers hit in HB:" << nlayHB << " HE:" << nlayHE << " HO:" << nlayHO - << " HF:" << nlayHF; + edm::LogVerbatim("MaterialBudget") << "Number of layers hit in HB:" << nlayHB << " HE:" << nlayHE << " HO:" << nlayHO << " HF:" << nlayHF; if (fillHistos) { fillHisto(maxSet - 1); fillLayer(); } if (printSum) { for (unsigned int ii = 0; ii < matList.size(); ii++) { - edm::LogInfo("MaterialBudget") << matList[ii] << "\t" << stepLength[ii] << "\t" << radLength[ii] << "\t" - << intLength[ii]; + edm::LogVerbatim("MaterialBudget") << matList[ii] << "\t" << stepLength[ii] << "\t" << radLength[ii] << "\t" << intLength[ii]; } } } @@ -243,10 +227,7 @@ void MaterialBudgetHcalHistos::book() { << "please add it to config file"; double maxPhi = pi; - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos: Booking user " - << "histos === with " << binEta << " bins " - << "in eta from " << -maxEta << " to " << maxEta << " and " << binPhi << " bins " - << "in phi from " << -maxPhi << " to " << maxPhi; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Booking user histos === with " << binEta << " bins in eta from " << -maxEta << " to " << maxEta << " and " << binPhi << " bins in phi from " << -maxPhi << " to " << maxPhi; std::string iter; // total X0 @@ -323,15 +304,13 @@ void MaterialBudgetHcalHistos::book() { maxEta); } - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos: Booking user " - << "histos done ==="; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Booking user histos done ==="; } void MaterialBudgetHcalHistos::fillHisto(int ii) { - LogDebug("MaterialBudget") << "MaterialBudgetHcalHistos:FillHisto called " - << "with index " << ii << " integrated step " << stepLen << " X0 " << radLen << " Lamda " - << intLen; - +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos:FillHisto called with index " << ii << " integrated step " << stepLen << " X0 " << radLen << " Lamda " << intLen; +#endif if (ii >= 0 && ii < maxSet) { me100[ii]->Fill(eta, radLen); me200[ii]->Fill(eta, intLen); @@ -394,8 +373,7 @@ void MaterialBudgetHcalHistos::fillLayer() { } void MaterialBudgetHcalHistos::hend() { - edm::LogInfo("MaterialBudget") << "MaterialBudgetHcalHistos: Save user " - << "histos ==="; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Save user histos ==="; } std::vector MaterialBudgetHcalHistos::getNames(DDFilteredView& fv) { @@ -416,24 +394,21 @@ std::vector MaterialBudgetHcalHistos::getNames(DDFilteredView& fv) } std::vector MaterialBudgetHcalHistos::getDDDArray(const std::string& str, const DDsvalues_type& sv) { - LogDebug("MaterialBudget") << "MaterialBudgetHcalHistos:getDDDArray called " - << "for " << str; +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos:getDDDArray called for " << str; +#endif DDValue value(str); if (DDfetch(&sv, value)) { LogDebug("MaterialBudget") << value; const std::vector& fvec = value.doubles(); int nval = fvec.size(); if (nval < 1) { - edm::LogError("MaterialBudget") << "MaterialBudgetHcalHistos : # of " << str << " bins " << nval - << " < 1 ==> illegal"; - throw cms::Exception("Unknown", "MaterialBudgetHcalHistos") << "nval < 1 for array " << str << "\n"; + throw cms::Exception("MaterialBudgetHcalHistos") << "nval = " << nval << " < 1 for array " << str << "\n"; } return fvec; } else { - edm::LogError("MaterialBudget") << "MaterialBudgetHcalHistos : cannot get " - << "array " << str; - throw cms::Exception("Unknown", "MaterialBudgetHcalHistos") << "cannot get array " << str << "\n"; + throw cms::Exception("MaterialBudgetHcalHistos") << "cannot get array " << str << "\n"; } } @@ -447,13 +422,13 @@ bool MaterialBudgetHcalHistos::isSensitive(std::string name) { } bool MaterialBudgetHcalHistos::isItHF(const G4VTouchable* touch) { - // std::vector::const_iterator it = hfNames.begin(); int levels = ((touch->GetHistoryDepth()) + 1); for (unsigned int it = 0; it < hfNames.size(); it++) { if (levels >= hfLevels[it]) { std::string name = touch->GetVolume(levels - hfLevels[it])->GetName(); - if (name == hfNames[it]) + if (name == hfNames[it]) { return true; + } } } return false; diff --git a/Validation/Geometry/test/MatBudgetHcal.C b/Validation/Geometry/test/MatBudgetHcal.C index 5d5a7c2a537bd..28b4b61d1a2f4 100644 --- a/Validation/Geometry/test/MatBudgetHcal.C +++ b/Validation/Geometry/test/MatBudgetHcal.C @@ -45,60 +45,70 @@ int colorLayer[25] = { 2, 7, 9, 30, 34, 38, 14, 40, 41, 42, void etaPhiPlot(TString fileName="matbdg_HCAL.root", TString plot="IntLen", int ifirst=0, int ilast=21, int drawLeg=1, bool ifEta=true, - double maxEta=-1, bool debug=true); + double maxEta=-1, std::string tag="Run", bool debug=false); void etaPhiPlotHO(TString fileName="matbdg_HCAL.root", TString plot="IntLen", - int drawLeg=1, bool ifEta=true, double maxEta=-1); + int drawLeg=1, bool ifEta=true, double maxEta=-1, + std::string tag="Run"); void etaPhiPlotEC(TString fileName="matbdg_HCAL.root", TString plot="IntLen", - int drawLeg=1, bool ifEta=true, double maxEta=-1); + int drawLeg=1, bool ifEta=true, double maxEta=-1, + std::string tag="Run"); void etaPhiPlotHC(TString fileName="matbdg_HCAL.root", TString plot="IntLen", - int drawLeg=1, bool ifEta=true, double maxEta=-1); + int drawLeg=1, bool ifEta=true, double maxEta=-1, + std::string tag="Run"); void etaPhi2DPlot(TString fileName="matbdg_HCAL.root", TString plot="IntLen", - int ifirst=0, int ilast=19, int drawLeg=1); + int ifirst=0, int ilast=19, int drawLeg=1, + std::string tag="Run"); void etaPhi2DPlot(int nslice, int kslice, TString fileName="matbdg_HCAL.root", TString plot="IntLen", int ifirst=0, int ilast=21, - int drawLeg=1); + int drawLeg=1, std::string tag="Run"); void printTable (TString fileName="matbdg_HCAL.root", TString outputFileName="hcal.txt", TString inputFileName="None"); -void plotDiff (TString fileName="matbdg_HCAL.root", TString plot="IntLen"); +void plotDiff (TString fileName="matbdg_HCAL.root", TString plot="IntLen", + std::string tag="Run"); void getDiff (TString fileName="matbdg_HCAL.root", TString plot="IntLen"); -void plotHE(int flag=0, int logy=0, int save=0); +void plotHE(int flag=0, int logy=0, std::string tag="Run", int save=0); void etaPhiCastorPlot(TString fileName="matbdg_Castor.root", TString plot="IntLen", TString type="All", bool etaPlus=true, int drawLeg=1, bool ifEta=true, - bool debug=true); + std::string tag="Run", bool debug=false); void efficiencyPlot(TString fileName="matbdg_HCAL.root", TString type="All", - bool ifEtaPhi=true, double maxEta=-1, bool debug=false); + bool ifEtaPhi=true, double maxEta=-1, std::string tag="Run", + bool debug=false); void etaPhiFwdPlot(TString fileName="matbdg_Fwd.root", TString plot="IntLen", - int first=0, int last=9, int drawLeg=1, bool debug=false); + int first=0, int last=9, int drawLeg=1, std::string tag="Run" + , bool debug=false); void setStyle (); -void standardPlot (TString fileName="matbdg_HCAL.root", +void standardPlot (TString fileName="matbdgHCAL_run3.root", + std::string tag="Run3", int maxEta=4.8, TString outputFileName="hcal.txt") { - etaPhiPlot (fileName, "IntLen", 0, 21, 1, true, 4.8); - etaPhiPlot (fileName, "IntLen", 0, 20, 1, false, -1.); - etaPhiPlot (fileName, "RadLen", 0, 21, 1, true, 4.8); - etaPhiPlot (fileName, "RadLen", 0, 20, 1, false, -1); - etaPhiPlot (fileName, "StepLen",0, 21, 1, true, -1); - etaPhiPlot (fileName, "StepLen",0, 20, 1, false, -1); - etaPhiPlotHO(fileName, "IntLen", 1, true, 2.5); - etaPhiPlotEC(fileName, "IntLen", 1, true, 2.5); - plotDiff (fileName, "IntLen"); - plotDiff (fileName, "RadLen"); + etaPhiPlot (fileName, "IntLen", 0, 21, 1, true, maxEta, tag); + etaPhiPlot (fileName, "IntLen", 0, 20, 1, false, -1.0, tag); + etaPhiPlot (fileName, "RadLen", 0, 21, 1, true, maxEta, tag); + etaPhiPlot (fileName, "RadLen", 0, 20, 1, false, -1.0, tag); + etaPhiPlot (fileName, "StepLen",0, 21, 1, true, -1.0, tag); + etaPhiPlot (fileName, "StepLen",0, 20, 1, false, -1.0, tag); + etaPhiPlotHO(fileName, "IntLen", 1, true, 2.5, tag); + etaPhiPlotEC(fileName, "IntLen", 1, true, 2.5, tag); + plotDiff (fileName, "IntLen", tag); + plotDiff (fileName, "RadLen", tag); printTable (fileName, outputFileName); - etaPhi2DPlot(fileName, "IntLen", 0, 21, 1); - etaPhi2DPlot(fileName, "RadLen", 0, 21, 1); + etaPhi2DPlot(fileName, "IntLen", 0, 21, 1, tag); + etaPhi2DPlot(fileName, "RadLen", 0, 21, 1, tag); } void etaPhiPlot(TString fileName, TString plot, int ifirst, int ilast, - int drawLeg, bool ifEta, double maxEta, bool debug) { + int drawLeg, bool ifEta, double maxEta, std::string tag, + bool debug) { TFile* hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); setStyle(); TString xtit = TString("#eta"); + TString ztit = TString("eta"); TString ytit = "none"; int ymin = 0, ymax = 20, istart = 200; double xh = 0.90; @@ -115,6 +125,7 @@ void etaPhiPlot(TString fileName, TString plot, int ifirst, int ilast, if (!ifEta) { istart += 400; xtit = TString("#phi"); + ztit = TString("phi"); } TLegend *leg = new TLegend(xh-0.13, 0.60, xh, 0.90); @@ -136,6 +147,8 @@ void etaPhiPlot(TString fileName, TString plot, int ifirst, int ilast, prof[nplots]->GetXaxis()->SetRangeUser(-maxEta,maxEta); if (xh < 0.8) prof[nplots]->GetYaxis()->SetTitleOffset(1.7); + else + prof[nplots]->GetYaxis()->SetTitleOffset(1.0); int lay = ii - 1; if (lay > 0 && lay < 20) { sprintf(title, "Layer %d", lay); @@ -163,7 +176,7 @@ void etaPhiPlot(TString fileName, TString plot, int ifirst, int ilast, } } - TString cname = "c_" + plot + xtit; + TString cname = "c_" + plot + ztit + tag; TCanvas *cc1 = new TCanvas(cname, cname, 700, 600); if (xh < 0.8) { cc1->SetLeftMargin(0.15); cc1->SetRightMargin(0.05); @@ -176,7 +189,7 @@ void etaPhiPlot(TString fileName, TString plot, int ifirst, int ilast, } void etaPhiPlotHO(TString fileName, TString plot, int drawLeg, bool ifEta, - double maxEta) { + double maxEta, std::string tag) { TFile* hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); @@ -184,6 +197,7 @@ void etaPhiPlotHO(TString fileName, TString plot, int drawLeg, bool ifEta, int ihid[3] = {2, 18, 20}; TString xtit = TString("#eta"); + TString ztit = TString("eta"); TString ytit = "none"; int ymin = 0, ymax = 20, istart = 200; double xh = 0.90; @@ -200,6 +214,7 @@ void etaPhiPlotHO(TString fileName, TString plot, int drawLeg, bool ifEta, if (!ifEta) { istart += 400; xtit = TString("#phi"); + ztit = TString("phi"); } TLegend *leg = new TLegend(xh-0.25, 0.60, xh, 0.90); @@ -252,7 +267,7 @@ void etaPhiPlotHO(TString fileName, TString plot, int drawLeg, bool ifEta, } } - TString cname = "c_HO" + plot + xtit; + TString cname = "c_HO" + plot + ztit + tag; new TCanvas(cname, cname, 700, 400); prof[0]->Draw("h"); @@ -262,7 +277,7 @@ void etaPhiPlotHO(TString fileName, TString plot, int drawLeg, bool ifEta, } void etaPhiPlotEC(TString fileName, TString plot, int drawLeg, bool ifEta, - double maxEta) { + double maxEta, std::string tag) { TFile* hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); @@ -270,6 +285,7 @@ void etaPhiPlotEC(TString fileName, TString plot, int drawLeg, bool ifEta, int ihid[3] = {0, 1, 2}; TString xtit = TString("#eta"); + TString ztit = TString("eta"); TString ytit = "none"; int ymin = 0, ymax = 7, istart = 200, ymax1 = 5; double xh = 0.90, xh1 = 0.90; @@ -286,6 +302,7 @@ void etaPhiPlotEC(TString fileName, TString plot, int drawLeg, bool ifEta, if (!ifEta) { istart += 400; xtit = TString("#phi"); + ztit = TString("phi"); } TLegend *leg = new TLegend(xh-0.25, 0.60, xh, 0.90); @@ -323,7 +340,7 @@ void etaPhiPlotEC(TString fileName, TString plot, int drawLeg, bool ifEta, nplots++; } - TString cname = "c_EC1" + plot + xtit; + TString cname = "c_EC1" + plot + ztit + tag; new TCanvas(cname, cname, 700, 400); prof[0]->Draw("h"); @@ -365,14 +382,14 @@ void etaPhiPlotEC(TString fileName, TString plot, int drawLeg, bool ifEta, mlg->SetBorderSize(1); mlg->SetFillColor(10); mlg->SetMargin(0.30); mlg->SetTextSize(0.04); mlg->AddEntry(prof1, title, "lf"); - cname = "c_EC2" + plot + xtit; + cname = "c_EC2" + plot + ztit + tag; new TCanvas(cname, cname, 700, 400); prof1->Draw(); if (drawLeg > 0) mlg->Draw("sames"); } void etaPhiPlotHC(TString fileName, TString plot, int drawLeg, bool ifEta, - double maxEta) { + double maxEta, std::string tag) { TFile* hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); @@ -380,6 +397,7 @@ void etaPhiPlotHC(TString fileName, TString plot, int drawLeg, bool ifEta, int ihid[20] = {2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21}; TString xtit = TString("#eta"); + TString ztit = TString("eta"); TString ytit = "none"; int ymin = 0, ymax = 20, istart = 200; double xh = 0.90; @@ -396,6 +414,7 @@ void etaPhiPlotHC(TString fileName, TString plot, int drawLeg, bool ifEta, if (!ifEta) { istart += 400; xtit = TString("#phi"); + ztit = TString("phi"); } TLegend *leg = new TLegend(xh-0.25, 0.70, xh, 0.90); @@ -438,7 +457,7 @@ void etaPhiPlotHC(TString fileName, TString plot, int drawLeg, bool ifEta, nplots++; } - TString cname = "c_HC" + plot + xtit; + TString cname = "c_HC" + plot + ztit + tag; new TCanvas(cname, cname, 700, 400); prof[0]->Draw("h"); @@ -448,7 +467,7 @@ void etaPhiPlotHC(TString fileName, TString plot, int drawLeg, bool ifEta, } void etaPhi2DPlot(TString fileName, TString plot, int ifirst, int ilast, - int drawLeg) { + int drawLeg, std::string tag) { TFile* hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); @@ -456,6 +475,7 @@ void etaPhi2DPlot(TString fileName, TString plot, int ifirst, int ilast, TString xtit = TString("#eta"); TString ytit = TString("#phi"); + TString xytit= TString("etaphi"); TString ztit = TString("HCal Material Budget (#lambda)"); int ymin = 0, ymax = 20, istart = 1000; double xh=0.95, yh=0.95; @@ -498,7 +518,7 @@ void etaPhi2DPlot(TString fileName, TString plot, int ifirst, int ilast, nplots++; } - TString cname = "c_" + plot + xtit + ytit; + TString cname = "c_" + plot + xytit + tag; TCanvas *cc1 = new TCanvas(cname, cname, 700, 600); cc1->SetLeftMargin(0.15); cc1->SetRightMargin(0.05); @@ -509,7 +529,7 @@ void etaPhi2DPlot(TString fileName, TString plot, int ifirst, int ilast, } void etaPhi2DPlot(int nslice, int kslice, TString fileName, TString plot, - int ifirst, int ilast, int drawLeg) { + int ifirst, int ilast, int drawLeg, std::string tag) { char hname[200], title[50]; @@ -617,7 +637,7 @@ void etaPhi2DPlot(int nslice, int kslice, TString fileName, TString plot, ismax = ismin+1; } for (int is=ismin; isSetLeftMargin(0.15); cc1[is]->SetRightMargin(0.05); hist[0][is]->Draw(); @@ -733,7 +753,7 @@ void printTable (TString fileName, TString outputFileName, << nbadR << " discrepancies for X0\n"; } -void plotDiff (TString fileName, TString plot) { +void plotDiff (TString fileName, TString plot, std::string tag) { setStyle(); gStyle->SetTitleOffset(1.0,"Y"); @@ -807,7 +827,7 @@ void plotDiff (TString fileName, TString plot) { mg->Add(gr_eta25, "pc"); leg_mg->AddEntry(gr_eta25,"HE #eta = 26"); - TString cname = "c_diff_" + plot; + TString cname = "c_diff_" + plot + tag; new TCanvas(cname, cname, 700, 400); mg->Draw("a"); mg->GetXaxis()->SetTitle(xtit); @@ -917,7 +937,7 @@ void getDiff (TString fileName, TString plot) { */ } -void plotHE(int flag, int logy, int save) { +void plotHE(int flag, int logy, std::string tag, int save) { double angle[31] = {-2.5,-2.25,-2.00,-1.75,-1.50,-1.25,-1.00,-0.75,-0.50, -0.25,-0.20,-0.15,-0.10,-0.05,-0.025,0,0.025,0.05,0.10, @@ -1054,7 +1074,8 @@ void plotHE(int flag, int logy, int save) { gr[i]->GetXaxis()->SetRangeUser(-3.0,+3.0); } - TCanvas *c1 = new TCanvas("c1", name, 800, 500); + TString cname = "c1" + tag; + TCanvas *c1 = new TCanvas(cname, name, 800, 500); if (logy != 0) gPad->SetLogy(1); gr[kfirst]->Draw("alp"); for (int i=0; i<4; i++) { @@ -1079,11 +1100,11 @@ void plotHE(int flag, int logy, int save) { c1->SaveAs(fname); } - } void etaPhiCastorPlot(TString fileName, TString plot, TString type, - bool etaPlus, int drawLeg, bool ifEta, bool debug) { + bool etaPlus, int drawLeg, bool ifEta, std::string tag, + bool debug) { TFile* hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); @@ -1091,6 +1112,7 @@ void etaPhiCastorPlot(TString fileName, TString plot, TString type, if (debug) std::cout << fileName << " is opened at " << hcalFile << std::endl; TString xtit = TString("#eta"); + TString ztit = TString("eta"); char ytit[80], ytpart[10]; int ymin = 0, ymax = 20, istart = 200, ifirst=0; double xh = 0.90; @@ -1120,6 +1142,7 @@ void etaPhiCastorPlot(TString fileName, TString plot, TString type, if (!ifEta) { istart += 400; xtit = TString("#phi"); + ztit = TString("phi"); } if (debug) std::cout << "Title (x) " << xtit << " (y) " << ytit << " First " << ifirst << ":" << istart << std::endl; @@ -1151,7 +1174,7 @@ void etaPhiCastorPlot(TString fileName, TString plot, TString type, nplots++; } - TString cname = "c_" + plot + xtit; + TString cname = "c_" + plot + ztit + tag; TCanvas *cc1 = new TCanvas(cname, cname, 700, 600); if (xh < 0.8) { cc1->SetLeftMargin(0.15); cc1->SetRightMargin(0.05); @@ -1164,7 +1187,7 @@ void etaPhiCastorPlot(TString fileName, TString plot, TString type, } void efficiencyPlot(TString fileName, TString type, bool ifEtaPhi, - double maxEta, bool debug) { + double maxEta, std::string tag, bool debug) { TFile* hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); @@ -1173,15 +1196,15 @@ void efficiencyPlot(TString fileName, TString type, bool ifEtaPhi, int id0=1300, idpl1=8, idpl2=0; char hname[100], title[100]; if (type.CompareTo("HB") == 0) { - idpl1 = 1; idpl2 =2; sprintf (title, "Efficiency for HB"); + idpl1 = 1; idpl2 = 2; sprintf (title, "Efficiency for HB (%s)", tag.c_str()); } else if (type.CompareTo("HE") == 0) { - idpl1 = 3; idpl2 = 4; sprintf (title, "Efficiency for HE"); + idpl1 = 3; idpl2 = 4; sprintf (title, "Efficiency for HE (%s)", tag.c_str()); } else if (type.CompareTo("HO") == 0) { - idpl1 = 5; idpl2 = 6; sprintf (title, "Efficiency for HO"); + idpl1 = 5; idpl2 = 6; sprintf (title, "Efficiency for HO (%s)", tag.c_str()); } else if (type.CompareTo("HF") == 0) { - idpl1 = 7; sprintf (title, "Efficiency for HF"); + idpl1 = 7; sprintf (title, "Efficiency for HF (%s)", tag.c_str()); } else { - sprintf (title, "Efficiency for HCAL"); + sprintf (title, "Efficiency for HCAL (%s)", tag.c_str()); } TLegend *leg = new TLegend(0.70, 0.82, 0.90, 0.90); leg->SetBorderSize(1); leg->SetFillColor(10); leg->SetMargin(0.25); @@ -1296,13 +1319,14 @@ void efficiencyPlot(TString fileName, TString type, bool ifEtaPhi, } void etaPhiFwdPlot(TString fileName, TString plot, int first, int last, - int drawLeg, bool debug) { + int drawLeg, std::string tag, bool debug) { TFile* hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); setStyle(); TString xtit = TString("#eta"); + TString ztit = TString("eta"); TString ytit = "none"; int ymin = 0, ymax = 20, istart = 200; double xh = 0.90, xl=0.1; @@ -1358,7 +1382,7 @@ void etaPhiFwdPlot(TString fileName, TString plot, int first, int last, nplots++; } - TString cname = "c_" + plot + xtit; + TString cname = "c_" + plot + ztit + tag; TCanvas *cc1 = new TCanvas(cname, cname, 700, 600); if (xh < 0.8) { cc1->SetLeftMargin(0.15); cc1->SetRightMargin(0.05); diff --git a/Validation/Geometry/test/runP_Calo_cfg.py b/Validation/Geometry/test/runP_Calo_cfg.py index fb898c060ff6e..16f16e950205f 100644 --- a/Validation/Geometry/test/runP_Calo_cfg.py +++ b/Validation/Geometry/test/runP_Calo_cfg.py @@ -50,7 +50,8 @@ process.p1 = cms.Path(process.g4SimHits) process.g4SimHits.UseMagneticField = False -#process.g4SimHits.Physics.type = 'SimG4Core/Physics/DummyPhysics' +process.g4SimHits.StackingAction.TrackNeutrino = True +process.g4SimHits.Physics.type = 'SimG4Core/Physics/DummyPhysics' process.g4SimHits.Physics.DummyEMPhysics = True process.g4SimHits.Physics.CutsPerRegion = False process.g4SimHits.Generator.ApplyEtaCuts = False @@ -71,5 +72,3 @@ ), type = cms.string('MaterialBudget') )) - - diff --git a/Validation/Geometry/test/runP_Castor_Debug_cfg.py b/Validation/Geometry/test/runP_Castor_Debug_cfg.py index 66387251ee2be..e006dde5346e7 100644 --- a/Validation/Geometry/test/runP_Castor_Debug_cfg.py +++ b/Validation/Geometry/test/runP_Castor_Debug_cfg.py @@ -2,6 +2,7 @@ process = cms.Process("PROD") process.load("SimGeneral.HepPDTESSource.pythiapdt_cfi") +process.load('FWCore.MessageService.MessageLogger_cfi') #Geometry # @@ -20,26 +21,8 @@ process.RandomNumberGeneratorService.generator.initialSeed = 456789 process.RandomNumberGeneratorService.g4SimHits.initialSeed = 9876 -process.MessageLogger = cms.Service("MessageLogger", - debugModules = cms.untracked.vstring('*'), - cout = cms.untracked.PSet( - threshold = cms.untracked.string('DEBUG'), - default = cms.untracked.PSet( - limit = cms.untracked.int32(0) - ), - G4cerr = cms.untracked.PSet( - limit = cms.untracked.int32(0) - ), - G4cout = cms.untracked.PSet( - limit = cms.untracked.int32(0) - ), - MaterialBudget = cms.untracked.PSet( - limit = cms.untracked.int32(-1) - ) - ), - categories = cms.untracked.vstring('G4cerr','G4cout','MaterialBudget'), - destinations = cms.untracked.vstring('cout') -) +if hasattr(process,'MessageLogger'): + process.MessageLogger.categories.append('MaterialBudget') process.source = cms.Source("EmptySource", firstRun = cms.untracked.uint32(1), @@ -70,6 +53,7 @@ process.p1 = cms.Path(process.generator*process.g4SimHits) process.g4SimHits.UseMagneticField = False +process.g4SimHits.StackingAction.TrackNeutrino = True process.g4SimHits.Physics.type = 'SimG4Core/Physics/DummyPhysics' process.g4SimHits.Physics.DummyEMPhysics = True process.g4SimHits.Physics.CutsPerRegion = False diff --git a/Validation/Geometry/test/runP_Castor_cfg.py b/Validation/Geometry/test/runP_Castor_cfg.py index 72b2815033828..2975488f10785 100644 --- a/Validation/Geometry/test/runP_Castor_cfg.py +++ b/Validation/Geometry/test/runP_Castor_cfg.py @@ -48,6 +48,7 @@ process.p1 = cms.Path(process.g4SimHits) process.g4SimHits.UseMagneticField = False +process.g4SimHits.StackingAction.TrackNeutrino = True process.g4SimHits.Physics.type = 'SimG4Core/Physics/DummyPhysics' process.g4SimHits.Physics.DummyEMPhysics = True process.g4SimHits.Physics.CutsPerRegion = False diff --git a/Validation/Geometry/test/runP_ECAL_cfg.py b/Validation/Geometry/test/runP_ECAL_cfg.py index 978692e299a8f..2b73094954c0e 100644 --- a/Validation/Geometry/test/runP_ECAL_cfg.py +++ b/Validation/Geometry/test/runP_ECAL_cfg.py @@ -48,6 +48,7 @@ process.p1 = cms.Path(process.g4SimHits) process.g4SimHits.UseMagneticField = False +process.g4SimHits.StackingAction.TrackNeutrino = True process.g4SimHits.Physics.type = 'SimG4Core/Physics/DummyPhysics' process.g4SimHits.Physics.DummyEMPhysics = True process.g4SimHits.Physics.CutsPerRegion = False diff --git a/Validation/Geometry/test/runP_Forward_cfg.py b/Validation/Geometry/test/runP_Forward_cfg.py index bb8b38fb159fc..365029fde46cd 100644 --- a/Validation/Geometry/test/runP_Forward_cfg.py +++ b/Validation/Geometry/test/runP_Forward_cfg.py @@ -65,6 +65,7 @@ process.p1 = cms.Path(process.generator*process.g4SimHits) process.g4SimHits.UseMagneticField = False +process.g4SimHits.StackingAction.TrackNeutrino = True process.g4SimHits.Physics.type = 'SimG4Core/Physics/DummyPhysics' process.g4SimHits.Physics.DummyEMPhysics = True process.g4SimHits.Physics.CutsPerRegion = False diff --git a/Validation/Geometry/test/runP_H2TB_Debug_cfg.py b/Validation/Geometry/test/runP_H2TB_Debug_cfg.py index d82f065f4616a..3b423a7cb1a61 100644 --- a/Validation/Geometry/test/runP_H2TB_Debug_cfg.py +++ b/Validation/Geometry/test/runP_H2TB_Debug_cfg.py @@ -84,6 +84,7 @@ process.p1 = cms.Path(process.generator*process.VtxSmeared*process.g4SimHits) process.g4SimHits.NonBeamEvent = True process.g4SimHits.UseMagneticField = False +process.g4SimHits.StackingAction.TrackNeutrino = True process.g4SimHits.Physics.type = 'SimG4Core/Physics/DummyPhysics' process.g4SimHits.Physics.DummyEMPhysics = True process.g4SimHits.Physics.CutsPerRegion = False diff --git a/Validation/Geometry/test/runP_HCAL_Debug_cfg.py b/Validation/Geometry/test/runP_HCAL_Debug_cfg.py index f165ae54a3736..fddae13f81e24 100644 --- a/Validation/Geometry/test/runP_HCAL_Debug_cfg.py +++ b/Validation/Geometry/test/runP_HCAL_Debug_cfg.py @@ -1,13 +1,14 @@ import FWCore.ParameterSet.Config as cms -process = cms.Process("PROD") +from Configuration.Eras.Era_Run2_2018_cff import Run2_2018 +process = cms.Process('PROD',Run2_2018) process.load("SimGeneral.HepPDTESSource.pythiapdt_cfi") +process.load('FWCore.MessageService.MessageLogger_cfi') #Geometry # -process.load("Geometry.CMSCommonData.cmsIdealGeometryXML_cfi") -process.load("Geometry.TrackerNumberingBuilder.trackerNumberingGeometry_cfi") +process.load("Configuration.Geometry.GeometryExtended2018_cff") #Magnetic Field # @@ -21,20 +22,8 @@ process.RandomNumberGeneratorService.generator.initialSeed = 456789 process.RandomNumberGeneratorService.g4SimHits.initialSeed = 9876 -process.MessageLogger = cms.Service("MessageLogger", - debugModules = cms.untracked.vstring('*'), - cout = cms.untracked.PSet( - threshold = cms.untracked.string('DEBUG'), - default = cms.untracked.PSet( - limit = cms.untracked.int32(0) - ), - MaterialBudget = cms.untracked.PSet( - limit = cms.untracked.int32(-1) - ) - ), - categories = cms.untracked.vstring('MaterialBudget'), - destinations = cms.untracked.vstring('cout') -) +if hasattr(process,'MessageLogger'): + process.MessageLogger.categories.append('MaterialBudget') process.source = cms.Source("EmptySource", firstRun = cms.untracked.uint32(1), @@ -66,6 +55,7 @@ process.p1 = cms.Path(process.generator*process.g4SimHits) process.g4SimHits.UseMagneticField = False process.g4SimHits.Physics.type = 'SimG4Core/Physics/DummyPhysics' +process.g4SimHits.StackingAction.TrackNeutrino = True process.g4SimHits.Physics.DummyEMPhysics = True process.g4SimHits.Physics.CutsPerRegion = False process.g4SimHits.Generator.ApplyEtaCuts = False @@ -84,5 +74,3 @@ ), type = cms.string('MaterialBudgetHcal') )) - - diff --git a/Validation/Geometry/test/runP_HCAL_cfg.py b/Validation/Geometry/test/runP_HCAL_cfg.py index ac755bd630baf..4a520b4d99ad4 100644 --- a/Validation/Geometry/test/runP_HCAL_cfg.py +++ b/Validation/Geometry/test/runP_HCAL_cfg.py @@ -1,18 +1,18 @@ import FWCore.ParameterSet.Config as cms -process = cms.Process("PROD") +from Configuration.Eras.Era_Run3_cff import Run3 +process = cms.Process('PROD',Run3) process.load("SimGeneral.HepPDTESSource.pythiapdt_cfi") +process.load('FWCore.MessageService.MessageLogger_cfi') #Geometry # -process.load("Geometry.CMSCommonData.cmsIdealGeometryXML_cfi") -process.load("Geometry.TrackerNumberingBuilder.trackerNumberingGeometry_cfi") +process.load("Configuration.Geometry.GeometryExtended2021_cff") #Magnetic Field # process.load("Configuration.StandardSequences.MagneticField_cff") - process.load("Configuration.EventContent.EventContent_cff") # Detector simulation (Geant4-based) @@ -22,20 +22,9 @@ process.load("IOMC.RandomEngine.IOMC_cff") process.RandomNumberGeneratorService.g4SimHits.initialSeed = 9876 -process.MessageLogger = cms.Service("MessageLogger", - destinations = cms.untracked.vstring('cout'), - categories = cms.untracked.vstring('MaterialBudget'), -# debugModules = cms.untracked.vstring('*'), - cout = cms.untracked.PSet( -# threshold = cms.untracked.string('DEBUG'), - default = cms.untracked.PSet( - limit = cms.untracked.int32(0) - ), - MaterialBudget = cms.untracked.PSet( - limit = cms.untracked.int32(0) - ) - ) -) +if hasattr(process,'MessageLogger'): + process.MessageLogger.categories.append('HCalGeom') +# process.MessageLogger.categories.append('MaterialBudget') process.source = cms.Source("PoolSource", noEventSort = cms.untracked.bool(True), @@ -48,12 +37,13 @@ ) process.TFileService = cms.Service("TFileService", - fileName = cms.string('matbdg_HCAL.root') + fileName = cms.string('matbdgHCAL_run3.root') ) process.p1 = cms.Path(process.g4SimHits) process.g4SimHits.UseMagneticField = False process.g4SimHits.Physics.type = 'SimG4Core/Physics/DummyPhysics' +process.g4SimHits.StackingAction.TrackNeutrino = True process.g4SimHits.Physics.DummyEMPhysics = True process.g4SimHits.Physics.CutsPerRegion = False process.g4SimHits.Watchers = cms.VPSet(cms.PSet( diff --git a/Validation/Geometry/test/runP_InFrontOfECAL_cfg.py b/Validation/Geometry/test/runP_InFrontOfECAL_cfg.py index a4bd456de20ac..615919bfcf116 100644 --- a/Validation/Geometry/test/runP_InFrontOfECAL_cfg.py +++ b/Validation/Geometry/test/runP_InFrontOfECAL_cfg.py @@ -51,6 +51,7 @@ process.p1 = cms.Path(process.g4SimHits) process.g4SimHits.UseMagneticField = False +process.g4SimHits.StackingAction.TrackNeutrino = True process.g4SimHits.Physics.type = 'SimG4Core/Physics/DummyPhysics' process.g4SimHits.Physics.DummyEMPhysics = True process.g4SimHits.Physics.CutsPerRegion = False From b8b51a8542fb7c7473af0be3f57a3f8f3614e5cf Mon Sep 17 00:00:00 2001 From: Sunanda Date: Wed, 20 May 2020 17:12:04 +0200 Subject: [PATCH 2/2] Code Check --- Validation/Geometry/src/MaterialBudgetHcal.cc | 3 +- .../Geometry/src/MaterialBudgetHcalHistos.cc | 66 +- Validation/Geometry/test/MatBudgetHcal.C | 1537 +++++++++-------- 3 files changed, 896 insertions(+), 710 deletions(-) diff --git a/Validation/Geometry/src/MaterialBudgetHcal.cc b/Validation/Geometry/src/MaterialBudgetHcal.cc index e7ca6573a07a0..e329bd072f06c 100644 --- a/Validation/Geometry/src/MaterialBudgetHcal.cc +++ b/Validation/Geometry/src/MaterialBudgetHcal.cc @@ -22,7 +22,8 @@ MaterialBudgetHcal::MaterialBudgetHcal(const edm::ParameterSet& p) : theHistoHca rMax = m_p.getUntrackedParameter("RMax", 4.5) * CLHEP::m; zMax = m_p.getUntrackedParameter("ZMax", 13.0) * CLHEP::m; bool doHcal = m_p.getUntrackedParameter("DoHCAL", true); - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcal initialized with rMax " << rMax << " mm and zMax " << zMax << " mm doHcal is set to " << doHcal; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcal initialized with rMax " << rMax << " mm and zMax " << zMax + << " mm doHcal is set to " << doHcal; if (doHcal) theHistoHcal = new MaterialBudgetHcalHistos(m_p); else diff --git a/Validation/Geometry/src/MaterialBudgetHcalHistos.cc b/Validation/Geometry/src/MaterialBudgetHcalHistos.cc index 0f1b8c3a5246a..f1985414df346 100644 --- a/Validation/Geometry/src/MaterialBudgetHcalHistos.cc +++ b/Validation/Geometry/src/MaterialBudgetHcalHistos.cc @@ -24,7 +24,10 @@ MaterialBudgetHcalHistos::MaterialBudgetHcalHistos(const edm::ParameterSet& p) { etaHigh = p.getUntrackedParameter("EtaHigh", 5.2); fillHistos = p.getUntrackedParameter("FillHisto", true); printSum = p.getUntrackedParameter("PrintSummary", false); - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: FillHisto : " << fillHistos << " PrintSummary " << printSum << " == Eta plot: NX " << binEta << " Range " << -maxEta << ":" << maxEta << " Phi plot: NX " << binPhi << " Range " << -pi << ":" << pi << " (Eta limit " << etaLow << ":" << etaHigh << ")"; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: FillHisto : " << fillHistos << " PrintSummary " + << printSum << " == Eta plot: NX " << binEta << " Range " << -maxEta << ":" + << maxEta << " Phi plot: NX " << binPhi << " Range " << -pi << ":" << pi + << " (Eta limit " << etaLow << ":" << etaHigh << ")"; if (fillHistos) book(); } @@ -36,7 +39,8 @@ void MaterialBudgetHcalHistos::fillBeginJob(const DDCompactView& cpv) { DDSpecificsMatchesValueFilter filter1{DDValue(attribute, value, 0)}; DDFilteredView fv1(cpv, filter1); sensitives = getNames(fv1); - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Names to be tested for " << attribute << " = " << value << " has " << sensitives.size() << " elements"; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Names to be tested for " << attribute << " = " + << value << " has " << sensitives.size() << " elements"; for (unsigned int i = 0; i < sensitives.size(); i++) edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: sensitives[" << i << "] = " << sensitives[i]; @@ -48,11 +52,13 @@ void MaterialBudgetHcalHistos::fillBeginJob(const DDCompactView& cpv) { fv2.firstChild(); DDsvalues_type sv(fv2.mergedSpecifics()); std::vector temp = getDDDArray("Levels", sv); - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Names to be tested for " << attribute << " = " << value << " has " << hfNames.size() << " elements"; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Names to be tested for " << attribute << " = " + << value << " has " << hfNames.size() << " elements"; for (unsigned int i = 0; i < hfNames.size(); i++) { int level = static_cast(temp[i]); hfLevels.push_back(level); - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: HF[" << i << "] = " << hfNames[i] << " at level " << hfLevels[i]; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: HF[" << i << "] = " << hfNames[i] + << " at level " << hfLevels[i]; } std::string ecalRO[2] = {"EcalHitsEB", "EcalHitsEE"}; @@ -62,7 +68,8 @@ void MaterialBudgetHcalHistos::fillBeginJob(const DDCompactView& cpv) { DDSpecificsMatchesValueFilter filter3{DDValue(attribute, value, 0)}; DDFilteredView fv3(cpv, filter3); std::vector senstmp = getNames(fv3); - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Names to be tested for " << attribute << " = " << value << " has " << senstmp.size() << " elements"; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Names to be tested for " << attribute << " = " + << value << " has " << senstmp.size() << " elements"; for (unsigned int i = 0; i < senstmp.size(); i++) sensitiveEC.push_back(senstmp[i]); } @@ -93,7 +100,9 @@ void MaterialBudgetHcalHistos::fillStartTrack(const G4Track* aTrack) { intLength.clear(); } - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Track " << aTrack->GetTrackID() << " Code " << theID << " Energy " << theEnergy / CLHEP::GeV << " GeV; Eta " << eta << " Phi " << phi / CLHEP::deg << " PT " << dir.perp() / CLHEP::GeV << " GeV *****"; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Track " << aTrack->GetTrackID() << " Code " << theID + << " Energy " << theEnergy / CLHEP::GeV << " GeV; Eta " << eta << " Phi " + << phi / CLHEP::deg << " PT " << dir.perp() / CLHEP::GeV << " GeV *****"; } void MaterialBudgetHcalHistos::fillPerStep(const G4Step* aStep) { @@ -124,9 +133,19 @@ void MaterialBudgetHcalHistos::fillPerStep(const G4Step* aStep) { radLength.push_back(step / radl); intLength.push_back(step / intl); } - edm::LogVerbatim("MaterialBudget") << "Volume " << name << " id " << id << ":" << idOld << " Step " << step << " Material " << matName << " Old Length " << stepLen << " X0 " << step / radl << ":" << radLen << " Lambda " << step / intl << ":" << intLen; + edm::LogVerbatim("MaterialBudget") << "Volume " << name << " id " << id << ":" << idOld << " Step " << step + << " Material " << matName << " Old Length " << stepLen << " X0 " << step / radl + << ":" << radLen << " Lambda " << step / intl << ":" << intLen; } else { - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Step at " << name << " id " << id << ":" << idOld << " Length " << step << " in " << matName << " of density " << density << " g/cc; Radiation Length " << radl << " mm; Interaction Length " << intl << " mm\n" << " Position " << aStep->GetPreStepPoint()->GetPosition() << " Cylindrical R " << aStep->GetPreStepPoint()->GetPosition().perp() << " Length (so far) " << stepLen << " L/X0 " << step / radl << "/" << radLen << " L/Lambda " << step / intl << "/" << intLen; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Step at " << name << " id " << id << ":" << idOld + << " Length " << step << " in " << matName << " of density " << density + << " g/cc; Radiation Length " << radl << " mm; Interaction Length " << intl + << " mm\n" + << " Position " + << aStep->GetPreStepPoint()->GetPosition() << " Cylindrical R " + << aStep->GetPreStepPoint()->GetPosition().perp() << " Length (so far) " + << stepLen << " L/X0 " << step / radl << "/" << radLen << " L/Lambda " + << step / intl << "/" << intLen; } int det = 0, lay = 0; @@ -166,7 +185,8 @@ void MaterialBudgetHcalHistos::fillPerStep(const G4Step* aStep) { } } #ifdef EDM_ML_DEBUG - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Det " << det << " Layer " << lay << " Eta " << eta << " Phi " << phi / CLHEP::deg; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Det " << det << " Layer " << lay << " Eta " + << eta << " Phi " << phi / CLHEP::deg; #endif } else if (layer == 1) { det = -1; @@ -182,7 +202,8 @@ void MaterialBudgetHcalHistos::fillPerStep(const G4Step* aStep) { if (id > idOld) { #ifdef EDM_ML_DEBUG - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Step at " << name << " calls filHisto with " << (id - 1); + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Step at " << name << " calls filHisto with " + << (id - 1); #endif fillHisto(id - 1); } @@ -195,25 +216,29 @@ void MaterialBudgetHcalHistos::fillPerStep(const G4Step* aStep) { if (id == 21) { if (!isItHF(aStep->GetPostStepPoint()->GetTouchable())) { #ifdef EDM_ML_DEBUG - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: After HF in " << name << ":" << aStep->GetPostStepPoint()->GetTouchable()->GetVolume(0)->GetName() << " calls fillHisto with " << id; + edm::LogVerbatim("MaterialBudget") + << "MaterialBudgetHcalHistos: After HF in " << name << ":" + << aStep->GetPostStepPoint()->GetTouchable()->GetVolume(0)->GetName() << " calls fillHisto with " << id; #endif - fillHisto(idOld); - id++; - layer = 0; + fillHisto(idOld); + id++; + layer = 0; } } } } void MaterialBudgetHcalHistos::fillEndTrack() { - edm::LogVerbatim("MaterialBudget") << "Number of layers hit in HB:" << nlayHB << " HE:" << nlayHE << " HO:" << nlayHO << " HF:" << nlayHF; + edm::LogVerbatim("MaterialBudget") << "Number of layers hit in HB:" << nlayHB << " HE:" << nlayHE << " HO:" << nlayHO + << " HF:" << nlayHF; if (fillHistos) { fillHisto(maxSet - 1); fillLayer(); } if (printSum) { for (unsigned int ii = 0; ii < matList.size(); ii++) { - edm::LogVerbatim("MaterialBudget") << matList[ii] << "\t" << stepLength[ii] << "\t" << radLength[ii] << "\t" << intLength[ii]; + edm::LogVerbatim("MaterialBudget") << matList[ii] << "\t" << stepLength[ii] << "\t" << radLength[ii] << "\t" + << intLength[ii]; } } } @@ -227,7 +252,9 @@ void MaterialBudgetHcalHistos::book() { << "please add it to config file"; double maxPhi = pi; - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Booking user histos === with " << binEta << " bins in eta from " << -maxEta << " to " << maxEta << " and " << binPhi << " bins in phi from " << -maxPhi << " to " << maxPhi; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos: Booking user histos === with " << binEta + << " bins in eta from " << -maxEta << " to " << maxEta << " and " << binPhi + << " bins in phi from " << -maxPhi << " to " << maxPhi; std::string iter; // total X0 @@ -309,7 +336,8 @@ void MaterialBudgetHcalHistos::book() { void MaterialBudgetHcalHistos::fillHisto(int ii) { #ifdef EDM_ML_DEBUG - edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos:FillHisto called with index " << ii << " integrated step " << stepLen << " X0 " << radLen << " Lamda " << intLen; + edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalHistos:FillHisto called with index " << ii + << " integrated step " << stepLen << " X0 " << radLen << " Lamda " << intLen; #endif if (ii >= 0 && ii < maxSet) { me100[ii]->Fill(eta, radLen); @@ -403,7 +431,7 @@ std::vector MaterialBudgetHcalHistos::getDDDArray(const std::string& str const std::vector& fvec = value.doubles(); int nval = fvec.size(); if (nval < 1) { - throw cms::Exception("MaterialBudgetHcalHistos") << "nval = " << nval << " < 1 for array " << str << "\n"; + throw cms::Exception("MaterialBudgetHcalHistos") << "nval = " << nval << " < 1 for array " << str << "\n"; } return fvec; diff --git a/Validation/Geometry/test/MatBudgetHcal.C b/Validation/Geometry/test/MatBudgetHcal.C index 28b4b61d1a2f4..c57abb8e4124f 100644 --- a/Validation/Geometry/test/MatBudgetHcal.C +++ b/Validation/Geometry/test/MatBudgetHcal.C @@ -18,92 +18,116 @@ const int nlaymax = 25; const int nbinmax = 41; -double mean[nlaymax][nbinmax], diff[nlaymax][nbinmax]; - -double towLow[41] = { 0.000, 0.087, 0.174, 0.261, 0.348, - 0.435, 0.522, 0.609, 0.696, 0.783, - 0.870, 0.957, 1.044, 1.131, 1.218, - 1.305, 1.392, 1.479, 1.566, 1.653, - 1.740, 1.830, 1.930, 2.043, 2.172, - 2.322, 2.500, 2.650, 2.853, 2.964, - 3.139, 3.314, 3.489, 3.664, 3.839, - 4.013, 4.191, 4.363, 4.538, 4.716, - 4.889}; -double towHigh[41] = { 0.087, 0.174, 0.261, 0.348, 0.435, - 0.522, 0.609, 0.696, 0.783, 0.870, - 0.957, 1.044, 1.131, 1.218, 1.305, - 1.392, 1.479, 1.566, 1.653, 1.740, - 1.830, 1.930, 2.043, 2.172, 2.322, - 2.500, 2.650, 3.000, 2.964, 3.139, - 3.314, 3.489, 3.664, 3.839, 4.013, - 4.191, 4.363, 4.538, 4.716, 4.889, - 5.191}; - -int colorLayer[25] = { 2, 7, 9, 30, 34, 38, 14, 40, 41, 42, - 45, 46, 8, 49, 37, 28, 4, 1, 48, 50, - 3, 6, 5, 156, 159}; - -void etaPhiPlot(TString fileName="matbdg_HCAL.root", TString plot="IntLen", - int ifirst=0, int ilast=21, int drawLeg=1, bool ifEta=true, - double maxEta=-1, std::string tag="Run", bool debug=false); -void etaPhiPlotHO(TString fileName="matbdg_HCAL.root", TString plot="IntLen", - int drawLeg=1, bool ifEta=true, double maxEta=-1, - std::string tag="Run"); -void etaPhiPlotEC(TString fileName="matbdg_HCAL.root", TString plot="IntLen", - int drawLeg=1, bool ifEta=true, double maxEta=-1, - std::string tag="Run"); -void etaPhiPlotHC(TString fileName="matbdg_HCAL.root", TString plot="IntLen", - int drawLeg=1, bool ifEta=true, double maxEta=-1, - std::string tag="Run"); -void etaPhi2DPlot(TString fileName="matbdg_HCAL.root", TString plot="IntLen", - int ifirst=0, int ilast=19, int drawLeg=1, - std::string tag="Run"); -void etaPhi2DPlot(int nslice, int kslice, TString fileName="matbdg_HCAL.root", - TString plot="IntLen", int ifirst=0, int ilast=21, - int drawLeg=1, std::string tag="Run"); -void printTable (TString fileName="matbdg_HCAL.root", - TString outputFileName="hcal.txt", - TString inputFileName="None"); -void plotDiff (TString fileName="matbdg_HCAL.root", TString plot="IntLen", - std::string tag="Run"); -void getDiff (TString fileName="matbdg_HCAL.root", TString plot="IntLen"); -void plotHE(int flag=0, int logy=0, std::string tag="Run", int save=0); -void etaPhiCastorPlot(TString fileName="matbdg_Castor.root", - TString plot="IntLen", TString type="All", - bool etaPlus=true, int drawLeg=1, bool ifEta=true, - std::string tag="Run", bool debug=false); -void efficiencyPlot(TString fileName="matbdg_HCAL.root", TString type="All", - bool ifEtaPhi=true, double maxEta=-1, std::string tag="Run", - bool debug=false); -void etaPhiFwdPlot(TString fileName="matbdg_Fwd.root", TString plot="IntLen", - int first=0, int last=9, int drawLeg=1, std::string tag="Run" - , bool debug=false); -void setStyle (); - -void standardPlot (TString fileName="matbdgHCAL_run3.root", - std::string tag="Run3", int maxEta=4.8, - TString outputFileName="hcal.txt") { - - etaPhiPlot (fileName, "IntLen", 0, 21, 1, true, maxEta, tag); - etaPhiPlot (fileName, "IntLen", 0, 20, 1, false, -1.0, tag); - etaPhiPlot (fileName, "RadLen", 0, 21, 1, true, maxEta, tag); - etaPhiPlot (fileName, "RadLen", 0, 20, 1, false, -1.0, tag); - etaPhiPlot (fileName, "StepLen",0, 21, 1, true, -1.0, tag); - etaPhiPlot (fileName, "StepLen",0, 20, 1, false, -1.0, tag); +double mean[nlaymax][nbinmax], diff[nlaymax][nbinmax]; + +double towLow[41] = {0.000, 0.087, 0.174, 0.261, 0.348, 0.435, 0.522, 0.609, 0.696, 0.783, 0.870, 0.957, 1.044, 1.131, + 1.218, 1.305, 1.392, 1.479, 1.566, 1.653, 1.740, 1.830, 1.930, 2.043, 2.172, 2.322, 2.500, 2.650, + 2.853, 2.964, 3.139, 3.314, 3.489, 3.664, 3.839, 4.013, 4.191, 4.363, 4.538, 4.716, 4.889}; +double towHigh[41] = {0.087, 0.174, 0.261, 0.348, 0.435, 0.522, 0.609, 0.696, 0.783, 0.870, 0.957, 1.044, 1.131, 1.218, + 1.305, 1.392, 1.479, 1.566, 1.653, 1.740, 1.830, 1.930, 2.043, 2.172, 2.322, 2.500, 2.650, 3.000, + 2.964, 3.139, 3.314, 3.489, 3.664, 3.839, 4.013, 4.191, 4.363, 4.538, 4.716, 4.889, 5.191}; + +int colorLayer[25] = {2, 7, 9, 30, 34, 38, 14, 40, 41, 42, 45, 46, 8, 49, 37, 28, 4, 1, 48, 50, 3, 6, 5, 156, 159}; + +void etaPhiPlot(TString fileName = "matbdg_HCAL.root", + TString plot = "IntLen", + int ifirst = 0, + int ilast = 21, + int drawLeg = 1, + bool ifEta = true, + double maxEta = -1, + std::string tag = "Run", + bool debug = false); +void etaPhiPlotHO(TString fileName = "matbdg_HCAL.root", + TString plot = "IntLen", + int drawLeg = 1, + bool ifEta = true, + double maxEta = -1, + std::string tag = "Run"); +void etaPhiPlotEC(TString fileName = "matbdg_HCAL.root", + TString plot = "IntLen", + int drawLeg = 1, + bool ifEta = true, + double maxEta = -1, + std::string tag = "Run"); +void etaPhiPlotHC(TString fileName = "matbdg_HCAL.root", + TString plot = "IntLen", + int drawLeg = 1, + bool ifEta = true, + double maxEta = -1, + std::string tag = "Run"); +void etaPhi2DPlot(TString fileName = "matbdg_HCAL.root", + TString plot = "IntLen", + int ifirst = 0, + int ilast = 19, + int drawLeg = 1, + std::string tag = "Run"); +void etaPhi2DPlot(int nslice, + int kslice, + TString fileName = "matbdg_HCAL.root", + TString plot = "IntLen", + int ifirst = 0, + int ilast = 21, + int drawLeg = 1, + std::string tag = "Run"); +void printTable(TString fileName = "matbdg_HCAL.root", + TString outputFileName = "hcal.txt", + TString inputFileName = "None"); +void plotDiff(TString fileName = "matbdg_HCAL.root", TString plot = "IntLen", std::string tag = "Run"); +void getDiff(TString fileName = "matbdg_HCAL.root", TString plot = "IntLen"); +void plotHE(int flag = 0, int logy = 0, std::string tag = "Run", int save = 0); +void etaPhiCastorPlot(TString fileName = "matbdg_Castor.root", + TString plot = "IntLen", + TString type = "All", + bool etaPlus = true, + int drawLeg = 1, + bool ifEta = true, + std::string tag = "Run", + bool debug = false); +void efficiencyPlot(TString fileName = "matbdg_HCAL.root", + TString type = "All", + bool ifEtaPhi = true, + double maxEta = -1, + std::string tag = "Run", + bool debug = false); +void etaPhiFwdPlot(TString fileName = "matbdg_Fwd.root", + TString plot = "IntLen", + int first = 0, + int last = 9, + int drawLeg = 1, + std::string tag = "Run", + bool debug = false); +void setStyle(); + +void standardPlot(TString fileName = "matbdgHCAL_run3.root", + std::string tag = "Run3", + int maxEta = 4.8, + TString outputFileName = "hcal.txt") { + etaPhiPlot(fileName, "IntLen", 0, 21, 1, true, maxEta, tag); + etaPhiPlot(fileName, "IntLen", 0, 20, 1, false, -1.0, tag); + etaPhiPlot(fileName, "RadLen", 0, 21, 1, true, maxEta, tag); + etaPhiPlot(fileName, "RadLen", 0, 20, 1, false, -1.0, tag); + etaPhiPlot(fileName, "StepLen", 0, 21, 1, true, -1.0, tag); + etaPhiPlot(fileName, "StepLen", 0, 20, 1, false, -1.0, tag); etaPhiPlotHO(fileName, "IntLen", 1, true, 2.5, tag); etaPhiPlotEC(fileName, "IntLen", 1, true, 2.5, tag); - plotDiff (fileName, "IntLen", tag); - plotDiff (fileName, "RadLen", tag); - printTable (fileName, outputFileName); + plotDiff(fileName, "IntLen", tag); + plotDiff(fileName, "RadLen", tag); + printTable(fileName, outputFileName); etaPhi2DPlot(fileName, "IntLen", 0, 21, 1, tag); etaPhi2DPlot(fileName, "RadLen", 0, 21, 1, tag); } -void etaPhiPlot(TString fileName, TString plot, int ifirst, int ilast, - int drawLeg, bool ifEta, double maxEta, std::string tag, - bool debug) { - - TFile* hcalFile = new TFile(fileName); +void etaPhiPlot(TString fileName, + TString plot, + int ifirst, + int ilast, + int drawLeg, + bool ifEta, + double maxEta, + std::string tag, + bool debug) { + TFile *hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); setStyle(); @@ -114,38 +138,47 @@ void etaPhiPlot(TString fileName, TString plot, int ifirst, int ilast, double xh = 0.90; if (plot.CompareTo("RadLen") == 0) { ytit = TString("HCal Material Budget (X_{0})"); - ymin = 0; ymax = 200; istart = 100; + ymin = 0; + ymax = 200; + istart = 100; } else if (plot.CompareTo("StepLen") == 0) { ytit = TString("HCal Material Budget (Step Length)"); - ymin = 0; ymax = 15000; istart = 300; xh = 0.70; + ymin = 0; + ymax = 15000; + istart = 300; + xh = 0.70; } else { ytit = TString("HCal Material Budget (#lambda)"); - ymin = 0; ymax = 20; istart = 200; + ymin = 0; + ymax = 20; + istart = 200; } if (!ifEta) { istart += 400; - xtit = TString("#phi"); - ztit = TString("phi"); + xtit = TString("#phi"); + ztit = TString("phi"); } - - TLegend *leg = new TLegend(xh-0.13, 0.60, xh, 0.90); - leg->SetBorderSize(1); leg->SetFillColor(10); leg->SetMargin(0.25); + + TLegend *leg = new TLegend(xh - 0.13, 0.60, xh, 0.90); + leg->SetBorderSize(1); + leg->SetFillColor(10); + leg->SetMargin(0.25); leg->SetTextSize(0.018); - int nplots=0; + int nplots = 0; TProfile *prof[nlaymax]; - for (int ii=ilast; ii>=ifirst; ii--) { + for (int ii = ilast; ii >= ifirst; ii--) { char hname[10], title[50]; - sprintf(hname, "%i", istart+ii); - gDirectory->GetObject(hname,prof[nplots]); + sprintf(hname, "%i", istart + ii); + gDirectory->GetObject(hname, prof[nplots]); prof[nplots]->GetXaxis()->SetTitle(xtit); prof[nplots]->GetYaxis()->SetTitle(ytit); prof[nplots]->GetYaxis()->SetRangeUser(ymin, ymax); prof[nplots]->SetLineColor(colorLayer[ii]); prof[nplots]->SetFillColor(colorLayer[ii]); - if (ifEta && maxEta > 0) - prof[nplots]->GetXaxis()->SetRangeUser(-maxEta,maxEta); - if (xh < 0.8) + if (ifEta && maxEta > 0) + prof[nplots]->GetXaxis()->SetRangeUser(-maxEta, maxEta); + if (xh < 0.8) prof[nplots]->GetYaxis()->SetTitleOffset(1.7); else prof[nplots]->GetYaxis()->SetTitleOffset(1.0); @@ -154,7 +187,7 @@ void etaPhiPlot(TString fileName, TString plot, int ifirst, int ilast, sprintf(title, "Layer %d", lay); } else if (lay == 0) { sprintf(title, "After Crystal"); - } else if (lay >= 20 ) { + } else if (lay >= 20) { sprintf(title, "After HF"); } else { sprintf(title, "Before Crystal"); @@ -162,15 +195,15 @@ void etaPhiPlot(TString fileName, TString plot, int ifirst, int ilast, leg->AddEntry(prof[nplots], title, "lf"); nplots++; if (ii == ilast && debug) { - int nbinX = prof[0]->GetNbinsX(); - double xmin = prof[0]->GetXaxis()->GetXmin(); - double xmax = prof[0]->GetXaxis()->GetXmax(); - double dx = (xmax - xmin)/nbinX; + int nbinX = prof[0]->GetNbinsX(); + double xmin = prof[0]->GetXaxis()->GetXmin(); + double xmax = prof[0]->GetXaxis()->GetXmax(); + double dx = (xmax - xmin) / nbinX; cout << "Hist " << ii; - for (int ibx=0; ibxGetBinContent(ibx+1); - cout << " | " << ibx << "(" << xx1 << ":" << (xx1+dx) << ") " << cont; + for (int ibx = 0; ibx < nbinX; ibx++) { + double xx1 = xmin + ibx * dx; + double cont = prof[0]->GetBinContent(ibx + 1); + cout << " | " << ibx << "(" << xx1 << ":" << (xx1 + dx) << ") " << cont; } cout << "\n"; } @@ -179,19 +212,19 @@ void etaPhiPlot(TString fileName, TString plot, int ifirst, int ilast, TString cname = "c_" + plot + ztit + tag; TCanvas *cc1 = new TCanvas(cname, cname, 700, 600); if (xh < 0.8) { - cc1->SetLeftMargin(0.15); cc1->SetRightMargin(0.05); + cc1->SetLeftMargin(0.15); + cc1->SetRightMargin(0.05); } prof[0]->Draw("h"); - for(int i=1; iDraw("h sames"); - if (drawLeg > 0) leg->Draw("sames"); + if (drawLeg > 0) + leg->Draw("sames"); } -void etaPhiPlotHO(TString fileName, TString plot, int drawLeg, bool ifEta, - double maxEta, std::string tag) { - - TFile* hcalFile = new TFile(fileName); +void etaPhiPlotHO(TString fileName, TString plot, int drawLeg, bool ifEta, double maxEta, std::string tag) { + TFile *hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); setStyle(); @@ -203,30 +236,39 @@ void etaPhiPlotHO(TString fileName, TString plot, int drawLeg, bool ifEta, double xh = 0.90; if (plot.CompareTo("RadLen") == 0) { ytit = TString("Material Budget (X_{0})"); - ymin = 0; ymax = 200; istart = 100; + ymin = 0; + ymax = 200; + istart = 100; } else if (plot.CompareTo("StepLen") == 0) { ytit = TString("Material Budget (Step Length)"); - ymin = 0; ymax = 15000; istart = 300; xh = 0.70; + ymin = 0; + ymax = 15000; + istart = 300; + xh = 0.70; } else { ytit = TString("Material Budget (#lambda)"); - ymin = 0; ymax = 20; istart = 200; + ymin = 0; + ymax = 20; + istart = 200; } if (!ifEta) { istart += 400; - xtit = TString("#phi"); - ztit = TString("phi"); + xtit = TString("#phi"); + ztit = TString("phi"); } - - TLegend *leg = new TLegend(xh-0.25, 0.60, xh, 0.90); - leg->SetBorderSize(1); leg->SetFillColor(10); leg->SetMargin(0.45); + + TLegend *leg = new TLegend(xh - 0.25, 0.60, xh, 0.90); + leg->SetBorderSize(1); + leg->SetFillColor(10); + leg->SetMargin(0.45); leg->SetTextSize(0.04); - int nplots=0; + int nplots = 0; TProfile *prof[nlaymax]; - for (int ii=2; ii>=0; ii--) { + for (int ii = 2; ii >= 0; ii--) { char hname[10], title[50]; - int idpl = istart+ihid[ii]; + int idpl = istart + ihid[ii]; sprintf(hname, "%i", idpl); - gDirectory->GetObject(hname,prof[nplots]); + gDirectory->GetObject(hname, prof[nplots]); prof[nplots]->GetXaxis()->SetTitle(xtit); prof[nplots]->GetYaxis()->SetTitle(ytit); prof[nplots]->GetXaxis()->SetLabelSize(0.05); @@ -237,10 +279,10 @@ void etaPhiPlotHO(TString fileName, TString plot, int drawLeg, bool ifEta, prof[nplots]->GetYaxis()->SetRangeUser(ymin, ymax); prof[nplots]->SetLineColor(colorLayer[ii]); prof[nplots]->SetFillColor(colorLayer[ii]); - if (xh < 0.8) + if (xh < 0.8) prof[nplots]->GetYaxis()->SetTitleOffset(1.7); - if (ifEta && maxEta > 0) - prof[nplots]->GetXaxis()->SetRangeUser(0.0,maxEta); + if (ifEta && maxEta > 0) + prof[nplots]->GetXaxis()->SetRangeUser(0.0, maxEta); if (ii >= 2) { sprintf(title, "With HO"); } else if (ii == 1) { @@ -251,18 +293,18 @@ void etaPhiPlotHO(TString fileName, TString plot, int drawLeg, bool ifEta, leg->AddEntry(prof[nplots], title, "lf"); nplots++; if (ii == 1) { - for (int kk=0; kk<7; kk++) { - idpl--; - sprintf(hname, "%i", idpl); - gDirectory->GetObject(hname,prof[nplots]); - prof[nplots]->GetXaxis()->SetTitle(xtit); - prof[nplots]->GetYaxis()->SetTitle(ytit); - prof[nplots]->GetYaxis()->SetRangeUser(ymin, ymax); - prof[nplots]->SetLineColor(colorLayer[ii]); - prof[nplots]->SetFillColor(colorLayer[ii]); - if (ifEta && maxEta > 0) - prof[nplots]->GetXaxis()->SetRangeUser(0.0,maxEta); - nplots++; + for (int kk = 0; kk < 7; kk++) { + idpl--; + sprintf(hname, "%i", idpl); + gDirectory->GetObject(hname, prof[nplots]); + prof[nplots]->GetXaxis()->SetTitle(xtit); + prof[nplots]->GetYaxis()->SetTitle(ytit); + prof[nplots]->GetYaxis()->SetRangeUser(ymin, ymax); + prof[nplots]->SetLineColor(colorLayer[ii]); + prof[nplots]->SetFillColor(colorLayer[ii]); + if (ifEta && maxEta > 0) + prof[nplots]->GetXaxis()->SetRangeUser(0.0, maxEta); + nplots++; } } } @@ -271,15 +313,14 @@ void etaPhiPlotHO(TString fileName, TString plot, int drawLeg, bool ifEta, new TCanvas(cname, cname, 700, 400); prof[0]->Draw("h"); - for(int i=1; iDraw("h sames"); - if (drawLeg > 0) leg->Draw("sames"); + if (drawLeg > 0) + leg->Draw("sames"); } -void etaPhiPlotEC(TString fileName, TString plot, int drawLeg, bool ifEta, - double maxEta, std::string tag) { - - TFile* hcalFile = new TFile(fileName); +void etaPhiPlotEC(TString fileName, TString plot, int drawLeg, bool ifEta, double maxEta, std::string tag) { + TFile *hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); setStyle(); @@ -291,30 +332,41 @@ void etaPhiPlotEC(TString fileName, TString plot, int drawLeg, bool ifEta, double xh = 0.90, xh1 = 0.90; if (plot.CompareTo("RadLen") == 0) { ytit = TString("Material Budget (X_{0})"); - ymin = 0; ymax = 70; istart = 100; ymax1 = 50; + ymin = 0; + ymax = 70; + istart = 100; + ymax1 = 50; } else if (plot.CompareTo("StepLen") == 0) { ytit = TString("Material Budget (Step Length)"); - ymin = 0; ymax = 6000; istart = 300; xh = 0.35; ymax1 = 2500; + ymin = 0; + ymax = 6000; + istart = 300; + xh = 0.35; + ymax1 = 2500; } else { ytit = TString("Material Budget (#lambda)"); - ymin = 0; ymax = 7; istart = 200; + ymin = 0; + ymax = 7; + istart = 200; } if (!ifEta) { istart += 400; - xtit = TString("#phi"); - ztit = TString("phi"); + xtit = TString("#phi"); + ztit = TString("phi"); } - - TLegend *leg = new TLegend(xh-0.25, 0.60, xh, 0.90); - leg->SetBorderSize(1); leg->SetFillColor(10); leg->SetMargin(0.30); + + TLegend *leg = new TLegend(xh - 0.25, 0.60, xh, 0.90); + leg->SetBorderSize(1); + leg->SetFillColor(10); + leg->SetMargin(0.30); leg->SetTextSize(0.04); - int nplots=0; + int nplots = 0; TProfile *prof[nlaymax]; - for (int ii=2; ii>=0; ii--) { + for (int ii = 2; ii >= 0; ii--) { char hname[10], title[50]; - int idpl = istart+ihid[ii]; + int idpl = istart + ihid[ii]; sprintf(hname, "%i", idpl); - gDirectory->GetObject(hname,prof[nplots]); + gDirectory->GetObject(hname, prof[nplots]); prof[nplots]->GetXaxis()->SetTitle(xtit); prof[nplots]->GetYaxis()->SetTitle(ytit); prof[nplots]->GetXaxis()->SetLabelSize(0.05); @@ -325,10 +377,10 @@ void etaPhiPlotEC(TString fileName, TString plot, int drawLeg, bool ifEta, prof[nplots]->GetYaxis()->SetRangeUser(ymin, ymax); prof[nplots]->SetLineColor(colorLayer[ii]); prof[nplots]->SetFillColor(colorLayer[ii]); - if (xh < 0.8) + if (xh < 0.8) prof[nplots]->GetYaxis()->SetTitleOffset(1.05); - if (ifEta && maxEta > 0) - prof[nplots]->GetXaxis()->SetRangeUser(0.0,maxEta); + if (ifEta && maxEta > 0) + prof[nplots]->GetXaxis()->SetRangeUser(0.0, maxEta); if (ii >= 2) { sprintf(title, "Front of HCAL"); } else if (ii == 1) { @@ -344,24 +396,26 @@ void etaPhiPlotEC(TString fileName, TString plot, int drawLeg, bool ifEta, new TCanvas(cname, cname, 700, 400); prof[0]->Draw("h"); - for(int i=1; iDraw("h sames"); - if (drawLeg > 0) leg->Draw("sames"); + if (drawLeg > 0) + leg->Draw("sames"); - double xmin = prof[2]->GetXaxis()->GetXmin(); - double xmax = prof[2]->GetXaxis()->GetXmax(); - int nbins = prof[2]->GetNbinsX(); + double xmin = prof[2]->GetXaxis()->GetXmin(); + double xmax = prof[2]->GetXaxis()->GetXmax(); + int nbins = prof[2]->GetNbinsX(); TH1D *prof1 = new TH1D("Temp01", "Temp01", nbins, xmin, xmax); - for (int ii=0; iiGetBinLowEdge(ii+1); - double x2 = prof[0]->GetBinWidth(ii+1); - double v0 = prof[0]->GetBinContent(ii+1); - double v1 = prof[1]->GetBinContent(ii+1); - double v2 = prof[2]->GetBinContent(ii+1); - double xx = x1+0.5*x2; + for (int ii = 0; ii < nbins; ii++) { + double x1 = prof[0]->GetBinLowEdge(ii + 1); + double x2 = prof[0]->GetBinWidth(ii + 1); + double v0 = prof[0]->GetBinContent(ii + 1); + double v1 = prof[1]->GetBinContent(ii + 1); + double v2 = prof[2]->GetBinContent(ii + 1); + double xx = x1 + 0.5 * x2; double cont = v0 - v1; - prof1->Fill(xx,cont); - std::cout << "Bin " << ii << " Eta/Phi " << std::setw(4) << xx << " Material " << std::setw(6) << cont << " " << std::setw(6) << (v0-v2) << "\n"; + prof1->Fill(xx, cont); + std::cout << "Bin " << ii << " Eta/Phi " << std::setw(4) << xx << " Material " << std::setw(6) << cont << " " + << std::setw(6) << (v0 - v2) << "\n"; } prof1->GetXaxis()->SetTitle(xtit); prof1->GetYaxis()->SetTitle(ytit); @@ -373,25 +427,29 @@ void etaPhiPlotEC(TString fileName, TString plot, int drawLeg, bool ifEta, prof1->GetYaxis()->SetRangeUser(ymin, ymax1); prof1->SetLineColor(colorLayer[3]); prof1->SetFillColor(colorLayer[3]); - if (ifEta && maxEta > 0) prof1->GetXaxis()->SetRangeUser(0.0,maxEta); - if (xh < 0.8) prof1->GetYaxis()->SetTitleOffset(1.05); + if (ifEta && maxEta > 0) + prof1->GetXaxis()->SetRangeUser(0.0, maxEta); + if (xh < 0.8) + prof1->GetYaxis()->SetTitleOffset(1.05); - TLegend *mlg = new TLegend(xh1-0.3, 0.80, xh1, 0.90); + TLegend *mlg = new TLegend(xh1 - 0.3, 0.80, xh1, 0.90); char title[100]; - sprintf (title, "End crystal to Layer 0"); - mlg->SetBorderSize(1); mlg->SetFillColor(10); mlg->SetMargin(0.30); - mlg->SetTextSize(0.04); mlg->AddEntry(prof1, title, "lf"); - - cname = "c_EC2" + plot + ztit + tag; + sprintf(title, "End crystal to Layer 0"); + mlg->SetBorderSize(1); + mlg->SetFillColor(10); + mlg->SetMargin(0.30); + mlg->SetTextSize(0.04); + mlg->AddEntry(prof1, title, "lf"); + + cname = "c_EC2" + plot + ztit + tag; new TCanvas(cname, cname, 700, 400); prof1->Draw(); - if (drawLeg > 0) mlg->Draw("sames"); + if (drawLeg > 0) + mlg->Draw("sames"); } -void etaPhiPlotHC(TString fileName, TString plot, int drawLeg, bool ifEta, - double maxEta, std::string tag) { - - TFile* hcalFile = new TFile(fileName); +void etaPhiPlotHC(TString fileName, TString plot, int drawLeg, bool ifEta, double maxEta, std::string tag) { + TFile *hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); setStyle(); @@ -403,33 +461,44 @@ void etaPhiPlotHC(TString fileName, TString plot, int drawLeg, bool ifEta, double xh = 0.90; if (plot.CompareTo("RadLen") == 0) { ytit = TString("Material Budget (X_{0})"); - ymin = 0; ymax = 200; istart = 100; + ymin = 0; + ymax = 200; + istart = 100; } else if (plot.CompareTo("StepLen") == 0) { ytit = TString("Material Budget (Step Length)"); - ymin = 0; ymax = 15000; istart = 300; xh = 0.35; + ymin = 0; + ymax = 15000; + istart = 300; + xh = 0.35; } else { ytit = TString("Material Budget (#lambda)"); - ymin = 0; ymax = 20; istart = 200; + ymin = 0; + ymax = 20; + istart = 200; } if (!ifEta) { istart += 400; - xtit = TString("#phi"); - ztit = TString("phi"); + xtit = TString("#phi"); + ztit = TString("phi"); } - - TLegend *leg = new TLegend(xh-0.25, 0.70, xh, 0.90); - leg->SetBorderSize(1); leg->SetFillColor(10); leg->SetMargin(0.30); + + TLegend *leg = new TLegend(xh - 0.25, 0.70, xh, 0.90); + leg->SetBorderSize(1); + leg->SetFillColor(10); + leg->SetMargin(0.30); leg->SetTextSize(0.04); - int nplots=0; + int nplots = 0; TProfile *prof[nlaymax]; - for (int ii=19; ii>=0; ii--) { + for (int ii = 19; ii >= 0; ii--) { char hname[10], title[50]; - int idpl = istart+ihid[ii]; + int idpl = istart + ihid[ii]; sprintf(hname, "%i", idpl); int icol = colorLayer[0]; - if (ii >= 18) icol = colorLayer[2]; - else if (ii > 0) icol = colorLayer[1]; - gDirectory->GetObject(hname,prof[nplots]); + if (ii >= 18) + icol = colorLayer[2]; + else if (ii > 0) + icol = colorLayer[1]; + gDirectory->GetObject(hname, prof[nplots]); prof[nplots]->GetXaxis()->SetTitle(xtit); prof[nplots]->GetYaxis()->SetTitle(ytit); prof[nplots]->GetXaxis()->SetLabelSize(0.05); @@ -440,10 +509,10 @@ void etaPhiPlotHC(TString fileName, TString plot, int drawLeg, bool ifEta, prof[nplots]->GetYaxis()->SetRangeUser(ymin, ymax); prof[nplots]->SetLineColor(icol); prof[nplots]->SetFillColor(icol); - if (xh < 0.8) + if (xh < 0.8) prof[nplots]->GetYaxis()->SetTitleOffset(1.05); - if (ifEta && maxEta > 0) - prof[nplots]->GetXaxis()->SetRangeUser(0.0,maxEta); + if (ifEta && maxEta > 0) + prof[nplots]->GetXaxis()->SetRangeUser(0.0, maxEta); if (ii == 19) { sprintf(title, "End of HF"); leg->AddEntry(prof[nplots], title, "lf"); @@ -461,42 +530,47 @@ void etaPhiPlotHC(TString fileName, TString plot, int drawLeg, bool ifEta, new TCanvas(cname, cname, 700, 400); prof[0]->Draw("h"); - for(int i=1; iDraw("h sames"); - if (drawLeg > 0) leg->Draw("sames"); + if (drawLeg > 0) + leg->Draw("sames"); } - -void etaPhi2DPlot(TString fileName, TString plot, int ifirst, int ilast, - int drawLeg, std::string tag) { - TFile* hcalFile = new TFile(fileName); +void etaPhi2DPlot(TString fileName, TString plot, int ifirst, int ilast, int drawLeg, std::string tag) { + TFile *hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); setStyle(); TString xtit = TString("#eta"); TString ytit = TString("#phi"); - TString xytit= TString("etaphi"); + TString xytit = TString("etaphi"); TString ztit = TString("HCal Material Budget (#lambda)"); int ymin = 0, ymax = 20, istart = 1000; - double xh=0.95, yh=0.95; + double xh = 0.95, yh = 0.95; if (plot.CompareTo("RadLen") == 0) { ztit = TString("HCal Material Budget (X_{0})"); - ymin = 0; ymax = 200; istart = 900; + ymin = 0; + ymax = 200; + istart = 900; } else if (plot.CompareTo("StepLen") == 0) { ztit = TString("HCal Material Budget (Step Length)"); - ymin = 0; ymax = 15000; istart = 1100; + ymin = 0; + ymax = 15000; + istart = 1100; } - - TLegend *leg = new TLegend(xh-0.13, yh-0.30, xh, yh); - leg->SetBorderSize(1); leg->SetFillColor(10); leg->SetMargin(0.25); + + TLegend *leg = new TLegend(xh - 0.13, yh - 0.30, xh, yh); + leg->SetBorderSize(1); + leg->SetFillColor(10); + leg->SetMargin(0.25); leg->SetTextSize(0.018); - int nplots=0; + int nplots = 0; TProfile2D *prof[nlaymax]; - for (int ii=ilast; ii>=ifirst; ii--) { + for (int ii = ilast; ii >= ifirst; ii--) { char hname[10], title[50]; - sprintf(hname, "%i", istart+ii); - gDirectory->GetObject(hname,prof[nplots]); + sprintf(hname, "%i", istart + ii); + gDirectory->GetObject(hname, prof[nplots]); prof[nplots]->GetXaxis()->SetTitle(xtit); prof[nplots]->GetYaxis()->SetTitle(ytit); prof[nplots]->GetZaxis()->SetTitle(ztit); @@ -509,7 +583,7 @@ void etaPhi2DPlot(TString fileName, TString plot, int ifirst, int ilast, sprintf(title, "Layer %d", lay); } else if (lay == 0) { sprintf(title, "After Crystal"); - } else if (lay >= 20 ) { + } else if (lay >= 20) { sprintf(title, "After HF"); } else { sprintf(title, "Before Crystal"); @@ -520,109 +594,119 @@ void etaPhi2DPlot(TString fileName, TString plot, int ifirst, int ilast, TString cname = "c_" + plot + xytit + tag; TCanvas *cc1 = new TCanvas(cname, cname, 700, 600); - cc1->SetLeftMargin(0.15); cc1->SetRightMargin(0.05); + cc1->SetLeftMargin(0.15); + cc1->SetRightMargin(0.05); prof[0]->Draw("lego fb bb"); - for(int i=1; iDraw("lego fb bb sames"); - if (drawLeg > 0) leg->Draw("sames"); + if (drawLeg > 0) + leg->Draw("sames"); } -void etaPhi2DPlot(int nslice, int kslice, TString fileName, TString plot, - int ifirst, int ilast, int drawLeg, std::string tag) { - +void etaPhi2DPlot( + int nslice, int kslice, TString fileName, TString plot, int ifirst, int ilast, int drawLeg, std::string tag) { char hname[200], title[50]; - TFile* hcalFile = new TFile(fileName); + TFile *hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); setStyle(); TString xtit = TString("#eta"); TString ytit; int ymin, ymax, istart; - double xh=0.95, yh=0.90; + double xh = 0.95, yh = 0.90; char type[10]; if (plot.CompareTo("RadLen") == 0) { ytit = TString("HCal Material Budget (X_{0})"); - ymin = 0; ymax = 200; istart = 900; - sprintf (type, "Radlen"); + ymin = 0; + ymax = 200; + istart = 900; + sprintf(type, "Radlen"); } else if (plot.CompareTo("StepLen") == 0) { ytit = TString("HCal Material Budget (Step Length)"); - ymin = 0; ymax = 15000; istart = 1100; - sprintf (type, "Steplen"); + ymin = 0; + ymax = 15000; + istart = 1100; + sprintf(type, "Steplen"); } else { ytit = TString("HCal Material Budget (#lambda)"); - ymin = 0; ymax = 20; istart = 1000; - sprintf (type, "Intlen"); + ymin = 0; + ymax = 20; + istart = 1000; + sprintf(type, "Intlen"); } - int nplots=0; + int nplots = 0; TProfile2D *prof[nlaymax]; - for (int ii=ilast; ii>=ifirst; ii--) { - sprintf(hname, "%i", istart+ii); - gDirectory->GetObject(hname,prof[nplots]); + for (int ii = ilast; ii >= ifirst; ii--) { + sprintf(hname, "%i", istart + ii); + gDirectory->GetObject(hname, prof[nplots]); nplots++; } - - double xmin = prof[0]->GetXaxis()->GetXmin(); - double xmax = prof[0]->GetXaxis()->GetXmax(); - int nbinX = prof[0]->GetNbinsX(); + + double xmin = prof[0]->GetXaxis()->GetXmin(); + double xmax = prof[0]->GetXaxis()->GetXmax(); + int nbinX = prof[0]->GetNbinsX(); double ymin1 = prof[0]->GetYaxis()->GetXmin(); double ymax1 = prof[0]->GetYaxis()->GetXmax(); - int nbinY = prof[0]->GetNbinsY(); - int ngroup= nbinY/nslice; - double dy = (ymax1-ymin1)/nbinY; - cout << "X " << nbinX << "/" << xmin << "/" << xmax << " Slice " << nbinY << "/" << nslice << "/" << ngroup << " " << nplots << " " << dy << "\n"; + int nbinY = prof[0]->GetNbinsY(); + int ngroup = nbinY / nslice; + double dy = (ymax1 - ymin1) / nbinY; + cout << "X " << nbinX << "/" << xmin << "/" << xmax << " Slice " << nbinY << "/" << nslice << "/" << ngroup << " " + << nplots << " " << dy << "\n"; - istart= 0; + istart = 0; TLegend *leg[360]; - TH1D *hist[nlaymax][360]; - for (int is=0; isSetBorderSize(1); leg[is]->SetFillColor(10); - leg[is]->SetMargin(0.25); leg[is]->SetTextSize(0.023); - double y1 = (ymin1 + istart*dy)*180./3.1415926; - double y2 = y1 + ngroup*dy*180./3.1415926; + TH1D *hist[nlaymax][360]; + for (int is = 0; is < nslice; is++) { + leg[is] = new TLegend(xh - 0.13, yh - 0.43, xh, yh); + leg[is]->SetBorderSize(1); + leg[is]->SetFillColor(10); + leg[is]->SetMargin(0.25); + leg[is]->SetTextSize(0.023); + double y1 = (ymin1 + istart * dy) * 180. / 3.1415926; + double y2 = y1 + ngroup * dy * 180. / 3.1415926; if (y1 < 0) { y1 += 360.; y2 += 360.; } - sprintf (title, "#phi = %6.1f :%6.1f", y1, y2); + sprintf(title, "#phi = %6.1f :%6.1f", y1, y2); leg[is]->SetHeader(title); - for (int ii=0; iiGetBinContent(ibx+1, ibin+1); - // cout << "/" << ibin << " " << cont; - contb += cont; - } - contb /= ngroup; - // cout << " " << contb; - hist[ii][is]->SetBinContent(ibx+1, contb); + for (int ibx = 0; ibx < nbinX; ibx++) { + double contb = 0; + // cout << " / " << ibx; + for (int iby = 0; iby < ngroup; iby++) { + int ibin = iby + istart; + double cont = prof[ii]->GetBinContent(ibx + 1, ibin + 1); + // cout << "/" << ibin << " " << cont; + contb += cont; + } + contb /= ngroup; + // cout << " " << contb; + hist[ii][is]->SetBinContent(ibx + 1, contb); } // cout << "\n"; hist[ii][is]->GetXaxis()->SetTitle(xtit); hist[ii][is]->GetYaxis()->SetTitle(ytit); hist[ii][is]->GetYaxis()->SetRangeUser(ymin, ymax); - hist[ii][is]->SetLineColor(colorLayer[ilast-ii]); - hist[ii][is]->SetFillColor(colorLayer[ilast-ii]); + hist[ii][is]->SetLineColor(colorLayer[ilast - ii]); + hist[ii][is]->SetFillColor(colorLayer[ilast - ii]); hist[ii][is]->GetYaxis()->SetTitleOffset(0.8); int lay = ilast - ii - 1; if (lay > 0 && lay < 20) { - sprintf(title, "Layer %d", lay); + sprintf(title, "Layer %d", lay); } else if (lay == 0) { - sprintf(title, "After Crystal"); - } else if (lay >= 20 ) { - sprintf(title, "After HF"); + sprintf(title, "After Crystal"); + } else if (lay >= 20) { + sprintf(title, "After HF"); } else { - sprintf(title, "Before Crystal"); + sprintf(title, "Before Crystal"); } leg[is]->AddEntry(hist[ii][is], title, "lf"); } @@ -631,53 +715,57 @@ void etaPhi2DPlot(int nslice, int kslice, TString fileName, TString plot, cout << "All histograms created now plot\n"; TCanvas *cc1[360]; - int ismin=0, ismax=nslice; - if (kslice >=0 && kslice <= nslice) { + int ismin = 0, ismax = nslice; + if (kslice >= 0 && kslice <= nslice) { ismin = kslice; - ismax = ismin+1; + ismax = ismin + 1; } - for (int is=ismin; isSetLeftMargin(0.15); cc1[is]->SetRightMargin(0.05); + cc1[is]->SetLeftMargin(0.15); + cc1[is]->SetRightMargin(0.05); hist[0][is]->Draw(); - for(int i=1; iDraw("sames"); - if (drawLeg > 0) leg[is]->Draw("sames"); + if (drawLeg > 0) + leg[is]->Draw("sames"); } } -void printTable (TString fileName, TString outputFileName, - TString inputFileName) { - - double radl[nlaymax][nbinmax], intl[nlaymax][nbinmax]; +void printTable(TString fileName, TString outputFileName, TString inputFileName) { + double radl[nlaymax][nbinmax], intl[nlaymax][nbinmax]; bool compare = false; if (inputFileName != "None") { ifstream inp(inputFileName, ios::in); cout << "Opens " << inputFileName << "\n"; if (inp) { TString line; - int tower; - double eta; - for (int i = 0; i < 23; i++) - inp >> line; - for (int itow=0; itow> tower >> eta; - int laymax=19; - if (itow > 29) laymax = 2; - else if (itow > 3) laymax = 18; - for (int ilay=0; ilay> intl[ilay][tower]; + int tower; + double eta; + for (int i = 0; i < 23; i++) + inp >> line; + for (int itow = 0; itow < nbinmax; itow++) { + inp >> tower >> eta; + int laymax = 19; + if (itow > 29) + laymax = 2; + else if (itow > 3) + laymax = 18; + for (int ilay = 0; ilay < laymax; ilay++) + inp >> intl[ilay][tower]; } - for (int i = 0; i < 23; i++) - inp >> line; - for (int itow=0; itow> tower >> eta; - int laymax=19; - if (itow > 29) laymax = 2; - else if (itow > 3) laymax = 18; - for (int ilay=0; ilay> radl[ilay][tower]; + for (int i = 0; i < 23; i++) + inp >> line; + for (int itow = 0; itow < nbinmax; itow++) { + inp >> tower >> eta; + int laymax = 19; + if (itow > 29) + laymax = 2; + else if (itow > 3) + laymax = 18; + for (int ilay = 0; ilay < laymax; ilay++) + inp >> radl[ilay][tower]; } compare = true; inp.close(); @@ -686,93 +774,99 @@ void printTable (TString fileName, TString outputFileName, std::ofstream os; os.open(outputFileName); - int nbadI=0; - getDiff (fileName, "IntLen"); - os << "Interaction Length\n" << "==================\n" + int nbadI = 0; + getDiff(fileName, "IntLen"); + os << "Interaction Length\n" + << "==================\n" << "Eta Tower/Layer 0 1 2 3 4 5 " << " 6 7 8 9 10 11 12 13 " << " 14 15 16 17\n"; - for (int itow=0; itow 29) {laymax = 2; ioff=0;} - else if (itow > 3) {laymax = 18; ioff=1;} - for (int ilay=0; ilay 29) { + laymax = 2; + ioff = 0; + } else if (itow > 3) { + laymax = 18; + ioff = 1; + } + for (int ilay = 0; ilay < laymax; ilay++) { + os << setw(8) << setprecision(4) << diff[ilay + ioff][itow]; if (compare) { - double num = (diff[ilay+ioff][itow] - intl[ilay][itow]); - double den = (diff[ilay+ioff][itow] + intl[ilay][itow]); - double dd = (den == 0.? 0. : 2.0*num/den); - if (dd > 0.01) { - nbadI++; - cout << "Lambda::Tower " << setw(3) << itow << " Layer " << setw(3) - << ilay << " Old" << setw(8) << setprecision(4) - << intl[ilay][itow] << " New" << setw(8) << setprecision(4) - << diff[ilay+ioff][itow] << " Diff"<< setw(8) << setprecision(4) - << dd << "\n"; - } + double num = (diff[ilay + ioff][itow] - intl[ilay][itow]); + double den = (diff[ilay + ioff][itow] + intl[ilay][itow]); + double dd = (den == 0. ? 0. : 2.0 * num / den); + if (dd > 0.01) { + nbadI++; + cout << "Lambda::Tower " << setw(3) << itow << " Layer " << setw(3) << ilay << " Old" << setw(8) + << setprecision(4) << intl[ilay][itow] << " New" << setw(8) << setprecision(4) << diff[ilay + ioff][itow] + << " Diff" << setw(8) << setprecision(4) << dd << "\n"; + } } } os << "\n"; } int nbadR = 0; - getDiff (fileName, "RadLen"); - os << "\n\nRadiation Length\n" << "================\n" + getDiff(fileName, "RadLen"); + os << "\n\nRadiation Length\n" + << "================\n" << "Eta Tower/Layer 0 1 2 3 4 5 " << " 6 7 8 9 10 11 12 13 " << " 14 15 16 17\n"; - for (int itow=0; itow 29) {laymax = 2; ioff=0;} - else if (itow > 3) {laymax = 18; ioff=1;} - for (int ilay=0; ilay 29) { + laymax = 2; + ioff = 0; + } else if (itow > 3) { + laymax = 18; + ioff = 1; + } + for (int ilay = 0; ilay < laymax; ilay++) { + os << setw(8) << setprecision(4) << diff[ilay + ioff][itow]; if (compare) { - double num = (diff[ilay+ioff][itow] - radl[ilay][itow]); - double den = (diff[ilay+ioff][itow] + radl[ilay][itow]); - double dd = (den == 0.? 0. : 2.0*num/den); - if (dd > 0.01) { - nbadR++; - cout << "X0::Tower " << setw(3) << itow << " Layer " << setw(3) - << ilay << " Old" << setw(8) << setprecision(4) - << radl[ilay][itow] << " New" << setw(8) << setprecision(4) - << diff[ilay+ioff][itow] << " Diff"<< setw(8) << setprecision(4) - << dd << "\n"; - } + double num = (diff[ilay + ioff][itow] - radl[ilay][itow]); + double den = (diff[ilay + ioff][itow] + radl[ilay][itow]); + double dd = (den == 0. ? 0. : 2.0 * num / den); + if (dd > 0.01) { + nbadR++; + cout << "X0::Tower " << setw(3) << itow << " Layer " << setw(3) << ilay << " Old" << setw(8) + << setprecision(4) << radl[ilay][itow] << " New" << setw(8) << setprecision(4) << diff[ilay + ioff][itow] + << " Diff" << setw(8) << setprecision(4) << dd << "\n"; + } } } os << "\n"; } os.close(); - cout << "Comparison Results " << nbadI << " discrepancies for Lambda and " - << nbadR << " discrepancies for X0\n"; + cout << "Comparison Results " << nbadI << " discrepancies for Lambda and " << nbadR << " discrepancies for X0\n"; } -void plotDiff (TString fileName, TString plot, std::string tag) { - +void plotDiff(TString fileName, TString plot, std::string tag) { setStyle(); - gStyle->SetTitleOffset(1.0,"Y"); - getDiff (fileName, plot); + gStyle->SetTitleOffset(1.0, "Y"); + getDiff(fileName, plot); TString xtit = TString("Layer Number"); TString ytit = TString("HCal Material Budget (#lambda)"); - if (plot.CompareTo("RadLen") == 0) + if (plot.CompareTo("RadLen") == 0) ytit = TString("HCal Material Budget (X_{0})"); TMultiGraph *mg = new TMultiGraph(); - TLegend *leg_mg = new TLegend(.5,.5,.75,.80); - leg_mg->SetFillColor(10); leg_mg->SetBorderSize(1); - leg_mg->SetMargin(0.25); leg_mg->SetTextSize(0.04); + TLegend *leg_mg = new TLegend(.5, .5, .75, .80); + leg_mg->SetFillColor(10); + leg_mg->SetBorderSize(1); + leg_mg->SetMargin(0.25); + leg_mg->SetTextSize(0.04); leg_mg->SetHeader(ytit); - double diff_lay[19], idx[19]; - for (int ilay=1; ilay<20; ilay++) { - diff_lay[ilay-1] = diff[ilay][0]; - idx[ilay-1] = ilay-1; + double diff_lay[19], idx[19]; + for (int ilay = 1; ilay < 20; ilay++) { + diff_lay[ilay - 1] = diff[ilay][0]; + idx[ilay - 1] = ilay - 1; } TGraph *gr_eta1 = new TGraph(19, idx, diff_lay); gr_eta1->SetMarkerStyle(20); @@ -780,22 +874,22 @@ void plotDiff (TString fileName, TString plot, std::string tag) { gr_eta1->SetLineColor(2); gr_eta1->GetXaxis()->SetTitle(xtit); gr_eta1->GetYaxis()->SetTitle(ytit); - mg->Add(gr_eta1, "pc"); + mg->Add(gr_eta1, "pc"); leg_mg->AddEntry(gr_eta1, "HB #eta = 1"); - for (int ilay=1; ilay<20; ilay++) - diff_lay[ilay-1] = diff[ilay][6]; + for (int ilay = 1; ilay < 20; ilay++) + diff_lay[ilay - 1] = diff[ilay][6]; TGraph *gr_eta7 = new TGraph(18, idx, diff_lay); gr_eta7->SetMarkerStyle(22); gr_eta7->SetMarkerColor(4); gr_eta7->SetLineColor(4); - mg->Add(gr_eta7, "pc"); + mg->Add(gr_eta7, "pc"); gr_eta7->GetXaxis()->SetTitle(xtit); gr_eta7->GetYaxis()->SetTitle(ytit); leg_mg->AddEntry(gr_eta7, "HB #eta = 7"); - - for (int ilay=1; ilay<20; ilay++) - diff_lay[ilay-1] = diff[ilay][12]; + + for (int ilay = 1; ilay < 20; ilay++) + diff_lay[ilay - 1] = diff[ilay][12]; TGraph *gr_eta13 = new TGraph(18, idx, diff_lay); gr_eta13->SetMarkerStyle(29); gr_eta13->SetMarkerColor(kGreen); @@ -803,10 +897,10 @@ void plotDiff (TString fileName, TString plot, std::string tag) { gr_eta13->GetXaxis()->SetTitle(xtit); gr_eta13->GetYaxis()->SetTitle(ytit); mg->Add(gr_eta13, "pc"); - leg_mg->AddEntry(gr_eta13,"HB #eta = 13"); + leg_mg->AddEntry(gr_eta13, "HB #eta = 13"); - for (int ilay=1; ilay<20; ilay++) - diff_lay[ilay-1] = diff[ilay][19]; + for (int ilay = 1; ilay < 20; ilay++) + diff_lay[ilay - 1] = diff[ilay][19]; TGraph *gr_eta19 = new TGraph(18, idx, diff_lay); gr_eta19->SetMarkerStyle(24); gr_eta19->SetMarkerColor(kCyan); @@ -814,10 +908,10 @@ void plotDiff (TString fileName, TString plot, std::string tag) { gr_eta19->GetXaxis()->SetTitle(xtit); gr_eta19->GetYaxis()->SetTitle(ytit); mg->Add(gr_eta19, "pc"); - leg_mg->AddEntry(gr_eta19,"HE #eta = 20"); + leg_mg->AddEntry(gr_eta19, "HE #eta = 20"); - for(int ilay=1; ilay<20; ilay++) - diff_lay[ilay-1] = diff[ilay][25]; + for (int ilay = 1; ilay < 20; ilay++) + diff_lay[ilay - 1] = diff[ilay][25]; TGraph *gr_eta25 = new TGraph(18, idx, diff_lay); gr_eta25->SetMarkerStyle(26); gr_eta25->SetMarkerColor(kCyan); @@ -825,7 +919,7 @@ void plotDiff (TString fileName, TString plot, std::string tag) { gr_eta25->GetXaxis()->SetTitle(xtit); gr_eta25->GetYaxis()->SetTitle(ytit); mg->Add(gr_eta25, "pc"); - leg_mg->AddEntry(gr_eta25,"HE #eta = 26"); + leg_mg->AddEntry(gr_eta25, "HE #eta = 26"); TString cname = "c_diff_" + plot + tag; new TCanvas(cname, cname, 700, 400); @@ -835,60 +929,69 @@ void plotDiff (TString fileName, TString plot, std::string tag) { leg_mg->Draw("same"); } -void getDiff (TString fileName, TString plot) { - - TFile* hcalFile = new TFile(fileName); +void getDiff(TString fileName, TString plot) { + TFile *hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); - int istart = 200; + int istart = 200; if (plot.CompareTo("RadLen") == 0) { istart = 100; } else if (plot.CompareTo("StepLen") == 0) { - istart = 300; + istart = 300; } - for (int ilay=0; ilay<22; ilay++) { + for (int ilay = 0; ilay < 22; ilay++) { char hname[10]; - sprintf(hname, "%i", istart+ilay+1); + sprintf(hname, "%i", istart + ilay + 1); TProfile *prof; - gDirectory->GetObject(hname,prof); - int nbins = prof->GetNbinsX(); - for (int itow=0; itowGetObject(hname, prof); + int nbins = prof->GetNbinsX(); + for (int itow = 0; itow < nbinmax; itow++) { double ent = 0, value = 0; - for (int ii=0; iiGetBinLowEdge(ii+1); - double xu = prof->GetBinWidth(ii+1); - if (xl >= 0) { xu += xl;} - else { double tmp = xu; xu =-xl; xl = xu-tmp;} - double cont = (prof->GetBinContent(ii+1)); - double dx = 1; - if (cont > 0) { - if (xl >= towLow[itow] && xu <= towHigh[itow]) { - ent += dx; value += cont; - } else if (xl < towLow[itow] && xu > towLow[itow]) { - dx = (xu-towLow[itow])/(xu-xl); - ent += dx; value += dx*cont; - } else if (xu > towHigh[itow] && xl < towHigh[itow]) { - dx = (towHigh[itow]-xl)/(xu-xl); - ent += dx; value += dx*cont; - } - } + for (int ii = 0; ii < nbins; ii++) { + double xl = prof->GetBinLowEdge(ii + 1); + double xu = prof->GetBinWidth(ii + 1); + if (xl >= 0) { + xu += xl; + } else { + double tmp = xu; + xu = -xl; + xl = xu - tmp; + } + double cont = (prof->GetBinContent(ii + 1)); + double dx = 1; + if (cont > 0) { + if (xl >= towLow[itow] && xu <= towHigh[itow]) { + ent += dx; + value += cont; + } else if (xl < towLow[itow] && xu > towLow[itow]) { + dx = (xu - towLow[itow]) / (xu - xl); + ent += dx; + value += dx * cont; + } else if (xu > towHigh[itow] && xl < towHigh[itow]) { + dx = (towHigh[itow] - xl) / (xu - xl); + ent += dx; + value += dx * cont; + } + } } - if (ent > 0) mean[ilay][itow] = value/ent; - else mean[ilay][itow] = 0.; + if (ent > 0) + mean[ilay][itow] = value / ent; + else + mean[ilay][itow] = 0.; } } - for (int itow=30; itow0; ilay--) { - if (mean[ilay-1][itow] <= 0) { - mean[ilay-1][itow] = mean[ilay][itow]; + for (int itow = 0; itow < 30; itow++) { + for (int ilay = 20; ilay > 0; ilay--) { + if (mean[ilay - 1][itow] <= 0) { + mean[ilay - 1][itow] = mean[ilay][itow]; } } } - for (int itow=0; itow 4 && itow < 26) mean[19][itow] = 0; + for (int itow = 0; itow < nbinmax; itow++) { + if (itow > 4 && itow < 26) + mean[19][itow] = 0; diff[0][itow] = mean[0][itow]; } - for (int itow=15; itow<17; itow++) { - for (int ilay=1; ilay<22; ilay++) { - if (mean[ilay][itow] > mean[ilay+1][itow]) { - for (int jlay=ilay+1; jlay<22; jlay++) - mean[jlay][itow] = 0; - break; + for (int itow = 15; itow < 17; itow++) { + for (int ilay = 1; ilay < 22; ilay++) { + if (mean[ilay][itow] > mean[ilay + 1][itow]) { + for (int jlay = ilay + 1; jlay < 22; jlay++) + mean[jlay][itow] = 0; + break; } } } - for (int ilay=1; ilay<22; ilay++) { - for (int itow=0; itowSetTitleOffset(1.2,"Y"); +void plotHE(int flag, int logy, std::string tag, int save) { + double angle[31] = {-2.5, -2.25, -2.00, -1.75, -1.50, -1.25, -1.00, -0.75, -0.50, -0.25, -0.20, + -0.15, -0.10, -0.05, -0.025, 0, 0.025, 0.05, 0.10, 0.15, 0.20, 0.25, + 0.50, 0.75, 1.00, 1.25, 1.50, 1.75, 2.00, 2.25, 2.50}; + double lAir[31] = {11440, 11440, 11440, 11440, 11440, 11440, 11440, 11440, 11440, 11440, 11440, + 11440, 11440, 11440, 11440, 11569, 11569, 11440, 11440, 11440, 11440, 11440, + 11440, 11440, 11440, 11440, 11440, 11440, 11440, 11440, 11440}; + double xAir[31] = {0.037941, 0.037941, 0.037941, 0.037941, 0.037941, 0.037941, 0.037941, 0.037941, + 0.037941, 0.037941, 0.037941, 0.037941, 0.037941, 0.037941, 0.037941, 0.038369, + 0.038369, 0.037941, 0.037941, 0.037941, 0.037941, 0.037941, 0.037941, 0.037941, + 0.037941, 0.037941, 0.037941, 0.037941, 0.037941, 0.037941, 0.037941}; + double iAir[31] = {0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481, + 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0164314, + 0.0164314, 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481, + 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481, 0.0162481}; + double lPol[31] = {60.8381, 60.8381, 60.8381, 60.8381, 60.8381, 60.8381, 60.8381, 60.8381, 60.8381, 60.8381, 60.8381, + 60.8381, 60.8381, 60.8381, 60.8381, 0.00000, 0.00000, 60.8381, 60.8381, 60.8381, 60.8381, 60.8381, + 60.8381, 60.8381, 60.8381, 60.8381, 60.8381, 60.8381, 60.8381, 60.8381, 60.8381}; + double xPol[31] = {0.129083, 0.129083, 0.129083, 0.129083, 0.129083, 0.129083, 0.129083, 0.129083, + 0.129083, 0.129083, 0.129083, 0.129083, 0.129083, 0.129083, 0.129083, 0.000000, + 0.000000, 0.129083, 0.129083, 0.129083, 0.129083, 0.129083, 0.129083, 0.129083, + 0.129083, 0.129083, 0.129083, 0.129083, 0.129083, 0.129083, 0.129083}; + double iPol[31] = {0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129, + 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0000000, + 0.0000000, 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129, + 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129, 0.0854129}; + double lScn[31] = {68.2125, 68.2125, 68.2125, 68.2125, 68.2125, 68.2125, 68.2125, 68.2125, 68.2125, 68.2125, 68.2125, + 68.2125, 68.2125, 68.2125, 68.2125, 0.00000, 0.00000, 68.2125, 68.2125, 68.2125, 68.2125, 68.2125, + 68.2125, 68.2125, 68.2125, 68.2125, 68.2125, 68.2125, 68.2125, 68.2125, 68.2125}; + double xScn[31] = {0.160352, 0.160352, 0.160352, 0.160352, 0.160352, 0.160352, 0.160352, 0.160352, + 0.160352, 0.160352, 0.160352, 0.160352, 0.160352, 0.160352, 0.160352, 0.000000, + 0.000000, 0.160352, 0.160352, 0.160352, 0.160352, 0.160352, 0.160352, 0.160352, + 0.160352, 0.160352, 0.160352, 0.160352, 0.160352, 0.160352, 0.160352}; + double iScn[31] = {0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967, + 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0000000, + 0.0000000, 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967, + 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967, 0.0973967}; + double lBra[31] = {1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, + 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, + 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41, 1444.41}; + double xBra[31] = {96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, + 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, + 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848, 96.7848}; + double iBra[31] = {8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, + 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, + 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637, 8.79637}; + std::string nameMat[4] = {"Air", "Polythene", "Scintillator", "Brass"}; + int colMat[4] = {1, 2, 6, 4}; + int symbMat[4] = {24, 29, 25, 27}; + + setStyle(); + gStyle->SetTitleOffset(1.2, "Y"); char name[30], title[60], gname[12]; TGraph *gr[4]; int kfirst = 3; - double ymi=0, ymx=100; + double ymi = 0, ymx = 100; if (flag < 0) { - sprintf (name, "Step Length"); - sprintf (title, "Step Length (mm)"); - sprintf (gname, "stepLength"); + sprintf(name, "Step Length"); + sprintf(title, "Step Length (mm)"); + sprintf(gname, "stepLength"); gr[0] = new TGraph(31, angle, lAir); gr[1] = new TGraph(31, angle, lPol); gr[2] = new TGraph(31, angle, lScn); @@ -1032,12 +1110,13 @@ void plotHE(int flag, int logy, std::string tag, int save) { if (logy == 0) { ymx = 12000; } else { - ymi = 10.0; ymx = 20000; + ymi = 10.0; + ymx = 20000; } - } else if (flag>0) { - sprintf (name, "# Interaction Length"); - sprintf (title, "Number of Interaction Length"); - sprintf (gname, "intLength"); + } else if (flag > 0) { + sprintf(name, "# Interaction Length"); + sprintf(title, "Number of Interaction Length"); + sprintf(gname, "intLength"); gr[0] = new TGraph(31, angle, iAir); gr[1] = new TGraph(31, angle, iPol); gr[2] = new TGraph(31, angle, iScn); @@ -1045,12 +1124,13 @@ void plotHE(int flag, int logy, std::string tag, int save) { if (logy == 0) { ymx = 10; } else { - ymi = 0.01; ymx = 20; + ymi = 0.01; + ymx = 20; } } else { - sprintf (name, "# Radiation Length"); - sprintf (title, "Number of Radiation Length"); - sprintf (gname, "radLength"); + sprintf(name, "# Radiation Length"); + sprintf(title, "Number of Radiation Length"); + sprintf(gname, "radLength"); gr[0] = new TGraph(31, angle, xAir); gr[1] = new TGraph(31, angle, xPol); gr[2] = new TGraph(31, angle, xScn); @@ -1058,112 +1138,137 @@ void plotHE(int flag, int logy, std::string tag, int save) { if (logy == 0) { ymx = 100; } else { - ymi = 0.01; ymx = 200; + ymi = 0.01; + ymx = 200; } } gr[kfirst]->GetXaxis()->SetTitle("#phi ( ^{o})"); gr[kfirst]->GetYaxis()->SetTitle(title); gr[kfirst]->SetTitle(""); - for (int i=0; i<4; i++) { + for (int i = 0; i < 4; i++) { int icol = colMat[i]; int type = symbMat[i]; - gr[i]->SetMarkerSize(1.2); gr[i]->SetLineColor(icol); - gr[i]->SetLineStyle(i+1); gr[i]->SetLineWidth(2); - gr[i]->SetMarkerColor(icol); gr[i]->SetMarkerStyle(type); - gr[i]->GetYaxis()->SetRangeUser(ymi,ymx); - gr[i]->GetXaxis()->SetRangeUser(-3.0,+3.0); + gr[i]->SetMarkerSize(1.2); + gr[i]->SetLineColor(icol); + gr[i]->SetLineStyle(i + 1); + gr[i]->SetLineWidth(2); + gr[i]->SetMarkerColor(icol); + gr[i]->SetMarkerStyle(type); + gr[i]->GetYaxis()->SetRangeUser(ymi, ymx); + gr[i]->GetXaxis()->SetRangeUser(-3.0, +3.0); } TString cname = "c1" + tag; TCanvas *c1 = new TCanvas(cname, name, 800, 500); - if (logy != 0) gPad->SetLogy(1); + if (logy != 0) + gPad->SetLogy(1); gr[kfirst]->Draw("alp"); - for (int i=0; i<4; i++) { - if (i != kfirst) gr[i]->Draw("lp"); + for (int i = 0; i < 4; i++) { + if (i != kfirst) + gr[i]->Draw("lp"); } double ylow = 0.4; char list[20]; - TLegend *leg1 = new TLegend(0.60,ylow,0.90,ylow+0.2); - for (int i=0; i<4; i++) { - sprintf (list, "%s", nameMat[i].c_str()); - leg1->AddEntry(gr[i],list,"LP"); + TLegend *leg1 = new TLegend(0.60, ylow, 0.90, ylow + 0.2); + for (int i = 0; i < 4; i++) { + sprintf(list, "%s", nameMat[i].c_str()); + leg1->AddEntry(gr[i], list, "LP"); } - leg1->SetHeader(name); leg1->SetFillColor(0); + leg1->SetHeader(name); + leg1->SetFillColor(0); leg1->SetTextSize(0.04); leg1->Draw(); if (save != 0) { char fname[20]; - if (save > 0) sprintf (fname, "%s.eps", gname); - else sprintf (fname, "%s.gif", gname); + if (save > 0) + sprintf(fname, "%s.eps", gname); + else + sprintf(fname, "%s.gif", gname); c1->SaveAs(fname); } - } -void etaPhiCastorPlot(TString fileName, TString plot, TString type, - bool etaPlus, int drawLeg, bool ifEta, std::string tag, - bool debug) { - - TFile* hcalFile = new TFile(fileName); +void etaPhiCastorPlot( + TString fileName, TString plot, TString type, bool etaPlus, int drawLeg, bool ifEta, std::string tag, bool debug) { + TFile *hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); setStyle(); - if (debug) std::cout << fileName << " is opened at " << hcalFile << std::endl; + if (debug) + std::cout << fileName << " is opened at " << hcalFile << std::endl; TString xtit = TString("#eta"); TString ztit = TString("eta"); char ytit[80], ytpart[10]; - int ymin = 0, ymax = 20, istart = 200, ifirst=0; + int ymin = 0, ymax = 20, istart = 200, ifirst = 0; double xh = 0.90; - if (!etaPlus) ifirst = 10; + if (!etaPlus) + ifirst = 10; if (type.CompareTo("EC") == 0) { - sprintf (ytpart, "(EC)"); ifirst += 2; + sprintf(ytpart, "(EC)"); + ifirst += 2; } else if (type.CompareTo("HC") == 0) { - sprintf (ytpart, "(HC)"); ifirst += 4; + sprintf(ytpart, "(HC)"); + ifirst += 4; } else if (type.CompareTo("ED") == 0) { - sprintf (ytpart, "(Dead EC)"); ifirst += 6; + sprintf(ytpart, "(Dead EC)"); + ifirst += 6; } else if (type.CompareTo("HD") == 0) { - sprintf (ytpart, "(Dead HC)"); ifirst += 8; + sprintf(ytpart, "(Dead HC)"); + ifirst += 8; } else { - sprintf (ytpart, "(All)"); + sprintf(ytpart, "(All)"); } - if (debug) std::cout << type << " Gives " << ifirst << " Title " << ytpart << std::endl; + if (debug) + std::cout << type << " Gives " << ifirst << " Title " << ytpart << std::endl; if (plot.CompareTo("RadLen") == 0) { sprintf(ytit, "Castor %s Material Budget (X_{0})", ytpart); - ymin = 0; ymax = 200; istart = 100; + ymin = 0; + ymax = 200; + istart = 100; } else if (plot.CompareTo("StepLen") == 0) { sprintf(ytit, "Castor %s Material Budget (Step Length)", ytpart); - ymin = 0; ymax = 15000; istart = 300; xh = 0.70; + ymin = 0; + ymax = 15000; + istart = 300; + xh = 0.70; } else { sprintf(ytit, "Castor %s Material Budget (#lambda)", ytpart); - ymin = 0; ymax = 20; istart = 200; + ymin = 0; + ymax = 20; + istart = 200; } if (!ifEta) { istart += 400; - xtit = TString("#phi"); - ztit = TString("phi"); + xtit = TString("#phi"); + ztit = TString("phi"); } - if (debug) std::cout << "Title (x) " << xtit << " (y) " << ytit << " First " << ifirst << ":" << istart << std::endl; - - TLegend *leg = new TLegend(xh-0.13, 0.80, xh, 0.90); - leg->SetBorderSize(1); leg->SetFillColor(10); leg->SetMargin(0.25); + if (debug) + std::cout << "Title (x) " << xtit << " (y) " << ytit << " First " << ifirst << ":" << istart << std::endl; + + TLegend *leg = new TLegend(xh - 0.13, 0.80, xh, 0.90); + leg->SetBorderSize(1); + leg->SetFillColor(10); + leg->SetMargin(0.25); leg->SetTextSize(0.025); - int nplots=0; + int nplots = 0; TProfile *prof[2]; - for (int ii=ifirst+1; ii>=ifirst; ii--) { + for (int ii = ifirst + 1; ii >= ifirst; ii--) { char hname[10], title[50]; - sprintf(hname, "%i", istart+ii); - if (debug) std::cout << "[" << nplots << "] " << ii << " " << hname << "\n"; - gDirectory->GetObject(hname,prof[nplots]); - if (debug) std::cout << "Histogram[" << nplots << "] : " << hname << " at " << prof[nplots] << std::endl; + sprintf(hname, "%i", istart + ii); + if (debug) + std::cout << "[" << nplots << "] " << ii << " " << hname << "\n"; + gDirectory->GetObject(hname, prof[nplots]); + if (debug) + std::cout << "Histogram[" << nplots << "] : " << hname << " at " << prof[nplots] << std::endl; prof[nplots]->GetXaxis()->SetTitle(xtit); prof[nplots]->GetYaxis()->SetTitle(ytit); prof[nplots]->GetYaxis()->SetRangeUser(ymin, ymax); prof[nplots]->SetLineColor(colorLayer[nplots]); prof[nplots]->SetFillColor(colorLayer[nplots]); - if (xh < 0.8) + if (xh < 0.8) prof[nplots]->GetYaxis()->SetTitleOffset(1.7); if (ii == ifirst) { sprintf(title, "Front"); @@ -1177,37 +1282,46 @@ void etaPhiCastorPlot(TString fileName, TString plot, TString type, TString cname = "c_" + plot + ztit + tag; TCanvas *cc1 = new TCanvas(cname, cname, 700, 600); if (xh < 0.8) { - cc1->SetLeftMargin(0.15); cc1->SetRightMargin(0.05); + cc1->SetLeftMargin(0.15); + cc1->SetRightMargin(0.05); } prof[0]->Draw("h"); - for(int i=1; iDraw("h sames"); - if (drawLeg > 0) leg->Draw("sames"); + if (drawLeg > 0) + leg->Draw("sames"); } -void efficiencyPlot(TString fileName, TString type, bool ifEtaPhi, - double maxEta, std::string tag, bool debug) { - - TFile* hcalFile = new TFile(fileName); +void efficiencyPlot(TString fileName, TString type, bool ifEtaPhi, double maxEta, std::string tag, bool debug) { + TFile *hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); setStyle(); - int id0=1300, idpl1=8, idpl2=0; + int id0 = 1300, idpl1 = 8, idpl2 = 0; char hname[100], title[100]; - if (type.CompareTo("HB") == 0) { - idpl1 = 1; idpl2 = 2; sprintf (title, "Efficiency for HB (%s)", tag.c_str()); + if (type.CompareTo("HB") == 0) { + idpl1 = 1; + idpl2 = 2; + sprintf(title, "Efficiency for HB (%s)", tag.c_str()); } else if (type.CompareTo("HE") == 0) { - idpl1 = 3; idpl2 = 4; sprintf (title, "Efficiency for HE (%s)", tag.c_str()); + idpl1 = 3; + idpl2 = 4; + sprintf(title, "Efficiency for HE (%s)", tag.c_str()); } else if (type.CompareTo("HO") == 0) { - idpl1 = 5; idpl2 = 6; sprintf (title, "Efficiency for HO (%s)", tag.c_str()); + idpl1 = 5; + idpl2 = 6; + sprintf(title, "Efficiency for HO (%s)", tag.c_str()); } else if (type.CompareTo("HF") == 0) { - idpl1 = 7; sprintf (title, "Efficiency for HF (%s)", tag.c_str()); + idpl1 = 7; + sprintf(title, "Efficiency for HF (%s)", tag.c_str()); } else { - sprintf (title, "Efficiency for HCAL (%s)", tag.c_str()); + sprintf(title, "Efficiency for HCAL (%s)", tag.c_str()); } TLegend *leg = new TLegend(0.70, 0.82, 0.90, 0.90); - leg->SetBorderSize(1); leg->SetFillColor(10); leg->SetMargin(0.25); + leg->SetBorderSize(1); + leg->SetFillColor(10); + leg->SetMargin(0.25); leg->SetTextSize(0.03); if (ifEtaPhi) { @@ -1215,52 +1329,65 @@ void efficiencyPlot(TString fileName, TString type, bool ifEtaPhi, TH2F *hist0, *hist1, *hist2; sprintf(hname, "%i", id0); gDirectory->GetObject(hname, hist0); - sprintf(hname, "%i", id0+idpl1); + sprintf(hname, "%i", id0 + idpl1); gDirectory->GetObject(hname, hist1); if (idpl2 > 0) { - sprintf(hname, "%i", id0+idpl2); + sprintf(hname, "%i", id0 + idpl2); gDirectory->GetObject(hname, hist2); } else { hist2 = 0; } - if (debug) std::cout << "Get Histos at " <GetXaxis()->GetXmin(); - double xmax = hist0->GetXaxis()->GetXmax(); - int nbinX = hist0->GetNbinsX(); - double ymin = hist0->GetYaxis()->GetXmin(); - double ymax = hist0->GetYaxis()->GetXmax(); - int nbinY = hist0->GetNbinsY(); - if (debug) std::cout <<"NbinX " <GetBinContent(ibx+1,iby+1); - double contD = hist0->GetBinContent(ibx+1,iby+1); - double cont = contN/std::max(contD,1.0); - hist->SetBinContent(ibx+1, iby+1, cont); - if (hist2) { - contN = hist2->GetBinContent(ibx+1,iby+1); - cont = contN/std::max(contD,1.0); - histe->SetBinContent(ibx+1, iby+1, cont); - } - } + double xmin = hist0->GetXaxis()->GetXmin(); + double xmax = hist0->GetXaxis()->GetXmax(); + int nbinX = hist0->GetNbinsX(); + double ymin = hist0->GetYaxis()->GetXmin(); + double ymax = hist0->GetYaxis()->GetXmax(); + int nbinY = hist0->GetNbinsY(); + if (debug) + std::cout << "NbinX " << nbinX << " range " << std::setprecision(5) << xmin << ":" << std::setprecision(5) + << xmax << " " + << "NbinY " << nbinY << " range " << std::setprecision(5) << ymin << ":" << std::setprecision(5) + << ymax << "\n"; + TH2D *hist = new TH2D("hist", title, nbinX, xmin, xmax, nbinY, ymin, ymax); + TH2D *histe = 0; + if (hist2) + histe = new TH2D("histe", title, nbinX, xmin, xmax, nbinY, ymin, ymax); + for (int ibx = 0; ibx < nbinX; ++ibx) { + for (int iby = 0; iby < nbinY; ++iby) { + double contN = hist1->GetBinContent(ibx + 1, iby + 1); + double contD = hist0->GetBinContent(ibx + 1, iby + 1); + double cont = contN / std::max(contD, 1.0); + hist->SetBinContent(ibx + 1, iby + 1, cont); + if (hist2) { + contN = hist2->GetBinContent(ibx + 1, iby + 1); + cont = contN / std::max(contD, 1.0); + histe->SetBinContent(ibx + 1, iby + 1, cont); + } + } } - hist->GetXaxis()->SetTitle("#eta");hist->GetYaxis()->SetTitle("#phi"); - hist->GetZaxis()->SetTitle(title);hist->GetZaxis()->SetTitleOffset(.8); + hist->GetXaxis()->SetTitle("#eta"); + hist->GetYaxis()->SetTitle("#phi"); + hist->GetZaxis()->SetTitle(title); + hist->GetZaxis()->SetTitleOffset(.8); new TCanvas(title, title, 700, 400); - hist->SetLineColor(2); hist->SetLineStyle(1); hist->SetLineWidth(1); - if (maxEta > 0) hist->GetXaxis()->SetRangeUser(-maxEta,maxEta); - hist->Draw("lego fb bb"); leg->AddEntry(hist, "At least 1 layer", "l"); + hist->SetLineColor(2); + hist->SetLineStyle(1); + hist->SetLineWidth(1); + if (maxEta > 0) + hist->GetXaxis()->SetRangeUser(-maxEta, maxEta); + hist->Draw("lego fb bb"); + leg->AddEntry(hist, "At least 1 layer", "l"); if (histe) { - histe->SetLineColor(4);histe->SetLineStyle(2);histe->SetLineWidth(1); - if (maxEta > 0) histe->GetXaxis()->SetRangeUser(-maxEta,maxEta); - histe->Draw("lego fb bb sames");leg->AddEntry(histe,"All layers","l"); + histe->SetLineColor(4); + histe->SetLineStyle(2); + histe->SetLineWidth(1); + if (maxEta > 0) + histe->GetXaxis()->SetRangeUser(-maxEta, maxEta); + histe->Draw("lego fb bb sames"); + leg->AddEntry(histe, "All layers", "l"); } leg->Draw("sames"); } @@ -1268,60 +1395,69 @@ void efficiencyPlot(TString fileName, TString type, bool ifEtaPhi, TH1F *hist0, *hist1, *hist2; sprintf(hname, "%i", id0); gDirectory->GetObject(hname, hist0); - sprintf(hname, "%i", id0+idpl1); + sprintf(hname, "%i", id0 + idpl1); gDirectory->GetObject(hname, hist1); if (idpl2 > 0) { - sprintf(hname, "%i", id0+idpl2); + sprintf(hname, "%i", id0 + idpl2); gDirectory->GetObject(hname, hist2); } else { hist2 = 0; } - if (debug) std::cout << "Get Histos at " <GetXaxis()->GetXmin(); - double xmax = hist0->GetXaxis()->GetXmax(); - int nbinX = hist0->GetNbinsX(); - if (debug) std::cout <<"Nbin " <GetXaxis()->GetXmin(); + double xmax = hist0->GetXaxis()->GetXmax(); + int nbinX = hist0->GetNbinsX(); + if (debug) + std::cout << "Nbin " << nbinX << " range " << std::setprecision(5) << xmin << ":" << std::setprecision(5) + << xmax << "\n"; TH1D *hist = new TH1D("hist", title, nbinX, xmin, xmax); - TH1D *histe= 0; - if (hist2) histe = new TH1D("histe", title, nbinX, xmin, xmax); - for (int ib=0; ibGetBinContent(ib+1); - double contD = hist0->GetBinContent(ib+1); - double cont = contN/std::max(contD,1.0); - hist->SetBinContent(ib+1, cont); - /* + TH1D *histe = 0; + if (hist2) + histe = new TH1D("histe", title, nbinX, xmin, xmax); + for (int ib = 0; ib < nbinX; ++ib) { + double contN = hist1->GetBinContent(ib + 1); + double contD = hist0->GetBinContent(ib + 1); + double cont = contN / std::max(contD, 1.0); + hist->SetBinContent(ib + 1, cont); + /* double eror = std::sqrt(contN)/std::max(contD,1.0); hist->SetBinError(ib+1, eror); */ - if (hist2) { - contN = hist2->GetBinContent(ib+1); - cont = contN/std::max(contD,1.0); - histe->SetBinContent(ib+1, cont); - } + if (hist2) { + contN = hist2->GetBinContent(ib + 1); + cont = contN / std::max(contD, 1.0); + histe->SetBinContent(ib + 1, cont); + } } hist->GetXaxis()->SetTitle("#eta"); hist->GetYaxis()->SetTitle(title); hist->GetYaxis()->SetTitleOffset(0.8); new TCanvas(title, title, 700, 400); - hist->SetLineColor(2); hist->SetLineStyle(1); hist->SetLineWidth(1); - if (maxEta > 0) hist->GetXaxis()->SetRangeUser(-maxEta,maxEta); - hist->Draw(); leg->AddEntry(hist, "At least 1 layer", "l"); + hist->SetLineColor(2); + hist->SetLineStyle(1); + hist->SetLineWidth(1); + if (maxEta > 0) + hist->GetXaxis()->SetRangeUser(-maxEta, maxEta); + hist->Draw(); + leg->AddEntry(hist, "At least 1 layer", "l"); if (histe) { - histe->SetLineColor(4); histe->SetLineStyle(2); histe->SetLineWidth(1); - if (maxEta > 0) histe->GetXaxis()->SetRangeUser(-maxEta,maxEta); - histe->Draw("sames"); leg->AddEntry(histe, "All layers", "l"); + histe->SetLineColor(4); + histe->SetLineStyle(2); + histe->SetLineWidth(1); + if (maxEta > 0) + histe->GetXaxis()->SetRangeUser(-maxEta, maxEta); + histe->Draw("sames"); + leg->AddEntry(histe, "All layers", "l"); } leg->Draw("sames"); } } -} - -void etaPhiFwdPlot(TString fileName, TString plot, int first, int last, - int drawLeg, std::string tag, bool debug) { +} - TFile* hcalFile = new TFile(fileName); +void etaPhiFwdPlot(TString fileName, TString plot, int first, int last, int drawLeg, std::string tag, bool debug) { + TFile *hcalFile = new TFile(fileName); hcalFile->cd("g4SimHits"); setStyle(); @@ -1329,53 +1465,68 @@ void etaPhiFwdPlot(TString fileName, TString plot, int first, int last, TString ztit = TString("eta"); TString ytit = "none"; int ymin = 0, ymax = 20, istart = 200; - double xh = 0.90, xl=0.1; + double xh = 0.90, xl = 0.1; if (plot.CompareTo("RadLen") == 0) { ytit = TString("Material Budget (X_{0})"); - ymin = 0; ymax = 400; istart = 100; + ymin = 0; + ymax = 400; + istart = 100; } else if (plot.CompareTo("StepLen") == 0) { ytit = TString("Material Budget (Step Length)"); - ymin = 0; ymax = 20000; istart = 300; xh = 0.70, xl=0.15; + ymin = 0; + ymax = 20000; + istart = 300; + xh = 0.70, xl = 0.15; } else { ytit = TString("Material Budget (#lambda)"); - ymin = 0; ymax = 30; istart = 200; + ymin = 0; + ymax = 30; + istart = 200; } int index[10] = {9, 0, 1, 2, 3, 8, 4, 7, 5, 6}; - std::string label[10] = {"Empty", "Beam Pipe", "Tracker", "EM Calorimeter", - "Hadron Calorimeter", "Muon System", - "Forward Hadron Calorimeter", "Shielding", "TOTEM", - "CASTOR"}; - - TLegend *leg = new TLegend(xl, 0.75, xl+0.3, 0.90); - leg->SetBorderSize(1); leg->SetFillColor(10); leg->SetMargin(0.25); + std::string label[10] = {"Empty", + "Beam Pipe", + "Tracker", + "EM Calorimeter", + "Hadron Calorimeter", + "Muon System", + "Forward Hadron Calorimeter", + "Shielding", + "TOTEM", + "CASTOR"}; + + TLegend *leg = new TLegend(xl, 0.75, xl + 0.3, 0.90); + leg->SetBorderSize(1); + leg->SetFillColor(10); + leg->SetMargin(0.25); leg->SetTextSize(0.018); - int nplots=0; + int nplots = 0; TProfile *prof[nlaymax]; - for (int ii=last; ii>=first; ii--) { + for (int ii = last; ii >= first; ii--) { char hname[10], title[50]; - sprintf(hname, "%i", istart+index[ii]); - gDirectory->GetObject(hname,prof[nplots]); + sprintf(hname, "%i", istart + index[ii]); + gDirectory->GetObject(hname, prof[nplots]); prof[nplots]->GetXaxis()->SetTitle(xtit); prof[nplots]->GetYaxis()->SetTitle(ytit); prof[nplots]->GetYaxis()->SetRangeUser(ymin, ymax); prof[nplots]->SetLineColor(colorLayer[ii]); prof[nplots]->SetFillColor(colorLayer[ii]); - if (xh < 0.8) + if (xh < 0.8) prof[nplots]->GetYaxis()->SetTitleOffset(1.7); sprintf(title, "%s", label[ii].c_str()); leg->AddEntry(prof[nplots], title, "lf"); if (debug) { - int nbinX = prof[nplots]->GetNbinsX(); - double xmin = prof[nplots]->GetXaxis()->GetXmin(); - double xmax = prof[nplots]->GetXaxis()->GetXmax(); - double dx = (xmax - xmin)/nbinX; + int nbinX = prof[nplots]->GetNbinsX(); + double xmin = prof[nplots]->GetXaxis()->GetXmin(); + double xmax = prof[nplots]->GetXaxis()->GetXmax(); + double dx = (xmax - xmin) / nbinX; cout << "Hist " << ii; - for (int ibx=0; ibxGetBinContent(ibx+1); - cout << " | " << ibx << "(" << xx1 << ":" << (xx1+dx) << ") " << cont; + for (int ibx = 0; ibx < nbinX; ibx++) { + double xx1 = xmin + ibx * dx; + double cont = prof[nplots]->GetBinContent(ibx + 1); + cout << " | " << ibx << "(" << xx1 << ":" << (xx1 + dx) << ") " << cont; } cout << "\n"; } @@ -1385,24 +1536,30 @@ void etaPhiFwdPlot(TString fileName, TString plot, int first, int last, TString cname = "c_" + plot + ztit + tag; TCanvas *cc1 = new TCanvas(cname, cname, 700, 600); if (xh < 0.8) { - cc1->SetLeftMargin(0.15); cc1->SetRightMargin(0.05); + cc1->SetLeftMargin(0.15); + cc1->SetRightMargin(0.05); } prof[0]->Draw("h"); - for(int i=1; iDraw("h sames"); - if (drawLeg > 0) leg->Draw("sames"); + if (drawLeg > 0) + leg->Draw("sames"); } -void setStyle () { - - gStyle->SetCanvasBorderMode(0); gStyle->SetCanvasColor(kWhite); - gStyle->SetPadColor(kWhite); gStyle->SetFrameBorderMode(0); - gStyle->SetFrameBorderSize(1); gStyle->SetFrameFillColor(0); - gStyle->SetFrameFillStyle(0); gStyle->SetFrameLineColor(1); - gStyle->SetFrameLineStyle(1); gStyle->SetFrameLineWidth(1); - gStyle->SetOptStat(0); gStyle->SetLegendBorderSize(1); - gStyle->SetOptTitle(0); gStyle->SetTitleOffset(2.5,"Y"); - +void setStyle() { + gStyle->SetCanvasBorderMode(0); + gStyle->SetCanvasColor(kWhite); + gStyle->SetPadColor(kWhite); + gStyle->SetFrameBorderMode(0); + gStyle->SetFrameBorderSize(1); + gStyle->SetFrameFillColor(0); + gStyle->SetFrameFillStyle(0); + gStyle->SetFrameLineColor(1); + gStyle->SetFrameLineStyle(1); + gStyle->SetFrameLineWidth(1); + gStyle->SetOptStat(0); + gStyle->SetLegendBorderSize(1); + gStyle->SetOptTitle(0); + gStyle->SetTitleOffset(2.5, "Y"); } -