CaloHitMakers/src/CaloHitMaker.cc

0001 #include "FastSimulation/CaloHitMakers/interface/CaloHitMaker.h"
0002 #include "FastSimulation/CaloGeometryTools/interface/CaloGeometryHelper.h"
0003 #include "FastSimulation/CalorimeterProperties/interface/CalorimeterProperties.h"
0004 #include "FastSimulation/CalorimeterProperties/interface/PreshowerProperties.h"
0005 #include "FastSimulation/CalorimeterProperties/interface/PreshowerLayer1Properties.h"
0006 #include "FastSimulation/CalorimeterProperties/interface/ECALProperties.h"
0007 #include "FastSimulation/CalorimeterProperties/interface/HCALProperties.h"
0008 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
0009
0010 typedef ROOT::Math::Plane3D::Point Point;
0011
0012 CaloHitMaker::CaloHitMaker(
0013     const CaloGeometryHelper* theCalo, DetId::Detector basedet, int subdetn, int cal, unsigned sht)
0014     : myCalorimeter(theCalo),
0015       theCaloProperties(nullptr),
0016       base_(basedet),
0017       subdetn_(subdetn),
0018       onCal_(cal),
0019       showerType_(sht) {
0020   //  std::cout << " FamosCalorimeter " << basedet << " " << cal << std::endl;
0021   EMSHOWER = (sht == 0);
0022   HADSHOWER = (sht == 1);
0023   MIP = (sht == 2);
0024   if (base_ == DetId::Ecal && (subdetn_ == EcalBarrel || subdetn == EcalEndcap) && onCal_)
0025     theCaloProperties = myCalorimeter->ecalProperties(onCal_);
0026   // is it really necessary to cast here ?
0027   if (base_ == DetId::Ecal && subdetn_ == EcalPreshower && onCal_)
0028     theCaloProperties = myCalorimeter->layer1Properties(onCal_);
0029   if (base_ == DetId::Hcal && cal)
0030     theCaloProperties = myCalorimeter->hcalProperties(onCal_);
0031
0032   if (theCaloProperties) {
0033     moliereRadius = theCaloProperties->moliereRadius();
0034     interactionLength = theCaloProperties->interactionLength();
0035   } else {
0036     moliereRadius = 999;
0037     interactionLength = 999;
0038   }
0039 }
0040
0041 CaloHitMaker::XYZPoint CaloHitMaker::intersect(
0042     const Plane3D& p, const XYZPoint& a, const XYZPoint& b, double& t, bool segment, bool debug) {
0043   t = -9999.;
0044   // En Attendant //
0045   XYZVector normal = p.Normal();
0046   double AAA = normal.X();
0047   double BBB = normal.Y();
0048   double CCC = normal.Z();
0049   //  double DDD = p.Distance(Point(0.,0.,0.));
0050   double DDD = p.HesseDistance();
0051   //  double denom = p.A()*(b.X()-a.X()) + p.B()*(b.Y()-a.Y()) + p.C()*(b.Z()-a.Z());
0052   double denom = AAA * (b.X() - a.X()) + BBB * (b.Y() - a.Y()) + CCC * (b.Z() - a.Z());
0053   if (denom != 0.) {
0054     // t=-(p.A()*a.X()+p.B()*a.Y()+p.C()*a.Z()+p.D());
0055     t = -(AAA * a.X() + BBB * a.Y() + CCC * a.Z() + DDD);
0056     t /= denom;
0057     if (debug)
0058       std::cout << " T = " << t << std::endl;
0059     if (segment) {
0060       if (t >= 0 && t <= 1)
0061         return XYZPoint(a.X() + (b.X() - a.X()) * t, a.Y() + (b.Y() - a.Y()) * t, a.Z() + (b.Z() - a.Z()) * t);
0062     } else {
0063       return XYZPoint(a.X() + (b.X() - a.X()) * t, a.Y() + (b.Y() - a.Y()) * t, a.Z() + (b.Z() - a.Z()) * t);
0064     }
0065   }
0066
0067   return XYZPoint(0., 0., 0.);
0068 }