00001 #include "FastSimulation/CaloHitMakers/interface/CaloHitMaker.h"
00002 #include "FastSimulation/CaloGeometryTools/interface/CaloGeometryHelper.h"
00003 #include "FastSimulation/CalorimeterProperties/interface/CalorimeterProperties.h"
00004 #include "FastSimulation/CalorimeterProperties/interface/PreshowerProperties.h"
00005 #include "FastSimulation/CalorimeterProperties/interface/HCALProperties.h"
00006 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
00007
00008 typedef ROOT::Math::Plane3D::Point Point;
00009
00010 CaloHitMaker::CaloHitMaker(const CaloGeometryHelper * theCalo,DetId::Detector basedet,int subdetn,int cal,unsigned sht)
00011 :myCalorimeter(theCalo),theCaloProperties(NULL),base_(basedet),subdetn_(subdetn),onCal_(cal),showerType_(sht)
00012 {
00013
00014 EMSHOWER=(sht==0);
00015 HADSHOWER=(sht==1);
00016 if(base_==DetId::Ecal&&(subdetn_==EcalBarrel||subdetn==EcalEndcap)&&onCal_)
00017 theCaloProperties = (CalorimeterProperties*)myCalorimeter->ecalProperties(onCal_);
00018
00019 if(base_==DetId::Ecal&&subdetn_==EcalPreshower&&onCal_)
00020 theCaloProperties = (PreshowerProperties*)myCalorimeter->layer1Properties(onCal_);
00021 if(base_==DetId::Hcal&&cal) theCaloProperties = myCalorimeter->hcalProperties(onCal_);
00022
00023 if(theCaloProperties)
00024 {
00025 moliereRadius=theCaloProperties->moliereRadius();
00026 interactionLength=theCaloProperties->interactionLength();
00027 }
00028 else
00029 {
00030 moliereRadius=999;
00031 interactionLength=999;
00032 }
00033 }
00034
00035
00036 CaloHitMaker::XYZPoint
00037 CaloHitMaker::intersect(const Plane3D& p,const XYZPoint& a,const XYZPoint& b,
00038 double& t,bool segment, bool debug)
00039 {
00040 t=-9999.;
00041
00042 XYZVector normal = p.Normal();
00043 double AAA = normal.X();
00044 double BBB = normal.Y();
00045 double CCC = normal.Z();
00046
00047 double DDD = p.HesseDistance();
00048
00049 double denom = AAA*(b.X()-a.X()) + BBB*(b.Y()-a.Y()) + CCC*(b.Z()-a.Z());
00050 if(denom!=0.)
00051 {
00052
00053 t=-(AAA*a.X()+BBB*a.Y()+CCC*a.Z()+DDD);
00054 t/=denom;
00055 if(debug) std::cout << " T = " << t <<std::endl;
00056 if(segment)
00057 {
00058 if(t>=0&&t<=1)
00059 return XYZPoint(a.X()+(b.X()-a.X())*t,
00060 a.Y()+(b.Y()-a.Y())*t,
00061 a.Z()+(b.Z()-a.Z())*t);
00062 }
00063 else
00064 {
00065 return XYZPoint(a.X()+(b.X()-a.X())*t,
00066 a.Y()+(b.Y()-a.Y())*t,
00067 a.Z()+(b.Z()-a.Z())*t);
00068 }
00069
00070
00071 }
00072
00073 return XYZPoint(0.,0.,0.);
00074 }