Skip to content

Commit 415c635

Browse files
author
Sunanda
committed
Change ShowerLibrary Producer by Lev
1 parent 67f6b89 commit 415c635

File tree

4 files changed

+873
-900
lines changed

4 files changed

+873
-900
lines changed

‎SimG4CMS/Calo/src/HFFibre.cc‎

Lines changed: 145 additions & 153 deletions
Original file line numberDiff line numberDiff line change
@@ -1,153 +1,145 @@
1-
///////////////////////////////////////////////////////////////////////////////
2-
// File: HFFibre.cc
3-
// Description: Loads the table for attenuation length and calculates it
4-
///////////////////////////////////////////////////////////////////////////////
5-
6-
#include "SimG4CMS/Calo/interface/HFFibre.h"
7-
#include "FWCore/Utilities/interface/Exception.h"
8-
9-
#include "CLHEP/Units/GlobalSystemOfUnits.h"
10-
#include "CLHEP/Units/GlobalPhysicalConstants.h"
11-
#include <iostream>
12-
#include <sstream>
13-
14-
//#define EDM_ML_DEBUG
15-
16-
HFFibre::HFFibre(const std::string& name,
17-
const HcalDDDSimConstants* hcons,
18-
const HcalSimulationParameters* hps,
19-
edm::ParameterSet const& p)
20-
: hcalConstant_(hcons), hcalsimpar_(hps) {
21-
edm::ParameterSet m_HF =
22-
(p.getParameter<edm::ParameterSet>("HFShower")).getParameter<edm::ParameterSet>("HFShowerBlock");
23-
cFibre = c_light * (m_HF.getParameter<double>("CFibre"));
24-
25-
edm::LogVerbatim("HFShower") << "HFFibre:: Speed of light in fibre " << cFibre << " m/ns";
26-
27-
// Attenuation length
28-
attL = hcalsimpar_->attenuationLength_;
29-
nBinAtt = static_cast<int>(attL.size());
30-
#ifdef EDM_ML_DEBUG
31-
std::stringstream ss1;
32-
for (int it = 0; it < nBinAtt; it++) {
33-
if (it / 10 * 10 == it) {
34-
ss1 << "\n";
35-
}
36-
ss1 << " " << attL[it] * CLHEP::cm;
37-
}
38-
edm::LogVerbatim("HFShower") << "HFFibre: " << nBinAtt << " attL(1/cm): " << ss1.str();
39-
#endif
40-
// Limits on Lambda
41-
std::vector<int> nvec = hcalsimpar_->lambdaLimits_;
42-
lambLim[0] = nvec[0];
43-
lambLim[1] = nvec[1];
44-
#ifdef EDM_ML_DEBUG
45-
edm::LogVerbatim("HFShower") << "HFFibre: Limits on lambda " << lambLim[0] << " and " << lambLim[1];
46-
#endif
47-
// Fibre Lengths
48-
longFL = hcalsimpar_->longFiberLength_;
49-
#ifdef EDM_ML_DEBUG
50-
std::stringstream ss2;
51-
for (unsigned int it = 0; it < longFL.size(); it++) {
52-
if (it / 10 * 10 == it) {
53-
ss2 << "\n";
54-
}
55-
ss2 << " " << longFL[it] / CLHEP::cm;
56-
}
57-
edm::LogVerbatim("HFShower") << "HFFibre: " << longFL.size() << " Long Fibre Length(cm):" << ss2.str();
58-
#endif
59-
shortFL = hcalsimpar_->shortFiberLength_;
60-
#ifdef EDM_ML_DEBUG
61-
std::stringstream ss3;
62-
for (unsigned int it = 0; it < shortFL.size(); it++) {
63-
if (it / 10 * 10 == it) {
64-
ss3 << "\n";
65-
}
66-
ss3 << " " << shortFL[it] / CLHEP::cm;
67-
}
68-
edm::LogVerbatim("HFShower") << "HFFibre: " << shortFL.size() << " Short Fibre Length(cm):" << ss3.str();
69-
#endif
70-
71-
// Now geometry parameters
72-
gpar = hcalConstant_->getGparHF();
73-
radius = hcalConstant_->getRTableHF();
74-
75-
nBinR = static_cast<int>(radius.size());
76-
#ifdef EDM_ML_DEBUG
77-
std::stringstream sss;
78-
for (int i = 0; i < nBinR; ++i) {
79-
if (i / 10 * 10 == i) {
80-
sss << "\n";
81-
}
82-
sss << " " << radius[i] / CLHEP::cm;
83-
}
84-
edm::LogVerbatim("HFShower") << "HFFibre: " << radius.size() << " rTable(cm):" << sss.str();
85-
#endif
86-
}
87-
88-
double HFFibre::attLength(double lambda) {
89-
int i = int(nBinAtt * (lambda - lambLim[0]) / (lambLim[1] - lambLim[0]));
90-
91-
int j = i;
92-
if (i >= nBinAtt)
93-
j = nBinAtt - 1;
94-
else if (i < 0)
95-
j = 0;
96-
double att = attL[j];
97-
#ifdef EDM_ML_DEBUG
98-
edm::LogVerbatim("HFShower") << "HFFibre::attLength for Lambda " << lambda << " index " << i << " " << j
99-
<< " Att. Length " << att;
100-
#endif
101-
return att;
102-
}
103-
104-
double HFFibre::tShift(const G4ThreeVector& point, int depth, int fromEndAbs) {
105-
double zFibre = zShift(point, depth, fromEndAbs);
106-
double time = zFibre / cFibre;
107-
#ifdef EDM_ML_DEBUG
108-
edm::LogVerbatim("HFShower") << "HFFibre::tShift for point " << point << " ( depth = " << depth
109-
<< ", traversed length = " << zFibre / CLHEP::cm << " cm) = " << time / CLHEP::ns
110-
<< " ns";
111-
#endif
112-
return time;
113-
}
114-
115-
double HFFibre::zShift(const G4ThreeVector& point, int depth, int fromEndAbs) { // point is z-local
116-
117-
double zFibre = 0;
118-
double hR = sqrt((point.x()) * (point.x()) + (point.y()) * (point.y()));
119-
int ieta = 0;
120-
double length = 250 * CLHEP::cm;
121-
if (fromEndAbs < 0) {
122-
zFibre = 0.5 * gpar[1] - point.z(); // Never, as fromEndAbs=0 (?)
123-
} else {
124-
// Defines the Radius bin by radial subdivision
125-
for (int i = nBinR - 1; i > 0; --i)
126-
if (hR < radius[i])
127-
ieta = nBinR - i - 1;
128-
// define the length of the fibre
129-
if (depth == 2) {
130-
if ((int)(shortFL.size()) > ieta)
131-
length = shortFL[ieta];
132-
} else {
133-
if ((int)(longFL.size()) > ieta)
134-
length = longFL[ieta];
135-
}
136-
zFibre = length;
137-
if (fromEndAbs > 0) {
138-
zFibre -= gpar[1]; // Never, as fromEndAbs=0 (M.K. ?)
139-
} else {
140-
double zz = 0.5 * gpar[1] + point.z();
141-
zFibre -= zz;
142-
}
143-
if (depth == 2)
144-
zFibre += gpar[0]; // here zFibre is reduced for Short
145-
}
146-
147-
#ifdef EDM_ML_DEBUG
148-
edm::LogVerbatim("HFShower") << "HFFibre::zShift for point " << point << " (R = " << hR / CLHEP::cm
149-
<< " cm, Index = " << ieta << ", depth = " << depth
150-
<< ", Fibre Length = " << length / CLHEP::cm << " cm = " << zFibre / CLHEP::cm << " cm)";
151-
#endif
152-
return zFibre;
153-
}
1+
///////////////////////////////////////////////////////////////////////////////
2+
// File: HFFibre.cc
3+
// Description: Loads the table for attenuation length and calculates it
4+
///////////////////////////////////////////////////////////////////////////////
5+
6+
#include "SimG4CMS/Calo/interface/HFFibre.h"
7+
#include "FWCore/Utilities/interface/Exception.h"
8+
9+
#include "CLHEP/Units/GlobalSystemOfUnits.h"
10+
#include "CLHEP/Units/GlobalPhysicalConstants.h"
11+
#include <iostream>
12+
#include <sstream>
13+
14+
//#define EDM_ML_DEBUG
15+
16+
HFFibre::HFFibre(const std::string& name,
17+
const HcalDDDSimConstants* hcons,
18+
const HcalSimulationParameters* hps,
19+
edm::ParameterSet const& p)
20+
: hcalConstant_(hcons), hcalsimpar_(hps) {
21+
edm::ParameterSet m_HF = (p.getParameter<edm::ParameterSet>("HFShower")).getParameter<edm::ParameterSet>("HFShowerBlock");
22+
cFibre = c_light * (m_HF.getParameter<double>("CFibre"));
23+
24+
edm::LogVerbatim("HFShower") << "HFFibre:: Speed of light in fibre " << cFibre << " m/ns";
25+
26+
// Attenuation length
27+
attL = hcalsimpar_->attenuationLength_;
28+
nBinAtt = static_cast<int>(attL.size());
29+
#ifdef EDM_ML_DEBUG
30+
std::stringstream ss1;
31+
for (int it = 0; it < nBinAtt; it++) {
32+
if (it / 10 * 10 == it) {
33+
ss1 << "\n";
34+
}
35+
ss1 << " " << attL[it] * CLHEP::cm;
36+
}
37+
edm::LogVerbatim("HFShower") << "HFFibre: " << nBinAtt << " attL(1/cm): " << ss1.str();
38+
#endif
39+
// Limits on Lambda
40+
std::vector<int> nvec = hcalsimpar_->lambdaLimits_;
41+
lambLim[0] = nvec[0];
42+
lambLim[1] = nvec[1];
43+
#ifdef EDM_ML_DEBUG
44+
edm::LogVerbatim("HFShower") << "HFFibre: Limits on lambda " << lambLim[0] << " and " << lambLim[1];
45+
#endif
46+
// Fibre Lengths
47+
longFL = hcalsimpar_->longFiberLength_;
48+
#ifdef EDM_ML_DEBUG
49+
std::stringstream ss2;
50+
for (unsigned int it = 0; it < longFL.size(); it++) {
51+
if (it / 10 * 10 == it) {
52+
ss2 << "\n";
53+
}
54+
ss2 << " " << longFL[it] / CLHEP::cm;
55+
}
56+
edm::LogVerbatim("HFShower") << "HFFibre: " << longFL.size() << " Long Fibre Length(cm):" << ss2.str();
57+
#endif
58+
shortFL = hcalsimpar_->shortFiberLength_;
59+
#ifdef EDM_ML_DEBUG
60+
std::stringstream ss3;
61+
for (unsigned int it = 0; it < shortFL.size(); it++) {
62+
if (it / 10 * 10 == it) {
63+
ss3 << "\n";
64+
}
65+
ss3 << " " << shortFL[it] / CLHEP::cm;
66+
}
67+
edm::LogVerbatim("HFShower") << "HFFibre: " << shortFL.size() << " Short Fibre Length(cm):" << ss3.str();
68+
#endif
69+
70+
// Now geometry parameters
71+
gpar = hcalConstant_->getGparHF();
72+
radius = hcalConstant_->getRTableHF();
73+
74+
nBinR = static_cast<int>(radius.size());
75+
#ifdef EDM_ML_DEBUG
76+
std::stringstream sss;
77+
for (int i = 0; i < nBinR; ++i) {
78+
if (i / 10 * 10 == i) {
79+
sss << "\n";
80+
}
81+
sss << " " << radius[i] / CLHEP::cm;
82+
}
83+
edm::LogVerbatim("HFShower") << "HFFibre: " << radius.size() << " rTable(cm):" << sss.str();
84+
#endif
85+
}
86+
87+
double HFFibre::attLength(double lambda) {
88+
int i = int(nBinAtt * (lambda - lambLim[0]) / (lambLim[1] - lambLim[0]));
89+
90+
int j = i;
91+
if (i >= nBinAtt)
92+
j = nBinAtt - 1;
93+
else if (i < 0)
94+
j = 0;
95+
double att = attL[j];
96+
#ifdef EDM_ML_DEBUG
97+
edm::LogVerbatim("HFShower") << "HFFibre::attLength for Lambda " << lambda << " index " << i << " " << j
98+
<< " Att. Length " << att;
99+
#endif
100+
return att;
101+
}
102+
103+
double HFFibre::tShift(const G4ThreeVector& point, int depth, int fromEndAbs) {
104+
double zFibre = zShift(point, depth, fromEndAbs);
105+
double time = zFibre / cFibre;
106+
#ifdef EDM_ML_DEBUG
107+
edm::LogVerbatim("HFShower") << "HFFibre::tShift for point " << point << " ( depth = " << depth
108+
<< ", traversed length = " << zFibre / CLHEP::cm << " cm) = " << time / CLHEP::ns
109+
<< " ns";
110+
#endif
111+
return time;
112+
}
113+
114+
double HFFibre::zShift(const G4ThreeVector& point, int depth, int fromEndAbs) {
115+
// point is z-local
116+
double zFibre = 0;
117+
double hR = sqrt((point.x()) * (point.x()) + (point.y()) * (point.y()));
118+
int ieta = 0;
119+
double length = 250 * CLHEP::cm;
120+
if (fromEndAbs < 0) {
121+
zFibre = 0.5 * gpar[1] - point.z(); //
122+
} else {
123+
// Defines the Radius bin by radial subdivision
124+
for (int i = nBinR - 1; i > 0; --i)
125+
if (hR < radius[i])
126+
ieta = nBinR - i - 1;
127+
// Defines the full length of the fibre (For onlyLong)
128+
if (static_cast<int>(longFL.size()) > ieta)
129+
length = longFL[ieta];
130+
zFibre = length; // from beginning of abs (full length)
131+
if (fromEndAbs > 0) {
132+
zFibre -= gpar[1]; // length from end of HF
133+
} else {
134+
double zz = 0.5 * gpar[1] + point.z(); // depth of point of photon emission (from beginning of HF)
135+
zFibre -= zz; // length of fiber from point of photon emission
136+
}
137+
}
138+
139+
#ifdef EDM_ML_DEBUG
140+
edm::LogVerbatim("HFShower") << "HFFibre::zShift for point " << point << " (R = " << hR / CLHEP::cm
141+
<< " cm, Index = " << ieta << ", depth = " << depth
142+
<< ", Fibre Length = " << length / CLHEP::cm << " cm = " << zFibre / CLHEP::cm << " cm)";
143+
#endif
144+
return zFibre;
145+
}

0 commit comments

Comments
 (0)