From 5bead954c91a87ad2c503d55ef61d0eeff55b6d0 Mon Sep 17 00:00:00 2001 From: wonyongc Date: Fri, 25 Oct 2024 15:47:49 -0400 Subject: [PATCH] Simpler dr crystal hits for SCEPCal --- detector/calorimeter/SCEPCalConstructor.cpp | 54 +++++++------- .../detectorSegmentations/DRCrystalHit.h | 54 ++++---------- ...le.h => SCEPCalSegmentationHandle_k4geo.h} | 26 +++---- ...entation.h => SCEPCalSegmentation_k4geo.h} | 14 ++-- detectorSegmentations/src/DRCrystalHit.cpp | 45 ++++-------- ...pp => SCEPCalSegmentationHandle_k4geo.cpp} | 4 +- ...tion.cpp => SCEPCalSegmentation_k4geo.cpp} | 40 +++++------ .../src/plugins/SegmentationFactories.cpp | 2 +- plugins/Geant4Output2EDM4hepDRCrystalHit.cpp | 15 ++-- plugins/SCEPCalSDActionDRCrystalHit.cpp | 72 +++++++------------ 10 files changed, 127 insertions(+), 199 deletions(-) rename detectorSegmentations/include/detectorSegmentations/{SCEPCalSegmentationHandle.h => SCEPCalSegmentationHandle_k4geo.h} (80%) rename detectorSegmentations/include/detectorSegmentations/{SCEPCalSegmentation.h => SCEPCalSegmentation_k4geo.h} (89%) rename detectorSegmentations/src/{SCEPCalSegmentationHandle.cpp => SCEPCalSegmentationHandle_k4geo.cpp} (53%) rename detectorSegmentations/src/{SCEPCalSegmentation.cpp => SCEPCalSegmentation_k4geo.cpp} (87%) diff --git a/detector/calorimeter/SCEPCalConstructor.cpp b/detector/calorimeter/SCEPCalConstructor.cpp index 97eee3fcc..0af2bfa3c 100644 --- a/detector/calorimeter/SCEPCalConstructor.cpp +++ b/detector/calorimeter/SCEPCalConstructor.cpp @@ -3,7 +3,7 @@ // Princeton University //=============================== -#include "detectorSegmentations/SCEPCalSegmentationHandle.h" +#include "detectorSegmentations/SCEPCalSegmentationHandle_k4geo.h" #include "DD4hep/DetFactoryHelper.h" #include "DD4hep/DetectorTools.h" #include "DD4hep/Printout.h" @@ -34,7 +34,7 @@ create_detector_SCEPCal(dd4hep::Detector &theDetector,xml_h xmlElement,dd4hep::S xml_comp_t timingLgXML =detectorXML.child(_Unicode(timingLayerLg)); xml_comp_t sipmLgXML =detectorXML.child(_Unicode(sipmLg)); xml_comp_t sipmTrXML =detectorXML.child(_Unicode(sipmTr)); - xml_comp_t instXML =detectorXML.child(_Unicode(inst)); + // xml_comp_t instXML =detectorXML.child(_Unicode(inst)); xml_comp_t timingAssemblyGlobalVisXML =detectorXML.child(_Unicode(timingAssemblyGlobalVis)); xml_comp_t barrelAssemblyGlobalVisXML =detectorXML.child(_Unicode(barrelAssemblyGlobalVis)); xml_comp_t endcapAssemblyGlobalVisXML =detectorXML.child(_Unicode(endcapAssemblyGlobalVis)); @@ -45,7 +45,7 @@ create_detector_SCEPCal(dd4hep::Detector &theDetector,xml_h xmlElement,dd4hep::S dd4hep::Material timingLgMat =theDetector.material(timingLgXML.materialStr()); dd4hep::Material sipmLgMat =theDetector.material(sipmLgXML.materialStr()); dd4hep::Material sipmTrMat =theDetector.material(sipmTrXML.materialStr()); - dd4hep::Material instMat =theDetector.material(instXML.materialStr()); + // dd4hep::Material instMat =theDetector.material(instXML.materialStr()); const double EBz =dimXML.attr(_Unicode(barrelHalfZ)); const double Rin =dimXML.attr(_Unicode(barrelInnerR)); const double nomfw =dimXML.attr(_Unicode(crystalFaceWidthNominal)); @@ -56,8 +56,8 @@ create_detector_SCEPCal(dd4hep::Detector &theDetector,xml_h xmlElement,dd4hep::S const int PHI_SEGMENTS =dimXML.attr(_Unicode(phiSegments)); const int N_PROJECTIVE_FILL =dimXML.attr(_Unicode(projectiveFill)); const bool CONSTRUCT_TIMING =timingXML.attr(_Unicode(construct)); - const int TIMING_PHI_START =timingXML.attr(_Unicode(phistart)); - const int TIMING_PHI_END =timingXML.attr(_Unicode(phiend)); + // const int TIMING_PHI_START =timingXML.attr(_Unicode(phistart)); + // const int TIMING_PHI_END =timingXML.attr(_Unicode(phiend)); const bool CONSTRUCT_BARREL =barrelXML.attr(_Unicode(construct)); const int BARREL_PHI_START =barrelXML.attr(_Unicode(phistart)); const int BARREL_PHI_END =barrelXML.attr(_Unicode(phiend)); @@ -93,30 +93,30 @@ create_detector_SCEPCal(dd4hep::Detector &theDetector,xml_h xmlElement,dd4hep::S int nCy =floor(lT/nomth); double actY =lT/nCy; double actX =2*wT/nCy; - double r2slice_end =r0slice_end+Fdz+Rdz; - double z2slice_end =r2slice_end*cos(thC_end)+PROJECTIVE_GAP; - double Rin2slice_end =r2slice_end*sin(thC_end); - double y2slice_end =r2slice_end*tan(D_THETA_BARREL/2.); - double slice_front_jut2 =y2slice_end*sin(M_PI/2-thC_end); - double slice_side_jut2 =y2slice_end*cos(M_PI/2-thC_end); + // double r2slice_end =r0slice_end+Fdz+Rdz; + // double z2slice_end =r2slice_end*cos(thC_end)+PROJECTIVE_GAP; + // double Rin2slice_end =r2slice_end*sin(thC_end); + // double y2slice_end =r2slice_end*tan(D_THETA_BARREL/2.); + // double slice_front_jut2 =y2slice_end*sin(M_PI/2-thC_end); + // double slice_side_jut2 =y2slice_end*cos(M_PI/2-thC_end); double barrelSlice_z1 =z0slice_end+slice_side_jut; double barrelSlice_z2 =y2slice; - double barrelSlice_rmin2 =Rin2slice_end-slice_front_jut2; + // double barrelSlice_rmin2 =Rin2slice_end-slice_front_jut2; double thCEnd =THETA_SIZE_ENDCAP-D_THETA_ENDCAP/2; double thCBeg =D_THETA_ENDCAP/2+ENDCAP_THETA_START*D_THETA_ENDCAP; double r0eEnd =EBz/cos(thCEnd); - double r2eEnd =r0eEnd+Fdz+Rdz; - double y2eEnd =r2eEnd*tan(D_THETA_ENDCAP/2.); + // double r2eEnd =r0eEnd+Fdz+Rdz; + // double y2eEnd =r2eEnd*tan(D_THETA_ENDCAP/2.); double r0eBeg =EBz/cos(thCBeg); double r2eBeg =r0eBeg+Fdz+Rdz; double y2eBeg =r2eBeg*tan(D_THETA_ENDCAP/2.); double aEnd =r0eEnd/cos(D_THETA_ENDCAP/2); - double bEnd =sqrt(r2eEnd*r2eEnd+y2eEnd*y2eEnd); + // double bEnd =sqrt(r2eEnd*r2eEnd+y2eEnd*y2eEnd); double z1End =aEnd*cos(thCEnd+D_THETA_ENDCAP/2); - double z2End =bEnd*cos(thCEnd-D_THETA_ENDCAP/2); - double aBeg =r0eBeg/cos(D_THETA_ENDCAP/2); + // double z2End =bEnd*cos(thCEnd-D_THETA_ENDCAP/2); + // double aBeg =r0eBeg/cos(D_THETA_ENDCAP/2); double bBeg =sqrt(r2eBeg*r2eBeg+y2eBeg*y2eBeg); - double z1Beg =aBeg*cos(thCBeg+D_THETA_ENDCAP/2); + // double z1Beg =aBeg*cos(thCBeg+D_THETA_ENDCAP/2); double z2Beg =bBeg*cos(thCBeg-D_THETA_ENDCAP/2); double z1rmaxE =z1End*tan(thCEnd+D_THETA_ENDCAP/2); double z2rmaxE =z2Beg*tan(thCEnd+D_THETA_ENDCAP/2); @@ -175,7 +175,7 @@ create_detector_SCEPCal(dd4hep::Detector &theDetector,xml_h xmlElement,dd4hep::S dd4hep::Readout readout=sens.readout(); dd4hep::Segmentation geomseg=readout.segmentation(); dd4hep::Segmentation* _geoSeg=&geomseg; - auto segmentation=dynamic_cast(_geoSeg->segmentation()); + auto segmentation=dynamic_cast(_geoSeg->segmentation()); segmentation->setGeomParams(Fdz,Rdz,nomfw,nomth,EBz,Rin,sipmth,PHI_SEGMENTS,N_PROJECTIVE_FILL); std::vector zTimingPolyhedra ={-barrelSlice_z1,barrelSlice_z1}; @@ -211,10 +211,10 @@ create_detector_SCEPCal(dd4hep::Detector &theDetector,xml_h xmlElement,dd4hep::S int endcapAssemblyVolId32=segmentation->getFirst32bits(endcapAssemblyVolId); auto endcap1AssemblyVolId =segmentation->setVolumeID(2,1,0,0); int endcap1AssemblyVolId32=segmentation->getFirst32bits(endcap1AssemblyVolId); - dd4hep::PlacedVolume timingPlacedVol =experimentalHall.placeVolume(timingAssemblyVol,timingAssemblyVolId32); + experimentalHall.placeVolume(timingAssemblyVol,timingAssemblyVolId32); dd4hep::PlacedVolume barrelPlacedVol =experimentalHall.placeVolume(barrelAssemblyVol,barrelAssemblyVolId32); - dd4hep::PlacedVolume endcapPlacedVol =experimentalHall.placeVolume(endcapAssemblyVol,endcapAssemblyVolId32); - dd4hep::PlacedVolume endcap1PlacedVol=experimentalHall.placeVolume(endcap1AssemblyVol,endcap1AssemblyVolId32); + experimentalHall.placeVolume(endcapAssemblyVol,endcapAssemblyVolId32); + experimentalHall.placeVolume(endcap1AssemblyVol,endcap1AssemblyVolId32); ScepcalDetElement.setPlacement(barrelPlacedVol); @@ -231,7 +231,7 @@ create_detector_SCEPCal(dd4hep::Detector &theDetector,xml_h xmlElement,dd4hep::S RotationZ rotZPhi(phiEnvBarrel); double rTimingAssembly=rT+nomth; Position dispTimingAssembly(rTimingAssembly*cos(phiEnvBarrel),rTimingAssembly*sin(phiEnvBarrel),0); - dd4hep::PlacedVolume timingPhiAssemblyPlacedVol=timingAssemblyVol.placeVolume(timingPhiAssemblyVolume,Transform3D(rotZPhi,dispTimingAssembly)); + timingAssemblyVol.placeVolume(timingPhiAssemblyVolume,Transform3D(rotZPhi,dispTimingAssembly)); dd4hep::Box timingCrystalLg(nomth/2,actX/2,lT/2-sipmth); dd4hep::Box timingCrystalTr(nomth/2,wT-sipmth,actY/2); @@ -255,7 +255,7 @@ create_detector_SCEPCal(dd4hep::Detector &theDetector,xml_h xmlElement,dd4hep::S dd4hep::Volume tileAssemblyVolume("tileAssembly",tileAssemblyShape,theDetector.material("Vacuum")); tileAssemblyVolume.setVisAttributes(theDetector,scepcalAssemblyXML.visStr()); Position dispTileAssembly(0,0,-y1slice+nTile*lT+lT/2); - dd4hep::PlacedVolume tileAssemblyPlacedVol=timingPhiAssemblyVolume.placeVolume(tileAssemblyVolume,dispTileAssembly); + timingPhiAssemblyVolume.placeVolume(tileAssemblyVolume,dispTileAssembly); for (int nC=0;nC nwavelen_cer; - std::array nwavelen_scint; - - std::array ntime_cer; - std::array ntime_scint; + int nCerenkovProd; + int nScintillationProd; + double tAvgC; + double tAvgS; public: DRCrystalHit(); @@ -63,13 +44,4 @@ namespace SCEPCal { }; }; -#if defined(__CINT__) || defined(__MAKECINT__) || defined(__CLING__) || defined(__ROOTCLING__) -#pragma link off all globals; -#pragma link off all classes; -#pragma link off all functions; -#pragma link C++ namespace dd4hep; -#pragma link C++ namespace dd4hep::sim; -#pragma link C++ namespace SCEPCal; -#pragma link C++ class SCEPCal::DRCrystalHit+; -#endif #endif diff --git a/detectorSegmentations/include/detectorSegmentations/SCEPCalSegmentationHandle.h b/detectorSegmentations/include/detectorSegmentations/SCEPCalSegmentationHandle_k4geo.h similarity index 80% rename from detectorSegmentations/include/detectorSegmentations/SCEPCalSegmentationHandle.h rename to detectorSegmentations/include/detectorSegmentations/SCEPCalSegmentationHandle_k4geo.h index a6c48d10b..31104d251 100644 --- a/detectorSegmentations/include/detectorSegmentations/SCEPCalSegmentationHandle.h +++ b/detectorSegmentations/include/detectorSegmentations/SCEPCalSegmentationHandle_k4geo.h @@ -2,9 +2,9 @@ // Author: Wonyong Chung // Princeton University //=============================== -#ifndef SCEPCalSegmentationHandle_h -#define SCEPCalSegmentationHandle_h 1 -#include "detectorSegmentations/SCEPCalSegmentation.h" +#ifndef SCEPCalSegmentationHandle_k4geo_h +#define SCEPCalSegmentationHandle_k4geo_h 1 +#include "detectorSegmentations/SCEPCalSegmentation_k4geo.h" #include "DD4hep/Segmentations.h" #include "DD4hep/detail/SegmentationsInterna.h" @@ -14,21 +14,21 @@ class Segmentation; template class SegmentationWrapper; -typedef Handle> SCEPCalSegmentationHandle; +typedef Handle> SCEPCalSegmentationHandle_k4geo; -class SCEPCalSegmentation : public SCEPCalSegmentationHandle { +class SCEPCalSegmentation_k4geo : public SCEPCalSegmentationHandle_k4geo { public: - typedef SCEPCalSegmentationHandle::Object Object; + typedef SCEPCalSegmentationHandle_k4geo::Object Object; public: - SCEPCalSegmentation() = default; - SCEPCalSegmentation(const SCEPCalSegmentation& e) = default; - SCEPCalSegmentation(const Segmentation& e) : Handle(e) {} - SCEPCalSegmentation(const Handle& e) : Handle(e) {} + SCEPCalSegmentation_k4geo() = default; + SCEPCalSegmentation_k4geo(const SCEPCalSegmentation_k4geo& e) = default; + SCEPCalSegmentation_k4geo(const Segmentation& e) : Handle(e) {} + SCEPCalSegmentation_k4geo(const Handle& e) : Handle(e) {} template - SCEPCalSegmentation(const Handle& e) : Handle(e) {} - SCEPCalSegmentation& operator=(const SCEPCalSegmentation& seg) = default; - bool operator==(const SCEPCalSegmentation& seg) const { return m_element == seg.m_element; } + SCEPCalSegmentation_k4geo(const Handle& e) : Handle(e) {} + SCEPCalSegmentation_k4geo& operator=(const SCEPCalSegmentation_k4geo& seg) = default; + bool operator==(const SCEPCalSegmentation_k4geo& seg) const { return m_element == seg.m_element; } inline Position position(const CellID& id) const { return Position(access()->implementation->position(id)); diff --git a/detectorSegmentations/include/detectorSegmentations/SCEPCalSegmentation.h b/detectorSegmentations/include/detectorSegmentations/SCEPCalSegmentation_k4geo.h similarity index 89% rename from detectorSegmentations/include/detectorSegmentations/SCEPCalSegmentation.h rename to detectorSegmentations/include/detectorSegmentations/SCEPCalSegmentation_k4geo.h index e25754498..26f184b29 100644 --- a/detectorSegmentations/include/detectorSegmentations/SCEPCalSegmentation.h +++ b/detectorSegmentations/include/detectorSegmentations/SCEPCalSegmentation_k4geo.h @@ -2,8 +2,8 @@ // Author: Wonyong Chung // Princeton University //=============================== -#ifndef SCEPCalSegmentation_h -#define SCEPCalSegmentation_h 1 +#ifndef SCEPCalSegmentation_k4geo_h +#define SCEPCalSegmentation_k4geo_h 1 #include "DDSegmentation/Segmentation.h" #include "TVector3.h" #include "DD4hep/DetFactoryHelper.h" @@ -13,15 +13,15 @@ namespace dd4hep { namespace DDSegmentation { -class SCEPCalSegmentation : public Segmentation { +class SCEPCalSegmentation_k4geo : public Segmentation { public: - SCEPCalSegmentation(const std::string& aCellEncoding); - SCEPCalSegmentation(const BitFieldCoder* decoder); - virtual ~SCEPCalSegmentation() override; + SCEPCalSegmentation_k4geo(const std::string& aCellEncoding); + SCEPCalSegmentation_k4geo(const BitFieldCoder* decoder); + virtual ~SCEPCalSegmentation_k4geo() override; virtual Vector3D position(const CellID& aCellID) const; - virtual Vector3D myPosition(const CellID& aCellID) ; + virtual Vector3D myPosition(const CellID& aCellID) const; virtual CellID cellID(const Vector3D& aLocalPosition, const Vector3D& aGlobalPosition, diff --git a/detectorSegmentations/src/DRCrystalHit.cpp b/detectorSegmentations/src/DRCrystalHit.cpp index 775d0b2a2..0cb35d027 100644 --- a/detectorSegmentations/src/DRCrystalHit.cpp +++ b/detectorSegmentations/src/DRCrystalHit.cpp @@ -5,45 +5,28 @@ #include #include #include -#include -#include -#include -#include -#include +// #include +// #include +// #include +// #include +// #include #include "detectorSegmentations/DRCrystalHit.h" -#include "G4Track.hh" +// #include "G4Track.hh" -using namespace dd4hep::sim; -using namespace dd4hep; -using namespace std; -using namespace SCEPCal; +SCEPCal::DRCrystalHit::DRCrystalHit() +: dd4hep::sim::Geant4HitData(), position(), truth(), energyDeposit(0), nCerenkovProd(0), nScintillationProd(0), tAvgC(0), tAvgS(0) { -DRCrystalHit::DRCrystalHit() -: Geant4HitData(), position(), truth(), energyDeposit(0), eta(0), phi(0), depth(0), system(0), ncerenkov(0), nscintillator(0) { + dd4hep::InstanceCount::increment(this); - InstanceCount::increment(this); - - for( int i=0; i #include #include @@ -10,8 +10,8 @@ namespace dd4hep { namespace DDSegmentation { -SCEPCalSegmentation::SCEPCalSegmentation(const std::string& cellEncoding) : Segmentation(cellEncoding) { - _type = "SCEPCalSegmentation"; +SCEPCalSegmentation_k4geo::SCEPCalSegmentation_k4geo(const std::string& cellEncoding) : Segmentation(cellEncoding) { + _type = "SCEPCalSegmentation_k4geo"; _description = "SCEPCal segmentation based on side/eta/phi/depth/S/C"; registerIdentifier("identifier_system", "Cell ID identifier for numSystem", fSystemId, "system"); registerIdentifier("identifier_eta", "Cell ID identifier for numEta", fEtaId, "eta"); @@ -19,8 +19,8 @@ SCEPCalSegmentation::SCEPCalSegmentation(const std::string& cellEncoding) : Segm registerIdentifier("identifier_depth", "Cell ID identifier for numDepth", fDepthId, "depth"); } -SCEPCalSegmentation::SCEPCalSegmentation(const BitFieldCoder* decoder) : Segmentation(decoder) { - _type = "SCEPCalSegmentation"; +SCEPCalSegmentation_k4geo::SCEPCalSegmentation_k4geo(const BitFieldCoder* decoder) : Segmentation(decoder) { + _type = "SCEPCalSegmentation_k4geo"; _description = "SCEPCal segmentation based on side/eta/phi/depth/S/C"; registerIdentifier("identifier_system", "Cell ID identifier for numSystem", fSystemId, "system"); registerIdentifier("identifier_eta", "Cell ID identifier for Eta", fEtaId, "eta"); @@ -28,13 +28,13 @@ SCEPCalSegmentation::SCEPCalSegmentation(const BitFieldCoder* decoder) : Segment registerIdentifier("identifier_depth", "Cell ID identifier for Depth", fDepthId, "depth"); } -SCEPCalSegmentation::~SCEPCalSegmentation() {} +SCEPCalSegmentation_k4geo::~SCEPCalSegmentation_k4geo() {} -Vector3D SCEPCalSegmentation::position(const CellID& cID) const { - return Vector3D(0,0,0); +Vector3D SCEPCalSegmentation_k4geo::position(const CellID& cID) const { + return myPosition(cID); }; -Vector3D SCEPCalSegmentation::myPosition(const CellID& cID) { +Vector3D SCEPCalSegmentation_k4geo::myPosition(const CellID& cID) const { int copyNum = (int)cID; @@ -54,7 +54,7 @@ Vector3D SCEPCalSegmentation::myPosition(const CellID& cID) { double D_PHI_GLOBAL =2*M_PI/PHI_SEGMENTS; double PROJECTIVE_GAP =(N_PROJECTIVE_FILL*nomfw)/2; - double THETA_SIZE_BARREL =atan(EBz/Rin); + // double THETA_SIZE_BARREL =atan(EBz/Rin); double THETA_SIZE_ENDCAP =atan(Rin/EBz); int N_THETA_BARREL =2*floor(EBz/nomfw); int N_THETA_ENDCAP =floor(Rin/nomfw); @@ -69,7 +69,7 @@ Vector3D SCEPCalSegmentation::myPosition(const CellID& cID) { double y0slice_end =r0slice_end*tan(D_THETA_BARREL/2.); double slice_front_jut =y0slice_end*sin(M_PI/2-thC_end); double z1slice =Rin -slice_front_jut; - double z2slice =Rin +Fdz +Rdz +slice_front_jut; + // double z2slice =Rin +Fdz +Rdz +slice_front_jut; double y1slice =z1slice*tan(M_PI/2-THETA_SIZE_ENDCAP) +PROJECTIVE_GAP; double phiTiming =nPhi_in*D_PHI_GLOBAL; double rT =z1slice -2*nomth; @@ -180,13 +180,13 @@ Vector3D SCEPCalSegmentation::myPosition(const CellID& cID) { return Vector3D(0,0,0); } -CellID SCEPCalSegmentation::cellID(const Vector3D& /*localPosition*/, +CellID SCEPCalSegmentation_k4geo::cellID(const Vector3D& /*localPosition*/, const Vector3D& /*globalPosition*/, const VolumeID& vID) const { return setCellID(System(vID), Eta(vID), Phi(vID), Depth(vID) ); } -VolumeID SCEPCalSegmentation::setVolumeID(int System, int Eta, int Phi, int Depth) const { +VolumeID SCEPCalSegmentation_k4geo::setVolumeID(int System, int Eta, int Phi, int Depth) const { VolumeID SystemId = static_cast(System); VolumeID EtaId = static_cast(Eta); VolumeID PhiId = static_cast(Phi); @@ -199,7 +199,7 @@ VolumeID SCEPCalSegmentation::setVolumeID(int System, int Eta, int Phi, int Dept return vID; } -CellID SCEPCalSegmentation::setCellID(int System, int Eta, int Phi, int Depth) const { +CellID SCEPCalSegmentation_k4geo::setCellID(int System, int Eta, int Phi, int Depth) const { VolumeID SystemId = static_cast(System); VolumeID EtaId = static_cast(Eta); VolumeID PhiId = static_cast(Phi); @@ -212,33 +212,33 @@ CellID SCEPCalSegmentation::setCellID(int System, int Eta, int Phi, int Depth) c return vID; } -int SCEPCalSegmentation::System(const CellID& aCellID) const { +int SCEPCalSegmentation_k4geo::System(const CellID& aCellID) const { VolumeID System = static_cast(_decoder->get(aCellID, fSystemId)); return static_cast(System); } -int SCEPCalSegmentation::Eta(const CellID& aCellID) const { +int SCEPCalSegmentation_k4geo::Eta(const CellID& aCellID) const { VolumeID Eta = static_cast(_decoder->get(aCellID, fEtaId)); return static_cast(Eta); } -int SCEPCalSegmentation::Phi(const CellID& aCellID) const { +int SCEPCalSegmentation_k4geo::Phi(const CellID& aCellID) const { VolumeID Phi = static_cast(_decoder->get(aCellID, fPhiId)); return static_cast(Phi); } -int SCEPCalSegmentation::Depth(const CellID& aCellID) const { +int SCEPCalSegmentation_k4geo::Depth(const CellID& aCellID) const { VolumeID Depth = static_cast(_decoder->get(aCellID, fDepthId)); return static_cast(Depth); } -int SCEPCalSegmentation::getLast32bits(const CellID& aCellID) const { +int SCEPCalSegmentation_k4geo::getLast32bits(const CellID& aCellID) const { CellID aId64 = aCellID >> sizeof(int)*CHAR_BIT; int aId32 = (int)aId64; return aId32; } -CellID SCEPCalSegmentation::convertLast32to64(const int aId32) const { +CellID SCEPCalSegmentation_k4geo::convertLast32to64(const int aId32) const { CellID aId64 = (CellID)aId32; aId64 <<= sizeof(int)*CHAR_BIT; return aId64; diff --git a/detectorSegmentations/src/plugins/SegmentationFactories.cpp b/detectorSegmentations/src/plugins/SegmentationFactories.cpp index ea2bc93a6..a7081097f 100644 --- a/detectorSegmentations/src/plugins/SegmentationFactories.cpp +++ b/detectorSegmentations/src/plugins/SegmentationFactories.cpp @@ -42,4 +42,4 @@ DECLARE_SEGMENTATION(FCCSWHCalPhiTheta_k4geo, create_segmentation) #include "detectorSegmentations/SCEPCalSegmentation_k4geo.h" -DECLARE_SEGMENTATION(SCEPCalSegmentation_k4geo, create_segmentation) \ No newline at end of file +DECLARE_SEGMENTATION(SCEPCalSegmentation_k4geo, create_segmentation) diff --git a/plugins/Geant4Output2EDM4hepDRCrystalHit.cpp b/plugins/Geant4Output2EDM4hepDRCrystalHit.cpp index 869dd512e..1840f0f3b 100644 --- a/plugins/Geant4Output2EDM4hepDRCrystalHit.cpp +++ b/plugins/Geant4Output2EDM4hepDRCrystalHit.cpp @@ -461,20 +461,13 @@ void Geant4Output2EDM4hepDRCrystalHit::saveCollection(OutputContext& /* sch.setCellID( hit->cellID ); sch.setPosition( hitpos ); - sch.setEta( hit->eta ); - sch.setPhi( hit->phi ); - sch.setDepth( hit->depth ); - sch.setSystem( hit->system ); sch.setEnergy( hit->energyDeposit/CLHEP::GeV ); - sch.setNcerenkov( hit->ncerenkov ); - sch.setNscintillator( hit->nscintillator ); + sch.setNCerenkovProd( hit->nCerenkovProd ); + sch.setNScintillationProd(hit->nScintillationProd ); - sch.setNwavelen_cer( hit->nwavelen_cer ); - sch.setNwavelen_scint( hit->nwavelen_scint ); - - sch.setNtime_cer( hit->ntime_cer ); - sch.setNtime_scint( hit->ntime_scint ); + sch.setTAvgS( hit->tAvgC ); + sch.setTAvgS( hit->tAvgS ); for(auto ci=hit->truth.begin(); ci != hit->truth.end(); ++ci){ diff --git a/plugins/SCEPCalSDActionDRCrystalHit.cpp b/plugins/SCEPCalSDActionDRCrystalHit.cpp index 0bd01ccfa..127c81298 100644 --- a/plugins/SCEPCalSDActionDRCrystalHit.cpp +++ b/plugins/SCEPCalSDActionDRCrystalHit.cpp @@ -3,14 +3,17 @@ // Princeton University //=============================== #include "detectorSegmentations/DRCrystalHit.h" -#include "detectorSegmentations/SCEPCalSegmentation.h" +#include "detectorSegmentations/SCEPCalSegmentation_k4geo.h" #include "DDG4/Geant4SensDetAction.inl" #include "DDG4/Factories.h" +#include "G4OpticalPhoton.hh" +#include "G4VProcess.hh" +#include namespace SCEPCal { G4double convertEvtoNm(G4double energy) { - return 1239.84187/energy*1000.; + return 1239.84187/energy*1000.; //GeV to nm } class SegmentedCrystalCalorimeterSD_DRHit { public: @@ -45,74 +48,51 @@ namespace dd4hep { Geant4HitCollection* coll = collection(m_collectionID); dd4hep::Segmentation* _geoSeg = &m_segmentation; - auto segmentation=dynamic_cast(_geoSeg->segmentation()); + auto segmentation=dynamic_cast(_geoSeg->segmentation()); auto copyNum64 = segmentation->convertFirst32to64(thePreStepTouchable->GetCopyNumber(0)); int cellID = (int)copyNum64; SegmentedCrystalCalorimeterSD_DRHit::Hit* hit = coll->findByKey(cellID); + if(!hit) { DDSegmentation::Vector3D pos = segmentation->myPosition(copyNum64); Position global(pos.x(),pos.y(),pos.z()); hit = new SegmentedCrystalCalorimeterSD_DRHit::Hit(global); hit->cellID = cellID; - hit->system = segmentation->System(copyNum64); - hit->eta = segmentation->Eta(copyNum64); - hit->phi = segmentation->Phi(copyNum64); - hit->depth = segmentation->Depth(copyNum64); coll->add(cellID, hit); } + G4Track * track = step->GetTrack(); + if(track->GetDefinition()==G4OpticalPhoton::OpticalPhotonDefinition()) { - float wavelength=convertEvtoNm(track->GetTotalEnergy()/eV); - int ibin=-1; - float binsize=(hit->wavelen_max-hit->wavelen_min)/hit->nbins; - ibin = (wavelength-hit->wavelen_min)/binsize; + float avgarrival=(pretime+posttime)/2.; - int jbin=-1; - float tbinsize=(hit->time_max-hit->time_min)/hit->nbins; - jbin = (avgarrival-hit->time_min)/tbinsize; + + // count 1st and kill + // apply scale factor and poisson smearing offline int phstep = track->GetCurrentStepNumber(); if (track->GetCreatorProcess()->G4VProcess::GetProcessName()=="CerenkovPhys") { - std::string amedia = ((track->GetMaterial())->GetName()); - if(amedia.find("Silicon")!=std::string::npos) - { - if(phstep>1) { - hit->ncerenkov+=1; - if(ibin>-1 && ibinnbins) ((hit->nwavelen_cer).at(ibin)) +=1; - if(jbin>-1 && jbinnbins) ((hit->ntime_cer).at(jbin)) +=1; - } - track->SetTrackStatus(fStopAndKill); - } - else { - if(phstep==1) { - hit->ncerenkov+=1; - if(ibin>-1 && ibinnbins) ((hit->nwavelen_cer).at(ibin)) +=1; - if(jbin>-1 && jbinnbins) ((hit->ntime_cer).at(jbin)) +=1; - } + if(phstep==1) { + float tAvgC_new = (((hit->tAvgC)*(hit->nCerenkovProd)) +avgarrival)/(hit->nCerenkovProd+1); + hit->nCerenkovProd+=1; + hit->tAvgC = tAvgC_new; } + track->SetTrackStatus(fStopAndKill); } + else if (track->GetCreatorProcess()->G4VProcess::GetProcessName()=="ScintillationPhys") { - std::string amedia = ((track->GetMaterial())->GetName()); - if(amedia.find("Silicon")!=std::string::npos) - { - if(phstep>1) { - hit->nscintillator+=1; - if((ibin>-1)&&(ibinnbins)) ((hit->nwavelen_scint).at(ibin))+=1; - if(jbin>-1&&jbinnbins) ((hit->ntime_scint).at(jbin))+=1; - } - track->SetTrackStatus(fStopAndKill); - } - else { - if((track->GetCurrentStepNumber()==1)) { - hit->nscintillator+=1; - if((ibin>-1)&&(ibinnbins)) ((hit->nwavelen_scint).at(ibin))+=1; - if(jbin>-1&&jbinnbins) ((hit->ntime_scint).at(jbin))+=1; - } + if(phstep==1) { + float tAvgS_new = (((hit->tAvgS)*(hit->nScintillationProd)) +avgarrival)/(hit->nScintillationProd+1); + hit->nScintillationProd+=1; + hit->tAvgS = tAvgS_new; } + track->SetTrackStatus(fStopAndKill); } + } + hit->truth.emplace_back(contrib); hit->energyDeposit+=edep; mark(h.track);