From a892427beca0b5de979224d5b33ae9c8e796a9d5 Mon Sep 17 00:00:00 2001 From: FU Chengdong <fucd@ihep.ac.cn> Date: Sun, 20 Oct 2024 05:31:38 +0000 Subject: [PATCH] REC: add one of OTK version as patch --- .../CRD_common_v01/OTKEndcap_v01_01.xml | 26 +- .../CRD_common_v02/FTD_SkewRing_v01_09.xml | 50 ++ .../TDR_o1_v01/TDR_o1_v01-patchOTK.xml | 78 +++ Detector/DetCRD/scripts/TDR_o1_v01/sim.py | 3 +- .../DetCRD/scripts/TDR_o1_v01/tracking.py | 8 +- .../Tracker/SiTracker_itkbarrel_v02_geo.cpp | 4 +- .../Tracker/SiTracker_otkendcap_v01_geo.cpp | 185 ++++-- Digitization/DigiSimple/src/PlanarDigiAlg.cpp | 11 +- Digitization/DigiSimple/src/SmearDigiTool.cpp | 11 +- Examples/options/tracking-patchOTK.py | 308 ++++++++++ .../RecSiTracking/src/SiliconTrackingAlg.cpp | 2 +- .../FullLDCTracking/FullLDCTrackingAlg.cpp | 126 ++-- .../src/FullLDCTracking/FullLDCTrackingAlg.h | 2 +- Service/GearSvc/src/GearSvc.cpp | 90 ++- Service/GearSvc/src/GearSvc.h | 1 + Service/TrackSystemSvc/src/MarlinKalTest.cc | 39 +- .../TrackSystemSvc/src/MarlinKalTestTrack.cc | 7 +- .../kaldet/CEPCOTKEndcapKalDetector.h} | 44 +- .../kaldet/CEPCSegmentedDiscMeasLayer.h} | 89 ++- .../KalDet/include/kaldet/ILDDiscMeasLayer.h | 16 +- .../KalDet/include/kaldet/ILDFTDKalDetector.h | 4 +- .../KalDet/include/kaldet/ILDVXDKalDetector.h | 3 +- .../ild/common/CEPCSegmentedDiscMeasLayer.cc | 485 ++++++++++++++++ .../KalDet/src/ild/common/ILDConeMeasLayer.cc | 2 +- .../KalDet/src/ild/common/ILDConeMeasLayer.h | 99 ---- .../KalDet/src/ild/common/ILDCylinderHit.cc | 4 +- .../KalDet/src/ild/common/ILDCylinderHit.h | 38 -- .../src/ild/common/ILDCylinderMeasLayer.cc | 29 +- .../src/ild/common/ILDCylinderMeasLayer.h | 93 --- .../KalDet/src/ild/common/ILDDiscMeasLayer.cc | 104 ++-- .../KalDet/src/ild/common/ILDDiscMeasLayer.h | 101 ---- .../ILDMeasurementSurfaceStoreFiller.cc | 2 +- .../common/ILDMeasurementSurfaceStoreFiller.h | 81 --- .../ild/common/ILDParallelPlanarMeasLayer.cc | 2 +- .../ild/common/ILDParallelPlanarMeasLayer.h | 80 --- .../common/ILDParallelPlanarStripMeasLayer.cc | 4 +- .../common/ILDParallelPlanarStripMeasLayer.h | 62 -- .../KalDet/src/ild/common/ILDPlanarHit.cc | 8 +- .../KalDet/src/ild/common/ILDPlanarHit.h | 41 -- .../src/ild/common/ILDPlanarMeasLayer.cc | 4 +- .../src/ild/common/ILDPlanarMeasLayer.h | 84 --- .../src/ild/common/ILDPlanarStripHit.cc | 8 +- .../KalDet/src/ild/common/ILDPlanarStripHit.h | 42 -- .../ild/common/ILDPolygonBarrelMeasLayer.cc | 2 +- .../ild/common/ILDPolygonBarrelMeasLayer.h | 138 ----- .../src/ild/common/ILDRotatedTrapMeaslayer.cc | 4 +- .../src/ild/common/ILDRotatedTrapMeaslayer.h | 106 ---- .../ild/common/ILDSegmentedDiscMeasLayer.cc | 4 +- .../common/ILDSegmentedDiscStripMeasLayer.cc | 5 +- .../common/ILDSegmentedDiscStripMeasLayer.h | 80 --- .../KalDet/src/ild/common/ILDVMeasLayer.cc | 4 +- .../KalDet/src/ild/common/ILDVMeasLayer.h | 89 --- .../KalDet/src/ild/common/ILDVTrackHit.h | 33 -- .../KalDet/src/ild/common/MaterialDataBase.cc | 15 +- .../KalDet/src/ild/common/MaterialDataBase.h | 76 --- .../src/ild/etd/CEPCOTKEndcapKalDetector.cc | 547 ++++++++++++++++++ .../src/ild/ftd/ILDFTDDiscBasedKalDetector.cc | 3 +- .../src/ild/ftd/ILDFTDDiscBasedKalDetector.h | 44 -- .../KalDet/src/ild/ftd/ILDFTDKalDetector.cc | 3 +- .../KalDet/src/ild/lctpc/LCTPCKalDetector.h | 39 -- .../KalDet/src/ild/set/ILDSETKalDetector.h | 58 -- .../KalDet/src/ild/sit/CEPCITKKalDetector.cc | 7 +- .../src/ild/sit/ILDSITCylinderKalDetector.h | 44 -- .../KalDet/src/ild/sit/ILDSITKalDetector.cc | 9 +- .../KalDet/src/ild/sit/ILDSITKalDetector.h | 60 -- .../src/ild/support/ILDSupportKalDetector.cc | 12 +- .../src/ild/support/ILDSupportKalDetector.h | 37 -- .../KalDet/src/ild/tpc/ILDTPCKalDetector.cc | 11 +- .../KalDet/src/ild/tpc/ILDTPCKalDetector.h | 30 - .../KalDet/src/ild/vxd/CEPCVTXKalDetector.cc | 9 +- .../KalDet/src/ild/vxd/ILDVXDKalDetector.cc | 12 +- .../KalDet/src/ild/vxd/ILDVXDKalDetector.h | 74 --- 72 files changed, 2036 insertions(+), 1999 deletions(-) create mode 100644 Detector/DetCRD/compact/CRD_common_v02/FTD_SkewRing_v01_09.xml create mode 100644 Detector/DetCRD/compact/TDR_o1_v01/TDR_o1_v01-patchOTK.xml create mode 100644 Examples/options/tracking-patchOTK.py rename Utilities/KalDet/{src/ild/ftd/ILDFTDKalDetector.h => include/kaldet/CEPCOTKEndcapKalDetector.h} (62%) rename Utilities/KalDet/{src/ild/common/ILDSegmentedDiscMeasLayer.h => include/kaldet/CEPCSegmentedDiscMeasLayer.h} (59%) create mode 100644 Utilities/KalDet/src/ild/common/CEPCSegmentedDiscMeasLayer.cc delete mode 100644 Utilities/KalDet/src/ild/common/ILDConeMeasLayer.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDCylinderHit.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDCylinderMeasLayer.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDDiscMeasLayer.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDMeasurementSurfaceStoreFiller.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDParallelPlanarMeasLayer.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDParallelPlanarStripMeasLayer.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDPlanarHit.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDPlanarMeasLayer.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDPlanarStripHit.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDPolygonBarrelMeasLayer.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDRotatedTrapMeaslayer.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDSegmentedDiscStripMeasLayer.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDVMeasLayer.h delete mode 100644 Utilities/KalDet/src/ild/common/ILDVTrackHit.h delete mode 100644 Utilities/KalDet/src/ild/common/MaterialDataBase.h create mode 100644 Utilities/KalDet/src/ild/etd/CEPCOTKEndcapKalDetector.cc delete mode 100644 Utilities/KalDet/src/ild/ftd/ILDFTDDiscBasedKalDetector.h delete mode 100644 Utilities/KalDet/src/ild/lctpc/LCTPCKalDetector.h delete mode 100644 Utilities/KalDet/src/ild/set/ILDSETKalDetector.h delete mode 100644 Utilities/KalDet/src/ild/sit/ILDSITCylinderKalDetector.h delete mode 100644 Utilities/KalDet/src/ild/sit/ILDSITKalDetector.h delete mode 100644 Utilities/KalDet/src/ild/support/ILDSupportKalDetector.h delete mode 100644 Utilities/KalDet/src/ild/tpc/ILDTPCKalDetector.h delete mode 100644 Utilities/KalDet/src/ild/vxd/ILDVXDKalDetector.h diff --git a/Detector/DetCRD/compact/CRD_common_v01/OTKEndcap_v01_01.xml b/Detector/DetCRD/compact/CRD_common_v01/OTKEndcap_v01_01.xml index 5f6b2f83..9142b728 100644 --- a/Detector/DetCRD/compact/CRD_common_v01/OTKEndcap_v01_01.xml +++ b/Detector/DetCRD/compact/CRD_common_v01/OTKEndcap_v01_01.xml @@ -21,13 +21,13 @@ <constant name="OTKEndcap_r8" value="1520*mm" /> <constant name="OTKEndcap_r9" value="1660*mm" /> <constant name="OTKEndcap_r10" value="1800*mm" /> + <constant name="OTKEndcap_piece_num" value="24" /> <constant name="OTKEndcap_inner_radius" value="OTKEndcap_r0" /> <constant name="OTKEndcap_outer_radius" value="OTKEndcap_r10" /> <constant name="OTKEndcap_half_length" value="2900*mm" /> <constant name="OTKEndcap_piece_deg" value="7.5" /><!-- variable OTKEndcap_piece_deg needs no dimension --> <constant name="OTKEndcap_dead_deg" value="0.5" /> <constant name="OTKEndcap_dead_thickness" value="1*mm" /> - <constant name="OTKEndcap_piece_num" value="24" /> <constant name="OTKEndcap_asic_num_0" value="5" /> <constant name="OTKEndcap_asic_num_1" value="7" /> <constant name="OTKEndcap_asic_num_2" value="9" /> @@ -48,9 +48,9 @@ <constant name="OTKEndcap_module_num_7" value="2" /> <constant name="OTKEndcap_module_num_8" value="2" /> <constant name="OTKEndcap_module_num_9" value="2" /> - <constant name="OTKEndcap_layer_thickness" value="25*mm" /> - <constant name="OTKEndcap_layer0_zpos" value="2850*mm" /> - <constant name="OTKEndcap_layer1_zpos" value="2875*mm" /> + <constant name="OTKEndcap_layer_thickness" value="15*mm" /> + <constant name="OTKEndcap_layer0_zpos" value="OTKEndcap_half_length+OTKEndcap_layer_thickness/2" /> + <constant name="OTKEndcap_layer1_zpos" value="OTKEndcap_layer0_zpos+OTKEndcap_layer_thickness" /> <constant name="OTKEndcap_support_thickness" value="1*mm" /> <constant name="OTKEndcap_sensor_thickness" value="500*um" /> <constant name="OTKEndcap_sensor_gap" value="0.5*mm" /> @@ -66,14 +66,16 @@ </define> <detectors> - <detector id="DetID_OTKEndcap" limits="otk_limits" name="OTKEndcap" type="SiTracker_otkendcap_v01" vis="OTKEndcapVis" readout="OTKEndcapCollection" insideTrackingVolume="true"> + <detector id="DetID_OTKEndcap" name="OTKEndcap" type="SiTracker_otkendcap_v01" vis="OTKEndcapVis" readout="OTKEndcapCollection" insideTrackingVolume="true" + limits="otk_limits" combineHits="true"> <envelope> - <shape type="BooleanShape" operation="Union" material="Air" > - <shape type="Tube" rmin="OTKEndcap_inner_radius" rmax="OTKEndcap_outer_radius" dz="OTKEndcap_half_length" /> - </shape> + <shape type="BooleanShape" operation="Subtraction" material="Air" > + <shape type="Tube" rmin="OTKEndcap_inner_radius" rmax="OTKEndcap_outer_radius" dz="Ecal_endcap_zmin"/> + <shape type="Tube" rmin="OTKEndcap_inner_radius-env_safety" rmax="OTKEndcap_outer_radius+env_safety" dz="OTKEndcap_half_length"/> + </shape> </envelope> - <display support="WhiteVis" sens_env="SeeThrough" sens="GrayVis" deadsensor="GreenVis" - pcb="GreenVis" asic="YellowVis" dead="BlackVis"/> + <display support="GrayVis" sens_env="SeeThrough" sens="MagentaVis" deadsensor="CyanVis" + pcb="RedVis" asic="OrangeVis" dead="BlackVis"/> <type_flags type="DetType_TRACKER + DetType_ENDCAP + DetType_STRIP "/> <global sensitive_mat="G4_Si" support_mat="G4_C" sensitive_threshold_KeV="64*keV"/> <support thickness="OTKEndcap_support_thickness" inner_radius="OTKEndcap_inner_radius" outer_radius="OTKEndcap_outer_radius" mat="CarbonFiber" dead_mat="epoxy"/> @@ -92,7 +94,7 @@ <readouts> <readout name="OTKEndcapCollection"> - <id>system:5,side:-2,layer:9,module:8,active:8,sensor:8</id> + <id>system:5,side:-2,layer:9,module:8,sensor:8,active:8</id> </readout> </readouts> -</lccdd> \ No newline at end of file +</lccdd> diff --git a/Detector/DetCRD/compact/CRD_common_v02/FTD_SkewRing_v01_09.xml b/Detector/DetCRD/compact/CRD_common_v02/FTD_SkewRing_v01_09.xml new file mode 100644 index 00000000..f7a31383 --- /dev/null +++ b/Detector/DetCRD/compact/CRD_common_v02/FTD_SkewRing_v01_09.xml @@ -0,0 +1,50 @@ +<lccdd> + <define> + <constant name="SiliconThickness" value="0.2*mm"/> + <constant name="SupportThickness" value="1.4925*mm"/> <!--equivalent from carbon fiber to carbon/--> + <constant name="ModuleZGap" value="2.0*mm"/> + <constant name="ModuleRPhiGap" value="-10*mm"/> + <constant name="FTDPetalNumber" value="16"/> + </define> + + <detectors> + <detector id="DetID_FTD" name="FTD" type="SiTrackerSkewRing_v01" vis="FTDVis" readout="FTDCollection" insideTrackingVolume="true" reflect="true"> + <envelope> + <shape type="Assembly"/> + </envelope> + + <type_flags type="DetType_TRACKER + DetType_ENDCAP + DetType_PIXEL "/> + + <reconstruction strip_width="0.05*mm" strip_length="92*mm" strip_pitch="0" strip_angle="0"/> + + <layer id="0" z="SiTracker_endcap_z1" dz="0.5*ModuleZGap" inner_r="SiTracker_endcap_inner_radius1" outer_r="SiTracker_endcap_outer_radius1" + skew="0" phi0="0" gap="ModuleRPhiGap" is_pixel="true" nmodules="FTDPetalNumber" vis="SeeThrough"> + <component material="G4_Si" thickness="SiliconThickness" vis="FTDSensitiveVis" sensitive="yes"/> + <component material="Carbon" thickness="SupportThickness" vis="FTDSupportVis"/> + </layer> + <layer id="1" z="SiTracker_endcap_z2" dz="0.5*ModuleZGap" inner_r="SiTracker_endcap_inner_radius2" outer_r="SiTracker_endcap_outer_radius2" + skew="0" phi0="0" gap="ModuleRPhiGap" is_pixel="true" nmodules="FTDPetalNumber" vis="SeeThrough"> + <component material="G4_Si" thickness="SiliconThickness" vis="FTDSensitiveVis" sensitive="yes"/> + <component material="Carbon" thickness="SupportThickness" vis="FTDSupportVis"/> + </layer> + <layer id="2" z="SiTracker_endcap_z3" dz="0.5*ModuleZGap" inner_r="SiTracker_endcap_inner_radius3" outer_r="SiTracker_endcap_outer_radius3" + skew="0" phi0="0" gap="ModuleRPhiGap" is_pixel="true" nmodules="FTDPetalNumber" vis="SeeThrough"> + <component material="G4_Si" thickness="SiliconThickness" vis="FTDSensitiveVis" sensitive="yes"/> + <component material="Carbon" thickness="SupportThickness" vis="FTDSupportVis"/> + </layer> + <layer id="3" z="SiTracker_endcap_z4" dz="0.5*ModuleZGap" inner_r="SiTracker_endcap_inner_radius4" outer_r="SiTracker_endcap_outer_radius4" + skew="0" phi0="0" gap="ModuleRPhiGap" is_pixel="true" nmodules="FTDPetalNumber" vis="SeeThrough"> + <component material="G4_Si" thickness="SiliconThickness" vis="FTDSensitiveVis" sensitive="yes"/> + <component material="Carbon" thickness="SupportThickness" vis="FTDSupportVis"/> + </layer> + + </detector> + </detectors> + + <readouts> + <readout name="FTDCollection"> + <id>system:5,side:-2,layer:9,module:8,sensor:8</id> + </readout> + </readouts> + +</lccdd> diff --git a/Detector/DetCRD/compact/TDR_o1_v01/TDR_o1_v01-patchOTK.xml b/Detector/DetCRD/compact/TDR_o1_v01/TDR_o1_v01-patchOTK.xml new file mode 100644 index 00000000..47402751 --- /dev/null +++ b/Detector/DetCRD/compact/TDR_o1_v01/TDR_o1_v01-patchOTK.xml @@ -0,0 +1,78 @@ +<?xml version="1.0" encoding="UTF-8"?> +<lccdd xmlns:compact="http://www.lcsim.org/schemas/compact/1.0" + xmlns:xs="http://www.w3.org/2001/XMLSchema" + xs:noNamespaceSchemaLocation="http://www.lcsim.org/schemas/compact/1.0/compact.xsd"> + <info name="TDR_o1_v01" + title="CepC reference detctor for TDR" + author="" + url="http://cepc.ihep.ac.cn" + status="developing" + version="v01"> + <comment>CepC reference detector simulation models used for TDR </comment> + </info> + + <includes> + <gdmlFile ref="${DD4hepINSTALL}/DDDetectors/compact/elements.xml"/> + <gdmlFile ref="../CRD_common_v02/materials.xml"/> + </includes> + + <define> + <constant name="world_size" value="10*m"/> + <constant name="world_x" value="world_size"/> + <constant name="world_y" value="world_size"/> + <constant name="world_z" value="world_size"/> + + <include ref="${DD4hepINSTALL}/DDDetectors/compact/detector_types.xml"/> + </define> + + <include ref="./TDR_Dimensions_v01_01.xml"/> + + <include ref="../CRD_common_v02/Beampipe_v01_04.xml"/> + <!--preliminary vertex and tracker, to update/--> + <!--include ref="../CRD_common_v02/VXD_StaggeredLadder_v02_01.xml"/--> + <include ref="../CRD_common_v02/VXD_Composite_v01_01.xml"/> + <include ref="../CRD_common_v02/FTD_SkewRing_v01_09.xml"/> + <!--include ref="../CRD_common_v02/SIT_SimplePixel_v01_04.xml"/--> + <include ref="../CRD_common_v02/SIT_StaggeredStave_v02.xml"/> + <!--include ref="../CRD_common_v01/TPC_Simple_v10_02.xml"/--> + <!--use 10 rows clustering version/--> + <include ref="../CRD_common_v02/TPC_ModularEndcap_o1_v02.xml"/> + <!--include ref="../CRD_common_v01/SET_SimplePixel_v01_01.xml"/--> + <include ref="../CRD_common_v01/OTKBarrel_v01_01.xml"/> + <include ref="../CRD_common_v01/OTKEndcap_v01_01.xml"/> + + <include ref="../CRD_common_v01/Ecal_Crystal_Barrel_v01_02.xml"/> + <include ref="../CRD_common_v01/Ecal_Crystal_Endcap_v01_01.xml"/> + <include ref="../CRD_common_v01/SHcalGlass_Barrel_v05.xml"/> + <include ref="../CRD_common_v01/SHcalGlass_Endcaps_v01.xml"/> + + <!--Lumical to update--> + <include ref="../CRD_common_v01/Lumical_o1_v01.xml"/> + <!--preliminary Magnet, to update/--> + <include ref="../CRD_common_v02/Coil_Simple_v01_02.xml"/> + <!--preliminary Muon, obselete/--> + <!--include ref="../CRD_common_v02/Yoke_Polyhedra_Barrel_v01_01.xml"/> + <include ref="../CRD_common_v02/Yoke_Polyhedra_Endcaps_v01_01.xml"/--> + + <!--muon detector--> + <include ref="../CRD_common_v01/Muon_Barrel_v01_01.xml"/> + <include ref="../CRD_common_v01/Muon_Endcap_v01_02.xml"/> + + <fields> + <field name="InnerSolenoid" type="solenoid" + inner_field="Field_nominal_value" + outer_field="0" + zmax="SolenoidCoil_half_length" + inner_radius="SolenoidCoil_center_radius" + outer_radius="Solenoid_outer_radius"> + </field> + <field name="OuterSolenoid" type="solenoid" + inner_field="0" + outer_field="Field_outer_nominal_value" + zmax="SolenoidCoil_half_length" + inner_radius="Solenoid_outer_radius" + outer_radius="Yoke_barrel_inner_radius"> + </field> + </fields> + +</lccdd> diff --git a/Detector/DetCRD/scripts/TDR_o1_v01/sim.py b/Detector/DetCRD/scripts/TDR_o1_v01/sim.py index 8e871ed1..575297e7 100644 --- a/Detector/DetCRD/scripts/TDR_o1_v01/sim.py +++ b/Detector/DetCRD/scripts/TDR_o1_v01/sim.py @@ -17,6 +17,7 @@ rndmgensvc.Engine = rndmengine.name() # option for standalone tracker study #geometry_option = "TDR_o1_v01/TDR_o1_v01-oldVersion.xml" +#geometry_option = "TDR_o1_v01/TDR_o1_v01-patchOTK.xml" geometry_option = "TDR_o1_v01/TDR_o1_v01.xml" if not os.getenv("DETCRDROOT"): @@ -97,7 +98,7 @@ from Configurables import ApplicationMgr mgr = ApplicationMgr( TopAlg = [genalg, detsimalg, out], EvtSel = 'NONE', - EvtMax = 5, + EvtMax = 50, ExtSvc = [rndmengine, rndmgensvc, dsvc, geosvc], HistogramPersistency = 'ROOT', OutputLevel = ERROR diff --git a/Detector/DetCRD/scripts/TDR_o1_v01/tracking.py b/Detector/DetCRD/scripts/TDR_o1_v01/tracking.py index 5e7ac32a..8c009a18 100644 --- a/Detector/DetCRD/scripts/TDR_o1_v01/tracking.py +++ b/Detector/DetCRD/scripts/TDR_o1_v01/tracking.py @@ -248,12 +248,14 @@ full.FTDRawHits = ftdhitname full.TPCTracks = "ClupatraTracks" # add standalone TPC track full.SiTracks = "SubsetTracks" full.OutputTracks = "CompleteTracks" # default name -full.VTXHitToTrackDistance = 5. +#full.VTXHitToTrackDistance = 5. full.FTDHitToTrackDistance = 5. full.SITHitToTrackDistance = 3. full.SETHitToTrackDistance = 5. full.MinChi2ProbForSiliconTracks = 0 -full.MaxChi2PerHit = 500 +full.MaxChi2PerHit = 200 +#full.ForceSiTPCMerging = True +#full.ForceTPCSegmentsMerging = True #full.OutputLevel = DEBUG from Configurables import TPCDndxAlg @@ -287,7 +289,7 @@ from Configurables import ApplicationMgr mgr = ApplicationMgr( TopAlg = [podioinput, digiVXD, digiSIT, digiSET, digiFTD, digiTPC, digiMuon, tracking, forward, subset, clupatra, full, tpr, tpc_dndx, tmt, out], EvtSel = 'NONE', - EvtMax = 5, + EvtMax = 50, ExtSvc = [rndmengine, rndmgensvc, dsvc, evtseeder, geosvc, gearsvc, tracksystemsvc, pidsvc], HistogramPersistency = 'ROOT', OutputLevel = ERROR diff --git a/Detector/DetCRD/src/Tracker/SiTracker_itkbarrel_v02_geo.cpp b/Detector/DetCRD/src/Tracker/SiTracker_itkbarrel_v02_geo.cpp index 51bd66e5..037fe844 100644 --- a/Detector/DetCRD/src/Tracker/SiTracker_itkbarrel_v02_geo.cpp +++ b/Detector/DetCRD/src/Tracker/SiTracker_itkbarrel_v02_geo.cpp @@ -734,13 +734,13 @@ static dd4hep::Ref_t create_element(dd4hep::Detector& theDetector, xml_h e, dd4h Layer.thicknessSupport = support_thickness / 2.0; Layer.offsetSupport = stave_radius*sin(stave_phi0);//-stave_offset; Layer.widthSupport = support_width; - Layer.zHalfSupport = support_half_length; + Layer.zHalfSupport = module_length*n_modules_per_stave/2.0;//support_half_length; Layer.distanceSensitive = stave_radius*cos(stave_phi0) + StaveSupportenv_start_height + support_thickness; //sensitive_radius + tube_outer_radius*2. + support_thickness; Layer.thicknessSensitive = sensor_thickness;//module_thickness; Layer.offsetSensitive = stave_radius*sin(stave_phi0);//-stave_offset;//stave_offset/2.0; Layer.widthSensitive = module_width; - Layer.zHalfSensitive = support_half_length; + Layer.zHalfSensitive = module_length*n_modules_per_stave/2.0;//support_half_length; zPlanarData->layers.push_back(Layer); } diff --git a/Detector/DetCRD/src/Tracker/SiTracker_otkendcap_v01_geo.cpp b/Detector/DetCRD/src/Tracker/SiTracker_otkendcap_v01_geo.cpp index e8d2be84..080f0dd9 100644 --- a/Detector/DetCRD/src/Tracker/SiTracker_otkendcap_v01_geo.cpp +++ b/Detector/DetCRD/src/Tracker/SiTracker_otkendcap_v01_geo.cpp @@ -4,6 +4,7 @@ #include "DDRec/Surface.h" #include "DDRec/DetectorData.h" #include "XML/Utilities.h" +#include "Math/AxisAngle.h" #include <cmath> using namespace std; @@ -19,9 +20,89 @@ using dd4hep::Rotation3D; using dd4hep::Volume; using dd4hep::_toString; using dd4hep::rec::volSurfaceList; -using dd4hep::rec::ZPlanarData; +using dd4hep::rec::ZDiskPetalsData; using dd4hep::mm; +namespace dd4hep { + namespace rec { + class VolPlectaneImpl : public VolSurfaceBase { + + public: + /// default c'tor + VolPlectaneImpl() : VolSurfaceBase() { } + + /** The standard constructor. The origin vector points to the origin of the coordinate system on the cylinder, + * its rho defining the radius of the cylinder. The measurement direction v is set to be (0,0,1), the normal is + * chosen to be parallel to the origin vector and u = n X v. + */ + VolPlectaneImpl(Volume vol, SurfaceType type, double thickness_inner, double thickness_outer, Vector3D n, Vector3D o) : + VolSurfaceBase(type, thickness_inner, thickness_outer, Vector3D(), Vector3D(), n, o, vol, 0) { + Vector3D o_rphi(o.x(), o.y(), 0.); + Vector3D v_val = o_rphi.unit(); + Vector3D u_val = v_val.cross(n); + + setU(u_val); + setV(v_val); + + _type.setProperty(SurfaceType::Plane , false); + _type.setProperty(SurfaceType::Cylinder , true); + _type.setProperty(SurfaceType::Cone , false); + _type.checkParallelToZ(*this); + _type.checkOrthogonalToZ(*this); + } + + /** First direction of measurement U - rotated to point projected onto the plectane. + * No check is done whether the point actually is on the surface + */ + virtual Vector3D u(const Vector3D& point = Vector3D()) const { + Vector3D o_rphi(point.x(), point.y(), 0.); + Vector3D v_val = o_rphi.unit(); + Vector3D u_val = v_val.cross(normal()); + return u_val; + } + + /** First direction of measurement V - rotated to point projected onto the plectane. + * No check is done whether the point actually is on the surface + */ + virtual Vector3D v(const Vector3D& point = Vector3D()) const { + Vector3D o_rphi(point.x(), point.y(), 0.); + Vector3D v_val = o_rphi.unit(); + return v_val; + } + + /** Distance to surface */ + virtual double distance(const Vector3D& point) const { + return (point - origin()) * normal(); + } + + /** Convert the global position to the local position (u,v) on the surface - v runs along the axis of radial, u is r*phi */ + virtual Vector2D globalToLocal(const Vector3D& point) const { + double phi = point.phi() - origin().phi() ; + while( phi < -M_PI ) phi += 2.*M_PI; + while( phi > M_PI ) phi -= 2.*M_PI; + + return Vector2D(point.rho()*phi, point.rho() - origin().rho()); + } + + /** Convert the local position (u,v) on the surface to the global position - v runs along the axis of radial, u is r*phi*/ + virtual Vector3D localToGlobal(const Vector2D& point) const { + double z = origin().z(); + double rho = point.v() + origin().rho(); + double phi = point.u()/rho + origin().phi(); + while( phi < -M_PI ) phi += 2.*M_PI; + while( phi > M_PI ) phi -= 2.*M_PI; + + return Vector3D(rho, phi, z, Vector3D::cylindrical); + } + }; + + class VolPlectane : public VolSurface { + public: + VolPlectane(Volume vol, SurfaceType typ_val, double thickness_inner, double thickness_outer, Vector3D n, Vector3D o) : + VolSurface(new VolPlectaneImpl(vol, typ_val, thickness_inner, thickness_outer, n, o)) {} + }; + } +} static dd4hep::Ref_t create_element(dd4hep::Detector& theDetector, xml_h e, dd4hep::SensitiveDetector sens) { @@ -45,13 +126,13 @@ static dd4hep::Ref_t create_element(dd4hep::Detector& theDetector, xml_h e, dd4h } std::cout << " ** building SiTracker_otkendcap ..." << std::endl ; - dd4hep::rec::ZPlanarData* zPlanarData = new dd4hep::rec::ZPlanarData; + dd4hep::rec::ZDiskPetalsData* zDiskPetalsData = new dd4hep::rec::ZDiskPetalsData; //*****************************************************************// // Reading parameters from the .xml file //*****************************************************************// - + bool reflect = x_det.reflect(true); //fetch the geometry parameters const double inner_radius = theDetector.constant<double>("OTKEndcap_inner_radius"); const double outer_radius = theDetector.constant<double>("OTKEndcap_outer_radius"); @@ -292,7 +373,8 @@ static dd4hep::Ref_t create_element(dd4hep::Detector& theDetector, xml_h e, dd4h //create and place utilities in each ring std::vector<dd4hep::PlacedVolume> sensor_pv; - std::vector<dd4hep::rec::VolPlane> sensor_surf; + //std::vector<dd4hep::rec::VolPlane> sensor_surf; + std::vector<dd4hep::rec::VolPlectane> sensor_surf; dd4hep::rec::Vector3D o(0., 0., 0.); dd4hep::rec::Vector3D u(0., 1., 0.); dd4hep::rec::Vector3D v(1., 0., 0.); @@ -309,6 +391,7 @@ static dd4hep::Ref_t create_element(dd4hep::Detector& theDetector, xml_h e, dd4h double ring_asic_angle; double ring_asic_interval; double asic_mid_angle; + int sensor_number_per_piece(0); for(int i=0;i<10;i++){ ring_inner_radius = ring_outer_radius + sensor_dead_gap * 2; ring_outer_radius = r[i+1] - sensor_dead_gap; @@ -364,11 +447,17 @@ static dd4hep::Ref_t create_element(dd4hep::Detector& theDetector, xml_h e, dd4h SensorLogical.setSensitiveDetector(sens); SensorLogical.setVisAttributes(theDetector.visAttributes(sensVis)); if (x_det.hasAttr(_U(limits))) SensorLogical.setLimitSet(theDetector, x_det.limitsStr()); - dd4hep::rec::VolPlane surfsens( SensorLogical, + /* dd4hep::rec::VolPlane surfsens( SensorLogical, dd4hep::rec::SurfaceType(dd4hep::rec::SurfaceType::Sensitive), inner_thick, outer_thick, - u,v,n,o); + u,v,n,o);*/ + dd4hep::rec::Vector3D ocyl((ring_inner_radius+ring_outer_radius)/2.0, 0., 0.); + dd4hep::rec::VolPlectane surfsens(SensorLogical, dd4hep::rec::SurfaceType(dd4hep::rec::SurfaceType::Sensitive), + 0.5*sensor_thickness, 0.5*sensor_thickness, n, ocyl); + //surfsens.setU(u); + //surfsens.setV(v); + //surfsens.setNormal(n); sensor_surf.push_back(surfsens); DeadSensorLogicalA.setVisAttributes(theDetector.visAttributes(deadsensVis)); DeadSensorLogicalB.setVisAttributes(theDetector.visAttributes(deadsensVis)); @@ -382,6 +471,8 @@ static dd4hep::Ref_t create_element(dd4hep::Detector& theDetector, xml_h e, dd4h pv = SensorEnvLogical.placeVolume(DeadSensorLogicalC, Position(0, 0, -(pcb_thickness + asic_thickness) / 2.0)); pv = SensorEnvLogical.placeVolume(DeadSensorLogicalD, Position(0, 0, -(pcb_thickness + asic_thickness) / 2.0)); + sensor_number_per_piece++; + if(ring_module_number==2){ dd4hep::Tube SensorSolid( ring_inner_radius, ring_outer_radius, @@ -411,11 +502,8 @@ static dd4hep::Ref_t create_element(dd4hep::Detector& theDetector, xml_h e, dd4h SensorLogical.setSensitiveDetector(sens); SensorLogical.setVisAttributes(theDetector.visAttributes(sensVis)); if (x_det.hasAttr(_U(limits))) SensorLogical.setLimitSet(theDetector, x_det.limitsStr()); - dd4hep::rec::VolPlane surfsens( SensorLogical , - dd4hep::rec::SurfaceType(dd4hep::rec::SurfaceType::Sensitive), - inner_thick, - outer_thick, - u,v,n,o); + dd4hep::rec::VolPlectane surfsens(SensorLogical, dd4hep::rec::SurfaceType(dd4hep::rec::SurfaceType::Sensitive), + 0.5*sensor_thickness, 0.5*sensor_thickness, n, ocyl); sensor_surf.push_back(surfsens); DeadSensorLogicalA.setVisAttributes(theDetector.visAttributes(sensVis)); DeadSensorLogicalB.setVisAttributes(theDetector.visAttributes(sensVis)); @@ -424,6 +512,8 @@ static dd4hep::Ref_t create_element(dd4hep::Detector& theDetector, xml_h e, dd4h sensor_pv.push_back(pv); pv = SensorEnvLogical.placeVolume(DeadSensorLogicalA, Position(0, 0, -(pcb_thickness + asic_thickness) / 2.0)); pv = SensorEnvLogical.placeVolume(DeadSensorLogicalB, Position(0, 0, -(pcb_thickness + asic_thickness) / 2.0)); + + sensor_number_per_piece++; } //create and place pcb logical volume @@ -484,9 +574,8 @@ static dd4hep::Ref_t create_element(dd4hep::Detector& theDetector, xml_h e, dd4h // Assembling //*****************************************************************// - + float rot = layer_id*0.5; for(int i=0;i<piece_number; i++){ - float rot = layer_id*0.5; std::stringstream piece_enum; piece_enum << "otkendcap_piece_" << layer_id << "_" << i; DetElement pieceDE(layerDE, piece_enum.str(), x_det.id()); @@ -497,14 +586,14 @@ static dd4hep::Ref_t create_element(dd4hep::Detector& theDetector, xml_h e, dd4h for(int iring=0;iring<10;iring++){ sensor_num = sensor_num + module_num[iring]; } - for(int isensor=0;isensor<sensor_num;++isensor){ + for(int isensor=0;isensor<sensor_number_per_piece;++isensor){ std::stringstream sensor_str; sensor_str << piece_enum.str() << "_" << isensor; DetElement sensorDE(pieceDE, sensor_str.str(), x_det.id()); sensorDE.setPlacement(sensor_pv[isensor]); volSurfaceList(sensorDE)->push_back(sensor_surf[isensor]); } - Transform3D trA ( RotationZYX(deg_interval*(i+rot), + Transform3D trA ( RotationZYX(-deg_interval*(i+rot)+180*dd4hep::degree-deg-dead_deg, 180*dd4hep::degree, 0.), Position(0., @@ -516,37 +605,53 @@ static dd4hep::Ref_t create_element(dd4hep::Detector& theDetector, xml_h e, dd4h Position(0., 0., -layer_zpos)); + pv = layer_assembly.placeVolume(PieceEnvLogical,trA); - pv.addPhysVolID("module", i*2); - pieceDE.setPlacement(pv); - pv = layer_assembly.placeVolume(PieceEnvLogical,trB); - pv.addPhysVolID("module", i*2+1); + pv.addPhysVolID("module", i).addPhysVolID("side", +1); pieceDE.setPlacement(pv); - std::cout << piece_enum.str() << " done." << std::endl; + + if(reflect){ + pv = layer_assembly.placeVolume(PieceEnvLogical,trB); + pv.addPhysVolID("module", i).addPhysVolID("side", -1); + DetElement r_pieceDE(layerDE, piece_enum.str()+"_neg", x_det.id()); + r_pieceDE.setPlacement(pv); + for(int isensor=0;isensor<sensor_number_per_piece;++isensor){ + std::stringstream sensor_str; + sensor_str << piece_enum.str() << "_" << isensor; + DetElement r_sensorDE(r_pieceDE, sensor_str.str(), x_det.id()); + r_sensorDE.setPlacement(sensor_pv[isensor]); + volSurfaceList(r_sensorDE)->push_back(sensor_surf[isensor]); + } + } } // package the reconstruction data - dd4hep::rec::ZPlanarData::LayerLayout otkendcapLayer; - - otkendcapLayer.ladderNumber = piece_number; - otkendcapLayer.phi0 = 0.; - otkendcapLayer.sensorsPerLadder = 15; - otkendcapLayer.lengthSensor = r1-r0-sensor_dead_gap*2; - otkendcapLayer.distanceSupport = support_thickness/2.0; - otkendcapLayer.thicknessSupport = support_thickness/2.0; - otkendcapLayer.offsetSupport = 0.; - otkendcapLayer.widthSupport = support_inner_radius; - otkendcapLayer.zHalfSupport = support_outer_radius; - otkendcapLayer.distanceSensitive = support_thickness; - otkendcapLayer.thicknessSensitive = sensor_thickness; - otkendcapLayer.offsetSensitive = 0.; - otkendcapLayer.widthSensitive = 0.; - otkendcapLayer.zHalfSensitive = 0.; - - zPlanarData->layers.push_back(otkendcapLayer); + dd4hep::rec::ZDiskPetalsData::LayerLayout thisLayer; + thisLayer.typeFlags[ dd4hep::rec::ZDiskPetalsData::SensorType::DoubleSided ] = false; + thisLayer.typeFlags[ dd4hep::rec::ZDiskPetalsData::SensorType::Pixel ] = true; + thisLayer.petalHalfAngle = (deg+dead_deg)/2; + thisLayer.alphaPetal = 0; + thisLayer.zPosition = layer_zpos + layer_thickness/2.0 - support_thickness - sensor_thickness/2.0; + thisLayer.petalNumber = piece_number; + thisLayer.sensorsPerPetal = sensor_number_per_piece; + thisLayer.phi0 = 0.5*(deg+dead_deg) + deg_interval*rot; + thisLayer.zOffsetSupport = 0; + thisLayer.distanceSupport = r0; + thisLayer.thicknessSupport = support_thickness; + thisLayer.widthInnerSupport = 0; + thisLayer.widthOuterSupport = 0; + thisLayer.lengthSupport = r10 - r0; + thisLayer.zOffsetSensitive = 0; + thisLayer.distanceSensitive = r0; + thisLayer.thicknessSensitive = sensor_thickness; + thisLayer.widthInnerSensitive = 0; + thisLayer.widthOuterSensitive = 0; + thisLayer.lengthSensitive = r10 - r0; + + zDiskPetalsData->layers.push_back(thisLayer); } - std::cout << (*zPlanarData) << endl; - otkendcap.addExtension< ZPlanarData >(zPlanarData); + std::cout << (*zDiskPetalsData) << endl; + otkendcap.addExtension< ZDiskPetalsData >(zDiskPetalsData); if ( x_det.hasAttr(_U(combineHits)) ) { otkendcap.setCombineHits(x_det.attr<bool>(_U(combineHits)),sens); } diff --git a/Digitization/DigiSimple/src/PlanarDigiAlg.cpp b/Digitization/DigiSimple/src/PlanarDigiAlg.cpp index 29083e3f..3632c210 100644 --- a/Digitization/DigiSimple/src/PlanarDigiAlg.cpp +++ b/Digitization/DigiSimple/src/PlanarDigiAlg.cpp @@ -155,10 +155,11 @@ StatusCode PlanarDigiAlg::execute() encoder.setValue(SimTHit.getCellID()) ; det_id = encoder[lcio::ILDCellID0::subdet] ; - if ( det_id == lcio::ILDDetID::VXD ){} - else if( det_id == lcio::ILDDetID::SIT ){} - else if( det_id == lcio::ILDDetID::SET ){} - else if( det_id == lcio::ILDDetID::FTD ){} + if ( det_id == CEPCConf::DetID::VXD ){} + else if( det_id == CEPCConf::DetID::SIT ){} + else if( det_id == CEPCConf::DetID::SET ){} + else if( det_id == CEPCConf::DetID::FTD ){} + else if( det_id == CEPCConf::DetID::ETD ){} else { fatal() << "unsupported detector ID = " << det_id << " CellID = " << SimTHit.getCellID() << ": file " << __FILE__ << " line " << __LINE__ @@ -381,7 +382,7 @@ StatusCode PlanarDigiAlg::execute() } if( _isStrip || (resU!=0&&resV==0) ){ - trkHit.setType( UTIL::set_bit( trkHit.getType() , UTIL::ILDTrkHitTypeBit::ONE_DIMENSIONAL ) ) ; + trkHit.setType( UTIL::set_bit( trkHit.getType() , CEPCConf::TrkHitTypeBit::ONE_DIMENSIONAL ) ) ; } trkHit.setEDep( SimTHit.getEDep() ); trkHit.setTime( SimTHit.getTime() + gsl_ran_gaussian(_rng, resT) ); diff --git a/Digitization/DigiSimple/src/SmearDigiTool.cpp b/Digitization/DigiSimple/src/SmearDigiTool.cpp index e3fa7c00..a99c6ac4 100644 --- a/Digitization/DigiSimple/src/SmearDigiTool.cpp +++ b/Digitization/DigiSimple/src/SmearDigiTool.cpp @@ -89,7 +89,7 @@ StatusCode SmearDigiTool::Call(edm4hep::SimTrackerHit simhit, edm4hep::TrackerHi auto cellId = simhit.getCellID(); int system = m_decoder->get(cellId, "system"); - int side = m_decoder->get(cellId, "side" ); + int side = m_decoder->get(cellId, "side" ); int layer = m_decoder->get(cellId, "layer" ); int module = m_decoder->get(cellId, "module"); int sensor = m_decoder->get(cellId, "sensor"); @@ -176,10 +176,11 @@ StatusCode SmearDigiTool::Call(edm4hep::SimTrackerHit simhit, edm4hep::TrackerHi dd4hep::rec::Vector2D localPoint = surface->globalToLocal(oldPos); dd4hep::rec::Vector3D local3D(localPoint.u(), localPoint.v(), 0); // A small check, if the hit is in the boundaries: - if (!surface->insideBounds(oldPos, 1e100)) { - error() << " global: (" << oldPos.x() << " " << oldPos.y() << " " << oldPos.z() - << ") local: (" << localPoint.u() << ", " << localPoint.v() << " )" - << " is not within boundaries. Hit is skipped." << endmsg; + if (!surface->insideBounds(oldPos)) { + warning() << " global: (" << oldPos.x() << " " << oldPos.y() << " " << oldPos.z() + << ") local: (" << localPoint.u() << ", " << localPoint.v() << " )" + << " step: " << simhit.getPathLength() + << " is not within boundaries. Hit is skipped." << endmsg; return StatusCode::SUCCESS; } dd4hep::rec::Vector3D globalPointSmeared;//CLHEP::Hep3Vector globalPoint(pos[0],pos[1],pos[2]); diff --git a/Examples/options/tracking-patchOTK.py b/Examples/options/tracking-patchOTK.py new file mode 100644 index 00000000..77b345c6 --- /dev/null +++ b/Examples/options/tracking-patchOTK.py @@ -0,0 +1,308 @@ +#!/usr/bin/env python +import os +from Gaudi.Configuration import * + +from Configurables import k4DataSvc +dsvc = k4DataSvc("EventDataSvc", input="sim00.root") + +from Configurables import RndmGenSvc, HepRndm__Engine_CLHEP__RanluxEngine_ +seed = [12340] +# rndmengine = HepRndm__Engine_CLHEP__RanluxEngine_() # The default engine in Gaudi +rndmengine = HepRndm__Engine_CLHEP__HepJamesRandom_("RndmGenSvc.Engine") # The default engine in Geant4 +rndmengine.SetSingleton = True +rndmengine.Seeds = seed + +rndmgensvc = RndmGenSvc("RndmGenSvc") +rndmgensvc.Engine = rndmengine.name() + +geometry_option = "TDR_o1_v01/TDR_o1_v01-patchOTK.xml" + +if not os.getenv("DETCRDROOT"): + print("Can't find the geometry. Please setup envvar DETCRDROOT." ) + sys.exit(-1) + +geometry_path = os.path.join(os.getenv("DETCRDROOT"), "compact", geometry_option) +if not os.path.exists(geometry_path): + print("Can't find the compact geometry file: %s"%geometry_path) + sys.exit(-1) + +from Configurables import DetGeomSvc +geosvc = DetGeomSvc("GeomSvc") +geosvc.compact = geometry_path + +from Configurables import MarlinEvtSeeder +evtseeder = MarlinEvtSeeder("EventSeeder") + +from Configurables import GearSvc +gearsvc = GearSvc("GearSvc") + +from Configurables import TrackSystemSvc +tracksystemsvc = TrackSystemSvc("TrackSystemSvc") + +from Configurables import SimplePIDSvc +pidsvc = SimplePIDSvc("SimplePIDSvc") +cepcswdatatop = "/cvmfs/cepcsw.ihep.ac.cn/prototype/releases/data/latest" +pidsvc.ParFile = os.path.join(cepcswdatatop, "CEPCSWData/offline-data/Service/SimplePIDSvc/data/dNdx_TPC.root") + +from Configurables import PodioInput +podioinput = PodioInput("PodioReader", collections=[ +# "EventHeader", + "MCParticle", + "VXDCollection", + "SITCollection", + "TPCCollection", +# "SETCollection", + "OTKBarrelCollection", + "FTDCollection", + "OTKEndcapCollection", + "MuonBarrelCollection" + ]) + + +################## +# Digitization +################## + +## Config ## +vxdhitname = "VXDTrackerHits" +sithitname = "SITTrackerHits" +gashitname = "TPCTrackerHits" +sethitname = "OTKBarrelTrackerHits" +ftdhitname = "FTDTrackerHits" +etdhitname = "OTKEndcapTrackerHits" +from Configurables import SmearDigiTool +vxdtool = SmearDigiTool("VXD") +vxdtool.ResolutionU = [0.004, 0.004, 0.004, 0.004, 0.004, 0.004] +vxdtool.ResolutionV = [0.004, 0.004, 0.004, 0.004, 0.004, 0.004] +#vxdtool.OutputLevel = DEBUG + +otketool = SmearDigiTool("OTKE") +otketool.DetName = "OTKEndcap" +otketool.Readout = "OTKEndcapCollection" +otketool.ResolutionU = [0.010] +otketool.ResolutionV = [1.000] +otketool.OutputLevel = DEBUG + +## VXD ## +from Configurables import SiTrackerDigiAlg +digiVXD = SiTrackerDigiAlg("VXDDigi") +digiVXD.SimTrackHitCollection = "VXDCollection" +digiVXD.TrackerHitCollection = vxdhitname +digiVXD.TrackerHitAssociationCollection = "VXDTrackerHitAssociation" +digiVXD.DigiTool = "SmearDigiTool/VXD" +#digiVXD.OutputLevel = DEBUG + +## SIT ## +from Configurables import PlanarDigiAlg +digiSIT = PlanarDigiAlg("SITDigi") +digiSIT.IsStrip = False +digiSIT.SimTrackHitCollection = "SITCollection" +digiSIT.TrackerHitCollection = sithitname +digiSIT.TrackerHitAssociationCollection = "SITTrackerHitAssociation" +digiSIT.ResolutionU = [0.0098] +digiSIT.ResolutionV = [0.0433] +digiSIT.UsePlanarTag = True +digiSIT.ParameterizeResolution = False +digiSIT.ParametersU = [2.29655e-03, 0.965899e-03, 0.584699e-03, 17.0856, 84.566, 12.4695e-03, -0.0643059, 0.168662, 1.87998e-03, 0.514452] +digiSIT.ParametersV = [1.44629e-02, 2.20108e-03, 1.03044e-02, 4.39195e+00, 3.29641e+00, 1.55167e+18, -5.41954e+01, 5.72986e+00, -6.80699e-03, 5.04095e-01] +#digiSIT.OutputLevel = DEBUG + +## SET ## +digiSET = PlanarDigiAlg("SETDigi") +digiSET.IsStrip = False +digiSET.SimTrackHitCollection = "OTKBarrelCollection" +digiSET.TrackerHitCollection = sethitname +digiSET.TrackerHitAssociationCollection = "OTKBarrelTrackerHitAssociation" +digiSET.ResolutionU = [0.010] +digiSET.ResolutionV = [1.000] +digiSET.UsePlanarTag = True +digiSET.ParameterizeResolution = False +digiSET.ParametersU = [2.29655e-03, 0.965899e-03, 0.584699e-03, 17.0856, 84.566, 12.4695e-03, -0.0643059, 0.168662, 1.87998e-03, 0.514452] +digiSET.ParametersV = [1.44629e-02, 2.20108e-03, 1.03044e-02, 4.39195e+00, 3.29641e+00, 1.55167e+18, -5.41954e+01, 5.72986e+00, -6.80699e-03, 5.04095e-01] +#digiSET.OutputLevel = DEBUG + + +## FTD ## +digiFTD = PlanarDigiAlg("FTDDigi") +digiFTD.IsStrip = False +digiFTD.SimTrackHitCollection = "FTDCollection" +digiFTD.TrackerHitCollection = ftdhitname +digiFTD.TrackerHitAssociationCollection = "FTDTrackerHitAssociation" +digiFTD.ResolutionU = [0.0072] +digiFTD.ResolutionV = [0.086] +digiFTD.UsePlanarTag = True +digiFTD.ParameterizeResolution = False +digiFTD.ParametersU = [2.29655e-03, 0.965899e-03, 0.584699e-03, 17.0856, 84.566, 12.4695e-03, -0.0643059, 0.168662, 1.87998e-03, 0.514452] +digiFTD.ParametersV = [1.44629e-02, 2.20108e-03, 1.03044e-02, 4.39195e+00, 3.29641e+00, 1.55167e+18, -5.41954e+01, 5.72986e+00, -6.80699e-03, 5.04095e-01] +#digiFTD.OutputLevel = DEBUG + +## OTKEndcap ## +digiOTKE = SiTrackerDigiAlg("OTKEDigi") +digiOTKE.SimTrackHitCollection = "OTKEndcapCollection" +digiOTKE.TrackerHitCollection = etdhitname +digiOTKE.TrackerHitAssociationCollection = "OTKEndcapTrackerHitAssociation" +digiOTKE.DigiTool = "SmearDigiTool/OTKE" +#digiOTKE.OutputLevel = DEBUG + +## TPC ## +from Configurables import TPCDigiAlg +digiTPC = TPCDigiAlg("TPCDigi") +digiTPC.TPCCollection = "TPCCollection" +digiTPC.TPCLowPtCollection = "TPCLowPtCollection" +digiTPC.TPCTrackerHitsCol = gashitname +#default value, modify them according to future Garfield simulation results +#digiTPC.PixelClustering = True +#digiTPC.PointResolutionRPhi = 0.144 +#digiTPC.DiffusionCoeffRPhi = 0.0323 +#digiTPC.PointResolutionZ = 0.4 +#digiTPC.DiffusionCoeffZ = 0.23 +#digiTPC.N_eff = 30 +#digiTPC.OutputLevel = DEBUG + + +## Muon Detector ## +from Configurables import MuonDigiAlg +digiMuon = MuonDigiAlg("MuonDigiAlg") +digiMuon.MuonBarrelHitsCollection = "MuonBarrelCollection" +digiMuon.MuonBarrelTrackerHits = "MuonBarrelTrackerHits" +digiMuon.WriteNtuple = 0 +digiMuon.OutFileName = "Digi_MUON.root" +######################################### + +################ +# Tracking +################ +from Configurables import KalTestTool +# Close multiple scattering and smooth, used by clupatra +kt010 = KalTestTool("KalTest010") +kt010.MSOn = False +kt010.Smooth = False +#kt010.OutputLevel = DEBUG + +# Open multiple scattering, energy loss and smooth (default) +kt111 = KalTestTool("KalTest111") +#kt111.OutputLevel = DEBUG + +# Close smooth +kt110 = KalTestTool("KalTest110") +kt110.Smooth = False +#kt110.OutputLevel = DEBUG + +from Configurables import SiliconTrackingAlg +tracking = SiliconTrackingAlg("SiliconTracking") +# for 3 layer ITK + 4s-2d bent-planar VTX, s single d double +tracking.LayerCombinations = [8,7,6, 8,7,5, 8,7,4, 8,6,5, 8,6,4, 8,6,3, 7,6,5, 7,6,4, 7,6,3, 7,5,3, 7,5,2, 7,4,3, 7,4,2, + 6,5,3, 6,5,2, 6,4,3, 6,4,2, 6,3,2, 6,3,1, 5,3,2, 5,3,1, 5,2,1, 5,2,0, 4,3,2, 4,3,1, 4,2,1, + 4,2,0, 3,2,1, 3,2,0, 3,1,0, 2,1,0] +tracking.LayerCombinationsFTD = [] +tracking.HeaderCol = "EventHeader" +tracking.VTXHitCollection = vxdhitname +tracking.SITHitCollection = sithitname +tracking.FTDPixelHitCollection = ftdhitname +tracking.FTDSpacePointCollection = "NULL" +tracking.SITRawHitCollection = "NotNeedForPixelSIT" +tracking.FTDRawHitCollection = "NotNeedForPixelFTD" +tracking.UseSIT = True +tracking.SmoothOn = False +tracking.NDivisionsInTheta = 10 +tracking.NDivisionsInPhi = 60 +tracking.NDivisionsInPhiFTD = 16 +tracking.MinDistCutAttach = 50 +# for p=1GeV, theta=10degree, Chi2FitCut = 1500, HelixMaxChi2 = 1000000, Chi2WZ = 0.02 +tracking.Chi2FitCut = 200 +tracking.MaxChi2PerHit = 200 +tracking.HelixMaxChi2 = 50000 +tracking.Chi2WZTriplet = 0.1 +tracking.Chi2WZQuartet = 0.1 +tracking.Chi2WZSeptet = 0.1 +#tracking.FitterTool = "KalTestTool/KalTest111" +tracking.OutputLevel = DEBUG + +from Configurables import ForwardTrackingAlg +forward = ForwardTrackingAlg("ForwardTracking") +forward.FTDPixelHitCollection = ftdhitname +forward.FTDSpacePointCollection = "NULL" +forward.FTDRawHitCollection = "NotNeedForPixelFTD" +forward.Chi2ProbCut = 0.0 +forward.HitsPerTrackMin = 3 +forward.BestSubsetFinder = "SubsetSimple" +forward.Criteria = ["Crit2_DeltaPhi","Crit2_StraightTrackRatio","Crit3_3DAngle","Crit3_ChangeRZRatio","Crit3_IPCircleDist","Crit4_3DAngleChange","Crit4_DistToExtrapolation", + "Crit2_DeltaRho","Crit2_RZRatio","Crit3_PT"] +forward.CriteriaMin = [0, 0.9, 0, 0.995, 0, 0.8, 0, 20, 1.002, 0.1, 0, 0.99, 0, 0.999, 0, 0.99, 0] +forward.CriteriaMax = [30, 1.02, 10, 1.015, 20, 1.3, 1.0, 150, 1.08, 99999999, 0.8, 1.01, 0.35, 1.001, 1.5, 1.01, 0.05] +#forward.FitterTool = "KalTestTool/KalTest110" +#forward.OutputLevel = DEBUG + +from Configurables import TrackSubsetAlg +subset = TrackSubsetAlg("TrackSubset") +subset.TrackInputCollections = ["ForwardTracks", "SiTracks"] +subset.RawTrackerHitCollections = [vxdhitname, sithitname, ftdhitname] +subset.TrackSubsetCollection = "SubsetTracks" +#subset.FitterTool = "KalTestTool/KalTest111" +subset.OutputLevel = DEBUG + +from Configurables import ClupatraAlg +clupatra = ClupatraAlg("Clupatra") +clupatra.TPCHitCollection = gashitname +#clupatra.OutputLevel = DEBUG + +from Configurables import FullLDCTrackingAlg +full = FullLDCTrackingAlg("FullTracking") +full.VTXTrackerHits = vxdhitname +full.SITTrackerHits = sithitname +full.TPCTrackerHits = gashitname +full.SETTrackerHits = sethitname +full.FTDPixelTrackerHits = ftdhitname +full.FTDSpacePoints = "NULL" +full.ETDTrackerHits = etdhitname +full.SITRawHits = "NotNeedForPixelSIT" +full.SETRawHits = "NotNeedForPixelSET" +full.FTDRawHits = "NotNeedForPixelFTD" +full.TPCTracks = "ClupatraTracks" # add standalone TPC track +full.SiTracks = "SiTracks" #"SubsetTracks" +full.OutputTracks = "CompleteTracks" # default name +#full.VTXHitToTrackDistance = 5. +full.FTDHitToTrackDistance = 5. +full.SITHitToTrackDistance = 3. +full.SETHitToTrackDistance = 5. +full.MinChi2ProbForSiliconTracks = 0 +full.MaxChi2PerHit = 200 +#full.ForceTPCSegmentsMerging = True +full.OutputLevel = DEBUG + +from Configurables import TPCDndxAlg +tpc_dndx = TPCDndxAlg("TPCDndxAlg") +tpc_dndx.Method = "Simple" + +from Configurables import TrackParticleRelationAlg +tpr = TrackParticleRelationAlg("Track2Particle") +tpr.MCParticleCollection = "MCParticle" +tpr.TrackList = ["CompleteTracks", "ClupatraTracks"] +tpr.TrackerAssociationList = ["VXDTrackerHitAssociation", "SITTrackerHitAssociation", "OTKBarrelTrackerHitAssociation", "FTDTrackerHitAssociation", "OTKEndcapTrackerHitAssociation", "TPCTrackerHitAss"] +#tpr.OutputLevel = DEBUG + +from Configurables import TrueMuonTagAlg +tmt = TrueMuonTagAlg("TrueMuonTag") +tmt.MCParticleCollection = "MCParticle" +tmt.TrackList = ["CompleteTracks"] +tmt.TrackerAssociationList = ["VXDTrackerHitAssociation", "SITTrackerHitAssociation", "OTKBarrelTrackerHitAssociation", "FTDTrackerHitAssociation", "TPCTrackerHitAss"] +tmt.MuonTagEfficiency = 0.95 # muon true tag efficiency, default is 1.0 (100%) +tmt.MuonDetTanTheta = 1.2 # muon det barrel/endcap separation tan(theta) +#tmt.OutputLevel = DEBUG + +# output +from Configurables import PodioOutput +out = PodioOutput("outputalg") +out.filename = "rec00.root" +out.outputCommands = ["keep *"] + +# ApplicationMgr +from Configurables import ApplicationMgr +mgr = ApplicationMgr( + TopAlg = [podioinput, digiVXD, digiSIT, digiSET, digiFTD, digiOTKE, digiTPC, digiMuon, tracking, forward, subset, clupatra, full, tpr, tpc_dndx, tmt, out], + EvtSel = 'NONE', + EvtMax = 50, + ExtSvc = [rndmengine, rndmgensvc, dsvc, evtseeder, geosvc, gearsvc, tracksystemsvc, pidsvc], + HistogramPersistency = 'ROOT', + OutputLevel = ERROR +) diff --git a/Reconstruction/RecSiTracking/src/SiliconTrackingAlg.cpp b/Reconstruction/RecSiTracking/src/SiliconTrackingAlg.cpp index d2b5e8b7..73f833e6 100644 --- a/Reconstruction/RecSiTracking/src/SiliconTrackingAlg.cpp +++ b/Reconstruction/RecSiTracking/src/SiliconTrackingAlg.cpp @@ -3013,7 +3013,7 @@ StatusCode SiliconTrackingAlg::setupGearGeom(){ try{ - debug() << " filling VXD parameters from gear::SITParameters " << endmsg ; + debug() << " filling VXD parameters from gear::VXDParameters " << endmsg ; pVXDDetMain = &gearMgr->getVXDParameters(); pVXDLayerLayout = &(pVXDDetMain->getVXDLayerLayout()); diff --git a/Reconstruction/RecTrkGlobal/src/FullLDCTracking/FullLDCTrackingAlg.cpp b/Reconstruction/RecTrkGlobal/src/FullLDCTracking/FullLDCTrackingAlg.cpp index 7b7ffd87..e1899a16 100755 --- a/Reconstruction/RecTrkGlobal/src/FullLDCTracking/FullLDCTrackingAlg.cpp +++ b/Reconstruction/RecTrkGlobal/src/FullLDCTracking/FullLDCTrackingAlg.cpp @@ -590,12 +590,14 @@ fitstart: int nhits_in_sit = track.getSubdetectorHitNumbers(2); int nhits_in_tpc = track.getSubdetectorHitNumbers(3); int nhits_in_set = track.getSubdetectorHitNumbers(4); + int nhits_in_etd = track.getSubdetectorHitNumbers(5); #else int nhits_in_vxd = track.getSubDetectorHitNumbers(0); int nhits_in_ftd = track.getSubDetectorHitNumbers(1); int nhits_in_sit = track.getSubDetectorHitNumbers(2); int nhits_in_tpc = track.getSubDetectorHitNumbers(3); int nhits_in_set = track.getSubDetectorHitNumbers(4); + int nhits_in_etd = track.getSubDetectorHitNumbers(5); #endif //int nhits_in_vxd = Track->subdetectorHitNumbers()[ 2 * lcio::ILDDetID::VXD - 2 ]; @@ -612,6 +614,7 @@ fitstart: << " sit hits = " << nhits_in_sit << " tpc hits = " << nhits_in_tpc << " set hits = " << nhits_in_set + << " etd hits = " << nhits_in_etd << endmsg; if (nhits_in_vxd > 0) track.setType( track.getType()| (1<<lcio::ILDDetID::VXD) ) ; @@ -619,6 +622,7 @@ fitstart: if (nhits_in_sit > 0) track.setType( track.getType()| (1<<lcio::ILDDetID::SIT) ) ; if (nhits_in_tpc > 0) track.setType( track.getType()| (1<<lcio::ILDDetID::TPC) ) ; if (nhits_in_set > 0) track.setType( track.getType()| (1<<lcio::ILDDetID::SET) ) ; + if (nhits_in_etd > 0) track.setType( track.getType()| (1<<lcio::ILDDetID::ETD) ) ; bool rejectTrack_onTPCHits = (nhits_in_tpc < _cutOnTPCHits) && (nHitsSi<=0); //bool rejectTrack_onITKHits = ( (nhits_in_tpc<=0) && (nhits_in_sit<1 && nhits_in_ftd<1) ); @@ -682,25 +686,27 @@ fitstart: debug() << " Add Track to final Collection: ID = " << track.id() << " for trkCand "<< trkCand << endmsg; } - float omega = trkState.omega; - float tanLambda = trkState.tanLambda; - float phi0 = trkState.phi; - float d0 = trkState.D0; - float z0 = trkState.Z0; - - HelixClass helix; - helix.Initialize_Canonical(phi0,d0,z0,omega,tanLambda,_bField); - - float trkPx = helix.getMomentum()[0]; - float trkPy = helix.getMomentum()[1]; - float trkPz = helix.getMomentum()[2]; - float trkP = sqrt(trkPx*trkPx+trkPy*trkPy+trkPz*trkPz); - - eTot += trkP; - pxTot += trkPx; - pyTot += trkPy; - pzTot += trkPz; - nTotTracks++; + if (trkState.location == edm4hep::TrackState::AtIP) { + float omega = trkState.omega; + float tanLambda = trkState.tanLambda; + float phi0 = trkState.phi; + float d0 = trkState.D0; + float z0 = trkState.Z0; + + HelixClass helix; + helix.Initialize_Canonical(phi0,d0,z0,omega,tanLambda,_bField); + + float trkPx = helix.getMomentum()[0]; + float trkPy = helix.getMomentum()[1]; + float trkPz = helix.getMomentum()[2]; + float trkP = sqrt(trkPx*trkPx+trkPy*trkPy+trkPz*trkPz); + + eTot += trkP; + pxTot += trkPx; + pyTot += trkPy; + pzTot += trkPz; + nTotTracks++; + } } if(m_tuple) m_timeKalman = stopwatch.RealTime()*1000; // streamlog_out(DEBUG5) << endmsg; @@ -834,8 +840,7 @@ void FullLDCTrackingAlg::prepareVectors() { pos[i] = hit.getPosition()[i]; } - debug() << " FTD Pixel Hit added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " layer = " << layer << endmsg; - + debug() << " FTD Pixel Hit added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " layer = " << layer << endmsg; } } @@ -910,9 +915,7 @@ void FullLDCTrackingAlg::prepareVectors() { pos[i] = hit.getPosition()[i]; } - debug() << " FTD SpacePoint Hit added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " layer = " << layer << endmsg; - - + debug() << " FTD SpacePoint Hit added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " layer = " << layer << endmsg; } } @@ -1046,7 +1049,7 @@ void FullLDCTrackingAlg::prepareVectors() { pos[i] = trkhit.getPosition()[i]; } - debug() << " SIT Hit " << trkhit.id() << " added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " layer = " << layer << endmsg; + debug() << " SIT Hit " << trkhit.id() << " added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " layer = " << layer << endmsg; } } @@ -1180,7 +1183,7 @@ void FullLDCTrackingAlg::prepareVectors() { pos[i] = trkhit.getPosition()[i]; } - debug() << " SET Hit " << trkhit.id() << " added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " layer = " << layer << endmsg; + debug() << " SET Hit " << trkhit.id() << " added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " layer = " << layer << endmsg; } } @@ -1213,17 +1216,6 @@ void FullLDCTrackingAlg::prepareVectors() { double dz(NAN); for(edm4hep::TrackerHit trkhit : *hitETDCol){ - // hit could be of the following type - // 1) TrackerHit, either ILDTrkHitTypeBit::COMPOSITE_SPACEPOINT or just standard TrackerHit - // 2) TrackerHitPlane, either 1D or 2D - // 3) TrackerHitZCylinder, if coming from a simple cylinder design as in the LOI - - // Establish which of these it is in the following order of likelyhood - // i) ILDTrkHitTypeBit::ONE_DIMENSIONAL (TrackerHitPlane) Should Never Happen: SpacePoints Must be Used Instead - // ii) ILDTrkHitTypeBit::COMPOSITE_SPACEPOINT (TrackerHit) - // iii) TrackerHitPlane (Two dimentional) - // iv) TrackerHitZCylinder - // v) Must be standard TrackerHit int layer = getLayerID(trkhit); if (layer < 0 || (unsigned)layer >= _nLayersETD) { @@ -1233,7 +1225,7 @@ void FullLDCTrackingAlg::prepareVectors() { // first check that we have not been given 1D hits by mistake, as they won't work here if ( UTIL::BitSet32( trkhit.getType() )[ UTIL::ILDTrkHitTypeBit::ONE_DIMENSIONAL ] ) { - fatal() << "SiliconTrackingAlg => fatal error in SIT : layer is outside allowed range : " << layer << endmsg; + fatal() << "SiliconTrackingAlg => fatal error in ETD : layer is outside allowed range : " << layer << endmsg; exit(1); } // most likely case: COMPOSITE_SPACEPOINT hits formed from stereo strip hits @@ -1256,13 +1248,13 @@ void FullLDCTrackingAlg::prepareVectors() { const float eps = 1.0e-07; // V must be the global z axis - if( fabs(1.0 - V.dot(Z)) > eps ) { - fatal() << "PIXEL ETD Hit measurment vectors V is not equal to the global Z axis. \n exit(1) called from file " << __FILE__ << " : " << __LINE__ << endmsg; + if (fabs(V.dot(Z)) > eps) { + fatal() << "PIXEL ETD Hit measurment vectors V is not int the global X-Y plane. \n exit(1) called from file " << __FILE__ << " : " << __LINE__ << endmsg; exit(1); } // U must be normal to the global z axis - if( fabs(U.dot(Z)) > eps ) { + if (fabs(U.dot(Z)) > eps) { fatal() << "PIXEL ETD Hit measurment vectors U is not in the global X-Y plane. \n exit(1) called from file " << __FILE__ << " : " << __LINE__ << endmsg; exit(1); } @@ -1270,20 +1262,9 @@ void FullLDCTrackingAlg::prepareVectors() { // FIXME should use the correct // drphi = trkhit_P->getdU(); // dz = trkhit_P->getdV(); - drphi = trkhit.getCovMatrix()[2]; - dz = trkhit.getCovMatrix()[5]; - } - // or a simple cylindrical design, as used in the LOI - /*FIXME, fucd - else if ( true ) { - trkhit_C = hitCollection->at( ielem ); - // FIXME - // drphi = trkhit_C->getdRPhi(); - // dz = trkhit_C->getdZ(); - drphi = 1.0; - dz = 1.0; + drphi = sqrt(trkhit.getCovMatrix()[2]*trkhit.getCovMatrix()[2]+trkhit.getCovMatrix()[5]*trkhit.getCovMatrix()[5]); + dz = 0.1; } - */ // this would be very unlikely, but who knows ... just an ordinary TrackerHit, which is not a COMPOSITE_SPACEPOINT else { // SJA:FIXME: fudge for now by a factor of two and ignore covariance @@ -1312,7 +1293,7 @@ void FullLDCTrackingAlg::prepareVectors() { pos[i] = trkhit.getPosition()[i]; } - debug() << " ETD Hit " << trkhit.id() << " added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " layer = " << layer << endmsg; + debug() << " ETD Hit " << trkhit.id() << " added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " layer = " << layer << endmsg; } } @@ -1849,7 +1830,7 @@ TrackExtended * FullLDCTrackingAlg::CombineTracks(TrackExtended * tpcTrack, Trac std::auto_ptr<MarlinTrk::IMarlinTrack> marlin_trk_autop(_trksystem->createTrack()); MarlinTrk::IMarlinTrack& marlin_trk = *marlin_trk_autop.get(); - + edm4hep::TrackState pre_fit ; int errorCode = IMarlinTrack::success; @@ -3786,16 +3767,14 @@ void FullLDCTrackingAlg::AssignOuterHitsToTracks(TrackerHitExtendedVec hitVec, f TrackHitPair * trkHitPair = pairs[iP]; TrackExtended * trkExt = trkHitPair->getTrackExtended(); - TrackerHitExtended * trkHitExt = - - trkHitPair->getTrackerHitExtended(); + TrackerHitExtended * trkHitExt = trkHitPair->getTrackerHitExtended(); // check if the track or hit is still free to be combined if (flagTrack[trkExt] && flagHit[trkHitExt]) { if (refit==0) { // just set the association trkExt->addTrackerHitExtended( trkHitExt ); - trkHitExt->setUsedInFit( false ); + trkHitExt->setUsedInFit( true ); trkHitExt->setTrackExtended( trkExt ); } @@ -3870,8 +3849,10 @@ void FullLDCTrackingAlg::AssignOuterHitsToTracks(TrackerHitExtendedVec hitVec, f } debug() << "AssignOuterHitsToTracks: Start Fitting: AddHits: number of hits to fit " << trkHits.size() << endmsg; - + MarlinTrk::IMarlinTrack* marlin_trk = _trksystem->createTrack(); + //std::auto_ptr<MarlinTrk::IMarlinTrack> marlin_trk_autop(_trksystem->createTrack()); + //MarlinTrk::IMarlinTrack* marlin_trk = marlin_trk_autop.get(); edm4hep::TrackState pre_fit ; pre_fit.D0 = trkExt->getD0(); @@ -3900,16 +3881,16 @@ void FullLDCTrackingAlg::AssignOuterHitsToTracks(TrackerHitExtendedVec hitVec, f covMatrix[14] = ( _initialTrackError_tanL ); //sigma_tanl^2 pre_fit.covMatrix = covMatrix; - - int error = MarlinTrk::createFit( trkHits, marlin_trk, &pre_fit, _bField, IMarlinTrack::backward , _maxChi2PerHit ); - + debug() << "AssignOuterHitsToTracks: before createFit" << endmsg; + int error = MarlinTrk::createFit( trkHits, marlin_trk, &pre_fit, _bField, IMarlinTrack::forward/*backward*/, _maxChi2PerHit ); + debug() << "AssignOuterHitsToTracks: after createFit" << endmsg; if ( error != IMarlinTrack::success ) { debug() << "FullLDCTrackingAlg::AssignOuterHitsToTracks: creation of fit fails with error " << error << endmsg; delete marlin_trk ; continue ; } - + std::vector<std::pair<edm4hep::TrackerHit , double> > outliers ; marlin_trk->getOutliers(outliers); @@ -3922,7 +3903,7 @@ void FullLDCTrackingAlg::AssignOuterHitsToTracks(TrackerHitExtendedVec hitVec, f delete marlin_trk ; continue ; } - + debug() << "AssignOuterHitsToTracks: before propagate" << endmsg; edm4hep::Vector3d point(0.,0.,0.); // nominal IP int return_code = 0; @@ -3930,7 +3911,7 @@ void FullLDCTrackingAlg::AssignOuterHitsToTracks(TrackerHitExtendedVec hitVec, f return_code = marlin_trk->propagate(point, trkState, chi2_D, ndf ) ; delete marlin_trk ; - + debug() << "AssignOuterHitsToTracks: after delete" << endmsg; if ( error != IMarlinTrack::success ) { debug() << "FullLDCTrackingAlg::AssignOuterHitsToTracks: propagate to IP fails with error " << error << endmsg; continue ; @@ -4162,9 +4143,10 @@ void FullLDCTrackingAlg::AssignTPCHitsToTracks(TrackerHitExtendedVec hitVec, if (tracksToAttach[iH]!=NULL) { tracksToAttach[iH]->addTrackerHitExtended(trkHitExt); trkHitExt->setTrackExtended( tracksToAttach[iH] ); - //by fucd + //by fucd: if set true, too much TPC hit added for those swirl tracks will cause silicon tracker hits lost + // only set true while no TPC track input and assign TPC hits into silicon tracks //trkHitExt->setUsedInFit( false ); - trkHitExt->setUsedInFit( true ); + trkHitExt->setUsedInFit(_assignTPCHits==2); } } @@ -4419,7 +4401,7 @@ void FullLDCTrackingAlg::AssignSiHitsToTracks(TrackerHitExtendedVec hitVec, delete marlin_trk ; continue ; } - + std::vector<std::pair<edm4hep::TrackerHit , double> > outliers ; marlin_trk->getOutliers(outliers); @@ -5178,8 +5160,8 @@ void FullLDCTrackingAlg::setupGearGeom(){ try { debug() << " filling ETD parameters from gear::ETDParameters " << endmsg; - const gear::GearParameters& pETDDet = gearMgr->getGearParameters("ETD"); - _nLayersETD = int(pETDDet.getDoubleVals("ETDLayerZ").size()); + const gear::GearParameters& pETDDet = gearMgr->getGearParameters("ETDParameters"); + _nLayersETD = int(pETDDet.getIntVals("ETDPetalNumber").size()); } catch (...) { debug() << " ### gear::GearParameters ETD Not Present in GEAR FILE" << endmsg; diff --git a/Reconstruction/RecTrkGlobal/src/FullLDCTracking/FullLDCTrackingAlg.h b/Reconstruction/RecTrkGlobal/src/FullLDCTracking/FullLDCTrackingAlg.h index 9864168f..a73866ca 100755 --- a/Reconstruction/RecTrkGlobal/src/FullLDCTracking/FullLDCTrackingAlg.h +++ b/Reconstruction/RecTrkGlobal/src/FullLDCTracking/FullLDCTrackingAlg.h @@ -516,7 +516,7 @@ protected: DataHandle<edm4hep::TrackerHitCollection> _inSITRawColHdl{"SITTrackerHits", Gaudi::DataHandle::Reader, this}; DataHandle<edm4hep::TrackerHitCollection> _inSETRawColHdl{"SETTrackerHits", Gaudi::DataHandle::Reader, this}; - DataHandle<edm4hep::TrackerHitCollection> _inETDRawColHdl{"SETTrackerHits", Gaudi::DataHandle::Reader, this}; + DataHandle<edm4hep::TrackerHitCollection> _inETDRawColHdl{"ETDTrackerHits", Gaudi::DataHandle::Reader, this}; DataHandle<edm4hep::TrackerHitCollection> _inFTDRawColHdl{"FTDStripTrackerHits", Gaudi::DataHandle::Reader, this}; //DataHandle<edm4hep::SimTrackerHitCollection> _inVTXRawColHdl{"VXDCollection", Gaudi::DataHandle::Reader, this}; diff --git a/Service/GearSvc/src/GearSvc.cpp b/Service/GearSvc/src/GearSvc.cpp index c8ddfc05..a29f4e99 100644 --- a/Service/GearSvc/src/GearSvc.cpp +++ b/Service/GearSvc/src/GearSvc.cpp @@ -6,6 +6,7 @@ #include "gearxml/GearXML.h" #include "gearimpl/GearMgrImpl.h" #include "gearimpl/ConstantBField.h" +#include "gearimpl/GearParametersImpl.h" #include "gearimpl/ZPlanarParametersImpl.h" #include "gearimpl/ZPlanarLayerLayoutImpl.h" #include "gearimpl/FTDParametersImpl.h" @@ -104,6 +105,9 @@ StatusCode GearSvc::initialize() else if(it->first=="SET" || it->first=="OTKBarrel"){ sc = convertSET(sub); } + else if(it->first=="ETD" || it->first=="OTKEndcap"){ + sc = convertETD(sub); + } else if(it->first=="EcalBarrel"||it->first=="EcalEndcap"||it->first=="EcalPlug"|| it->first=="HcalBarrel"||it->first=="HcalEndcap"||it->first=="HcalRing"|| it->first=="YokeBarrel"||it->first=="YokeEndcap"||it->first=="YokePlug"|| @@ -323,6 +327,15 @@ StatusCode GearSvc::convertComposite(dd4hep::DetElement& vtx){ l.zHalfSensitive/dd4hep::mm, l.widthSensitive/dd4hep::mm, 0.); } + { + const dd4hep::rec::ZPlanarData::LayerLayout& l = vtxData->layersPlanar[0] ; + double offset = l.offsetSupport; + dd4hep::rec::Vector3D a( l.distanceSensitive + l.thicknessSensitive, offset, 2.*dd4hep::mm); + dd4hep::rec::Vector3D b( l.distanceSupport + l.thicknessSupport, offset, 2.*dd4hep::mm); + gear::SimpleMaterialImpl* VXDSupportMaterial = CreateGearMaterial(a, b, "VXDSupportMaterial"); + m_gearMgr->registerSimpleMaterial(VXDSupportMaterial); + } + std::vector<int> ids; std::vector<double> zhalfs, rsens, tsens, rsups, tsups, phi0s, rgaps, dphis; @@ -352,11 +365,13 @@ StatusCode GearSvc::convertComposite(dd4hep::DetElement& vtx){ m_gearMgr->setVXDParameters(gearVTX); - const dd4hep::rec::CylindricalData::LayerLayout& l = vtxData->layersBent[0] ; - dd4hep::rec::Vector3D a(l.radiusSupport, l.phi0, 0., dd4hep::rec::Vector3D::cylindrical); - dd4hep::rec::Vector3D b(l.radiusSupport + l.thicknessSupport, l.phi0, 0., dd4hep::rec::Vector3D::cylindrical); - gear::SimpleMaterialImpl* VXDSupportMaterial = CreateGearMaterial(a, b, "VXDSupportMaterial"); - m_gearMgr->registerSimpleMaterial(VXDSupportMaterial); + { + const dd4hep::rec::CylindricalData::LayerLayout& l = vtxData->layersBent[0] ; + dd4hep::rec::Vector3D a(l.radiusSupport, l.phi0, 0., dd4hep::rec::Vector3D::cylindrical); + dd4hep::rec::Vector3D b(l.radiusSupport + l.thicknessSupport, l.phi0, 0., dd4hep::rec::Vector3D::cylindrical); + gear::SimpleMaterialImpl* VXDSupportMaterial = CreateGearMaterial(a, b, "VXDBentSupportMaterial"); + m_gearMgr->registerSimpleMaterial(VXDSupportMaterial); + } if (vtxData->rOuterShell>vtxData->rInnerShell) { dd4hep::rec::Vector3D a1( vtxData->rInnerShell, 0, 2.*dd4hep::mm); @@ -423,6 +438,71 @@ StatusCode GearSvc::convertFTD(dd4hep::DetElement& ftd){ return StatusCode::SUCCESS; } +StatusCode GearSvc::convertETD(dd4hep::DetElement& etd){ + dd4hep::rec::ZDiskPetalsData* etdData = nullptr; + try{ + etdData = etd.extension<dd4hep::rec::ZDiskPetalsData>(); + } + catch(std::runtime_error& e){ + warning() << e.what() << " " << etdData << endmsg; + return StatusCode::FAILURE; + } + + std::vector<dd4hep::rec::ZDiskPetalsData::LayerLayout>& etdlayers = etdData->layers; + int nLayers = etdlayers.size(); + + gear::GearParametersImpl* etdParam = new gear::GearParametersImpl(); + etdParam->setDoubleVal("strip_width_mm", etdData->widthStrip*CLHEP::cm); + etdParam->setDoubleVal("strip_length_mm", etdData->lengthStrip*CLHEP::cm); + etdParam->setDoubleVal("strip_pitch_mm", etdData->pitchStrip*CLHEP::cm); + etdParam->setDoubleVal("strip_angle_deg", etdData->angleStrip*rad_to_deg); + + std::vector<int> nPetals, nSensors; + std::vector<double> petalangles, phi0s, alphas, zpositions, zoffsets, supRinners, supThicknesss, supHeights, senRinners, senThicknesss, senHeights; + for(int layer = 0; layer < nLayers; layer++){ + dd4hep::rec::ZDiskPetalsData::LayerLayout& etdlayer = etdlayers[layer]; + nPetals.push_back(etdlayer.petalNumber); + petalangles.push_back(etdlayer.petalHalfAngle*2); + phi0s.push_back(etdlayer.phi0); + alphas.push_back(etdlayer.alphaPetal); + zpositions.push_back(etdlayer.zPosition*CLHEP::cm); + zoffsets.push_back(etdlayer.zOffsetSupport*CLHEP::cm); + + supRinners.push_back(etdlayer.distanceSupport*CLHEP::cm); + supThicknesss.push_back(etdlayer.thicknessSupport*CLHEP::cm); + supHeights.push_back(etdlayer.lengthSupport*CLHEP::cm); + senRinners.push_back(etdlayer.distanceSensitive*CLHEP::cm); + senThicknesss.push_back(etdlayer.thicknessSensitive*CLHEP::cm); + senHeights.push_back(etdlayer.lengthSensitive*CLHEP::cm); + + nSensors.push_back(etdlayer.sensorsPerPetal); + } + etdParam->setIntVals("ETDPetalNumber", nPetals); + etdParam->setIntVals("ETDSensorNumber", nSensors); + etdParam->setDoubleVals("ETDPetalAngle", petalangles); + etdParam->setDoubleVals("ETDDiskPhi0", phi0s); + etdParam->setDoubleVals("ETDDiskAlpha", alphas); + etdParam->setDoubleVals("ETDDiskPosition", zpositions); + etdParam->setDoubleVals("ETDDiskOffset", zoffsets); + etdParam->setDoubleVals("ETDSupportRmin", supRinners); + etdParam->setDoubleVals("ETDSupportThickness", supThicknesss); + etdParam->setDoubleVals("ETDSupportHeight", supHeights); + etdParam->setDoubleVals("ETDSensitiveRmin", senRinners); + etdParam->setDoubleVals("ETDSensitiveThickness", senThicknesss); + etdParam->setDoubleVals("ETDSensitiveHeight", senHeights); + + const dd4hep::rec::ZDiskPetalsData::LayerLayout& l = etdData->layers[0]; + double z = l.zPosition; + dd4hep::rec::Vector3D a(l.distanceSupport+0.5*l.lengthSupport, l.phi0, z-l.thicknessSupport, dd4hep::rec::Vector3D::cylindrical); + dd4hep::rec::Vector3D b(l.distanceSupport+0.5*l.lengthSupport, l.phi0, z, dd4hep::rec::Vector3D::cylindrical); + gear::SimpleMaterialImpl* ETDSupportMaterial = CreateGearMaterial(a, b, "OTKEndcapSupportMaterial"); + m_gearMgr->registerSimpleMaterial(ETDSupportMaterial); + + m_gearMgr->setGearParameters("ETDParameters", etdParam); + info() << "nftd = " << nLayers << endmsg; + return StatusCode::SUCCESS; +} + StatusCode GearSvc::convertSIT(dd4hep::DetElement& sit){ dd4hep::rec::ZPlanarData* sitData = nullptr; try{ diff --git a/Service/GearSvc/src/GearSvc.h b/Service/GearSvc/src/GearSvc.h index 38e831de..9b7cd597 100644 --- a/Service/GearSvc/src/GearSvc.h +++ b/Service/GearSvc/src/GearSvc.h @@ -29,6 +29,7 @@ class GearSvc : public extends<Service, IGearSvc> StatusCode convertDC (dd4hep::DetElement& dc); StatusCode convertSET(dd4hep::DetElement& set); StatusCode convertFTD(dd4hep::DetElement& ftd); + StatusCode convertETD(dd4hep::DetElement& etd); StatusCode convertCal(dd4hep::DetElement& cal); TGeoNode* FindNode(TGeoNode* mother, char* name); gear::SimpleMaterialImpl* CreateGearMaterial(const dd4hep::rec::Vector3D& a, const dd4hep::rec::Vector3D& b, const std::string name); diff --git a/Service/TrackSystemSvc/src/MarlinKalTest.cc b/Service/TrackSystemSvc/src/MarlinKalTest.cc index bd7fade4..8ff54501 100644 --- a/Service/TrackSystemSvc/src/MarlinKalTest.cc +++ b/Service/TrackSystemSvc/src/MarlinKalTest.cc @@ -15,6 +15,7 @@ #include "kaldet/ILDSITCylinderKalDetector.h" #include "kaldet/ILDSETKalDetector.h" #include "kaldet/CEPCOTKKalDetector.h" +#include "kaldet/CEPCOTKEndcapKalDetector.h" #include "kaldet/ILDFTDKalDetector.h" #include "kaldet/ILDFTDDiscBasedKalDetector.h" #include "kaldet/ILDTPCKalDetector.h" @@ -186,6 +187,16 @@ namespace MarlinTrk{ std::cout << "Warning: " << " MarlinKalTest - Simple Disc Based FTD missing in gear file: Simple Disc Based FTD Not Built " << std::endl ; } } + + try{ + CEPCOTKEndcapKalDetector* etddet = new CEPCOTKEndcapKalDetector(*_gearMgr, _geoSvc); + // store the measurement layer id's for the active layers + this->storeActiveMeasurementModuleIDs(etddet); + _det->Install(*etddet); + } + catch( gear::UnknownParameterException& e){ + std::cout << "Warning: " << " MarlinKalTest - ETD missing in gear file: Petal Based ETD Not Built " << std::endl; + } try{ ILDTPCKalDetector* tpcdet = new ILDTPCKalDetector( *_gearMgr, _geoSvc ) ; @@ -333,31 +344,23 @@ namespace MarlinTrk{ int sensitive_element_id = *it; this->_active_measurement_modules.insert(std::pair<int,const ILDVMeasLayer*>( sensitive_element_id, ml )); ++it; - + } int subdet_layer_id = ml->getLayerID() ; this->_active_measurement_modules_by_layer.insert(std::pair<int ,const ILDVMeasLayer*>(subdet_layer_id,ml)); - - //streamlog_out(DEBUG0) << "MarlinKalTest::storeActiveMeasurementLayerIDs added active layer with " - //<< " LayerID = " << subdet_layer_id << " and DetElementIDs " ; - - for (it = ml->getCellIDs().begin(); it!=ml->getCellIDs().end(); ++it) { - - //streamlog_out(DEBUG0) << " : " << *it ; - - } - - //streamlog_out(DEBUG0) << std::endl; - - - - +//#define DEBUG_CELLID 1 +#ifdef DEBUG_CELLID + std::cout << "MarlinKalTest::storeActiveMeasurementLayerIDs added active layer with " + << " LayerID = " << subdet_layer_id << " and DetElementIDs " ; + for (it = ml->getCellIDs().begin(); it!=ml->getCellIDs().end(); ++it) { + std::cout << " : " << *it ; + } + std::cout << std::endl; +#endif } - } - } const ILDVMeasLayer* MarlinKalTest::getLastMeasLayer(THelicalTrack const& hel, TVector3 const& point) const { diff --git a/Service/TrackSystemSvc/src/MarlinKalTestTrack.cc b/Service/TrackSystemSvc/src/MarlinKalTestTrack.cc index 839da062..4f3bf458 100644 --- a/Service/TrackSystemSvc/src/MarlinKalTestTrack.cc +++ b/Service/TrackSystemSvc/src/MarlinKalTestTrack.cc @@ -116,19 +116,20 @@ namespace MarlinTrk { int MarlinKalTestTrack::addHit( edm4hep::TrackerHit& trkhit) { - + //std::cout << "MarlinKalTestTrack::addHit: trkhit = " << trkhit.id() << " cell: " << trkhit.getCellID() << std::endl; return this->addHit( trkhit, _ktest->findMeasLayer( trkhit )) ; } int MarlinKalTestTrack::addHit( edm4hep::TrackerHit& trkhit, const ILDVMeasLayer* ml) { - //std::cout << "MarlinKalTestTrack::addHit: trkhit = " << trkhit.id() << " addr: " << trkhit << " ml = " << ml << std::endl ; + //std::cout << "MarlinKalTestTrack::addHit: trkhit = " << trkhit.id() << " cell: " << trkhit.getCellID() << " ml = " << ml << std::endl ; if( trkhit.isAvailable() && ml ) { //if(ml){ + //std::cout << "in ILDVMeasLayer " << ml->GetName() << std::endl; return this->addHit( trkhit, ml->ConvertLCIOTrkHit(trkhit), ml) ; } else { - //std::cout << "MarlinKalTestTrack::addHit: - bad inputs trkhit = " << trkhit.id() << " addr: " << trkhit << " ml = " << ml << std::endl ; + //std::cout << "MarlinKalTestTrack::addHit: - bad inputs trkhit = " << trkhit.id() << " addr: " << trkhit.getCellID() << " ml = " << ml << std::endl ; return bad_intputs ; } return bad_intputs ; diff --git a/Utilities/KalDet/src/ild/ftd/ILDFTDKalDetector.h b/Utilities/KalDet/include/kaldet/CEPCOTKEndcapKalDetector.h similarity index 62% rename from Utilities/KalDet/src/ild/ftd/ILDFTDKalDetector.h rename to Utilities/KalDet/include/kaldet/CEPCOTKEndcapKalDetector.h index b3f66522..7a02cd20 100644 --- a/Utilities/KalDet/src/ild/ftd/ILDFTDKalDetector.h +++ b/Utilities/KalDet/include/kaldet/CEPCOTKEndcapKalDetector.h @@ -1,10 +1,10 @@ -#ifndef __ILDFTDDETECTOR__ -#define __ILDFTDDETECTOR__ +#ifndef __CEPCOTKEndcapDETECTOR__ +#define __CEPCOTKEndcapDETECTOR__ -/** Petal based FTD to be used for ILD DBD studies +/** Petal based OTKEndcap to be used for CEPC TDR studies * WARNING: Still very experimental * - * @author S.Aplin DESY, Robin Glattauer HEPHY + * @author */ #include "kaltest/TVKalDetector.h" @@ -17,17 +17,14 @@ namespace gear{ class GearMgr ; } -class ILDFTDKalDetector : public TVKalDetector { +class CEPCOTKEndcapKalDetector : public TVKalDetector { public: - - /** Initialize the FTD from GEAR */ - ILDFTDKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* geoSvc ); - - + /** Initialize the OTKEndcap from GEAR */ + CEPCOTKEndcapKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* geoSvc ); + private: - struct FTD_Petal { - + struct OTKEndcap_Petal { int ipetal; double phi; double alpha; @@ -39,11 +36,9 @@ private: double supThickness; double senZPos; bool faces_ip; - }; - - struct FTD_Disk { + struct OTKEndcap_Disk { int nPetals; double phi0; double dphi; @@ -65,29 +60,22 @@ private: bool isStripReadout; int nSensors; - - }; - void build_staggered_design(); - //void create_petal(TVector3 measurement_plane_centre, FTD_Petal petal, int CellID); + //void create_petal(TVector3 measurement_plane_centre, OTKEndcap_Petal petal, int CellID); /** * @param zpos the z position of the front measurement surface (middle of front sensitive) */ void create_segmented_disk_layers(int idisk, int nsegments, bool even_petals, double phi0, double zpos ); + void setupGearGeom(const gear::GearMgr& gearMgr); + void setupGearGeom(IGeomSvc* geoSvc); - void setupGearGeom( const gear::GearMgr& gearMgr ) ; - void setupGearGeom( IGeomSvc* geoSvc ); - - int _nDisks ; - double _bZ ; - + int _nDisks; + double _bZ; - std::vector<FTD_Disk> _FTDgeo; - + std::vector<OTKEndcap_Disk> _OTKEndcapgeo; }; - #endif diff --git a/Utilities/KalDet/src/ild/common/ILDSegmentedDiscMeasLayer.h b/Utilities/KalDet/include/kaldet/CEPCSegmentedDiscMeasLayer.h similarity index 59% rename from Utilities/KalDet/src/ild/common/ILDSegmentedDiscMeasLayer.h rename to Utilities/KalDet/include/kaldet/CEPCSegmentedDiscMeasLayer.h index 68bf1c99..28eeef8c 100644 --- a/Utilities/KalDet/src/ild/common/ILDSegmentedDiscMeasLayer.h +++ b/Utilities/KalDet/include/kaldet/CEPCSegmentedDiscMeasLayer.h @@ -1,10 +1,10 @@ -#ifndef __ILDSEGMENTEDDISCMEASLAYER_H__ -#define __ILDSEGMENTEDDISCMEASLAYER_H__ +#ifndef __CEPCSEGMENTEDDISCMEASLAYER_H__ +#define __CEPCSEGMENTEDDISCMEASLAYER_H__ -/** ILDSegmentedDiscMeasLayer: User defined Segemented Disk Planar KalTest measurement layer class used with ILDPLanarTrackHit. Segments are isosolese trapezoids whose axis of symmetry points to the origin +/** CEPCSegmentedDiscMeasLayer: User defined Segemented Disk Planar KalTest measurement layer class used with ILDPLanarTrackHit. Following ILDSegmentedDiscMeasLayer * WARNING: ONLY IMPLEMENTED FOR X AND Y COORDINATES AT FIXED Z * - * @author S.Aplin DESY + * @author C.Fu CEPC */ #include "TVector3.h" @@ -22,38 +22,37 @@ class TVTrackHit; -class ILDSegmentedDiscMeasLayer : public ILDVMeasLayer, public TPlane { +class CEPCSegmentedDiscMeasLayer : public ILDVMeasLayer, public TPlane { public: // Ctors and Dtor - ILDSegmentedDiscMeasLayer(TMaterial &min, - TMaterial &mout, - double Bz, - double SortingPolicy, - int nsegments, - double zpos, - double phi0, // defined by the axis of symmerty of the first petal - double trap_rmin, - double trap_height, - double trap_innerBaseLength, - double trap_outerBaseLength, - bool is_active, - std::vector<int> CellIDs, - const Char_t *name = "ILDDiscMeasL"); - - ILDSegmentedDiscMeasLayer(TMaterial &min, - TMaterial &mout, - double Bz, - double SortingPolicy, - int nsegments, - double zpos, - double phi0, // defined by the axis of symmerty of the first petal - double trap_rmin, - double trap_height, - double trap_innerBaseLength, - double trap_outerBaseLength, - bool is_active, - const Char_t *name = "ILDDiscMeasL"); + CEPCSegmentedDiscMeasLayer(TMaterial &min, + TMaterial &mout, + double Bz, + double SortingPolicy, + int nsegments, + double zpos, + double phi0, // defined by the axis of symmerty of the first petal + double rmin, + double rmax, + double halfPetal, + std::vector<int> nsensors, + bool is_active, + std::vector<int> CellIDs, + const Char_t *name = "CEPCDiscMeasL"); + + CEPCSegmentedDiscMeasLayer(TMaterial &min, + TMaterial &mout, + double Bz, + double SortingPolicy, + int nsegments, + double zpos, + double phi0, // defined by the axis of symmerty of the first petal + double rmin, + double rmax, + double halfPetal, + bool is_active, + const Char_t *name = "CEPCDiscMeasL"); // Parrent's pure virtuals that must be implemented @@ -63,7 +62,6 @@ public: const TVector3 &xv) const { return this->XvToMv(xv); } - /** Global to Local coordinates */ virtual TKalMatrix XvToMv (const TVector3 &xv) const; @@ -93,11 +91,8 @@ public: Double_t eps = 1.e-8) const{ return CalcXingPointWith(hel,xx,phi,0,eps); - } - - /** Get the intersection and the CellID, needed for multilayers */ virtual int getIntersectionAndCellID(const TVTrack &hel, TVector3 &xx, @@ -105,8 +100,6 @@ public: Int_t &CellID, Int_t mode, Double_t eps = 1.e-8) const ; - - /** Check if global point is on surface */ virtual Bool_t IsOnSurface (const TVector3 &xx) const; @@ -118,6 +111,7 @@ protected: double angular_range_2PI( double phi ) const; unsigned int get_segment_index(double phi) const; + unsigned int get_sensor_index(double r, double phi) const; double get_segment_phi(unsigned int index) const; TVector3 get_segment_centre(unsigned int index) const; @@ -125,20 +119,13 @@ private: double _sortingPolicy; int _nsegments; - double _trap_rmin; - double _trap_height; - double _trap_inner_base_length; - double _trap_outer_base_length; - double _trap_tan_beta; // tan of the openning angle of the petal - + double _rmin; double _rmax; + double _halfPetal; + std::vector<int> _nsensors; + double _start_phi; // trailing edge of the first sector double _segment_dphi; - - - + int _nrow; }; - - - #endif diff --git a/Utilities/KalDet/include/kaldet/ILDDiscMeasLayer.h b/Utilities/KalDet/include/kaldet/ILDDiscMeasLayer.h index 2f6f6167..9f489bdd 100644 --- a/Utilities/KalDet/include/kaldet/ILDDiscMeasLayer.h +++ b/Utilities/KalDet/include/kaldet/ILDDiscMeasLayer.h @@ -41,7 +41,21 @@ public: _sortingPolicy(SortingPolicy), _rMin(rMin), _rMax(rMax) { /* no op */ } - + ILDDiscMeasLayer(TMaterial &min, + TMaterial &mout, + const TVector3 ¢er, + const TVector3 &normal, + double Bz, + double SortingPolicy, + double rMin, + double rMax, + const std::vector<int>& CellIDs, + Bool_t is_active, + const Char_t *name = "ILDDiscMeasL") + : ILDVMeasLayer(min, mout, Bz, CellIDs, is_active, name), + TPlane(center, normal), + _sortingPolicy(SortingPolicy), _rMin(rMin), _rMax(rMax) + { /* no op */ } // Parrent's pure virtuals that must be implemented diff --git a/Utilities/KalDet/include/kaldet/ILDFTDKalDetector.h b/Utilities/KalDet/include/kaldet/ILDFTDKalDetector.h index c6ff6fa0..b3f66522 100644 --- a/Utilities/KalDet/include/kaldet/ILDFTDKalDetector.h +++ b/Utilities/KalDet/include/kaldet/ILDFTDKalDetector.h @@ -17,12 +17,11 @@ namespace gear{ class GearMgr ; } - class ILDFTDKalDetector : public TVKalDetector { public: /** Initialize the FTD from GEAR */ - ILDFTDKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* geoSvc ); + ILDFTDKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* geoSvc ); private: @@ -81,6 +80,7 @@ private: void setupGearGeom( const gear::GearMgr& gearMgr ) ; + void setupGearGeom( IGeomSvc* geoSvc ); int _nDisks ; double _bZ ; diff --git a/Utilities/KalDet/include/kaldet/ILDVXDKalDetector.h b/Utilities/KalDet/include/kaldet/ILDVXDKalDetector.h index 4feebacd..d1bbaaa2 100644 --- a/Utilities/KalDet/include/kaldet/ILDVXDKalDetector.h +++ b/Utilities/KalDet/include/kaldet/ILDVXDKalDetector.h @@ -27,7 +27,8 @@ public: private: - void setupGearGeom( const gear::GearMgr& gearMgr, IGeomSvc* geoSvc) ; + void setupGearGeom( const gear::GearMgr& gearMgr ); + void setupGearGeom( IGeomSvc* geoSvc) ; int _nLayers ; double _bZ ; diff --git a/Utilities/KalDet/src/ild/common/CEPCSegmentedDiscMeasLayer.cc b/Utilities/KalDet/src/ild/common/CEPCSegmentedDiscMeasLayer.cc new file mode 100644 index 00000000..0fd31662 --- /dev/null +++ b/Utilities/KalDet/src/ild/common/CEPCSegmentedDiscMeasLayer.cc @@ -0,0 +1,485 @@ +#include "kaldet/CEPCSegmentedDiscMeasLayer.h" +#include "kaldet/ILDPlanarHit.h" + +#include <UTIL/BitField64.h> +#include <UTIL/ILDConf.h> +#include "DetIdentifier/CEPCConf.h" + +#include "kaltest/TVTrack.h" +#include "TVector3.h" +#include "TMath.h" +#include "TRotMatrix.h" +#include "TBRIK.h" +#include "TNode.h" +#include "TString.h" + +//#include <EVENT/TrackerHitPlane.h> +#include <bitset> +#include <math.h> +#include <assert.h> +#include <algorithm> + +// #include "streamlog/streamlog.h" +//#define DEBUGPRINT + +CEPCSegmentedDiscMeasLayer::CEPCSegmentedDiscMeasLayer(TMaterial &min, + TMaterial &mout, + double Bz, + double sortingPolicy, + int nsegments, + double zpos, + double phi0, // defined by the axis of symmerty of the first petal + double rmin, + double rmax, + double halfPetal, + std::vector<int> nsensors, + bool is_active, + std::vector<int> CellIDs, + const Char_t *name) : + ILDVMeasLayer(min, mout, Bz, CellIDs, is_active, name), + TPlane(TVector3(0.,0.,zpos), TVector3(0.,0.,zpos)), + _sortingPolicy(sortingPolicy),_nsegments(nsegments),_rmin(rmin),_rmax(rmax),_halfPetal(halfPetal),_nsensors(nsensors) { + + _segment_dphi = 2.0*M_PI / _nsegments; + + phi0 = angular_range_2PI(phi0); + + _start_phi = phi0 - _halfPetal; + + _start_phi = angular_range_2PI(_start_phi); + + _nrow = _nsensors.size(); + // now check for constistency + if (_halfPetal*2 > _segment_dphi ) { + std::cout << "CEPCSegmentedDiscMeasLayer::CEPCSegmentedDiscMeasLayer overlaps: exit(1) called from " << __FILE__ << " line " << __LINE__ << std::endl; + exit(1); + } +} + +CEPCSegmentedDiscMeasLayer::CEPCSegmentedDiscMeasLayer(TMaterial &min, + TMaterial &mout, + double Bz, + double sortingPolicy, + int nsegments, + double zpos, + double phi0, // defined by the axis of symmerty of the first petal + double rmin, + double rmax, + double halfPetal, + bool is_active, + const Char_t *name) : + ILDVMeasLayer(min, mout, Bz, is_active, -1, name), + TPlane(TVector3(0.,0.,zpos), TVector3(0.,0.,zpos)), + _sortingPolicy(sortingPolicy),_nsegments(nsegments),_rmin(rmin),_rmax(rmax),_halfPetal(halfPetal) { + + _segment_dphi = 2.0*M_PI / _nsegments; + + phi0 = angular_range_2PI(phi0); + + _start_phi = phi0 - _halfPetal; + + _start_phi = angular_range_2PI(_start_phi); + + // now check for constistency + if (_halfPetal*2 > _segment_dphi ) { + std::cout << "delta phi between two segments: " << _segment_dphi << " half petal: " << _halfPetal << std::endl; + std::cout << "CEPCSegmentedDiscMeasLayer::CEPCSegmentedDiscMeasLayer overlaps: exit(1) called from " << __FILE__ << " line " << __LINE__ << std::endl; + exit(1); + } +} + + +TKalMatrix CEPCSegmentedDiscMeasLayer::XvToMv(const TVector3 &xv) const { + // Calculate measurement vector (hit coordinates) from global coordinates: + // coordinate matrix to return + TKalMatrix mv(ILDPlanarHit_DIM,1); + + int segmentIndex = get_segment_index(xv.Phi()); + + TVector3 XC = this->get_segment_centre(segmentIndex); + double rc = XC.Perp(); + double phic = XC.Phi(); + + double u = rc*(xv.Phi() - phic); + double v = xv.Perp() - rc; + + mv(0,0) = u; + mv(1,0) = v; + +#ifdef DEBUGPRINT + std::cout << "CEPCSegmentedDiscMeasLayer::XvToMv: phic = " << phic << " phi = " << xv.Phi() << " rc = " << rc << " r = " << xv.Perp() << std::endl; + std::cout << "CEPCSegmentedDiscMeasLayer::XvToMv: " + << " mv(0,0) = " << mv(0,0) + << " mv(1,0) = " << mv(1,0) + << std::endl; +#endif + return mv; +} + + +TVector3 CEPCSegmentedDiscMeasLayer::HitToXv(const TVTrackHit &vht) const { + // std::cout << "CEPCSegmentedDiscMeasLayer::HitToXv: " + // << " vht(0,0) = " << vht(0,0) << " vht(1,0) = " << vht(1,0) << std::endl; + + const ILDPlanarHit &mv = dynamic_cast<const ILDPlanarHit &>(vht); + +// double x = mv(0,0) ; +// double y = mv(1,0) ; +// +// double z = this->GetXc().Z() ; + + UTIL::BitField64 encoder(lcio::ILDCellID0::encoder_string); + edm4hep::TrackerHit hit = mv.getLCIOTrackerHit(); + encoder.setValue(hit.getCellID()); + int segmentIndex = encoder[lcio::ILDCellID0::module]; + + TVector3 XC = this->get_segment_centre(segmentIndex); + double rc = XC.Perp(); + double phic = XC.Phi(); + double zc = XC.Z(); + + double u = mv(0,0); + double v = mv(1,0); + + double phi = u/rc + phic; + double r = v + rc; + + double x = r*cos(phi); + double y = r*sin(phi); + double z = zc; + + // std::cout << "CEPCSegmentedDiscMeasLayer::HitToXv: " + // << " x = " << x + // << " y = " << y + // << " z = " << z + // << std::endl; + + return TVector3(x,y,z); +} + +void CEPCSegmentedDiscMeasLayer::CalcDhDa(const TVTrackHit &vht, + const TVector3 &xxv, + const TKalMatrix &dxphiada, + TKalMatrix &H) const { + // Calculate + // H = (@h/@a) = (@phi/@a, @z/@a)^t + // where + // h(a) = (phi, z)^t: expected meas vector + // a = (drho, phi0, kappa, dz, tanl, t0) + // + Int_t sdim = H.GetNcols(); + Int_t hdim = TMath::Max(5,sdim-1); + + // assume cylinder center at (0,0) + Double_t phiv = xxv.Phi(); + Double_t xv = xxv.X(); + Double_t yv = xxv.Y(); + Double_t xxyy = xv * xv + yv * yv; + + for (Int_t i = 0; i < hdim; i++) { + H(0, i) = - (yv / xxyy) * dxphiada(0, i) + (xv / xxyy) * dxphiada(1, i); + H(0, i) *= xxv.Perp(); + + H(1, i) = cos(phiv) * dxphiada(0, i) + sin(phiv) * dxphiada(1, i); + } + + if (sdim == 6) { + H(0,sdim-1) = 0.0; + H(1,sdim-1) = 0.; + } +} + +Int_t CEPCSegmentedDiscMeasLayer::CalcXingPointWith(const TVTrack &hel, + TVector3 &xx, + Double_t &phi, + Int_t mode, + Double_t eps) const{ + + // streamlog_out(DEBUG0) << "CEPCSegmentedDiscMeasLayer::CalcXingPointWith" << std::endl; + + phi = 0.0; + + xx.SetX(0.0); + xx.SetY(0.0); + xx.SetZ(0.0); + + // check that direction has one of the correct values + if( !( mode == 0 || mode == 1 || mode == -1) ) return -1 ; + +// + + // get helix parameters + Double_t dr = hel.GetDrho(); + Double_t phi0 = hel.GetPhi0(); // + Double_t kappa = hel.GetKappa(); + Double_t rho = hel.GetRho(); + Double_t omega = 1.0 / rho; + Double_t z0 = hel.GetDz(); + Double_t tanl = hel.GetTanLambda(); + + TVector3 ref_point = hel.GetPivot(); + + // + // Check if charge is nonzero. + // + + Int_t chg = (Int_t)TMath::Sign(1.1,kappa); + if (!chg) { + // streamlog_out(ERROR) << ">>>> Error >>>> CEPCSegmentedDiscMeasLayer::CalcXingPointWith" << std::endl + // << " Kappa = 0 is invalid for a helix " << std::endl; + return -1; + } + + const double sin_phi0 = sin(phi0); + const double cos_phi0 = cos(phi0); + + const double x_pca = ref_point.x() + dr * cos_phi0 ; + const double y_pca = ref_point.y() + dr * sin_phi0 ; + const double z_pca = ref_point.z() + z0 ; + + const double z = this->GetXc().Z() ; + // get path length to crossing point + + const double s = ( z - z_pca ) / tanl ; + + +// streamlog_out(DEBUG0) << "CEPCSegmentedDiscMeasLayer::CalcXingPointWith " +// << " ref_point.z() = " << ref_point.z() +// << " z = " << z +// << " z0 = " << z0 +// << " z_pca = " << z_pca +// << " tanl = " << tanl +// << " z - z_pca = " << z - z_pca +// << std::endl; +// +// TVector3 xx_n; +// int cuts = TVSurface::CalcXingPointWith(hel, xx_n, phi, mode, eps); +// streamlog_out(DEBUG0) << "CEPCSegmentedDiscMeasLayer::CalcXingPointWith from Newton: cuts = " << cuts << " x = " << xx_n.x() << " y = "<< xx_n.y() << " z = " << xx_n.z() << " r = " << xx_n.Perp() << " phi = " << xx_n.Phi() << " dphi = " << phi << std::endl; + + + phi = -omega * s; + + const double delta_phi_half = -phi/2.0 ; + + double x; + double y; + + if( fabs(s) > FLT_MIN ){ // protect against starting on the plane + + x = x_pca - s * ( sin(delta_phi_half) / delta_phi_half ) * sin( phi0 - delta_phi_half ) ; + + y = y_pca + s * ( sin(delta_phi_half) / delta_phi_half ) * cos( phi0 - delta_phi_half ) ; + + } + else{ + // streamlog_out(DEBUG0) << "CEPCSegmentedDiscMeasLayer::CalcXingPointWith Using PCA values " << std::endl; + x = x_pca; + y = y_pca; + phi = 0.0; + } + + + xx.SetXYZ(x, y, z); + + + // streamlog_out(DEBUG0) << "CEPCSegmentedDiscMeasLayer::CalcXingPointWith : cuts = " << (IsOnSurface(xx) && (chg*phi*mode)<0) << " x = " << xx.x() << " y = "<< xx.y() << " z = " << xx.z() << " r = " << xx.Perp() << " phi = " << xx.Phi() << " dphi = " << phi << " s = " << s << " " << this->TVMeasLayer::GetName() << std::endl; +// +// streamlog_out(DEBUG0) << "CEPCSegmentedDiscMeasLayer::CalcXingPointWith : xdiff = " << xx.x() - xx_n.x() << " ydiff = "<< xx.y() - xx_n.y() << " zdiff = " << xx.z() - xx_n.z() << std::endl; + + // check if intersection with plane is within boundaries + + if( mode!=0 && fabs(phi)>1.e-10 ){ // (+1,-1) = (fwd,bwd) + if( chg*phi*mode > 0){ + return 0; + } + } + + + return (IsOnSurface(xx) ? 1 : 0); + + +} + +Bool_t CEPCSegmentedDiscMeasLayer::IsOnSurface(const TVector3 &xx) const { + bool onSurface = false ; + + if (TMath::Abs(xx.Z()-GetXc().Z()) < 1e-4) { + //std::cout << "CEPCSegmentedDiscMeasLayer::IsOnSurface z passed " << std::endl; + double r2 = xx.Perp2(); + + // quick check to see weather the hit lies inside the min max r + if (r2 <= _rmax*_rmax && r2 >= _rmin*_rmin) { + //std::cout << "CEPCSegmentedDiscMeasLayer::IsOnSurface r2 passed " << std::endl; + + double phi_point = angular_range_2PI(xx.Phi()); + + // get the angle in the local system + double gamma = angular_range_2PI(phi_point - _start_phi); + + // the angle local to the sector + double local_phi = fmod(gamma, _segment_dphi); + + if (local_phi < 2*_halfPetal) { + //std::cout << "CEPCSegmentedDiscMeasLayer::IsOnSurface dphi passed " << std::endl; + onSurface = true ; + } + } + } +#ifdef DEBUGPRINT + if (!onSurface) { + std::cout << "IsOnSurface: zc = " << GetXc().Z() << " z = " << xx.Z() << " r = " << xx.Perp() << " phi = " << xx.Phi() << " phi0 = " << _start_phi + << " dphi = " << _segment_dphi << " local = " << fmod(xx.Phi()-_start_phi, _segment_dphi) << " half = " << _halfPetal << std::endl; + } +#endif + return onSurface; +} + + +ILDVTrackHit* CEPCSegmentedDiscMeasLayer::ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit) const { + std::bitset<32> type(trkhit.getType()); + // remember here the "position" of the hit in fact defines the origin of the plane it defines so u and v are per definition 0. + const edm4hep::Vector3d& pos=trkhit.getPosition(); + const TVector3 hit(pos.x, pos.y, pos.z) ; + + // convert to layer coordinates + TKalMatrix h(ILDPlanarHit_DIM,1); + + h = this->XvToMv(hit); + + double x[ILDPlanarHit_DIM] ; + double dx[ILDPlanarHit_DIM] ; + + x[0] = h(0, 0); + x[1] = h(1, 0); + + if (type[CEPCConf::TrkHitTypeBit::PLANAR]) { + dx[0] = trkhit.getCovMatrix(2); + dx[1] = trkhit.getCovMatrix(5); + } + else if (type[CEPCConf::TrkHitTypeBit::CYLINDER]) { + dx[0] = sqrt(trkhit.getCovMatrix(0)+trkhit.getCovMatrix(2)); + dx[1] = sqrt(trkhit.getCovMatrix(5)); + } + else { + dx[0] = sqrt(trkhit.getCovMatrix(0)+trkhit.getCovMatrix(2)); + dx[1] = sqrt(trkhit.getCovMatrix(5)); + } + + bool hit_on_surface = IsOnSurface(hit); + +#ifdef DEBUGPRINT + std::cout << "CEPCSegmentedDiscMeasLayer::ConvertLCIOTrkHit: ILDPlanarHit created" + << " for CellID " << trkhit.getCellID() + << " u = " << x[0] + << " v = " << x[1] + << " du = " << dx[0] + << " dv = " << dx[1] + << " x = " << pos.x + << " y = " << pos.y + << " z = " << pos.z + << " onSurface = " << hit_on_surface + << std::endl; +#endif + //ILDPlanarHit hh( *this , x, dx, this->GetBz(),trkhit); + + //this->HitToXv(hh); + + return hit_on_surface ? new ILDPlanarHit( *this , x, dx, this->GetBz(),trkhit) : NULL; + +} + + +/** Get the intersection and the CellID, needed for multilayers */ +int CEPCSegmentedDiscMeasLayer::getIntersectionAndCellID(const TVTrack &hel, + TVector3 &xx, + Double_t &phi, + Int_t &CellID, + Int_t mode, + Double_t eps) const { + int crosses = this->CalcXingPointWith(hel, xx, phi, mode, eps); + + if ( crosses != 0 ) { + unsigned int segment = this->get_segment_index(xx.Phi()); + + const std::vector<int>& cellIds = this->getCellIDs(); + + lcio::BitField64 bf(UTIL::ILDCellID0::encoder_string); + bf.setValue(this->getCellIDs()[0]); // get the first cell_ID, module = 0, sensor = 0 + + double phic = get_segment_phi(segment); + double local_phi = xx.Phi()-phic; + if (local_phi>M_PI) local_phi -= 2.0 * M_PI; + if (local_phi<-M_PI) local_phi += 2.0 * M_PI; + unsigned int isensor = get_sensor_index(xx.Perp(), local_phi); + bf[lcio::ILDCellID0::module] = segment;//cellIds.at(segment); + bf[lcio::ILDCellID0::sensor] = isensor; + CellID = bf.lowWord(); + } + + return crosses; +} + +unsigned int CEPCSegmentedDiscMeasLayer::get_segment_index(double phi) const { + phi = angular_range_2PI(phi-_start_phi); + return unsigned(floor(phi/_segment_dphi)); +} + +unsigned int CEPCSegmentedDiscMeasLayer::get_sensor_index(double r, double phi) const { + unsigned int irow = floor((r-_rmin)/(_rmax-_rmin)*_nrow); + if (irow >= _nrow) { + std::cout << "CEPCSegmentedDiscMeasLayer::get_sensor_index wrong row range: " << irow << std::endl; + exit(1); + } + unsigned int isensor = irow; + int nphi = _nsensors[irow]; + if (nphi>1) { + int iphi = floor((phi+_halfPetal)/(2*_halfPetal)*nphi); + for (int i=1; i<iphi; i++) { + for (int j=0; j<_nrow; j++) { + if (_nsensors[j]>i) isensor++; + } + } + } + + return isensor; +} + +double CEPCSegmentedDiscMeasLayer::get_segment_phi(unsigned int index) const{ + return angular_range_2PI(_start_phi + 0.5*_segment_dphi + index * _segment_dphi); +} + +TVector3 CEPCSegmentedDiscMeasLayer::get_segment_centre(unsigned int index) const{ + +// streamlog_out(DEBUG0) << "CEPCSegmentedDiscMeasLayer::get_segment_centre index = " << index << std::endl; + + double phi = this->get_segment_phi(index); + +// streamlog_out(DEBUG0) << "CEPCSegmentedDiscMeasLayer::get_segment_centre phi = " << phi << std::endl; + + + double rc = 0.5*(_rmin + _rmax); + +// streamlog_out(DEBUG0) << "CEPCSegmentedDiscMeasLayer::get_segment_centre rc = " << rc << std::endl; + + double xc = rc * cos(phi); + double yc = rc * sin(phi);; + + double zc = this->GetXc().Z(); + + TVector3 XC(xc,yc,zc); + + return XC; +} + +double CEPCSegmentedDiscMeasLayer::angular_range_2PI(double phi) const { + //bring phi_point into range 0 < phi < +2PI + while (phi < 0) { + phi += 2.0 * M_PI; + } + while (phi >= 2.0*M_PI) { + phi -= 2.0 * M_PI; + } + + return phi; +} + diff --git a/Utilities/KalDet/src/ild/common/ILDConeMeasLayer.cc b/Utilities/KalDet/src/ild/common/ILDConeMeasLayer.cc index 7ed8a826..f870e7a5 100644 --- a/Utilities/KalDet/src/ild/common/ILDConeMeasLayer.cc +++ b/Utilities/KalDet/src/ild/common/ILDConeMeasLayer.cc @@ -10,7 +10,7 @@ //************************************************************************* // -#include "ILDConeMeasLayer.h" +#include "kaldet/ILDConeMeasLayer.h" #include "TMath.h" diff --git a/Utilities/KalDet/src/ild/common/ILDConeMeasLayer.h b/Utilities/KalDet/src/ild/common/ILDConeMeasLayer.h deleted file mode 100644 index 19d617be..00000000 --- a/Utilities/KalDet/src/ild/common/ILDConeMeasLayer.h +++ /dev/null @@ -1,99 +0,0 @@ -#ifndef ILDCONEMEASLAYER_H -#define ILDCONEMEASLAYER_H -//************************************************************************* -//* =================== -//* ILDConeMeasLayer Class -//* =================== -//* -//* (Update Recored) -//* 2012/01/19 K.Fujii Original version. (EXBPConeMeasLayer) -//* 2012/01/24 R.Glattauer Adapted to ILD common in KalDet -//* -//************************************************************************* -// -#include "TVector3.h" - -#include "kaltest/TKalMatrix.h" -#include "kaltest/TCutCone.h" -#include "kaltest/KalTrackDim.h" - -#include "ILDVMeasLayer.h" -#include "iostream" -/* #include "streamlog/streamlog.h" */ -#include "UTIL/ILDConf.h" -#include "edm4hep/TrackerHit.h" - -class ILDConeMeasLayer : public ILDVMeasLayer, public TCutCone { -public: - // Ctors and Dtor - /** Constructor Taking inner and outer materials, z and radius at start and end, B-Field, whether the layer is sensitive, Cell ID, and an optional name */ - ILDConeMeasLayer(TMaterial &min, - TMaterial &mout, - Double_t z1, - Double_t r1, - Double_t z2, - Double_t r2, - Double_t Bz, - Double_t SortingPolicy, - Bool_t is_active, - Int_t CellID = -1, - const Char_t *name = "BPCONEML"); - virtual ~ILDConeMeasLayer(); - - // Parrent's pure virtuals that must be implemented - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVTrackHit &ht, - const TVector3 &xv) const; - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVector3 &xv) const; - - /** Local to Global coordinates */ - virtual TVector3 HitToXv (const TVTrackHit &ht) const; - - /** Calculate Projector Matrix */ - virtual void CalcDhDa (const TVTrackHit &ht, - const TVector3 &xv, - const TKalMatrix &dxphiada, - TKalMatrix &H) const; - - Bool_t IsOnSurface(const TVector3 &xx) const; - - /** Convert LCIO Tracker Hit to an ILDCylinderHit */ - virtual ILDVTrackHit* ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit) const { - - /* streamlog_out( ERROR ) << "Don't use this, it's not implemented!"; */ - return NULL; - } - - /** Get the intersection and the CellID, needed for multilayers */ - virtual int getIntersectionAndCellID(const TVTrack &hel, - TVector3 &xx, - Double_t &phi, - Int_t &CellID, - Int_t mode, - Double_t eps = 1.e-8) const { - - CellID = this->getCellIDs()[0]; // not multilayer - return this->CalcXingPointWith(hel,xx,phi,0,eps); - - - } - - /** Get sorting policy for this plane */ - virtual double GetSortingPolicy() const { return fsortingPolicy; } - - - - -private: - Double_t fZ1; // z of front face - Double_t fR1; // r of front face - Double_t fZ2; // z of back end - Double_t fR2; // r of back end - Double_t fsortingPolicy; // used for sorting the layers in to out - - -}; - -#endif diff --git a/Utilities/KalDet/src/ild/common/ILDCylinderHit.cc b/Utilities/KalDet/src/ild/common/ILDCylinderHit.cc index 27f1fc6f..d4a6832b 100644 --- a/Utilities/KalDet/src/ild/common/ILDCylinderHit.cc +++ b/Utilities/KalDet/src/ild/common/ILDCylinderHit.cc @@ -1,6 +1,6 @@ -#include "ILDCylinderHit.h" -#include "ILDCylinderMeasLayer.h" +#include "kaldet/ILDCylinderHit.h" +#include "kaldet/ILDCylinderMeasLayer.h" #include "TMath.h" #include <iostream> diff --git a/Utilities/KalDet/src/ild/common/ILDCylinderHit.h b/Utilities/KalDet/src/ild/common/ILDCylinderHit.h deleted file mode 100644 index 3fb9f1fc..00000000 --- a/Utilities/KalDet/src/ild/common/ILDCylinderHit.h +++ /dev/null @@ -1,38 +0,0 @@ -#ifndef ILDCYLINDERHIT_H -#define ILDCYLINDERHIT_H - -/** ILDCylinderHit: User defined KalTest hit class using R and Rphi coordinates, which provides coordinate vector as defined by the MeasLayer - * - * @author S.Aplin DESY - */ - -#include "kaltest/KalTrackDim.h" -#include "ILDVTrackHit.h" - - -class ILDCylinderHit : public ILDVTrackHit { - -public: - - - /** Constructor Taking R and Rphi coordinates and associated measurement layer, with bfield */ - ILDCylinderHit(const TVMeasLayer &ms, Double_t *x, Double_t *dx, - Double_t bfield, edm4hep::TrackerHit trkhit ) - : ILDVTrackHit(ms, x, dx, bfield, 2, trkhit) - { /* no op */ } - - - // TVTrackHit's pure virtuals that must be implemented - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv(const TVector3 &xv, Double_t t0) const; - - /** Print Debug information */ - virtual void DebugPrint(Option_t *opt = "") const; - - -private: - - -}; -#endif diff --git a/Utilities/KalDet/src/ild/common/ILDCylinderMeasLayer.cc b/Utilities/KalDet/src/ild/common/ILDCylinderMeasLayer.cc index ebf70e2f..5360f3e6 100644 --- a/Utilities/KalDet/src/ild/common/ILDCylinderMeasLayer.cc +++ b/Utilities/KalDet/src/ild/common/ILDCylinderMeasLayer.cc @@ -6,8 +6,8 @@ #include "kaltest/TKalTrack.h" -#include "ILDCylinderMeasLayer.h" -#include "ILDCylinderHit.h" +#include "kaldet/ILDCylinderMeasLayer.h" +#include "kaldet/ILDCylinderHit.h" #include <lcio.h> #include <edm4hep/TrackerHit.h> @@ -146,20 +146,19 @@ ILDVTrackHit* ILDCylinderMeasLayer::ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit bool hit_on_surface = IsOnSurface(hit); - //debug() << "ILDCylinderMeasLayer::ConvertLCIOTrkHit ILDCylinderHit created" + //std::cout << "ILDCylinderMeasLayer::ConvertLCIOTrkHit ILDCylinderHit created" // streamlog_out(DEBUG1) << "ILDCylinderMeasLayer::ConvertLCIOTrkHit ILDCylinderHit created" - // << " R = " << hit.Perp() - // << " Layer R = " << this->GetR() - // << " RPhi = " << x[0] - // << " Z = " << x[1] - // << " dRPhi = " << dx[0] - // << " dZ = " << dx[1] - // << " x = " << pos.x - // << " y = " << pos.y - // << " z = " << pos.z - // << " onSurface = " << hit_on_surface - // << std::endl ; - //<<endmsg; + // << " R = " << hit.Perp() + // << " Layer R = " << this->GetR() + // << " RPhi = " << x[0] + // << " Z = " << x[1] + // << " dRPhi = " << dx[0] + // << " dZ = " << dx[1] + // << " x = " << pos.x + // << " y = " << pos.y + // << " z = " << pos.z + // << " onSurface = " << hit_on_surface + // << std::endl; return hit_on_surface ? new ILDCylinderHit( *this , x, dx, this->GetBz(), trkhit) : NULL; diff --git a/Utilities/KalDet/src/ild/common/ILDCylinderMeasLayer.h b/Utilities/KalDet/src/ild/common/ILDCylinderMeasLayer.h deleted file mode 100644 index 00914f65..00000000 --- a/Utilities/KalDet/src/ild/common/ILDCylinderMeasLayer.h +++ /dev/null @@ -1,93 +0,0 @@ -#ifndef ILDCYLINDERMEASLAYER_H -#define ILDCYLINDERMEASLAYER_H - -/** ILDCylinderMeasLayer: User defined KalTest measurement layer class - * - * @author S.Aplin DESY - */ - - -#include "ILDVMeasLayer.h" -#include <iostream> -#include <cmath> -/* #include "streamlog/streamlog.h" */ - - -class ILDCylinderMeasLayer : public ILDVMeasLayer, public TCylinder { - -public: - - /** Constructor Taking inner and outer materials, radius and half length, B-Field, whether the layer is sensitive, Cell ID, and an optional name */ - ILDCylinderMeasLayer(TMaterial &min, - TMaterial &mout, - Double_t r0, - Double_t lhalf, - Double_t x0, - Double_t y0, - Double_t z0, - Double_t Bz, - Bool_t is_active, - Int_t CellID = -1, - const Char_t *name = "ILDCylinderMeasL") - : ILDVMeasLayer(min, mout, Bz, is_active, CellID, name), - TCylinder(r0, lhalf,x0,y0,z0) - { /* no op */ } - - - Bool_t IsOnSurface(const TVector3 &xx) const { - - bool z = (xx.Z() >= GetZmin() && xx.Z() <= GetZmax()); - bool r = std::fabs( (xx-this->GetXc()).Perp() - this->GetR() ) < 1.e-3; // for very short, very stiff tracks this can be poorly defined, so we relax this here a bit to 1 micron - -// streamlog_out(DEBUG0) << "ILDCylinderMeasLayer IsOnSurface for " << this->TVMeasLayer::GetName() << " R = " << this->GetR() << " GetZmin() = " << GetZmin() << " GetZmax() = " << GetZmax() -// << " dr = " << std::fabs( (xx-this->GetXc()).Perp() - this->GetR() ) << " r = " << r << " z = " << z -// << std::endl; - - return r && z; - } - - - - // Parent's pure virtuals that must be implemented - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVector3 &xv) const; - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVTrackHit &ht, - const TVector3 &xv) const - - { return this->XvToMv(xv); } - - - /** Local to Global coordinates */ - virtual TVector3 HitToXv (const TVTrackHit &ht) const; - - - /** Calculate Projector Matrix */ - virtual void CalcDhDa (const TVTrackHit &ht, - const TVector3 &xv, - const TKalMatrix &dxphiada, - TKalMatrix &H) const; - - /** Convert LCIO Tracker Hit to an ILDCylinderHit */ - virtual ILDVTrackHit* ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit) const ; - - /** Get the intersection and the CellID, needed for multilayers */ - virtual int getIntersectionAndCellID(const TVTrack &hel, - TVector3 &xx, - Double_t &phi, - Int_t &CellID, - Int_t mode, - Double_t eps = 1.e-8) const { - - CellID = this->getCellIDs()[0]; // not multilayer - return this->CalcXingPointWith(hel,xx,phi,0,eps); - - } - - -private: - -}; -#endif diff --git a/Utilities/KalDet/src/ild/common/ILDDiscMeasLayer.cc b/Utilities/KalDet/src/ild/common/ILDDiscMeasLayer.cc index cd0600e2..2bd0503e 100644 --- a/Utilities/KalDet/src/ild/common/ILDDiscMeasLayer.cc +++ b/Utilities/KalDet/src/ild/common/ILDDiscMeasLayer.cc @@ -1,8 +1,8 @@ #include <iostream> -#include "ILDDiscMeasLayer.h" -#include "ILDPlanarHit.h" +#include "kaldet/ILDDiscMeasLayer.h" +#include "kaldet/ILDPlanarHit.h" #include "kaltest/TVTrack.h" #include "TVector3.h" @@ -16,6 +16,8 @@ #include "gearimpl/Vector3D.h" +#include "DetIdentifier/CEPCConf.h" +#include <bitset> // #include "streamlog/streamlog.h" @@ -191,15 +193,20 @@ Bool_t ILDDiscMeasLayer::IsOnSurface(const TVector3 &xx) const TKalMatrix mv = XvToMv(xx); // check whether the hit lies in the same plane as the surface - if( TMath::Abs((xx.X()-GetXc().X())*GetNormal().X() + (xx.Y()-GetXc().Y())*GetNormal().Y() + (xx.Z()-GetXc().Z())*GetNormal().Z()) < 1e-4){ + if (TMath::Abs((xx.X()-GetXc().X())*GetNormal().X() + (xx.Y()-GetXc().Y())*GetNormal().Y() + (xx.Z()-GetXc().Z())*GetNormal().Z()) < 1e-4) { // check whether the hit lies within the boundary of the surface double r2 = mv(0,0) * mv(0,0) + mv(1,0) * mv(1,0) ; - if( r2 <= _rMax*_rMax && r2 >= _rMin*_rMin ) - { - onSurface = true ; - } + if (r2 <= _rMax*_rMax && r2 >= _rMin*_rMin) { + onSurface = true ; + } + else { + //std::cout << "r2: " << r2 << " r2min: " << _rMin*_rMin << " r2max: " << _rMax*_rMax << std::endl; + } + } + else { + //std::cout << "Xc: " << GetXc().X() << " " << GetXc().Y() << " " << GetXc().Z() << " Normal: " << GetNormal().X() << " " << GetNormal().Y() << " " << GetNormal().Z() << std::endl; } return onSurface; @@ -211,31 +218,46 @@ ILDVTrackHit* ILDDiscMeasLayer::ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit) co //edm4hep::TrackerHitPlane* plane_hit = dynamic_cast<EVENT::TrackerHitPlane*>( trkhit ) ; //edm4hep::TrackerHitPlane* plane_hit = trkhit; - if((trkhit.getType()&8)!=8) return NULL; + std::cout << "ILDDiscMeasLayer::ConvertLCIOTrkHit type = " << trkhit.getType() << std::endl; + std::bitset<32> type(trkhit.getType()); + //if (!type[CEPCConf::TrkHitTypeBit::PLANAR]) return NULL; //edm4hep::TrackerHit plane_hit = trkhit; //if( plane_hit == NULL ) return NULL; // SJA:FIXME: should be replaced with an exception //gear::Vector3D U(1.0,plane_hit.getU()[1],plane_hit.getU()[0],gear::Vector3D::spherical); //gear::Vector3D V(1.0,plane_hit.getV()[1],plane_hit.getV()[0],gear::Vector3D::spherical); - gear::Vector3D U(1.0,trkhit.getCovMatrix(1),trkhit.getCovMatrix(0),gear::Vector3D::spherical); - gear::Vector3D V(1.0,trkhit.getCovMatrix(5),trkhit.getCovMatrix(4),gear::Vector3D::spherical); - gear::Vector3D X(1.0,0.0,0.0); - gear::Vector3D Y(0.0,1.0,0.0); - - const float eps = 1.0e-07; + if (type[CEPCConf::TrkHitTypeBit::PLANAR]) { + gear::Vector3D U(1.0,trkhit.getCovMatrix(1),trkhit.getCovMatrix(0),gear::Vector3D::spherical); + gear::Vector3D V(1.0,trkhit.getCovMatrix(4),trkhit.getCovMatrix(3),gear::Vector3D::spherical); + gear::Vector3D X(1.0,0.0,0.0); + gear::Vector3D Y(0.0,1.0,0.0); + gear::Vector3D Z(0.0,0.0,1.0); + + const float eps = 1.0e-07; + // only require vertical to Z axis + if (U.dot(Z) > eps) { + std::cout << "ILDDiscMeasLayer: TrackerHit measurment vectors U is not vertical to the global Z axis. \n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << std::endl; + exit(1); + } + if (V.dot(Z) > eps) { + std::cout << "ILDDiscMeasLayer: TrackerHit measurment vectors V is not vertical to the global Z axis. \n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << std::endl; + exit(1); + } + } + /* // U must be the global X axis - if( fabs(1.0 - U.dot(X)) > eps ) { - // streamlog_out(ERROR) << "ILDDiscMeasLayer: TrackerHitPlane measurment vectors U is not equal to the global X axis. \n\n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << std::endl; + if( fabs(1.0 - U.dot(Y)) > eps ) { + std::cout << "ILDDiscMeasLayer: TrackerHitPlane measurment vectors U is not equal to the global Y axis. \n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << std::endl; exit(1); } // V must be the global X axis - if( fabs(1.0 - V.dot(Y)) > eps ) { - // streamlog_out(ERROR) << "ILDDiscMeasLayer: TrackerHitPlane measurment vectors V is not equal to the global Y axis. \n\n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << std::endl; + if( fabs(1.0 - V.dot(X)) > eps ) { + std::cout << "ILDDiscMeasLayer: TrackerHitPlane measurment vectors V is not equal to the global X axis. \n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << std::endl; exit(1); } - + */ const edm4hep::Vector3d& pos=trkhit.getPosition(); const TVector3 hit(pos.x, pos.y, pos.z); @@ -247,25 +269,33 @@ ILDVTrackHit* ILDDiscMeasLayer::ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit) co x[0] = h(0, 0); x[1] = h(1, 0); - - //dx[0] = plane_hit.getdU() ; - //dx[1] = plane_hit.getdV() ; - dx[0] = trkhit.getCovMatrix(2); - dx[1] = trkhit.getCovMatrix(5); + + if (type[CEPCConf::TrkHitTypeBit::PLANAR]) { + dx[0] = trkhit.getCovMatrix(2); + dx[1] = trkhit.getCovMatrix(5); + } + else if (type[CEPCConf::TrkHitTypeBit::CYLINDER]) { + dx[0] = sqrt(trkhit.getCovMatrix(0)+trkhit.getCovMatrix(2)); + dx[1] = sqrt(trkhit.getCovMatrix(5)); + } + else { + dx[0] = sqrt(trkhit.getCovMatrix(0)+trkhit.getCovMatrix(2)); + dx[1] = sqrt(trkhit.getCovMatrix(5)); + } bool hit_on_surface = IsOnSurface(hit); - - // streamlog_out(DEBUG1) << "ILDDiscMeasLayer::ConvertLCIOTrkHit ILDPlanarHit created" - // << " u = " << x[0] - // << " v = " << x[1] - // << " du = " << dx[0] - // << " dv = " << dx[1] - // << " x = " << pos.x - // << " y = " << pos.y - // << " z = " << pos.z - // << " onSurface = " << hit_on_surface - // << std::endl ; - +//#define DEBUG_CONVERT 1 +#ifdef DEBUG_CONVERT + std::cout << "ILDDiscMeasLayer::ConvertLCIOTrkHit ILDPlanarHit created" + << " u = " << x[0] + << " v = " << x[1] + << " du = " << dx[0] + << " dv = " << dx[1] + << " x = " << pos.x + << " y = " << pos.y + << " z = " << pos.z + << " onSurface = " << hit_on_surface + << std::endl ; +#endif return hit_on_surface ? new ILDPlanarHit( *this , x, dx, this->GetBz(), trkhit) : NULL; - } diff --git a/Utilities/KalDet/src/ild/common/ILDDiscMeasLayer.h b/Utilities/KalDet/src/ild/common/ILDDiscMeasLayer.h deleted file mode 100644 index 2f6f6167..00000000 --- a/Utilities/KalDet/src/ild/common/ILDDiscMeasLayer.h +++ /dev/null @@ -1,101 +0,0 @@ -#ifndef __ILDDISCMEASLAYER__ -#define __ILDDISCMEASLAYER__ - -/** ILDDiscMeasLayer: User defined KalTest Disc measurement layer class used with ILDPLanarTrackHit. - * - * @author S.Aplin DESY - */ - -#include "TVector3.h" - -#include "kaltest/TKalMatrix.h" -#include "kaltest/TPlane.h" -#include "kaltest/KalTrackDim.h" -#include "ILDVMeasLayer.h" - -#include "TMath.h" -#include <sstream> - -class TVTrackHit; - - -class ILDDiscMeasLayer : public ILDVMeasLayer, public TPlane { - -public: - - /** Constructor Taking inner and outer materials, center and normal to the plane, B-Field, sorting policy, min and max r, whether the layer is sensitive, Cell ID, and an optional name */ - - ILDDiscMeasLayer(TMaterial &min, - TMaterial &mout, - const TVector3 ¢er, - const TVector3 &normal, - double Bz, - double SortingPolicy, - double rMin, - double rMax, - Bool_t is_active, - Int_t CellID = -1, - const Char_t *name = "ILDDiscMeasL") - : ILDVMeasLayer(min, mout, Bz, is_active, CellID, name), - TPlane(center, normal), - _sortingPolicy(SortingPolicy), _rMin(rMin), _rMax(rMax) - { /* no op */ } - - - - // Parrent's pure virtuals that must be implemented - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVTrackHit &ht, - const TVector3 &xv) const - { return this->XvToMv(xv); } - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVector3 &xv) const; - - /** Local to Global coordinates */ - virtual TVector3 HitToXv (const TVTrackHit &ht) const; - - /** Calculate Projector Matrix */ - virtual void CalcDhDa (const TVTrackHit &ht, - const TVector3 &xv, - const TKalMatrix &dxphiada, - TKalMatrix &H) const; - - virtual Int_t CalcXingPointWith(const TVTrack &hel, - TVector3 &xx, - Double_t &phi, - Int_t mode, - Double_t eps) const; - - - /** Convert LCIO Tracker Hit to an ILDPLanarTrackHit */ - virtual ILDVTrackHit* ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit) const ; - - /** Check if global point is on surface */ - inline virtual Bool_t IsOnSurface (const TVector3 &xx) const; - - /** Get sorting policy for this plane */ - double GetSortingPolicy() const { return _sortingPolicy; } - - /** Get the intersection and the CellID, needed for multilayers */ - virtual int getIntersectionAndCellID(const TVTrack &hel, - TVector3 &xx, - Double_t &phi, - Int_t &CellID, - Int_t mode, - Double_t eps = 1.e-8) const { - - CellID = this->getCellIDs()[0]; // not multilayer - return this->CalcXingPointWith(hel,xx,phi,0,eps); - - } - -private: - double _sortingPolicy; - double _rMin; - double _rMax; - -}; - -#endif diff --git a/Utilities/KalDet/src/ild/common/ILDMeasurementSurfaceStoreFiller.cc b/Utilities/KalDet/src/ild/common/ILDMeasurementSurfaceStoreFiller.cc index 16fe6fe1..b5cf96ed 100644 --- a/Utilities/KalDet/src/ild/common/ILDMeasurementSurfaceStoreFiller.cc +++ b/Utilities/KalDet/src/ild/common/ILDMeasurementSurfaceStoreFiller.cc @@ -1,5 +1,5 @@ -#include "ILDMeasurementSurfaceStoreFiller.h" +#include "kaldet/ILDMeasurementSurfaceStoreFiller.h" #include "UTIL/ILDConf.h" diff --git a/Utilities/KalDet/src/ild/common/ILDMeasurementSurfaceStoreFiller.h b/Utilities/KalDet/src/ild/common/ILDMeasurementSurfaceStoreFiller.h deleted file mode 100644 index 79ff8e1e..00000000 --- a/Utilities/KalDet/src/ild/common/ILDMeasurementSurfaceStoreFiller.h +++ /dev/null @@ -1,81 +0,0 @@ -#ifndef ILDMEASUREMENTSURFACESTOREFILLER_H -#define ILDMEASUREMENTSURFACESTOREFILLER_H - -#include "gearsurf/MeasurementSurfaceStore.h" - -#include <vector> - -namespace gear{ - class GearMgr; - class ZPlanarParameters; - class FTDParameters; -} - -using namespace gear; - -class ILDMeasurementSurfaceStoreFiller : public MeasurementSurfaceStoreFiller{ - - public: - - - ILDMeasurementSurfaceStoreFiller(const gear::GearMgr& gear_mgr) : - _nVTXLayers(0), - _nSITLayers(0), - _nFTDLayers(0), - _nSETLayers(0) { - - this->get_gear_parameters(gear_mgr); - - } - - ~ILDMeasurementSurfaceStoreFiller() { /* no op */ } - - void getMeasurementSurfaces( std::vector<MeasurementSurface*>& surface_list ) const; - - std::string getName() const { return "ILDMeasurementSurfaceStoreFiller" ; } ; - - private: - - /** adds MeasurementSufaces to the store - * @param param: the ZPlanarParameters pointer of the detector, of which the measurement surfaces shall be added - * - * @param det_id: the detector id (as in ILDConf) - */ - void storeZPlanar( const gear::ZPlanarParameters* param , int det_id, std::vector<MeasurementSurface*>& surface_list ) const; - - void storeFTD( const gear::FTDParameters* param, std::vector<MeasurementSurface*>& surface_list ) const; - - - void get_gear_parameters(const gear::GearMgr& gear_mgr); - - /** the strip angles for every layer */ - std::vector< double > _VTXStripAngles; - std::vector< double > _SITStripAngles; - std::vector< double > _SETStripAngles; - - /** the number of sensors for every layer */ - std::vector< int > _VTXNSensors; - std::vector< int > _SITNSensors; - std::vector< int > _SETNSensors; - - - /** the strip angles for every layer and sensor */ - std::vector< std::vector< double > > _FTDStripAngles; - - unsigned _nVTXLayers; - unsigned _nSITLayers; - unsigned _nFTDLayers; - unsigned _nSETLayers; - - - const gear::ZPlanarParameters* _paramVXD; - const gear::ZPlanarParameters* _paramSIT; - const gear::ZPlanarParameters* _paramSET; - const gear::FTDParameters* _paramFTD; - - - -}; - -#endif - diff --git a/Utilities/KalDet/src/ild/common/ILDParallelPlanarMeasLayer.cc b/Utilities/KalDet/src/ild/common/ILDParallelPlanarMeasLayer.cc index 40eb8aeb..d03fd705 100644 --- a/Utilities/KalDet/src/ild/common/ILDParallelPlanarMeasLayer.cc +++ b/Utilities/KalDet/src/ild/common/ILDParallelPlanarMeasLayer.cc @@ -1,5 +1,5 @@ -#include "ILDParallelPlanarMeasLayer.h" +#include "kaldet/ILDParallelPlanarMeasLayer.h" #include "kaltest/TVTrack.h" // #include "streamlog/streamlog.h" diff --git a/Utilities/KalDet/src/ild/common/ILDParallelPlanarMeasLayer.h b/Utilities/KalDet/src/ild/common/ILDParallelPlanarMeasLayer.h deleted file mode 100644 index 7262a3e8..00000000 --- a/Utilities/KalDet/src/ild/common/ILDParallelPlanarMeasLayer.h +++ /dev/null @@ -1,80 +0,0 @@ -#ifndef __ILDParallelPlanarMeasLayer__ -#define __ILDParallelPlanarMeasLayer__ - -/** ILDParallelPlanarMeasLayer: User defined KalTest measurement layer class - * - * @author S.Aplin DESY - */ - - -#include "ILDPlanarMeasLayer.h" - -class ILDParallelPlanarMeasLayer : public ILDPlanarMeasLayer { - -public: - - /** Constructor Taking inner and outer materials, distance and phi of plane pca to origin, B-Field, Sorting policy, plane transverse witdth and offset of centre, longitudinal width, whether the layer is sensitive, Cell ID, and an optional name */ - ILDParallelPlanarMeasLayer(TMaterial &min, - TMaterial &mout, - Double_t r, - Double_t phi, - Double_t Bz, - Double_t SortingPolicy, - Double_t xiwidth, - Double_t zetawidth, - Double_t xioffset, - Double_t zoffset, - Double_t UOrigin, - Bool_t is_active, - Int_t CellID = -1, - const Char_t *name = "ILDParallelPlanarMeasLayer") - : - ILDPlanarMeasLayer(min,mout,TVector3(r*cos(phi),r*sin(phi),zoffset),TVector3(cos(phi),sin(phi),0.0),Bz,SortingPolicy,xiwidth,zetawidth,xioffset,UOrigin,is_active,CellID,name), _r(r),_phi(phi),_cos_phi(cos(_phi)),_sin_phi(sin(_phi)) - { /* no op */ } - - - // Parent's pure virtuals that must be implemented - - /** overloaded version of CalcXingPointWith using closed solution */ - virtual Int_t CalcXingPointWith(const TVTrack &hel, - TVector3 &xx, - Double_t &phi, - Int_t mode, - Double_t eps = 1.e-8) const; - - /** overloaded version of CalcXingPointWith using closed solution */ - virtual Int_t CalcXingPointWith(const TVTrack &hel, - TVector3 &xx, - Double_t &phi, - Double_t eps = 1.e-8) const{ - - return CalcXingPointWith(hel,xx,phi,0,eps); - - } - - /** Get the intersection and the CellID, needed for multilayers */ - virtual int getIntersectionAndCellID(const TVTrack &hel, - TVector3 &xx, - Double_t &phi, - Int_t &CellID, - Int_t mode, - Double_t eps = 1.e-8) const { - - CellID = this->getCellIDs()[0]; // not multilayer - return CalcXingPointWith(hel,xx,phi,0,eps); - - } - - -protected: - - Double_t _r; - Double_t _phi; - Double_t _cos_phi; - Double_t _sin_phi; - - -}; - -#endif - diff --git a/Utilities/KalDet/src/ild/common/ILDParallelPlanarStripMeasLayer.cc b/Utilities/KalDet/src/ild/common/ILDParallelPlanarStripMeasLayer.cc index d63cfe81..316a81b3 100644 --- a/Utilities/KalDet/src/ild/common/ILDParallelPlanarStripMeasLayer.cc +++ b/Utilities/KalDet/src/ild/common/ILDParallelPlanarStripMeasLayer.cc @@ -1,6 +1,6 @@ -#include "ILDParallelPlanarStripMeasLayer.h" -#include "ILDPlanarStripHit.h" +#include "kaldet/ILDParallelPlanarStripMeasLayer.h" +#include "kaldet/ILDPlanarStripHit.h" #include "kaltest/TVTrack.h" #include "kaltest/TVTrackHit.h" diff --git a/Utilities/KalDet/src/ild/common/ILDParallelPlanarStripMeasLayer.h b/Utilities/KalDet/src/ild/common/ILDParallelPlanarStripMeasLayer.h deleted file mode 100644 index f1bd711c..00000000 --- a/Utilities/KalDet/src/ild/common/ILDParallelPlanarStripMeasLayer.h +++ /dev/null @@ -1,62 +0,0 @@ -#ifndef __ILDParallelPlanarStripMeasLayer__ -#define __ILDParallelPlanarStripMeasLayer__ - -/** ILDParallelPlanarStripMeasLayer: User defined KalTest measurement layer class - * - * @author S.Aplin DESY - */ - -//#include "TKalMatrix.h" -//#include "TVector3.h" -//#include "TVTrackHit.h" - - -#include "ILDParallelPlanarMeasLayer.h" - -class ILDParallelPlanarStripMeasLayer : public ILDParallelPlanarMeasLayer { - -public: - - /** Constructor Taking inner and outer materials, distance and phi of plane pca to origin, B-Field, Sorting policy, plane transverse witdth and offset of centre, longitudinal width, whether the layer is sensitive, Cell ID, and an optional name */ - ILDParallelPlanarStripMeasLayer(TMaterial &min, - TMaterial &mout, - Double_t r, - Double_t phi, - Double_t Bz, - Double_t SortingPolicy, - Double_t xiwidth, - Double_t zetawidth, - Double_t xioffset, - Double_t zoffset, - Double_t UOrigin, - Double_t stripAngle, - Int_t CellID = -1, - const Char_t *name = "ILDParallelPlanarStripMeasLayer") - : - ILDParallelPlanarMeasLayer(min,mout,r,phi,Bz,SortingPolicy,xiwidth,zetawidth,xioffset,zoffset,UOrigin,true,CellID,name), _stripAngle(stripAngle) - - { /* no op */ } - - - // Parent's pure virtuals that must be implemented - - TKalMatrix XvToMv(const TVector3 &xv) const; - - TKalMatrix XvToMv(const TVTrackHit &, const TVector3 &xv) const { - return XvToMv(xv); - } - - TVector3 HitToXv(const TVTrackHit &vht) const ; - - void CalcDhDa(const TVTrackHit &vht, const TVector3 &xxv, const TKalMatrix &dxphiada, TKalMatrix &H) const; - - ILDVTrackHit* ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit) const; - -private: - - double _stripAngle; - -}; - -#endif - diff --git a/Utilities/KalDet/src/ild/common/ILDPlanarHit.cc b/Utilities/KalDet/src/ild/common/ILDPlanarHit.cc index c4fc4ae1..27fc2ac7 100644 --- a/Utilities/KalDet/src/ild/common/ILDPlanarHit.cc +++ b/Utilities/KalDet/src/ild/common/ILDPlanarHit.cc @@ -1,8 +1,8 @@ -#include "ILDPlanarHit.h" -#include "ILDPlanarMeasLayer.h" -#include "ILDSegmentedDiscMeasLayer.h" -#include "ILDDiscMeasLayer.h" +#include "kaldet/ILDPlanarHit.h" +#include "kaldet/ILDPlanarMeasLayer.h" +#include "kaldet/ILDSegmentedDiscMeasLayer.h" +#include "kaldet/ILDDiscMeasLayer.h" #include "TMath.h" #include <iostream> diff --git a/Utilities/KalDet/src/ild/common/ILDPlanarHit.h b/Utilities/KalDet/src/ild/common/ILDPlanarHit.h deleted file mode 100644 index f3b0bcce..00000000 --- a/Utilities/KalDet/src/ild/common/ILDPlanarHit.h +++ /dev/null @@ -1,41 +0,0 @@ -#ifndef ILDPLANARHIT_H -#define ILDPLANARHIT_H - -/** ILDPlanarHit: User defined KalTest hit class using u and v coordinates, which provides coordinate vector as defined by the MeasLayer - * - * @author S.Aplin DESY - */ - -#include "kaltest/KalTrackDim.h" - -#include "ILDVTrackHit.h" - -#define ILDPlanarHit_DIM 2 - -class ILDPlanarHit : public ILDVTrackHit { - -public: - - /** Constructor Taking u and v coordinates and associated measurement layer, with bfield */ - ILDPlanarHit(const TVMeasLayer &ms, - Double_t *x, - Double_t *dx, - Double_t bfield, - edm4hep::TrackerHit trkhit) - : ILDVTrackHit(ms, x, dx, bfield, ILDPlanarHit_DIM,trkhit) - { /* no op */ } - - // TVTrackHit's pure virtuals that must be implemented - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVector3 &xv, Double_t t0) const; - - /** Print Debug information */ - virtual void DebugPrint(Option_t *opt = "") const; - - -private: - - -}; -#endif diff --git a/Utilities/KalDet/src/ild/common/ILDPlanarMeasLayer.cc b/Utilities/KalDet/src/ild/common/ILDPlanarMeasLayer.cc index 35c4034d..c9984dfe 100644 --- a/Utilities/KalDet/src/ild/common/ILDPlanarMeasLayer.cc +++ b/Utilities/KalDet/src/ild/common/ILDPlanarMeasLayer.cc @@ -15,8 +15,8 @@ #include <iostream> #include <cmath> -#include "ILDPlanarMeasLayer.h" -#include "ILDPlanarHit.h" +#include "kaldet/ILDPlanarMeasLayer.h" +#include "kaldet/ILDPlanarHit.h" #include "kaltest/TVTrack.h" #include "TVector3.h" diff --git a/Utilities/KalDet/src/ild/common/ILDPlanarMeasLayer.h b/Utilities/KalDet/src/ild/common/ILDPlanarMeasLayer.h deleted file mode 100644 index f3bab1b7..00000000 --- a/Utilities/KalDet/src/ild/common/ILDPlanarMeasLayer.h +++ /dev/null @@ -1,84 +0,0 @@ -#ifndef __ILDPLANARMEASLAYER__ -#define __ILDPLANARMEASLAYER__ -//************************************************************************* -//* =================== -//* ILDPlanarMeasLayer Class -//* =================== -//* -//* (Description) -//* Planar measurement layer class used with ILDPLanarTrackHit. -//* (Requires) -//* ILDVMeasLayer -//* (Provides) -//* class ILDPlanarMeasLayer -//* (Update Recored) -//* 2003/09/30 Y.Nakashima Original version. -//* -//* 2011/06/17 D.Kamai Modified to handle ladder structure. -//************************************************************************* -// -#include "TVector3.h" - -#include "kaltest/TKalMatrix.h" -#include "kaltest/TPlane.h" -#include "kaltest/KalTrackDim.h" -#include "ILDVMeasLayer.h" - -#include "TMath.h" -#include <sstream> - -class TVTrackHit; - -class ILDPlanarMeasLayer : public ILDVMeasLayer, public TPlane { -public: - // Ctors and Dtor - - ILDPlanarMeasLayer(TMaterial &min, - TMaterial &mout, - const TVector3 ¢er, - const TVector3 &normal, - Double_t Bz, - Double_t SortingPolicy, - Double_t xiwidth, - Double_t zetawidth, - Double_t xioffset, - Double_t fUOrigin, - Bool_t is_active, - Int_t CellID = -1, - const Char_t *name = "ILDPlanarMeasL"); - - virtual ~ILDPlanarMeasLayer(); - - // Parrent's pure virtuals that must be implemented - - virtual TKalMatrix XvToMv (const TVTrackHit &ht, - const TVector3 &xv) const; - - virtual TKalMatrix XvToMv (const TVector3 &xv) const; - - virtual TVector3 HitToXv (const TVTrackHit &ht) const; - - virtual void CalcDhDa (const TVTrackHit &ht, - const TVector3 &xv, - const TKalMatrix &dxphiada, - TKalMatrix &H) const; - - virtual ILDVTrackHit* ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit) const ; - - virtual Bool_t IsOnSurface (const TVector3 &xx) const; - - Double_t GetSortingPolicy() const { return fSortingPolicy; } - Double_t GetXiwidth() const { return fXiwidth; } - Double_t GetZetawidth() const { return fZetawidth; } - Double_t GetXioffset() const { return fXioffset; } - -protected: - Double_t fSortingPolicy; - Double_t fXiwidth; - Double_t fZetawidth; - Double_t fXioffset; //determines how far the centre of the plane is translated in the direction positive rotation - Double_t fUOrigin; //determines origin of the transverse coordinate - -}; - -#endif diff --git a/Utilities/KalDet/src/ild/common/ILDPlanarStripHit.cc b/Utilities/KalDet/src/ild/common/ILDPlanarStripHit.cc index 613e34df..504633c4 100644 --- a/Utilities/KalDet/src/ild/common/ILDPlanarStripHit.cc +++ b/Utilities/KalDet/src/ild/common/ILDPlanarStripHit.cc @@ -1,8 +1,8 @@ -#include "ILDPlanarStripHit.h" -#include "ILDPlanarMeasLayer.h" -#include "ILDSegmentedDiscMeasLayer.h" -#include "ILDDiscMeasLayer.h" +#include "kaldet/ILDPlanarStripHit.h" +#include "kaldet/ILDPlanarMeasLayer.h" +#include "kaldet/ILDSegmentedDiscMeasLayer.h" +#include "kaldet/ILDDiscMeasLayer.h" #include "TMath.h" #include <iostream> diff --git a/Utilities/KalDet/src/ild/common/ILDPlanarStripHit.h b/Utilities/KalDet/src/ild/common/ILDPlanarStripHit.h deleted file mode 100644 index fe44a2dc..00000000 --- a/Utilities/KalDet/src/ild/common/ILDPlanarStripHit.h +++ /dev/null @@ -1,42 +0,0 @@ -#ifndef ILDPLANARSTRIPHIT_H -#define ILDPLANARSTRIPHIT_H - -/** ILDPlanarStripHit: User defined KalTest hit class using u coordinate, which provides coordinate vector as defined by the MeasLayer - * - * @author S.Aplin DESY - */ - -#include "kaltest/KalTrackDim.h" - -#include "ILDVTrackHit.h" - - -#define ILDPlanarStripHit_DIM 1 // set to 2 if one want to debug strip hits by using the 2nd dimention - -class ILDPlanarStripHit : public ILDVTrackHit { - -public: - - /** Constructor Taking a single coordinate and associated measurement layer, with bfield */ - ILDPlanarStripHit(const TVMeasLayer &ms, - Double_t *x, - Double_t *dx, - Double_t bfield, - edm4hep::TrackerHit trkhit) - : ILDVTrackHit(ms, x, dx, bfield, ILDPlanarStripHit_DIM,trkhit) - { /* no op */ } - - // TVTrackHit's pure virtuals that must be implemented - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVector3 &xv, Double_t t0) const; - - /** Print Debug information */ - virtual void DebugPrint(Option_t *opt = "") const; - - -private: - - -}; -#endif diff --git a/Utilities/KalDet/src/ild/common/ILDPolygonBarrelMeasLayer.cc b/Utilities/KalDet/src/ild/common/ILDPolygonBarrelMeasLayer.cc index d2668725..3c6676b2 100644 --- a/Utilities/KalDet/src/ild/common/ILDPolygonBarrelMeasLayer.cc +++ b/Utilities/KalDet/src/ild/common/ILDPolygonBarrelMeasLayer.cc @@ -1,5 +1,5 @@ -#include "ILDPolygonBarrelMeasLayer.h" +#include "kaldet/ILDPolygonBarrelMeasLayer.h" #include <UTIL/BitField64.h> #include <UTIL/ILDConf.h> diff --git a/Utilities/KalDet/src/ild/common/ILDPolygonBarrelMeasLayer.h b/Utilities/KalDet/src/ild/common/ILDPolygonBarrelMeasLayer.h deleted file mode 100644 index 791d3dbf..00000000 --- a/Utilities/KalDet/src/ild/common/ILDPolygonBarrelMeasLayer.h +++ /dev/null @@ -1,138 +0,0 @@ -#ifndef __ILDSEGMENTEDDISCMEASLAYER_H__ -#define __ILDSEGMENTEDDISCMEASLAYER_H__ - -/** ILDPolygonBarrelMeasLayer: User defined Polygonal Barrel KalTest measurement layer class to be used only for dead material. Segments are planes parallel to the z axis - * - * NOTE: ALL METHODS INVOLVING HITS ARE DISABLED AND CALL EXIT(1) - * THIS CLASS IS ONLY MEANT FOR DEAD MATERIAL - * - * @author S.Aplin DESY - */ - -#include "TVector3.h" - -#include "kaltest/TKalMatrix.h" -#include "kaltest/TVSurface.h" -#include "kaltest/KalTrackDim.h" -#include "ILDVMeasLayer.h" - -#include "TMath.h" -#include <sstream> - -#include "ILDParallelPlanarMeasLayer.h" - -#include <vector> - -class TVTrackHit; - - -class ILDPolygonBarrelMeasLayer : public ILDVMeasLayer, public TVSurface { -public: - // Ctors and Dtor - - - - ILDPolygonBarrelMeasLayer(TMaterial &min, - TMaterial &mout, - double Bz, - double SortingPolicy, - double r0, // min distance to the z-axis - double lhalf, // half length - int nsides, - double zpos, // z of the centre - double phi0, // phi of the first normal following the xaxis positive rotation - std::vector<int> CellIDs, - bool is_active, - const Char_t *name = "ILDPolygonBarrelMeasL"); - - - ~ILDPolygonBarrelMeasLayer(){ delete _enclosing_cylinder;} - - - // Parrent's pure virtuals that must be implemented - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVTrackHit &ht, - const TVector3 &xv) const - { return this->XvToMv(xv); } - - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVector3 &xv) const; - - /** Local to Global coordinates */ - virtual TVector3 HitToXv (const TVTrackHit &ht) const; - - /** Calculate Projector Matrix */ - virtual void CalcDhDa (const TVTrackHit &ht, - const TVector3 &xv, - const TKalMatrix &dxphiada, - TKalMatrix &H) const; - - /** Convert LCIO Tracker Hit to an ILDPLanarTrackHit */ - virtual ILDVTrackHit* ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit) const ; - - /** overloaded version of CalcXingPointWith using closed solution*/ - virtual Int_t CalcXingPointWith(const TVTrack &hel, - TVector3 &xx, - Double_t &phi, - Int_t mode, - Double_t eps = 1.e-8) const; - - /** overloaded version of CalcXingPointWith using closed solution*/ - virtual Int_t CalcXingPointWith(const TVTrack &hel, - TVector3 &xx, - Double_t &phi, - Double_t eps = 1.e-8) const{ - - return CalcXingPointWith(hel,xx,phi,0,eps); - - } - - - /** Get the intersection and the CellID, needed for multilayers */ - virtual int getIntersectionAndCellID(const TVTrack &hel, - TVector3 &xx, - Double_t &phi, - Int_t &CellID, - Int_t mode, - Double_t eps = 1.e-8) const ; - - - bool IsOutside(const TVector3 &xx) const; - - double CalcS(const TVector3 &xx) const; - - TMatrixD CalcDSDx(const TVector3 &xx) const; - - /** Check if global point is on surface */ - inline virtual Bool_t IsOnSurface (const TVector3 &xx) const; - - /** Get sorting policy for this plane */ - double GetSortingPolicy() const { return _sortingPolicy; } - -private: - - double angular_range_2PI( double phi ) const; - - unsigned int get_plane_index(double phi) const; - - double _sortingPolicy; - - double _r0; // min distance to the z-axis - double _lhalf; // half length - int _nsides; - double _zpos; // z of the centre - double _phi0; // phi of the first normal following the xaxis positive rotation - - double _segment_dphi; - double _start_phi; - double _rmax; - - std::vector<ILDParallelPlanarMeasLayer> _planes; - TCylinder* _enclosing_cylinder; -}; - - - -#endif diff --git a/Utilities/KalDet/src/ild/common/ILDRotatedTrapMeaslayer.cc b/Utilities/KalDet/src/ild/common/ILDRotatedTrapMeaslayer.cc index 7935ff81..2100d413 100644 --- a/Utilities/KalDet/src/ild/common/ILDRotatedTrapMeaslayer.cc +++ b/Utilities/KalDet/src/ild/common/ILDRotatedTrapMeaslayer.cc @@ -1,7 +1,7 @@ #include <iostream> -#include "ILDRotatedTrapMeaslayer.h" -#include "ILDPlanarHit.h" +#include "kaldet/ILDRotatedTrapMeaslayer.h" +#include "kaldet/ILDPlanarHit.h" #include "kaltest/TVTrack.h" #include "TVector3.h" diff --git a/Utilities/KalDet/src/ild/common/ILDRotatedTrapMeaslayer.h b/Utilities/KalDet/src/ild/common/ILDRotatedTrapMeaslayer.h deleted file mode 100644 index b375467d..00000000 --- a/Utilities/KalDet/src/ild/common/ILDRotatedTrapMeaslayer.h +++ /dev/null @@ -1,106 +0,0 @@ -#ifndef __ILDPLANARMEASLAYER__ -#define __ILDPLANARMEASLAYER__ -/** ILDRotatedTrapMeaslayer: User defined Rotated Trapezoid Planar KalTest measurement layer class used with ILDPLanarTrackHit. - * - * @author S.Aplin DESY - */ - -#include "TVector3.h" - -#include "kaltest/TKalMatrix.h" -#include "kaltest/TPlane.h" -#include "kaltest/KalTrackDim.h" -#include "ILDVMeasLayer.h" - -#include "TMath.h" -#include <sstream> -class TVTrackHit; - -class ILDRotatedTrapMeaslayer : public ILDVMeasLayer, public TPlane { - -public: - - /** Constructor Taking inner and outer materials, centre and normal of the plane, B-Field, sorting policy, height, inner base length, outer base length, tilt angle around axis of symmetry, represents full or half petal, whether the layer is sensitive, Cell ID, and an optional name */ - ILDRotatedTrapMeaslayer(TMaterial &min, - TMaterial &mout, - const TVector3 ¢er, - const TVector3 &normal, - Double_t Bz, - Double_t SortingPolicy, - Double_t height, - Double_t innerBaseLength, - Double_t outerBaseLength, - Double_t alpha, - Int_t half_petal, - Bool_t is_active, - Int_t CellID = -1, - const Char_t *name = "ILDRotatedTrapMeasL"); - - - // Parrent's pure virtuals that must be implemented - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVTrackHit &ht, - const TVector3 &xv) const - { return this->XvToMv(xv); } - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVector3 &xv) const; - - /** Local to Global coordinates */ - virtual TVector3 HitToXv (const TVTrackHit &ht) const; - - /** Calculate Projector Matrix */ - virtual void CalcDhDa (const TVTrackHit &ht, - const TVector3 &xv, - const TKalMatrix &dxphiada, - TKalMatrix &H) const; - - /** Convert LCIO Tracker Hit to an ILDPLanarTrackHit */ - virtual ILDVTrackHit* ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit) const ; - - /** Check if global point is on surface */ - inline virtual Bool_t IsOnSurface (const TVector3 &xx) const; - - /** Get sorting policy for this plane */ - Double_t GetSortingPolicy() const { return _sortingPolicy; } - - /** Get the intersection and the CellID, needed for multilayers */ - virtual int getIntersectionAndCellID(const TVTrack &hel, - TVector3 &xx, - Double_t &phi, - Int_t &CellID, - Int_t mode, - Double_t eps = 1.e-8) const { - - CellID = this->getCellIDs()[0]; // not multilayer - return this->CalcXingPointWith(hel,xx,phi,0,eps); - - } - - - -private: - Double_t _sortingPolicy ; - Double_t _signZ ; - Double_t _innerR ; - Double_t _outerR ; - Double_t _innerBaseLength ; - Double_t _outerBaseLength ; - Double_t _cosPhi ; //** cos of the azimuthal angle of the petal - Double_t _sinPhi ; //** sin of the azimuthal angle of the petal - Double_t _cosAlpha ; //** cos of the tilt angle of the petal - Double_t _sinAlpha ; //** sin of the tilt angle of the petal - Double_t _tanBeta ; //** tan of the openning angle of the petal - - // meaning of _halfPetal: - // 0 complete trapezoid - // +1 positive half only, i.e. the side of the petal in which the transverse coordinate, mv(0,0), is positive - // -1 negative half only, i.e. the side of the petal in which the transverse coordinate, mv(0,0), is negative - Int_t _halfPetal ; - - - -}; - -#endif diff --git a/Utilities/KalDet/src/ild/common/ILDSegmentedDiscMeasLayer.cc b/Utilities/KalDet/src/ild/common/ILDSegmentedDiscMeasLayer.cc index 08c5072a..47417c72 100644 --- a/Utilities/KalDet/src/ild/common/ILDSegmentedDiscMeasLayer.cc +++ b/Utilities/KalDet/src/ild/common/ILDSegmentedDiscMeasLayer.cc @@ -1,6 +1,6 @@ -#include "ILDSegmentedDiscMeasLayer.h" -#include "ILDPlanarHit.h" +#include "kaldet/ILDSegmentedDiscMeasLayer.h" +#include "kaldet/ILDPlanarHit.h" #include <UTIL/BitField64.h> #include <UTIL/ILDConf.h> diff --git a/Utilities/KalDet/src/ild/common/ILDSegmentedDiscStripMeasLayer.cc b/Utilities/KalDet/src/ild/common/ILDSegmentedDiscStripMeasLayer.cc index 60a3499b..24660abd 100644 --- a/Utilities/KalDet/src/ild/common/ILDSegmentedDiscStripMeasLayer.cc +++ b/Utilities/KalDet/src/ild/common/ILDSegmentedDiscStripMeasLayer.cc @@ -1,7 +1,6 @@ -#include "ILDSegmentedDiscStripMeasLayer.h" - -#include "ILDPlanarStripHit.h" +#include "kaldet/ILDSegmentedDiscStripMeasLayer.h" +#include "kaldet/ILDPlanarStripHit.h" #include <UTIL/BitField64.h> #include <UTIL/ILDConf.h> diff --git a/Utilities/KalDet/src/ild/common/ILDSegmentedDiscStripMeasLayer.h b/Utilities/KalDet/src/ild/common/ILDSegmentedDiscStripMeasLayer.h deleted file mode 100644 index 9031026b..00000000 --- a/Utilities/KalDet/src/ild/common/ILDSegmentedDiscStripMeasLayer.h +++ /dev/null @@ -1,80 +0,0 @@ -#ifndef __ILDSEGMENTEDDISCSTRIPMEASLAYER_H__ -#define __ILDSEGMENTEDDISCSTRIPMEASLAYER_H__ - -/** ILDSegmentedDiscStripMeasLayer: User defined Segemented Disk Planar KalTest measurement layer class used with ILDPLanarTrackHit. Segments are isosolese trapezoids whose axis of symmetry points to the origin - * WARNING: ONLY IMPLEMENTED FOR X AND Y COORDINATES AT FIXED Z - * - * @author S.Aplin DESY - */ -#include "TVector3.h" - -#include "kaltest/TKalMatrix.h" -#include "kaltest/TPlane.h" -#include "kaltest/KalTrackDim.h" -#include "ILDSegmentedDiscMeasLayer.h" - -#include "TMath.h" -#include <sstream> -#include <iostream> - -#include <vector> - -class TVTrackHit; - -class ILDSegmentedDiscStripMeasLayer : public ILDSegmentedDiscMeasLayer { - -public: - // Ctors and Dtor - - ILDSegmentedDiscStripMeasLayer(TMaterial &min, - TMaterial &mout, - double Bz, - double SortingPolicy, - int nsegments, - double zpos, - double phi0, // defined by the axis of symmerty of the first petal - double trap_rmin, - double trap_height, - double trap_innerBaseLength, - double trap_outerBaseLength, - double stripAngle, - bool is_active, - std::vector<int> CellIDs, - const Char_t *name = "ILDDiscMeasL") - : ILDSegmentedDiscMeasLayer(min,mout,Bz,SortingPolicy,nsegments,zpos,phi0,trap_rmin,trap_height,trap_innerBaseLength,trap_outerBaseLength,is_active,CellIDs,name), - _stripAngle(stripAngle) - { /* no op */ } - - - // Parrent's pure virtuals that must be implemented - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVTrackHit &ht, - const TVector3 &xv) const - { return this->XvToMv(xv); } - - - /** Global to Local coordinates */ - virtual TKalMatrix XvToMv (const TVector3 &xv) const; - - /** Local to Global coordinates */ - virtual TVector3 HitToXv (const TVTrackHit &ht) const; - - /** Calculate Projector Matrix */ - virtual void CalcDhDa (const TVTrackHit &ht, - const TVector3 &xv, - const TKalMatrix &dxphiada, - TKalMatrix &H) const; - - /** Convert LCIO Tracker Hit to an ILDPLanarTrackHit */ - virtual ILDVTrackHit* ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit) const ; - -private: - - double _stripAngle; - -}; - - - -#endif diff --git a/Utilities/KalDet/src/ild/common/ILDVMeasLayer.cc b/Utilities/KalDet/src/ild/common/ILDVMeasLayer.cc index 67705851..de36c3d6 100644 --- a/Utilities/KalDet/src/ild/common/ILDVMeasLayer.cc +++ b/Utilities/KalDet/src/ild/common/ILDVMeasLayer.cc @@ -1,5 +1,5 @@ -#include "ILDVMeasLayer.h" +#include "kaldet/ILDVMeasLayer.h" #include <UTIL/BitField64.h> #include <UTIL/ILDConf.h> @@ -46,7 +46,7 @@ _isMultiLayer(true) { if (cellIDs.size() == 0 ) { - // streamlog_out(ERROR) << __FILE__ << " line " << __LINE__ << " size of cellIDs == 0" << std::endl; + std::cout << __FILE__ << " line " << __LINE__ << " size of cellIDs == 0" << std::endl; } UTIL::BitField64 encoder( lcio::ILDCellID0::encoder_string ) ; diff --git a/Utilities/KalDet/src/ild/common/ILDVMeasLayer.h b/Utilities/KalDet/src/ild/common/ILDVMeasLayer.h deleted file mode 100644 index edd81744..00000000 --- a/Utilities/KalDet/src/ild/common/ILDVMeasLayer.h +++ /dev/null @@ -1,89 +0,0 @@ -#ifndef __ILDVMEASLAYER__ -#define __ILDVMEASLAYER__ - -/** ILDVMeasLayer: Virtual measurement layer class used by ILD[X]MeasLayer Classes. - * - * @author S.Aplin DESY - */ - -#include "TVector3.h" -#include "kaltest/TKalMatrix.h" -#include "kaltest/TCylinder.h" -#include "kaltest/TVMeasLayer.h" -#include "kaltest/TAttDrawable.h" -#include "kaltest/KalTrackDim.h" -#include "TString.h" -#include "edm4hep/TrackerHit.h" - -#include <vector> - -class TVTrackHit; -class TNode; -class ILDVTrackHit; - -namespace edm4hep{ - class TrackerHit; -} - -class ILDVMeasLayer : public TVMeasLayer { -public: - - static Bool_t kActive; - static Bool_t kDummy; - - /** Get the layer ID */ - inline int getLayerID() const { return _layerID ; } - - /** Get the Cell ID associated with this measurement layer */ - inline const std::vector<int>& getCellIDs() const { return _cellIDs ; } - - /** Get the number of Cell ID associated with this measurement layer */ - inline unsigned int getNCellIDs() const { return _cellIDs.size() ; } - - /** Get the Magnetic field at the measurement surface */ - inline Double_t GetBz() const { return _Bz; } - - /** Convert LCIO Tracker Hit to an ILDPLanarTrackHit */ - virtual ILDVTrackHit* ConvertLCIOTrkHit(edm4hep::TrackerHit trkhit) const = 0 ; - - /** Check whether the measurement layer represents a series of detector elements */ - bool isMultilayer() const { return _isMultiLayer; } - - /** Get the intersection and the CellID, needed for multilayers */ - virtual int getIntersectionAndCellID(const TVTrack &hel, - TVector3 &xx, - Double_t &phi, - Int_t &CellID, - Int_t mode, - Double_t eps = 1.e-8) const = 0 ; - -protected: - - ILDVMeasLayer(TMaterial &min, - TMaterial &mout, - Double_t Bz, - Bool_t is_active = ILDVMeasLayer::kActive, - int CellID = -1 , - const Char_t *name = "ILDMeasL"); - - ILDVMeasLayer(TMaterial &min, - TMaterial &mout, - Double_t Bz, - const std::vector<int>& cellIDs, - Bool_t is_active = ILDVMeasLayer::kActive, - const Char_t *name = "ILDMeasL"); - - - - Double_t _Bz ; // Magnitude of B-Field in Z - int _layerID ; - std::vector<int> _cellIDs ; - - bool _isMultiLayer; - -private: - - -}; - -#endif diff --git a/Utilities/KalDet/src/ild/common/ILDVTrackHit.h b/Utilities/KalDet/src/ild/common/ILDVTrackHit.h deleted file mode 100644 index a2c22400..00000000 --- a/Utilities/KalDet/src/ild/common/ILDVTrackHit.h +++ /dev/null @@ -1,33 +0,0 @@ -#ifndef ILDVTrackHIT_H -#define ILDVTrackHIT_H - -/** ILDVMeasLayer: Virtual hit class used by ILD[X]Hit Classes, which should provide coordinate vector as defined by the MeasLayer - * - * @author S.Aplin DESY - */ - - -#include "kaltest/TVTrackHit.h" - -#include "ILDVMeasLayer.h" - -#include "edm4hep/TrackerHit.h" - -class ILDVTrackHit : public TVTrackHit { - -public: - - /** Constructor Taking coordinates and associated measurement layer, with bfield and number of measurement dimentions*/ - ILDVTrackHit(const TVMeasLayer &ms, Double_t *x, Double_t *dx, - Double_t bfield , Int_t dim, edm4hep::TrackerHit trkhit) - : TVTrackHit(ms, x, dx, bfield, dim), _trkhit(trkhit) - { /* no op */ } - - edm4hep::TrackerHit getLCIOTrackerHit() const { return _trkhit; } - -private: - - edm4hep::TrackerHit _trkhit; - -}; -#endif diff --git a/Utilities/KalDet/src/ild/common/MaterialDataBase.cc b/Utilities/KalDet/src/ild/common/MaterialDataBase.cc index abf8e4e4..50a989e0 100644 --- a/Utilities/KalDet/src/ild/common/MaterialDataBase.cc +++ b/Utilities/KalDet/src/ild/common/MaterialDataBase.cc @@ -1,5 +1,5 @@ -#include "MaterialDataBase.h" +#include "kaldet/MaterialDataBase.h" #include <stdexcept> #include <vector> @@ -238,6 +238,19 @@ void MaterialDataBase::createMaterials(const gear::GearMgr& gearMgr, IGeomSvc* g catch( gear::UnknownParameterException& e){ std::cout << "Warning! cannot get material VXDShellMaterial from GeomSvc! GearMgr=" << &gearMgr << std::endl; } + try { + const gear::SimpleMaterial& etd_sup_mat = gearMgr.getSimpleMaterial("OTKEndcapSupportMaterial"); + A = etd_sup_mat.getA(); + Z = etd_sup_mat.getZ(); + density = etd_sup_mat.getDensity() * (1000.0/ 1000000.0); // kg/m^3 -> g/cm^3 + radlen = etd_sup_mat.getRadLength() / 10.0 ; // mm -> cm + name = etd_sup_mat.getName() ; + TMaterial &etdsupport = *new TMaterial(name.c_str(), "", A, Z, density, radlen, 0.); + this->addMaterial(&etdsupport, name); + } + catch( gear::UnknownParameterException& e){ + std::cout << "Warning! cannot get material OTKEndcapSupportMaterial from GeomSvc! GearMgr=" << &gearMgr << std::endl; + } try { const gear::SimpleMaterial& tpc_readout_mat = gearMgr.getSimpleMaterial("TPCReadoutMaterial"); A = tpc_readout_mat.getA(); diff --git a/Utilities/KalDet/src/ild/common/MaterialDataBase.h b/Utilities/KalDet/src/ild/common/MaterialDataBase.h deleted file mode 100644 index 398c0252..00000000 --- a/Utilities/KalDet/src/ild/common/MaterialDataBase.h +++ /dev/null @@ -1,76 +0,0 @@ -#ifndef MaterialDataBase_h -#define MaterialDataBase_h - -/** MaterialDataBase: Class to hold and manage collection of materials - * - * @author S.Aplin DESY - */ - -#include <string> -#include <map> -#include <exception> - -#include "lcio.h" -#include "Exceptions.h" - -class TMaterial; - -namespace gear{ - class GearMgr ; -} -class IGeomSvc; -// fg: define the MaterialDataBaseException as an lcio Exception to allow for -// messages to be printed in what() -typedef lcio::Exception MaterialDataBaseException ; - -class MaterialDataBase { - -public: - - /** Accessor Method */ - static MaterialDataBase& Instance() { - - static MaterialDataBase singleton; - - return singleton; - - } - - // Other non-static member functions - -public: - - /** Destructor */ - ~MaterialDataBase(); - - /** Get Material via name */ - TMaterial* getMaterial(std::string mat_name) ; - - void registerForService(const gear::GearMgr& gearMgr, IGeomSvc* geoSvc=0) ; - - -private: - - void initialise(const gear::GearMgr& gearMgr, IGeomSvc* geoSvc) ; - - MaterialDataBase() { _material_map.clear(); _isInitialised = false ; _gearMgr = 0; } // Private constructor - - - MaterialDataBase(const MaterialDataBase&) ; // Prevent copy-construction - MaterialDataBase& operator=(const MaterialDataBase&) ; // Prevent assignment - - void addMaterial(TMaterial* mat, std::string name); - void createMaterials(const gear::GearMgr& gearMgr, IGeomSvc* geoSvc); - - // private member variables - std::map<std::string,TMaterial* > _material_map; - - bool _isInitialised; - - const gear::GearMgr* _gearMgr; - -}; - - - -#endif diff --git a/Utilities/KalDet/src/ild/etd/CEPCOTKEndcapKalDetector.cc b/Utilities/KalDet/src/ild/etd/CEPCOTKEndcapKalDetector.cc new file mode 100644 index 00000000..07515374 --- /dev/null +++ b/Utilities/KalDet/src/ild/etd/CEPCOTKEndcapKalDetector.cc @@ -0,0 +1,547 @@ + +#include "kaldet/CEPCOTKEndcapKalDetector.h" + +#include "kaldet/MaterialDataBase.h" + +#include <sstream> + +#include "DetInterface/IGeomSvc.h" +#include "DD4hep/Detector.h" +#include "DDRec/DetectorData.h" +#include "DD4hep/DD4hepUnits.h" + +#include "gear/GEAR.h" +#include "gear/BField.h" +#include "gearimpl/Util.h" + +#include "kaldet/ILDRotatedTrapMeaslayer.h" +//#include "kaldet/ILDSegmentedDiscMeasLayer.h" +#include "kaldet/ILDSegmentedDiscStripMeasLayer.h" +#include "kaldet/CEPCSegmentedDiscMeasLayer.h" + +#include <UTIL/BitField64.h> +#include <UTIL/ILDConf.h> +//#include "Identifier/CEPCConf.h" + +#include "CLHEP/Units/SystemOfUnits.h" +#include "TVector3.h" + +CEPCOTKEndcapKalDetector::CEPCOTKEndcapKalDetector(const gear::GearMgr& gearMgr, IGeomSvc* geoSvc) : +TVKalDetector(300), _nDisks(0) // SJA:FIXME initial size, 300 looks reasonable for ILD, though this would be better stored as a const somewhere +{ + //std::cout << "CEPCOTKEndcapKalDetector building OTKEndcap detector using GEAR " << std::endl ; + + MaterialDataBase::Instance().registerForService(gearMgr, geoSvc); + if(geoSvc){ + setupGearGeom(geoSvc); + } + else{ + setupGearGeom(gearMgr); + } + + this->build_staggered_design(); + + SetOwner(); + +} + +void CEPCOTKEndcapKalDetector::build_staggered_design() { + //std::cout << "CEPCOTKEndcapKalDetector::build_staggered_design " << std::endl; + std::string name = "OTKEndcap"; + + UTIL::BitField64 encoder(lcio::ILDCellID0::encoder_string); + + for (int idisk = 0; idisk < _nDisks; ++idisk) { + //std::cout << "CEPCOTKEndcapKalDetector::build_staggered_design build disk " << idisk << std::endl; + int npetals = _OTKEndcapgeo[idisk].nPetals; + double phi0 = _OTKEndcapgeo[idisk].phi0; + + double senZPos_even_front = _OTKEndcapgeo[idisk].senZPos_even_front; + double senZPos_odd_front = _OTKEndcapgeo[idisk].senZPos_odd_front; + + // check that the number of petals is divisible by 2 + int nsegments = npetals/2; + + // even segments forward + this->create_segmented_disk_layers(idisk, npetals, false, phi0, senZPos_odd_front); + + // even segments backwards + this->create_segmented_disk_layers(idisk, npetals, false, phi0, -senZPos_odd_front); + + // odd segements + // update phi0 by the angular distance of one petal + phi0 += 2.0 * M_PI / npetals; + + // odd segments forward + //this->create_segmented_disk_layers(idisk, nsegments, true, phi0, senZPos_even_front); + + // odd segments backward + //this->create_segmented_disk_layers(idisk, nsegments, true, phi0, -senZPos_even_front); + /* + // make the air disks + TMaterial & air = *MaterialDataBase::Instance().getMaterial("air") ; + + Bool_t dummy = false; + + // place air discs to help transport during track extrapolation + if (idisk < _nDisks-1) { + // place the disc half way between the two discs + double z = (_OTKEndcapgeo[idisk].senZPos_even_front + _OTKEndcapgeo[idisk+1].senZPos_even_front) * 0.5; + + TVector3 xc_fwd(0.0, 0.0, z); + TVector3 normal_fwd(xc_fwd); + normal_fwd.SetMag(1.0); + + double eps1 = 1.0e-04; // offset for disk number + double eps2 = 1.0e-05; // odd or even + double eps4 = 1.0e-08; // offset for forwards and backwards + + double height = _OTKEndcapgeo[idisk].height; + double rInner = _OTKEndcapgeo[idisk].rInner; + + double sort_policy = rInner+height + eps1 * idisk + eps2 * 2; // we need to be after the even and odd + + //std::cout << "CEPCOTKEndcapKalDetector::create air support disk at " << xc_fwd.z() << " sort_policy = " << sort_policy << std::endl; + + Add(new ILDDiscMeasLayer(air, air, xc_fwd, normal_fwd, _bZ, sort_policy, + rInner, rInner+height, dummy, -1, "OTKEndcapAirSupportDiscFront")); + + TVector3 xc_bwd(0.0, 0.0, -z); + TVector3 normal_bwd(xc_bwd); + normal_bwd.SetMag(1.0); + + // offset needed for rear disks + sort_policy += eps4; + + //std::cout << "CEPCOTKEndcapKalDetector::create air support disk at " << xc_bwd.z() << " sort_policy = " << sort_policy << std::endl; + Add(new ILDDiscMeasLayer(air, air, xc_bwd, normal_bwd, _bZ, sort_policy, + rInner, rInner+height, dummy, -1, "OTKEndcapAirSupportDiscRear")); + } + */ + } +} + +void CEPCOTKEndcapKalDetector::create_segmented_disk_layers(int idisk, int nsegments, bool even_petals, double phi0, double zpos) { + + Bool_t active = true; + Bool_t dummy = false; + + TMaterial & air = *MaterialDataBase::Instance().getMaterial("air"); + TMaterial & silicon = *MaterialDataBase::Instance().getMaterial("silicon"); + TMaterial & carbon = *MaterialDataBase::Instance().getMaterial("carbon"); + TMaterial & stripsupport = *MaterialDataBase::Instance().getMaterial("OTKEndcapSupportMaterial"); + + double halfPetal = _OTKEndcapgeo[idisk].dphi/2.0; + double senThickness = _OTKEndcapgeo[idisk].senThickness; + double supThickness = _OTKEndcapgeo[idisk].supThickness; + double innerBaseLength = _OTKEndcapgeo[idisk].innerBaseLength; + double outerBaseLength = _OTKEndcapgeo[idisk].outerBaseLength; + double height = _OTKEndcapgeo[idisk].height; + double rInner = _OTKEndcapgeo[idisk].rInner; + bool isDoubleSided = _OTKEndcapgeo[idisk].isDoubleSided; + int nSensors = _OTKEndcapgeo[idisk].nSensors; + int zsign = zpos > 0 ? +1 : -1; + double stripAngle = _OTKEndcapgeo[idisk].stripAngle; + bool isStripReadout = _OTKEndcapgeo[idisk].isStripReadout; + + //SJA:FIXME: due to the space frame design of the strip layers there is far too much support so just leave it out for now ... + TMaterial & support = isStripReadout == false ? carbon : stripsupport; + + UTIL::BitField64 encoder(lcio::ILDCellID0::encoder_string); + encoder.reset(); // reset to 0 + + encoder[lcio::ILDCellID0::subdet] = lcio::ILDDetID::ETD; + encoder[lcio::ILDCellID0::side] = zsign; + encoder[lcio::ILDCellID0::layer] = idisk; + + int start_index = even_petals ? 0 : 1 ; + std::vector<int> sensors_front; + std::vector<int> sensors_back; + std::vector<int> module_ids_front; + std::vector<int> module_ids_back; + std::vector<int> nsensors(10,0); // FIXME: assume 10 row + + if (isDoubleSided) { // sensors on front and back double supZPos_odd_front = _OTKEndcapgeo[idisk].supZPos_odd; + // first half is on the front, second half on the back, sensors start with sensor number 1 + for (int iSensor=0; iSensor < nSensors/2; iSensor++) { + sensors_front.push_back(iSensor); + int irow = std::fmod(iSensor, 10); + nsensors[10-1-irow]++; + } + for (int iSensor=nSensors/2; iSensor < nSensors; iSensor++) sensors_back.push_back(iSensor); + } + else { // only sensors on the front + for (int iSensor=0; iSensor < nSensors; iSensor++) { + sensors_front.push_back(iSensor); + int irow = std::fmod(iSensor, 10); + nsensors[10-1-irow]++; + } + } + + for (int i=0; i<nsegments; ++i) { + encoder[lcio::ILDCellID0::module] = i;//start_index + i*2 ; + + for (unsigned j=0; j<sensors_front.size(); j++) { + encoder[lcio::ILDCellID0::sensor] = sensors_front[j]; + module_ids_front.push_back(encoder.lowWord()); + //std::cout << "will add front cell id: " << encoder.lowWord() << std::endl; + } + + for (unsigned j=0; j<sensors_back.size(); j++) { + encoder[lcio::ILDCellID0::sensor] = sensors_back[j]; + module_ids_back.push_back(encoder.lowWord()); + //std::cout << "will add back cell id: " << encoder.lowWord() << std::endl; + } + } + //std::cout << "total " << nsegments*nSensors << " cell IDs added, module: 0-" << nsegments-1 << " sensor: 0-" << nSensors-1 << std::endl; + + // create segmented disk + // strip: SegmentedDisc, trapezoid + // pixel: Disc, tube + + // front face of sensitive + double z = zpos - zsign*0.5*(senThickness); + // double sort_policy = fabs(z) ; + double eps1 = 1.0e-04 ; // disk + double eps2 = 1.0e-05 ; // odd or even + double eps3 = 1.0e-06 ; // layer in disk + double eps4 = 1.0e-08 ; // forward or backwards + + double sort_policy = rInner+height + eps1 * idisk + eps3 * 1; + if (!even_petals) sort_policy += eps2; + + // if this is the negative z disk add epsilon to the policy + if (z < 0) sort_policy += eps4 ; + const char *name1 = z > 0 ? "OTKEndcapSenFrontPositiveZ" : "OTKEndcapSenFrontNegativeZ"; + + //std::cout << "CEPCOTKEndcapKalDetector::create_segmented_disk_layers add front face of sensitive at " << z << " sort_policy = " << sort_policy << std::endl; + if (isStripReadout) { + Add(new ILDSegmentedDiscMeasLayer(air, silicon, _bZ, sort_policy, nsegments, z, phi0, rInner, height, innerBaseLength, outerBaseLength, dummy, name1)); + } + else { + TVector3 xc(0.0, 0.0, z); + TVector3 normal(xc); + normal.SetMag(1.0); + Add(new CEPCSegmentedDiscMeasLayer(air, silicon, _bZ, sort_policy, nsegments, z, phi0, rInner, rInner+height, halfPetal, dummy, name1)); + } + + // measurement plane + z += zsign*0.5*senThickness; + //sort_policy = fabs(z) ; + sort_policy = rInner+height + eps1 * idisk + eps3 * 2 ; + if (z < 0) sort_policy += eps4 ; + + if (!even_petals) sort_policy += eps2; + if (!even_petals) { + const char *name2 = z > 0 ? "OTKEndcapMeasLayerFrontPositiveZOdd" : "OTKEndcapMeasLayerFrontNegativeZOdd"; + //std::cout << "CEPCOTKEndcapKalDetector::create_segmented_disk_layers add measurement plane at " << z << " sort_policy = " << sort_policy << " Strip Readout = " << isStripReadout << " number of module_ids = " << module_ids_front.size(); + if (isStripReadout) { + Add(new ILDSegmentedDiscStripMeasLayer(silicon, silicon, _bZ, sort_policy, nsegments, z, phi0, rInner, height, innerBaseLength, outerBaseLength, stripAngle, active, module_ids_front, name2)); + //std::cout << " stripAngle = " << stripAngle ; + } + else { + TVector3 xc(0.0, 0.0, z); + TVector3 normal(xc); + normal.SetMag(1.0); + Add(new CEPCSegmentedDiscMeasLayer(silicon, silicon, _bZ, sort_policy, nsegments, z, phi0, rInner, rInner+height, halfPetal, nsensors, active, module_ids_front, name2)); + } + //std::cout << std::endl;*/ + } + else { + const char *name2 = z > 0 ? "OTKEndcapMeasLayerFrontPositiveZEven" : "OTKEndcapMeasLayerFrontNegativeZEven"; + std::cout << "CEPCOTKEndcapKalDetector::create_segmented_disk_layers add measurement plane at " << z << " sort_policy = " << sort_policy << " Strip Readout = " << isStripReadout << " number of module_ids = " << module_ids_front.size(); + if (isStripReadout) { + Add(new ILDSegmentedDiscStripMeasLayer(silicon, silicon, _bZ, sort_policy, nsegments, z, phi0, rInner, height, innerBaseLength, outerBaseLength, stripAngle, active, module_ids_front, name2)); + //std::cout << " stripAngle = " << stripAngle; + } + else { + TVector3 xc(0.0, 0.0, z); + TVector3 normal(xc); + normal.SetMag(1.0); + Add(new CEPCSegmentedDiscMeasLayer(silicon, silicon, _bZ, sort_policy, nsegments, z, phi0, rInner, rInner+height, halfPetal, nsensors, active, module_ids_front, name2)); + } + //std::cout << std::endl; + } + + // interface between sensitive and support + z += zsign*0.5*senThickness; + // sort_policy = fabs(z) ; + sort_policy = rInner+height + eps1 * idisk + eps3 * 3; + if (z < 0) sort_policy += eps4; + if (!even_petals) sort_policy += eps2; + + const char *name3 = z > 0 ? "OTKEndcapSenSupportIntfPositiveZ" : "OTKEndcapSenSupportIntfNegativeZ"; + + //std::cout << "CEPCOTKEndcapKalDetector::create_segmented_disk_layers add interface between sensitive and support at " << z << " sort_policy = " << sort_policy << std::endl; + if (isStripReadout) { + Add(new ILDSegmentedDiscMeasLayer(silicon, support, _bZ, sort_policy, nsegments, z, phi0, rInner, height, innerBaseLength, outerBaseLength, dummy, name3)); + } + else { + TVector3 xc(0.0, 0.0, z); + TVector3 normal(xc); + normal.SetMag(1.0); + Add(new CEPCSegmentedDiscMeasLayer(silicon, support, _bZ, sort_policy, nsegments, z, phi0, rInner, rInner+height, halfPetal, dummy, name3)); + //Add(new CEPCSegmentedDiscMeasLayer(silicon, support, xc, normal, _bZ, sort_policy, rInner, rInner+height, dummy, -1, name3)); + //std::cout << "add ILDDiscMeasLayer at z = " << xc[2] << std::endl; + } + + if (isDoubleSided) { + // interface between support and sensitive + z += zsign*supThickness; + // sort_policy = fabs(z) ; + sort_policy = rInner+height + eps1 * idisk + eps3 * 4; + if (z < 0) sort_policy += eps4; + if (!even_petals) sort_policy += eps2; + + const char *name4 = z > 0 ? "OTKEndcapSupportSenIntfPositiveZ" : "OTKEndcapSupportSenIntfNegativeZ"; + + //std::cout << "CEPCOTKEndcapKalDetector::create_segmented_disk_layers add interface between support and sensitive at " << z << " sort_policy = " << sort_policy << std::endl; + if (isStripReadout) { + Add(new ILDSegmentedDiscMeasLayer(support, silicon , _bZ, sort_policy, nsegments, z, phi0, rInner, height, innerBaseLength, outerBaseLength, dummy, name4)); + } + else { + TVector3 xc(0.0, 0.0, z); + TVector3 normal(xc); + normal.SetMag(1.0); + Add(new CEPCSegmentedDiscMeasLayer(support, silicon, _bZ, sort_policy, nsegments, z, phi0, rInner, rInner+height, halfPetal, dummy, name4)); + //Add(new CEPCSegmentedDiscMeasLayer(support, silicon, xc, normal, _bZ, sort_policy, rInner, rInner+height, dummy, -1, name4)); + } + // measurement plane at the back + z += zsign*0.5*senThickness; + // sort_policy = fabs(z) ; + sort_policy = rInner+height + eps1 * idisk + eps3 * 5; + if (z < 0) sort_policy += eps4; + if (!even_petals) { + sort_policy += eps2; + + const char *name5 = z > 0 ? "OTKEndcapMeasLayerBackPositiveZOdd" : "OTKEndcapMeasLayerBackNegativeZOdd"; + //std::cout << "CEPCOTKEndcapKalDetector::create_segmented_disk_layers add measurement plane at " << z << " sort_policy = " << sort_policy << " Strip Readout = " << isStripReadout << " number of module_ids = " << module_ids_back.size() ; + if (isStripReadout) { + Add(new ILDSegmentedDiscStripMeasLayer(silicon, silicon, _bZ, sort_policy, nsegments, z, phi0, rInner, height, innerBaseLength, outerBaseLength, -stripAngle, active, module_ids_back, name5)); + //std::cout << " stripAngle = " << -stripAngle; + } + else { + TVector3 xc(0.0, 0.0, z); + TVector3 normal(xc); + normal.SetMag(1.0); + Add(new CEPCSegmentedDiscMeasLayer(silicon, silicon, _bZ, sort_policy, nsegments, z, phi0, rInner, rInner+height, halfPetal, nsensors, active, module_ids_back, name5)); + //Add(new CEPCSegmentedDiscMeasLayer(silicon, silicon, xc, normal, _bZ, sort_policy, rInner, rInner+height, module_ids_back, active, name5)); + } + //std::cout << std::endl; + } + else{ + const char *name5 = z > 0 ? "OTKEndcapMeasLayerBackPositiveZEven" : "OTKEndcapMeasLayerBackNegativeZEven"; + //std::cout << "CEPCOTKEndcapKalDetector::create_segmented_disk_layers add measurement plane at " << z << " sort_policy = " << sort_policy<< " Strip Readout = " << isStripReadout << " number of module_ids = " << module_ids_back.size(); + if (isStripReadout) { + Add(new ILDSegmentedDiscStripMeasLayer(silicon, silicon, _bZ, sort_policy, nsegments, z, phi0, rInner, height, innerBaseLength, outerBaseLength, -stripAngle, active, module_ids_back, name5)); + //std::cout << " stripAngle = " << -stripAngle ; + } + else { + TVector3 xc(0.0, 0.0, z); + TVector3 normal(xc); + normal.SetMag(1.0); + Add(new CEPCSegmentedDiscMeasLayer(silicon, silicon, _bZ, sort_policy, nsegments, z, phi0, rInner, rInner+height, halfPetal, nsensors, active, module_ids_back, name5)); + //Add(new CEPCSegmentedDiscMeasLayer(silicon, silicon, xc, normal, _bZ, sort_policy, rInner, rInner+height, module_ids_back, active, name5)); + } + //std::cout << std::endl; + } + + // rear face of sensitive + z += zsign*0.5*senThickness; + // sort_policy = fabs(z) ; + sort_policy = rInner+height + eps1 * idisk + eps3 * 6; + if (z < 0) sort_policy += eps4; + if (!even_petals) sort_policy += eps2; + + const char *name6 = z > 0 ? "OTKEndcapSenRearPositiveZ" : "OTKEndcapSenRearNegativeZ"; + + //std::cout << "CEPCOTKEndcapKalDetector::create_segmented_disk_layers add rear face of sensitive at " << z << " sort_policy = " << sort_policy << std::endl; + if (isStripReadout) { + Add(new ILDSegmentedDiscMeasLayer(silicon, air, _bZ, sort_policy, nsegments, z, phi0, rInner, height, innerBaseLength, outerBaseLength, dummy, name6)); + } + else { + TVector3 xc(0.0, 0.0, z); + TVector3 normal(xc); + normal.SetMag(1.0); + Add(new CEPCSegmentedDiscMeasLayer(silicon, air, _bZ, sort_policy, nsegments, z, phi0, rInner, rInner+height, halfPetal, dummy, name6)); + //Add(new CEPCSegmentedDiscMeasLayer(silicon, air, xc, normal, _bZ, sort_policy, rInner, rInner+height, dummy, -1, name6)); + } + } + else{ + // rear face of support + z += zsign*supThickness; + // sort_policy = fabs(z) ; + sort_policy = rInner+height + eps1 * idisk + eps3 * 4; + if (z < 0) sort_policy += eps4; + if (!even_petals) sort_policy += eps2; + + const char *name4 = z > 0 ? "OTKEndcapSupRearPositiveZ" : "OTKEndcapSupRearNegativeZ"; + + //std::cout << "CEPCOTKEndcapKalDetector::create_segmented_disk_layers add rear face of support at " << z << " sort_policy = " << sort_policy << std::endl; + if (isStripReadout) { + Add(new ILDSegmentedDiscMeasLayer(support, air, _bZ, sort_policy, nsegments, z, phi0, rInner, height, innerBaseLength, outerBaseLength, dummy, name4)); + } + else { + TVector3 xc(0.0, 0.0, z); + TVector3 normal(xc); + normal.SetMag(1.0); + Add(new CEPCSegmentedDiscMeasLayer(support, air, _bZ, sort_policy, nsegments, z, phi0, rInner, rInner+height, halfPetal, dummy, name4)); + //Add(new CEPCSegmentedDiscMeasLayer(support, air, xc, normal, _bZ, sort_policy, rInner, rInner+height, dummy, -1, name4)); + } + } +} + +void CEPCOTKEndcapKalDetector::setupGearGeom(const gear::GearMgr& gearMgr) { + const gear::GearParameters& etdParams = gearMgr.getGearParameters("ETDParameters"); + //std::cout << etdParams << std::endl; + + _bZ = gearMgr.getBField().at(gear::Vector3D(0.,0.,0.)).z(); + + double strip_angle_deg = 0.0; + bool strip_angle_present = true; + try { + strip_angle_deg = etdParams.getDoubleVal("strip_angle_deg"); + } + catch (gear::UnknownParameterException& e) { + strip_angle_present = false; + } + //std::cout << "=============OTKEndcap strip angle: " << strip_angle_deg << "==============" << std::endl; + try { + std::vector<int> nPetals = etdParams.getIntVals("ETDPetalNumber"); + std::vector<int> nSensors = etdParams.getIntVals("ETDSensorNumber"); + std::vector<double> petalAngles = etdParams.getDoubleVals("ETDPetalAngle"); + std::vector<double> phi0s = etdParams.getDoubleVals("ETDDiskPhi0"); + std::vector<double> alphas = etdParams.getDoubleVals("ETDDiskAlpha"); + std::vector<double> zpositions = etdParams.getDoubleVals("ETDDiskPosition"); + std::vector<double> zoffsets = etdParams.getDoubleVals("ETDDiskOffset"); + std::vector<double> supRinners = etdParams.getDoubleVals("ETDSupportRmin"); + std::vector<double> supThicknesss = etdParams.getDoubleVals("ETDSupportThickness"); + std::vector<double> supHeights = etdParams.getDoubleVals("ETDSupportHeight"); + std::vector<double> senRinners = etdParams.getDoubleVals("ETDSensitiveRmin"); + std::vector<double> senThicknesss = etdParams.getDoubleVals("ETDSensitiveThickness"); + std::vector<double> senHeights = etdParams.getDoubleVals("ETDSensitiveHeight"); + + _nDisks = nPetals.size(); + _OTKEndcapgeo.resize(_nDisks); + + double eps = 1.0e-08; + + for(int disk=0; disk< _nDisks; ++disk){ + _OTKEndcapgeo[disk].nPetals = nPetals [disk]; + _OTKEndcapgeo[disk].nSensors = nSensors [disk]; + _OTKEndcapgeo[disk].dphi = petalAngles [disk]; + _OTKEndcapgeo[disk].phi0 = phi0s [disk]; + _OTKEndcapgeo[disk].alpha = alphas [disk]; + _OTKEndcapgeo[disk].rInner = senRinners [disk]; + _OTKEndcapgeo[disk].height = senHeights [disk]; + _OTKEndcapgeo[disk].innerBaseLength = 0; + _OTKEndcapgeo[disk].outerBaseLength = 0; + _OTKEndcapgeo[disk].senThickness = senThicknesss[disk]; + _OTKEndcapgeo[disk].supThickness = supThicknesss[disk]; + + _OTKEndcapgeo[disk].senZPos_even_front = zpositions[disk] - zoffsets[disk]; + _OTKEndcapgeo[disk].senZPos_odd_front = zpositions[disk] + zoffsets[disk]; + + // fixed single pixel layer (single long strip dealed as pixel), if design changed, to update + _OTKEndcapgeo[disk].isDoubleSided = false; + _OTKEndcapgeo[disk].isStripReadout = false; + + if (strip_angle_present) { + _OTKEndcapgeo[disk].stripAngle = strip_angle_deg * M_PI/180 ; + } + else { + _OTKEndcapgeo[disk].stripAngle = 0.0 ; + } + + //std::cout << _OTKEndcapgeo[disk].nPetals << " " << _OTKEndcapgeo[disk].dphi << " " << _OTKEndcapgeo[disk].phi0 << " " << _OTKEndcapgeo[disk].alpha << " " + // << _OTKEndcapgeo[disk].rInner << " " << _OTKEndcapgeo[disk].height << " " << _OTKEndcapgeo[disk].innerBaseLength << " " << _OTKEndcapgeo[disk].outerBaseLength << " " + // << _OTKEndcapgeo[disk].senThickness << " " << _OTKEndcapgeo[disk].supThickness << " " << _OTKEndcapgeo[disk].senZPos_even_front << " " << _OTKEndcapgeo[disk].senZPos_odd_front << " " + // << _OTKEndcapgeo[disk].isDoubleSided << " " << _OTKEndcapgeo[disk].isStripReadout << " " << _OTKEndcapgeo[disk].nSensors << " " << _OTKEndcapgeo[disk].stripAngle << std::endl; + //////////////////////////////////////////////////////////////////////////////////////////////////////// + // Assertions //////////////////////////////////////////////////////////////////////////////////// + + // there should be an even number of petals per disk if offset not zero + if (zoffsets[disk]!=0) assert(_OTKEndcapgeo[disk].nPetals%2 == 0); + + // support and sensitive should have the same size an position in xy + assert(fabs(supRinners[disk] - senRinners[disk]) < eps); + assert(fabs(supHeights[disk] - senHeights[disk]) < eps); + + // double sided petals should have as many sensors on the front as on the back + if(_OTKEndcapgeo[disk].isDoubleSided) assert(_OTKEndcapgeo[disk].nSensors%2 == 0); + + // petals should not be rotated around their symmetry axis + assert(fabs(alphas[disk]) < eps); // not still support alpha + } + } + catch (gear::UnknownParameterException& e) { + std::cout << e.what() << std::endl; + _nDisks = 0; + } +} + +void CEPCOTKEndcapKalDetector::setupGearGeom(IGeomSvc* geoSvc) { + dd4hep::DetElement world = geoSvc->getDD4HepGeo(); + dd4hep::DetElement ftd; + const std::map<std::string, dd4hep::DetElement>& subs = world.children(); + for(std::map<std::string, dd4hep::DetElement>::const_iterator it=subs.begin();it!=subs.end();it++){ + if(it->first!="OTKEndcap") continue; + ftd = it->second; + } + dd4hep::rec::ZDiskPetalsData* ftdData = nullptr; + try{ + ftdData = ftd.extension<dd4hep::rec::ZDiskPetalsData>(); + } + catch(std::runtime_error& e){ + std::cout << e.what() << " " << ftdData << std::endl; + throw GaudiException(e.what(), "FATAL", StatusCode::FAILURE); + } + + const dd4hep::Direction& field = geoSvc->lcdd()->field().magneticField(dd4hep::Position(0,0,0)); + _bZ = field.z()/dd4hep::tesla; + + double strip_angle_deg = ftdData->angleStrip/CLHEP::degree; + bool strip_angle_present = true; + + std::vector<dd4hep::rec::ZDiskPetalsData::LayerLayout>& ftdlayers = ftdData->layers; + _nDisks = ftdlayers.size() ; + + _OTKEndcapgeo.resize(_nDisks); + + double eps = 1.0e-08; + //std::cout << "=============OTKEndcap strip angle: " << strip_angle_deg << "==============" << std::endl; + for(int disk=0; disk< _nDisks; ++disk){ + dd4hep::rec::ZDiskPetalsData::LayerLayout& ftdlayer = ftdlayers[disk]; + _OTKEndcapgeo[disk].nPetals = ftdlayer.petalNumber; + _OTKEndcapgeo[disk].dphi = ftdlayer.petalHalfAngle*2; + _OTKEndcapgeo[disk].phi0 = ftdlayer.phi0; + _OTKEndcapgeo[disk].alpha = ftdlayer.alphaPetal; + _OTKEndcapgeo[disk].rInner = ftdlayer.distanceSensitive*CLHEP::cm; + _OTKEndcapgeo[disk].height = ftdlayer.lengthSensitive*CLHEP::cm; + _OTKEndcapgeo[disk].innerBaseLength = ftdlayer.widthInnerSensitive*CLHEP::cm; + _OTKEndcapgeo[disk].outerBaseLength = ftdlayer.widthOuterSensitive*CLHEP::cm; + _OTKEndcapgeo[disk].senThickness = ftdlayer.thicknessSensitive*CLHEP::cm; + _OTKEndcapgeo[disk].supThickness = ftdlayer.thicknessSupport*CLHEP::cm; + + _OTKEndcapgeo[disk].senZPos_even_front = ftdlayer.zPosition*CLHEP::cm - ftdlayer.zOffsetSensitive*CLHEP::cm; + _OTKEndcapgeo[disk].senZPos_odd_front = ftdlayer.zPosition*CLHEP::cm - ftdlayer.zOffsetSensitive*CLHEP::cm - 2*ftdlayer.zOffsetSupport*CLHEP::cm; + + _OTKEndcapgeo[disk].isDoubleSided = ftdlayer.typeFlags[dd4hep::rec::ZDiskPetalsData::SensorType::DoubleSided]; + _OTKEndcapgeo[disk].isStripReadout = !((bool)ftdlayer.typeFlags[dd4hep::rec::ZDiskPetalsData::SensorType::Pixel]); + _OTKEndcapgeo[disk].nSensors = ftdlayer.sensorsPerPetal; + + if (strip_angle_present) { + _OTKEndcapgeo[disk].stripAngle = strip_angle_deg * M_PI/180 ; + } else { + _OTKEndcapgeo[disk].stripAngle = 0.0 ; + } + + assert(_OTKEndcapgeo[disk].nPetals%2 == 0); + assert(fabs(ftdlayer.widthInnerSupport - ftdlayer.widthInnerSensitive) < eps); + assert(fabs(ftdlayer.widthOuterSupport - ftdlayer.widthOuterSensitive) < eps); + assert(fabs(ftdlayer.lengthSupport - ftdlayer.lengthSensitive) < eps); + assert(fabs(ftdlayer.distanceSupport - ftdlayer.distanceSensitive) < eps); + if (_OTKEndcapgeo[disk].isDoubleSided) assert(_OTKEndcapgeo[disk].nSensors%2 == 0); + assert(fabs(ftdlayer.alphaPetal) < eps); // not still support alpha + } +} diff --git a/Utilities/KalDet/src/ild/ftd/ILDFTDDiscBasedKalDetector.cc b/Utilities/KalDet/src/ild/ftd/ILDFTDDiscBasedKalDetector.cc index db2b578d..15fea781 100644 --- a/Utilities/KalDet/src/ild/ftd/ILDFTDDiscBasedKalDetector.cc +++ b/Utilities/KalDet/src/ild/ftd/ILDFTDDiscBasedKalDetector.cc @@ -1,6 +1,5 @@ -#include "ILDFTDDiscBasedKalDetector.h" - +#include "kaldet/ILDFTDDiscBasedKalDetector.h" #include "kaldet/MaterialDataBase.h" #include <sstream> diff --git a/Utilities/KalDet/src/ild/ftd/ILDFTDDiscBasedKalDetector.h b/Utilities/KalDet/src/ild/ftd/ILDFTDDiscBasedKalDetector.h deleted file mode 100644 index c8536c64..00000000 --- a/Utilities/KalDet/src/ild/ftd/ILDFTDDiscBasedKalDetector.h +++ /dev/null @@ -1,44 +0,0 @@ -#ifndef __ILDFTDDISCBASEDDETECTOR__ -#define __ILDFTDDISCBASEDDETECTOR__ - -/** Disk based version of the FTD alla LOI -* -* @author S.Aplin DESY -*/ - -#include "kaltest/TVKalDetector.h" - -class TNode; - -namespace gear{ - class GearMgr ; -} - - -class ILDFTDDiscBasedKalDetector : public TVKalDetector { -public: - - /** Initialize the FTD from GEAR */ - ILDFTDDiscBasedKalDetector( const gear::GearMgr& gearMgr ); - - -private: - - void setupGearGeom( const gear::GearMgr& gearMgr ) ; - - int _nDisks ; - double _bZ ; - - struct FTD_Disk { - double rInner; - double rOuter; - double senThickness; - double supThickness; - double zPos; - - }; - std::vector<FTD_Disk> _FTDgeo; - -}; - -#endif diff --git a/Utilities/KalDet/src/ild/ftd/ILDFTDKalDetector.cc b/Utilities/KalDet/src/ild/ftd/ILDFTDKalDetector.cc index b87a9d0f..59279be5 100644 --- a/Utilities/KalDet/src/ild/ftd/ILDFTDKalDetector.cc +++ b/Utilities/KalDet/src/ild/ftd/ILDFTDKalDetector.cc @@ -1,6 +1,5 @@ -#include "ILDFTDKalDetector.h" - +#include "kaldet/ILDFTDKalDetector.h" #include "kaldet/MaterialDataBase.h" #include <sstream> diff --git a/Utilities/KalDet/src/ild/lctpc/LCTPCKalDetector.h b/Utilities/KalDet/src/ild/lctpc/LCTPCKalDetector.h deleted file mode 100644 index 452aaaad..00000000 --- a/Utilities/KalDet/src/ild/lctpc/LCTPCKalDetector.h +++ /dev/null @@ -1,39 +0,0 @@ -#ifndef LCTPCKALDETECTOR_H -#define LCTPCKALDETECTOR_H - -#include "kaltest/TVKalDetector.h" - -#include "kaldet/ILDVMeasLayer.h" - -namespace gear{ - class GearMgr ; -} - -namespace kaldet{ - - /** - * The LCTPC implementation for a TPC which is completely instantiated from GEAR. - * - */ -class LCTPCKalDetector : public TVKalDetector { - -public: - - LCTPCKalDetector() {}; - - /** - * The constructor. All information to initialise the TPC is taken from GEAR. - * - * The class has been copied from GearTPCKalDetector class and adopted for the use of MarlinTrk - * You can find comments and necessary information in the original class - * - */ - LCTPCKalDetector(const gear::GearMgr& gearMgr); - - /// The destructor. - virtual ~LCTPCKalDetector(); - -}; - -}// namespace kaldet -#endif //LCTPCKALDETECTOR_H diff --git a/Utilities/KalDet/src/ild/set/ILDSETKalDetector.h b/Utilities/KalDet/src/ild/set/ILDSETKalDetector.h deleted file mode 100644 index 9327aa8a..00000000 --- a/Utilities/KalDet/src/ild/set/ILDSETKalDetector.h +++ /dev/null @@ -1,58 +0,0 @@ -#ifndef __ILDSETKALDETECTOR__ -#define __ILDSETKALDETECTOR__ - -/** Ladder based SET to be used for ILD DBD studies - * - * @author S.Aplin DESY - */ - -#include "kaltest/TVKalDetector.h" - -#include "TMath.h" - -class TNode; - -namespace gear{ - class GearMgr ; -} - - -class ILDSETKalDetector : public TVKalDetector { - -public: - - /** Initialize the SET from GEAR */ - ILDSETKalDetector( const gear::GearMgr& gearMgr ); - - -private: - - void setupGearGeom( const gear::GearMgr& gearMgr ) ; - - int _nLayers ; - double _bZ ; - - bool _isStripDetector; - - struct SET_Layer { - int nLadders; - int nSensorsPerLadder; - double phi0; - double dphi; - double senRMin; - double supRMin; - double length; - double width; - double offset; - double senThickness; - double supThickness; - double sensorLength; - double stripAngle; - }; - std::vector<SET_Layer> _SETgeo; - -}; - - - -#endif diff --git a/Utilities/KalDet/src/ild/sit/CEPCITKKalDetector.cc b/Utilities/KalDet/src/ild/sit/CEPCITKKalDetector.cc index ae897cc8..e4a6b179 100644 --- a/Utilities/KalDet/src/ild/sit/CEPCITKKalDetector.cc +++ b/Utilities/KalDet/src/ild/sit/CEPCITKKalDetector.cc @@ -1,10 +1,7 @@ #include "kaldet/CEPCITKKalDetector.h" - -#include "MaterialDataBase.h" - -#include "ILDParallelPlanarStripMeasLayer.h" - +#include "kaldet/MaterialDataBase.h" +#include "kaldet/ILDParallelPlanarStripMeasLayer.h" #include <UTIL/BitField64.h> #include <UTIL/ILDConf.h> diff --git a/Utilities/KalDet/src/ild/sit/ILDSITCylinderKalDetector.h b/Utilities/KalDet/src/ild/sit/ILDSITCylinderKalDetector.h deleted file mode 100644 index 29f13f6f..00000000 --- a/Utilities/KalDet/src/ild/sit/ILDSITCylinderKalDetector.h +++ /dev/null @@ -1,44 +0,0 @@ -#ifndef __ILDSITCYLINDERKALDETECTOR__ -#define __ILDSITCYLINDERKALDETECTOR__ - -/** SIT Cylinder based detector to be used for ILD DBD studies when using the old LOI base SIT - * - * @author S.Aplin DESY - */ - -#include "kaltest/TVKalDetector.h" - -class TNode; - -namespace gear{ - class GearMgr ; -} - - -class ILDSITCylinderKalDetector : public TVKalDetector { -public: - - /** Initialize the TPC from GEAR */ - ILDSITCylinderKalDetector( const gear::GearMgr& gearMgr ); - - -private: - - void setupGearGeom( const gear::GearMgr& gearMgr ) ; - - unsigned int _nLayers ; - double _bZ ; - - struct SIT_Layer { - double radius; - double half_length; - double senThickness; - double supThickness; - - }; - std::vector<SIT_Layer> _SITgeo; - - -}; - -#endif diff --git a/Utilities/KalDet/src/ild/sit/ILDSITKalDetector.cc b/Utilities/KalDet/src/ild/sit/ILDSITKalDetector.cc index 5d10c55f..a14b763e 100644 --- a/Utilities/KalDet/src/ild/sit/ILDSITKalDetector.cc +++ b/Utilities/KalDet/src/ild/sit/ILDSITKalDetector.cc @@ -1,10 +1,7 @@ -#include "ILDSITKalDetector.h" - -#include "MaterialDataBase.h" - -#include "ILDParallelPlanarStripMeasLayer.h" - +#include "kaldet/ILDSITKalDetector.h" +#include "kaldet/MaterialDataBase.h" +#include "kaldet/ILDParallelPlanarStripMeasLayer.h" #include <UTIL/BitField64.h> #include <UTIL/ILDConf.h> diff --git a/Utilities/KalDet/src/ild/sit/ILDSITKalDetector.h b/Utilities/KalDet/src/ild/sit/ILDSITKalDetector.h deleted file mode 100644 index 3f5a8636..00000000 --- a/Utilities/KalDet/src/ild/sit/ILDSITKalDetector.h +++ /dev/null @@ -1,60 +0,0 @@ -#ifndef __ILDSITKALDETECTOR__ -#define __ILDSITKALDETECTOR__ - -/** Ladder based SIT to be used for ILD DBD studies - * - * @author S.Aplin DESY - */ - -#include "kaltest/TVKalDetector.h" - -#include "TMath.h" - -class TNode; - -namespace gear{ - class GearMgr ; -} - -class IGeomSvc; - -class ILDSITKalDetector : public TVKalDetector { - -public: - - /** Initialize the SIT from GEAR */ - ILDSITKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* geoSvc ); - - -private: - - void setupGearGeom( const gear::GearMgr& gearMgr ) ; - void setupGearGeom( IGeomSvc* geoSvc ); - - int _nLayers ; - double _bZ ; - - bool _isStripDetector; - - struct SIT_Layer { - int nLadders; - int nSensorsPerLadder; - double phi0; - double dphi; - double senRMin; - double supRMin; - double length; - double width; - double offset; - double senThickness; - double supThickness; - double sensorLength; - double stripAngle; - }; - std::vector<SIT_Layer> _SITgeo; - -}; - - - -#endif diff --git a/Utilities/KalDet/src/ild/support/ILDSupportKalDetector.cc b/Utilities/KalDet/src/ild/support/ILDSupportKalDetector.cc index 9db65dd2..44484c7d 100644 --- a/Utilities/KalDet/src/ild/support/ILDSupportKalDetector.cc +++ b/Utilities/KalDet/src/ild/support/ILDSupportKalDetector.cc @@ -1,9 +1,9 @@ -#include "ILDSupportKalDetector.h" -#include "ILDCylinderMeasLayer.h" -#include "ILDConeMeasLayer.h" -#include "ILDPolygonBarrelMeasLayer.h" -#include "ILDDiscMeasLayer.h" +#include "kaldet/ILDSupportKalDetector.h" +#include "kaldet/ILDCylinderMeasLayer.h" +#include "kaldet/ILDConeMeasLayer.h" +#include "kaldet/ILDPolygonBarrelMeasLayer.h" +#include "kaldet/ILDDiscMeasLayer.h" #include "TMath.h" #include "TTUBE.h" @@ -11,7 +11,7 @@ #include <UTIL/BitField64.h> #include <UTIL/ILDConf.h> -#include "MaterialDataBase.h" +#include "kaldet/MaterialDataBase.h" #include <sstream> #include <cmath> diff --git a/Utilities/KalDet/src/ild/support/ILDSupportKalDetector.h b/Utilities/KalDet/src/ild/support/ILDSupportKalDetector.h deleted file mode 100644 index 373fe51b..00000000 --- a/Utilities/KalDet/src/ild/support/ILDSupportKalDetector.h +++ /dev/null @@ -1,37 +0,0 @@ -#ifndef __ILDSUPPORTDETECTOR__ -#define __ILDSUPPORTDETECTOR__ - -/** Support Material to be used for ILD DBD studies - * - * @author S.Aplin DESY - */ - -#include "kaltest/TVKalDetector.h" - -class TNode; - -namespace gear{ - class GearMgr ; -} -class IGeomSvc; - -class ILDCylinderMeasLayer; - -class ILDSupportKalDetector : public TVKalDetector { -public: - - /** Initialize the support structures from GEAR */ - ILDSupportKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* geoSvc ); - - /** Returns the special layer inside the Beam Pipe used for propagation to the IP */ - ILDCylinderMeasLayer* getIPLayer() { return _ipLayer; } - -private: - - ILDCylinderMeasLayer* _ipLayer; - - std::vector<std::string> _surface_names; - -}; - -#endif diff --git a/Utilities/KalDet/src/ild/tpc/ILDTPCKalDetector.cc b/Utilities/KalDet/src/ild/tpc/ILDTPCKalDetector.cc index 3344de26..301fb9fd 100644 --- a/Utilities/KalDet/src/ild/tpc/ILDTPCKalDetector.cc +++ b/Utilities/KalDet/src/ild/tpc/ILDTPCKalDetector.cc @@ -1,13 +1,13 @@ -#include "ILDTPCKalDetector.h" -#include "ILDCylinderMeasLayer.h" -#include "ILDCylinderHit.h" -#include "ILDDiscMeasLayer.h" +#include "kaldet/ILDTPCKalDetector.h" +#include "kaldet/ILDCylinderMeasLayer.h" +#include "kaldet/ILDCylinderHit.h" +#include "kaldet/ILDDiscMeasLayer.h" #include "TMath.h" #include "TTUBE.h" -#include "MaterialDataBase.h" +#include "kaldet/MaterialDataBase.h" #include <sstream> @@ -29,7 +29,6 @@ // #include "streamlog/streamlog.h" - ILDTPCKalDetector::ILDTPCKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* geoSvc ) : TVKalDetector(250) // SJA:FIXME initial size, 250 looks reasonable for ILD, though this would be better stored as a const somewhere { diff --git a/Utilities/KalDet/src/ild/tpc/ILDTPCKalDetector.h b/Utilities/KalDet/src/ild/tpc/ILDTPCKalDetector.h deleted file mode 100644 index 48012991..00000000 --- a/Utilities/KalDet/src/ild/tpc/ILDTPCKalDetector.h +++ /dev/null @@ -1,30 +0,0 @@ -#ifndef __ILDTPCDETECTOR__ -#define __ILDTPCDETECTOR__ - -/** TPC to be used for ILD DBD studies - * - * @author S.Aplin DESY - */ - -#include "kaltest/TVKalDetector.h" - -class TNode; - -namespace gear{ - class GearMgr ; -} - -class IGeomSvc; - -class ILDTPCKalDetector : public TVKalDetector { -public: - - /** Initialize the TPC from GEAR */ - ILDTPCKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* geoSvc=0 ); - - -private: - -}; - -#endif diff --git a/Utilities/KalDet/src/ild/vxd/CEPCVTXKalDetector.cc b/Utilities/KalDet/src/ild/vxd/CEPCVTXKalDetector.cc index 093b6bda..8f9a3304 100644 --- a/Utilities/KalDet/src/ild/vxd/CEPCVTXKalDetector.cc +++ b/Utilities/KalDet/src/ild/vxd/CEPCVTXKalDetector.cc @@ -1,11 +1,10 @@ #include "kaldet/CEPCVTXKalDetector.h" +#include "kaldet/MaterialDataBase.h" -#include "MaterialDataBase.h" - -#include "ILDParallelPlanarMeasLayer.h" -#include "ILDCylinderMeasLayer.h" -#include "ILDDiscMeasLayer.h" +#include "kaldet/ILDParallelPlanarMeasLayer.h" +#include "kaldet/ILDCylinderMeasLayer.h" +#include "kaldet/ILDDiscMeasLayer.h" #include <UTIL/BitField64.h> #include <UTIL/ILDConf.h> diff --git a/Utilities/KalDet/src/ild/vxd/ILDVXDKalDetector.cc b/Utilities/KalDet/src/ild/vxd/ILDVXDKalDetector.cc index 3cf66696..c1ee0c97 100644 --- a/Utilities/KalDet/src/ild/vxd/ILDVXDKalDetector.cc +++ b/Utilities/KalDet/src/ild/vxd/ILDVXDKalDetector.cc @@ -1,11 +1,9 @@ -#include "ILDVXDKalDetector.h" - -#include "MaterialDataBase.h" - -#include "ILDParallelPlanarMeasLayer.h" -#include "ILDCylinderMeasLayer.h" -#include "ILDDiscMeasLayer.h" +#include "kaldet/ILDVXDKalDetector.h" +#include "kaldet/ILDParallelPlanarMeasLayer.h" +#include "kaldet/ILDCylinderMeasLayer.h" +#include "kaldet/ILDDiscMeasLayer.h" +#include "kaldet/MaterialDataBase.h" #include <UTIL/BitField64.h> #include <UTIL/ILDConf.h> diff --git a/Utilities/KalDet/src/ild/vxd/ILDVXDKalDetector.h b/Utilities/KalDet/src/ild/vxd/ILDVXDKalDetector.h deleted file mode 100644 index d1bbaaa2..00000000 --- a/Utilities/KalDet/src/ild/vxd/ILDVXDKalDetector.h +++ /dev/null @@ -1,74 +0,0 @@ -#ifndef __ILDVXDKALDETECTOR__ -#define __ILDVXDKALDETECTOR__ - -/** Ladder based VXD to be used for ILD DBD studies - * - * @author S.Aplin DESY - */ - -#include "kaltest/TVKalDetector.h" - -#include "TMath.h" - -class TNode; -class IGeomSvc; - -namespace gear{ - class GearMgr ; -} - - -class ILDVXDKalDetector : public TVKalDetector { - -public: - - ILDVXDKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* geoSvc); - - -private: - - void setupGearGeom( const gear::GearMgr& gearMgr ); - void setupGearGeom( IGeomSvc* geoSvc) ; - - int _nLayers ; - double _bZ ; - - double _relative_position_of_measurement_surface ; - - struct VXD_Layer { - int nLadders; - double phi0; - double dphi; - double senRMin; - double supRMin; - double length; - double width; - double offset; - double senThickness; - double supThickness; - }; - std::vector<VXD_Layer> _VXDgeo; - - - struct VXD_Cryostat { - double alRadius; - double alThickness; - double alInnerR; - double alZEndCap; - double alHalfZ; - double shellRadius; - double shellThickness; - double shellInnerR; - double shellZEndCap; - double shelllHalfZ; - - bool exists; - }; - - VXD_Cryostat _vxd_Cryostat; - -}; - - - -#endif -- GitLab