File indexing completed on 2024-04-06 12:11:12
0001
0002 #include "FastSimulation/CaloGeometryTools/interface/CaloSegment.h"
0003 #include "FastSimulation/CaloGeometryTools/interface/CaloGeometryHelper.h"
0004 #include "FastSimulation/CalorimeterProperties/interface/PreshowerLayer1Properties.h"
0005 #include "FastSimulation/CalorimeterProperties/interface/PreshowerLayer2Properties.h"
0006 #include "FastSimulation/CalorimeterProperties/interface/HCALProperties.h"
0007 #include "FastSimulation/CalorimeterProperties/interface/ECALProperties.h"
0008
0009 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
0010
0011 CaloSegment::CaloSegment(const CaloPoint& in,
0012 const CaloPoint& out,
0013 double si,
0014 double siX0,
0015 double siL0,
0016 Material mat,
0017 const CaloGeometryHelper* myCalorimeter)
0018 : entrance_(in),
0019 exit_(out),
0020 sentrance_(si),
0021 sX0entrance_(siX0),
0022 sL0entrance_(siL0),
0023 material_(mat)
0024
0025 {
0026 sexit_ = sentrance_ + std::sqrt((exit_ - entrance_).mag2());
0027
0028 double radLenIncm = 999999;
0029 double intLenIncm = 999999;
0030 detector_ = in.whichDetector();
0031 if (detector_ != out.whichDetector() && mat != CRACK && mat != GAP && mat != ECALHCALGAP) {
0032 std::cout << " Problem in the segments " << detector_ << " " << out.whichDetector() << std::endl;
0033 }
0034 switch (mat) {
0035 case PbWO4: {
0036 int det = 0;
0037 if (in.whichSubDetector() == EcalBarrel)
0038 det = 1;
0039 if (in.whichSubDetector() == EcalEndcap)
0040 det = 2;
0041
0042 radLenIncm = myCalorimeter->ecalProperties(det)->radLenIncm();
0043 intLenIncm = myCalorimeter->ecalProperties(det)->interactionLength();
0044 } break;
0045 case CRACK: {
0046 radLenIncm = 8.9;
0047 intLenIncm = 35.4;
0048 } break;
0049 case PS: {
0050 radLenIncm = myCalorimeter->layer1Properties(1)->radLenIncm();
0051 intLenIncm = myCalorimeter->layer1Properties(1)->interactionLength();
0052 } break;
0053 case HCAL: {
0054 radLenIncm = myCalorimeter->hcalProperties(1)->radLenIncm();
0055 intLenIncm = myCalorimeter->hcalProperties(1)->interactionLength();
0056 } break;
0057 case ECALHCALGAP: {
0058
0059 radLenIncm = 22.3;
0060 intLenIncm = 140;
0061 } break;
0062 case PSEEGAP: {
0063
0064
0065
0066 radLenIncm = myCalorimeter->layer2Properties(1)->pseeRadLenIncm();
0067 intLenIncm = myCalorimeter->layer2Properties(1)->pseeIntLenIncm();
0068 } break;
0069 default:
0070 radLenIncm = 999999;
0071 }
0072 sX0exit_ = sX0entrance_ + (sexit_ - sentrance_) / radLenIncm;
0073 sL0exit_ = sL0entrance_ + (sexit_ - sentrance_) / intLenIncm;
0074 if (mat == GAP) {
0075 sX0exit_ = sX0entrance_;
0076 sL0exit_ = sL0entrance_;
0077 }
0078 length_ = sexit_ - sentrance_;
0079 X0length_ = sX0exit_ - sX0entrance_;
0080 L0length_ = sL0exit_ - sL0entrance_;
0081 }
0082
0083 CaloSegment::XYZPoint CaloSegment::positionAtDepthincm(double depth) const {
0084 if (depth < sentrance_ || depth > sexit_)
0085 return XYZPoint();
0086 return XYZPoint(entrance_ + ((depth - sentrance_) / (sexit_ - sentrance_) * (exit_ - entrance_)));
0087 }
0088
0089 CaloSegment::XYZPoint CaloSegment::positionAtDepthinX0(double depth) const {
0090 if (depth < sX0entrance_ || depth > sX0exit_)
0091 return XYZPoint();
0092 return XYZPoint(entrance_ + ((depth - sX0entrance_) / (sX0exit_ - sX0entrance_) * (exit_ - entrance_)));
0093 }
0094
0095 CaloSegment::XYZPoint CaloSegment::positionAtDepthinL0(double depth) const {
0096 if (depth < sL0entrance_ || depth > sL0exit_)
0097 return XYZPoint();
0098 return XYZPoint(entrance_ + ((depth - sL0entrance_) / (sL0exit_ - sL0entrance_) * (exit_ - entrance_)));
0099 }
0100
0101 double CaloSegment::x0FromCm(double cm) const { return sX0entrance_ + cm / length_ * X0length_; }
0102
0103 std::ostream& operator<<(std::ostream& ost, const CaloSegment& seg) {
0104 ost << " DetId ";
0105 if (!seg.entrance().getDetId().null())
0106 ost << seg.entrance().getDetId()();
0107 else {
0108 ost << seg.entrance().whichDetector();
0109
0110 ost << " Point " << (math::XYZVector)seg.entrance() << std::endl;
0111 }
0112 ost << "DetId ";
0113 if (!seg.exit().getDetId().null())
0114 ost << seg.exit().getDetId()();
0115 else
0116 ost << seg.exit().whichDetector();
0117
0118
0119 ost << " Point " << (math::XYZVector)seg.exit() << " " << seg.length() << " cm " << seg.X0length() << " X0 "
0120 << seg.L0length() << " Lambda0 ";
0121 switch (seg.material()) {
0122 case CaloSegment::PbWO4:
0123 ost << "PbWO4 ";
0124 break;
0125 case CaloSegment::CRACK:
0126 ost << "CRACK ";
0127 break;
0128 case CaloSegment::PS:
0129 ost << "PS ";
0130 break;
0131 case CaloSegment::HCAL:
0132 ost << "HCAL ";
0133 break;
0134 case CaloSegment::ECALHCALGAP:
0135 ost << "ECAL-HCAL GAP ";
0136 break;
0137 case CaloSegment::PSEEGAP:
0138 ost << "PS-ECAL GAP";
0139 break;
0140 default:
0141 ost << "GAP ";
0142 }
0143 return ost;
0144 }