CaloGeometryTools/src/CaloSegment.cc

0001 //FAMOS headers
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   // Change this. CaloProperties from FamosShower should be used instead
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;  //cracks : Al
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       // From Olga's & Patrick's talk PRS/JetMET 21 Sept 2004
0059       radLenIncm = 22.3;
0060       intLenIncm = 140;
0061     } break;
0062     case PSEEGAP: {
0063       // according to Sunanda 0.19 X0 (0.08X0 of polyethylene), support (0.06X0 of aluminium)  + other stuff
0064       // in the geometry 12 cm between layer and entrance of EE. Polyethylene is rather 48 and Al 8.9 (PDG)
0065       // for the inLen, just rescale according to PDG (85cm)
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     //  ost<< " Entrance side " << seg.entrance().getSide()
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   //  ost << " Exit side " << seg.exit().getSide()
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 }