diff --git a/Detector/DetDriftChamber/CMakeLists.txt b/Detector/DetDriftChamber/CMakeLists.txt index b680276955849c32df5af4e92216e9cab21accf9..0c1cc85b2d26105a245427a369270ad4027769a5 100644 --- a/Detector/DetDriftChamber/CMakeLists.txt +++ b/Detector/DetDriftChamber/CMakeLists.txt @@ -27,9 +27,9 @@ gaudi_add_module(DetDriftChamber ${DetDriftChamber_src} INCLUDE_DIRS # DD4hep ROOT Geant4 src/include - LINK_LIBRARIES + LINK_LIBRARIES GaudiKernel DD4hep ROOT DetSegmentation # GaudiKernel - DD4hep ${DD4hep_COMPONENT_LIBRARIES} + #DD4hep ${DD4hep_COMPONENT_LIBRARIES} # ROOT Geant4 ) diff --git a/Detector/DetDriftChamber/compact/det.xml b/Detector/DetDriftChamber/compact/det.xml index 84a1ec7b05d00220ac3fda86c50f4d0111b550b7..b791565dba6598a8c2edc8a178b30afd2c30e25a 100644 --- a/Detector/DetDriftChamber/compact/det.xml +++ b/Detector/DetDriftChamber/compact/det.xml @@ -38,6 +38,7 @@ <constant name="SDT_inner_chamber_layer_number" value="67"/> <constant name="SDT_outer_chamber_layer_number" value="63"/> <constant name="SDT_chamber_layer_width" value="10*mm"/> + <constant name="Epsilon" value="0*deg"/> </define> @@ -52,14 +53,15 @@ <detectors> <detector id="1" name="DriftChamber" type="DriftChamber" readout="DriftChamberHitsCollection" vis="VisibleBlue" sensitive="true"> <!-- Use cm as unit if you want to use Pandora for reconstruction --> + <sensitive type="SimpleDriftChamber"/> </detector> </detectors> <readouts> <readout name="DriftChamberHitsCollection"> - <segmentation type="GridDriftChamber" delta_phi="8*deg" identifier_phi="phi"/> + <segmentation type="GridDriftChamber" cell_size="10*mm" offset_phi="0." epsilon0="Epsilon" detector_length="SDT_half_length" identifier_phi="cellID" /> - <id>system:8,chamber:1,layer:8,phi:16</id> + <id>system:8,chamber:1,layer:8,cellID:16</id> </readout> </readouts> diff --git a/Detector/DetDriftChamber/src/driftchamber/DriftChamber.cpp b/Detector/DetDriftChamber/src/driftchamber/DriftChamber.cpp index c71012156662c8b053567bf2667d535c097016a6..07b95b983df469f822a61b22d7ba2827e27b50eb 100644 --- a/Detector/DetDriftChamber/src/driftchamber/DriftChamber.cpp +++ b/Detector/DetDriftChamber/src/driftchamber/DriftChamber.cpp @@ -8,8 +8,10 @@ //==================================================================== #include "DD4hep/DetFactoryHelper.h" +#include "DD4hep/Printout.h" #include "XML/Layering.h" #include "XML/Utilities.h" +#include "XML/XMLElements.h" #include "DDRec/DetectorData.h" #include "DDSegmentation/Segmentation.h" #include "DetSegmentation/GridDriftChamber.h" @@ -24,6 +26,9 @@ using namespace dd4hep::rec ; static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector, xml_h e, dd4hep::SensitiveDetector sens) { + // ------- Lambda functions ---- // + auto delta_a_func = [](auto x, auto y) { return 0.5 * ( x + y ); }; + // ======================================================================= // Parameter Definition // ======================================================================= @@ -33,6 +38,8 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector, std::string det_name = x_det.nameStr(); std::string det_type = x_det.typeStr(); + dd4hep::SensitiveDetector sd = sens; + // - global double chamber_radius_min = theDetector.constant<double>("SDT_radius_min"); double chamber_radius_max = theDetector.constant<double>("SDT_radius_max"); @@ -53,6 +60,8 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector, int outer_chamber_layer_number = theDetector.constant<int>("SDT_outer_chamber_layer_number"); double chamber_layer_width = theDetector.constant<double>("SDT_chamber_layer_width"); + double epsilon = theDetector.constant<double>("Epsilon"); + // ======================================================================= // Detector Construction // ======================================================================= @@ -75,36 +84,55 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector, dd4hep::Tube det_outer_chamber_solid(outer_chamber_radius_min, outer_chamber_radius_max, outer_chamber_length*0.5); dd4hep::Volume det_outer_chamber_vol(det_name+"_outer_chamber_vol", det_outer_chamber_solid, det_mat); + //Initialize the segmentation + dd4hep::Readout readout = sd.readout(); + dd4hep::Segmentation geomseg = readout.segmentation(); + dd4hep::Segmentation* _geoSeg = &geomseg; + + auto DCHseg = dynamic_cast<dd4hep::DDSegmentation::GridDriftChamber*>(_geoSeg->segmentation()); // - layer - for(int layer_id=0; layer_id<(inner_chamber_layer_number+outer_chamber_layer_number-1);layer_id++) { - double rmin,rmax; - std::string layer_name; - dd4hep::Volume* current_vol_ptr = nullptr; - - if(layer_id<inner_chamber_layer_number){ - current_vol_ptr = &det_inner_chamber_vol; - rmin = inner_chamber_radius_min+(layer_id*chamber_layer_width); - rmax = rmin+chamber_layer_width; - layer_name = det_name+"_inner_chamber_vol"+_toString(layer_id,"_layer%d"); - }else{ - current_vol_ptr = &det_outer_chamber_vol; - rmin = outer_chamber_radius_min+((layer_id-inner_chamber_layer_number)*chamber_layer_width); - rmax = rmin+chamber_layer_width; - layer_name = det_name+"_outer_chamber_vol"+_toString(layer_id,"_layer%d"); + for(int layer_id = 0; layer_id < (inner_chamber_layer_number+outer_chamber_layer_number); layer_id++) { + double rmin,rmax,offset; + std::string layer_name; + dd4hep::Volume* current_vol_ptr = nullptr; + + if( layer_id < inner_chamber_layer_number ) { + current_vol_ptr = &det_inner_chamber_vol; + rmin = inner_chamber_radius_min+(layer_id*chamber_layer_width); + rmax = rmin+chamber_layer_width; + layer_name = det_name+"_inner_chamber_vol"+_toString(layer_id,"_layer%d"); + } + else { + current_vol_ptr = &det_outer_chamber_vol; + rmin = outer_chamber_radius_min+((layer_id-inner_chamber_layer_number)*chamber_layer_width); + rmax = rmin+chamber_layer_width; + layer_name = det_name+"_outer_chamber_vol"+_toString(layer_id,"_layer%d"); + } + + //Construction of drift chamber layers + double rmid = delta_a_func(rmin,rmax); + double ilayer_cir = 2 * M_PI * rmid; + double ncell = ilayer_cir / chamber_layer_width; + int ncell_layer = ceil(ncell); + int numWire = ncell_layer; + double layer_Phi = 2*M_PI / ncell_layer; + if(layer_id %2 ==0){ offset = 0.; } + else { offset = 0.5 * layer_Phi; } + + DCHseg->setGeomParams(layer_id, layer_Phi, rmid, epsilon, offset); + DCHseg->setWiresInLayer(layer_id, numWire); + + dd4hep::Tube layer_solid(rmin,rmax,chamber_length*0.5); + dd4hep::Volume layer_vol(layer_name,layer_solid,det_mat); + dd4hep::Transform3D transform_layer(dd4hep::Rotation3D(),dd4hep::Position(0.,0.,0.)); + dd4hep::PlacedVolume layer_phy = (*current_vol_ptr).placeVolume(layer_vol, transform_layer); + layer_phy.addPhysVolID("layer",layer_id); + + //Set drift chamber layers to sensitive detector + layer_vol.setSensitiveDetector(sens); + sd.setType("tracker"); } - /// Construction of drift chamber layers - dd4hep::Tube layer_solid(rmin,rmax,chamber_length*0.5); - dd4hep::Volume layer_vol(layer_name,layer_solid,det_mat); - dd4hep::Transform3D transform_layer(dd4hep::Rotation3D(),dd4hep::Position(0.,0.,0.)); - dd4hep::PlacedVolume layer_phy = (*current_vol_ptr).placeVolume(layer_vol, transform_layer); - layer_phy.addPhysVolID("layer",layer_id); - - /// Set drift chamber layers to sensitive detector - dd4hep::SensitiveDetector sd = sens; - layer_vol.setSensitiveDetector(sens); - sd.setType("tracker"); - } // - place in det // inner diff --git a/Detector/DetSegmentation/DetSegmentation/GridDriftChamber.h b/Detector/DetSegmentation/DetSegmentation/GridDriftChamber.h index 379ad4a19fdbce33a63f009577970534f3f995c3..fc01bf114d8d1780690ecbc9db2c23a82cc8cbe5 100644 --- a/Detector/DetSegmentation/DetSegmentation/GridDriftChamber.h +++ b/Detector/DetSegmentation/DetSegmentation/GridDriftChamber.h @@ -6,6 +6,7 @@ #include "TVector3.h" #include <cmath> #include <iostream> +#include <map> /** GridDriftChamber Detector/DetSegmentation/DetSegmentation/GridDriftChamber.h GridDriftChamber.h * @@ -14,6 +15,21 @@ * @author nalipour */ + +typedef struct Layer + { + double layerphi; + double R; + double eps; + double offset; + Layer(){}; + Layer(double x, double y, double z, double k):layerphi(x),R(y),eps(z),offset(k){}; + bool operator < (const Layer &a) const + { + return layerphi < a.layerphi; + } + } LAYER; + namespace dd4hep { namespace DDSegmentation { class GridDriftChamber : public Segmentation { @@ -28,62 +44,102 @@ public: virtual Vector3D position(const CellID& aCellID) const; virtual CellID cellID(const Vector3D& aLocalPosition, const Vector3D& aGlobalPosition, const VolumeID& aVolumeID) const; + virtual double distanceTrackWire(const CellID& cID, const TVector3& hit_start, const TVector3& hit_end) const; -// inline double innerRadius() const { return m_innerRadius; } -// inline double detectorLength() const { return m_detectorLength; } - inline double offsetPhi() const { return m_offsetPhi; } - inline double delta_phi() const{ return m_delta_phi; } + inline double cell_Size() const { return m_cellSize; } + inline double offset_phi() const { return m_offsetPhi; } + inline double epsilon0() const { return m_epsilon0; } + inline double detectorLength() const { return m_detectorLength; } inline const std::string& fieldNamePhi() const { return m_phiID; } // Setters -// inline void setGeomParams(int layer, double sizePhi) { -// layer_params[layer] = {sizePhi}; -// } - -// void updateParams(int layer) const { -// auto it_end = layer_params.cend(); -// --it_end; -// double size = it_end->second[0]; -// double radius = it_end->second[1]; -// double eps = it_end->second[2]; -// -// auto map_it = layer_params.find(layer); -// if (map_it != layer_params.cend()) { -// size = map_it->second[0]; -// radius = map_it->second[1]; -// eps = map_it->second[2]; -// } -// -// _currentGridSizePhi = size; -// _currentRadius = radius; -// m_epsilon = eps; -// } - inline double phiFromXY(const Vector3D& aposition) const { return std::atan2(aposition.Y, aposition.X) + M_PI ; } -// inline int returnLayer(double x, double y) const { -// // Hit R position -// double R = std::sqrt(x * x + y * y); -// // Layer -// int layer = int((R - m_innerRadius) / m_cellSize); -// return layer; -// } + inline void setGeomParams(int layer, double layerphi, double R, double eps, double offset ) { + // layer_params[layer] = {layerphi,R,eps}; + layer_params.insert(std::pair<int,LAYER>(layer,LAYER(layerphi,R,eps,offset))); + } + + inline void setWiresInLayer(int layer, int numWires) + { + double phi0; + updateParams(layer); + for (int i = 0; i<numWires; ++i) { + +// if(layer % 2 == 0) { phi0 = 0.; } +// else { phi0 = 0.5 * _currentLayerphi; } + double phi0 = m_offset; + + auto phi_start = _currentLayerphi * i + phi0; + if(phi_start > 2 * M_PI) { phi_start = phi_start - 2 * M_PI; } + auto phi_end = phi_start + _currentLayerphi; + + TVector3 Wstart = returnWirePosition(phi_start, 1); + TVector3 Wend = returnWirePosition(phi_end, -1); + + TVector3 Wmid = (Wstart+Wend)*(1/2.0); + TVector3 Wdirection = (Wend - Wstart); + + m_wiresPositions[layer].push_back(std::make_pair(Wmid, Wdirection)); + } + } + + inline auto returnAllWires() const { return m_wiresPositions; } + + inline TVector3 returnWirePosition(double angle, int sign) const { + TVector3 w(0, 0, 0); + w.SetX(_currentRadius * std::cos(angle)); + w.SetY(_currentRadius * std::sin(angle)); + w.SetZ(sign * m_detectorLength / 2.0); + return w; + } + + void updateParams(int layer) const{ + auto it_end = layer_params.cend(); + --it_end; + double layerphi = it_end->second.layerphi; + double radius = it_end->second.R; + double eps = it_end->second.eps; + double offset = it_end->second.offset; + + auto map_it = layer_params.find(layer); + if (map_it != layer_params.cend()) { + layerphi = map_it->second.layerphi; + radius = map_it->second.R; + eps = map_it->second.eps; + offset = map_it->second.offset; + } + _currentLayerphi = layerphi; + _currentRadius = radius; + m_epsilon = eps; + m_offset = offset; + } + + inline double returnAlpha() const { + double alpha = 2 * std::asin(m_detectorLength * std::tan(m_epsilon0)/(2 * _currentRadius)); + return alpha; + } protected: /* *** nalipour *** */ double phi(const CellID& cID) const; + std::map<int,LAYER> layer_params; // <layer, {layerphi, R, eps, offset}> + std::map<int, std::vector<std::pair<TVector3, TVector3> >> m_wiresPositions; // < layer, vec<WireMidpoint, WireDirection> > - + double m_cellSize; double m_offsetPhi; - double m_delta_phi; + double m_epsilon0; + double m_detectorLength; std::string m_phiID; // Current parameters of the layer: sizePhi -// mutable double _currentGridSizePhi; // current size Phi -// mutable double _currentRadius; // current size radius -// mutable double m_epsilon; + mutable double _currentLayerphi; + mutable double _currentRadius; + mutable double m_epsilon; + mutable double m_offset; + }; } } diff --git a/Detector/DetSegmentation/src/GridDriftChamber.cpp b/Detector/DetSegmentation/src/GridDriftChamber.cpp index 4b1ed216fb0825e1b41a210871ca2451954ee22e..57a6d7ebacd9fb221dbb7e59dfedd2a14a52894b 100644 --- a/Detector/DetSegmentation/src/GridDriftChamber.cpp +++ b/Detector/DetSegmentation/src/GridDriftChamber.cpp @@ -1,4 +1,5 @@ #include "DetSegmentation/GridDriftChamber.h" +#include <map> namespace dd4hep { namespace DDSegmentation { @@ -9,17 +10,22 @@ GridDriftChamber::GridDriftChamber(const std::string& cellEncoding) : Segmentati _type = "GridDriftChamber"; _description = "Drift chamber segmentation in the global coordinates"; - registerIdentifier("identifier_phi", "Cell ID identifier for phi", m_phiID, "phi"); - registerParameter("delta_phi", "delta phi", m_delta_phi, 0., SegmentationParameter::LengthUnit); + registerParameter("cell_size", "cell size", m_cellSize, 0., SegmentationParameter::LengthUnit); + registerParameter("offset_phi", "offset in phi", m_offsetPhi, 0., SegmentationParameter::LengthUnit, true); + registerParameter("detector_length", "Length of the wire", m_detectorLength, 1., SegmentationParameter::LengthUnit); + registerIdentifier("identifier_phi", "Cell ID identifier for phi", m_phiID, "cellID"); } GridDriftChamber::GridDriftChamber(const BitFieldCoder* decoder) : Segmentation(decoder) { - // define type and description + _type = "GridDriftChamber"; _description = "Drift chamber segmentation in the global coordinates"; - registerIdentifier("identifier_phi", "Cell ID identifier for phi", m_phiID, "phi"); - registerParameter("delta_phi", "delta phi", m_delta_phi, 0., SegmentationParameter::LengthUnit); + registerParameter("cell_size", "cell size", m_cellSize, 1., SegmentationParameter::LengthUnit); + registerParameter("offset_phi", "offset in phi", m_offsetPhi, 0., SegmentationParameter::LengthUnit, true); + registerParameter("epsilon0", "epsilon", m_epsilon0, 0., SegmentationParameter::AngleUnit, true); + registerParameter("detector_length", "Length of the wire", m_detectorLength, 1., SegmentationParameter::LengthUnit); + registerIdentifier("identifier_phi", "Cell ID identifier for phi", m_phiID, "cellID"); } Vector3D GridDriftChamber::position(const CellID& /*cID*/) const { @@ -27,28 +33,78 @@ Vector3D GridDriftChamber::position(const CellID& /*cID*/) const { return cellPosition; } + CellID GridDriftChamber::cellID(const Vector3D& /*localPosition*/, const Vector3D& globalPosition, const VolumeID& vID) const { CellID cID = vID; + unsigned int layerID = _decoder->get(vID, "layer"); + updateParams(layerID); double phi_hit = phiFromXY(globalPosition); double posx = globalPosition.X; double posy = globalPosition.Y; - - int lphi = (int) (phi_hit/m_delta_phi); + double offsetphi= m_offset; + int _lphi; +// if(layerID % 2 == 0) { +// offsetphi = 0.; +// _lphi = (int) (phi_hit / _currentLayerphi); +// } +// else { +// offsetphi = _currentLayerphi / 2.; + if(phi_hit >= offsetphi) { + _lphi = (int) ((phi_hit - offsetphi)/ _currentLayerphi); + } + else { + _lphi = (int) ((phi_hit - offsetphi + 2 * M_PI)/ _currentLayerphi); + } + int lphi = _lphi; _decoder->set(cID, m_phiID, lphi); -// std::cout << " myliu: " -// << " x: " << posx -// << " y: " << posy -//// << " pre: " << phi_pre -// << " phi_hit: " << phi_hit -// << " lphi: " << lphi -// << std::endl; +//std::cout << "#######################################: " +// << " offset : " << m_offset +// << " offsetphi: " << offsetphi +// << " layerID: " << layerID +// << " r: " << _currentRadius +// << " layerphi: " << _currentLayerphi +// << std::endl; + return cID; } +double GridDriftChamber::phi(const CellID& cID) const { + CellID phiValue = _decoder->get(cID, m_phiID); + return binToPosition(phiValue, _currentLayerphi, m_offsetPhi); +} + +double GridDriftChamber::distanceTrackWire(const CellID& cID, const TVector3& hit_start, + const TVector3& hit_end) const { + + auto layerIndex = _decoder->get(cID, "layer"); + updateParams(layerIndex); + + double phi_start = phi(cID); + double phi_end = phi_start + returnAlpha(); + + TVector3 Wstart = returnWirePosition(phi_start, 1); + TVector3 Wend = returnWirePosition(phi_end, -1); + + TVector3 a = hit_end - hit_start; + TVector3 b = Wend - Wstart; + TVector3 c = Wstart - hit_start; + + double num = std::abs(c.Dot(a.Cross(b))); + double denum = (a.Cross(b)).Mag(); + + double DCA = 0; + + if (denum) { + DCA = num / denum; + } + + return DCA; +} + REGISTER_SEGMENTATION(GridDriftChamber) } diff --git a/Digitisers/DCHDigi/CMakeLists.txt b/Digitisers/DCHDigi/CMakeLists.txt index 33fa374effccc19a7fabd6fb841a3e34f368c067..a2d07c1b835859cf4c3e4deff65fcb055cad2a09 100644 --- a/Digitisers/DCHDigi/CMakeLists.txt +++ b/Digitisers/DCHDigi/CMakeLists.txt @@ -1,12 +1,12 @@ gaudi_subdir(DCHDigi v0r0) +find_package(CLHEP REQUIRED;CONFIG) find_package(DD4hep COMPONENTS DDG4 REQUIRED) find_package(EDM4HEP REQUIRED ) message("EDM4HEP_INCLUDE_DIRS: ${EDM4HEP_INCLUDE_DIR}") message("EDM4HEP_LIB: ${EDM4HEP_LIBRARIES}") include_directories(${EDM4HEP_INCLUDE_DIR}) -find_package(CLHEP REQUIRED) find_package(podio REQUIRED ) set(srcs @@ -18,8 +18,8 @@ gaudi_depends_on_subdirs( ) ## Modules gaudi_add_module(DCHDigi ${srcs} - INCLUDE_DIRS FWCore GaudiKernel GaudiAlgLib CLHEP DD4hep - LINK_LIBRARIES FWCore GaudiKernel GaudiAlgLib CLHEP DD4hep ${DD4hep_COMPONENT_LIBRARIES} DDRec + INCLUDE_DIRS FWCore GaudiKernel GaudiAlgLib ${CLHEP_INCLUDE_DIR} DD4hep + LINK_LIBRARIES FWCore GaudiKernel GaudiAlgLib ${CLHEP_LIBRARIES} DD4hep ${DD4hep_COMPONENT_LIBRARIES} DDRec -Wl,--no-as-needed EDM4HEP::edm4hep EDM4HEP::edm4hepDict ) diff --git a/Digitisers/G2CDArbor/CMakeLists.txt b/Digitisers/G2CDArbor/CMakeLists.txt index 8b35957aebdcf428770521dd13169739ae081115..d1895d1c6a2d45af08bec4cf6ea4854985974d2d 100644 --- a/Digitisers/G2CDArbor/CMakeLists.txt +++ b/Digitisers/G2CDArbor/CMakeLists.txt @@ -1,5 +1,6 @@ gaudi_subdir(G2CDArbor v0r0) +find_package(CLHEP REQUIRED;CONFIG) find_package(DD4hep COMPONENTS DDG4 REQUIRED) find_package(EDM4HEP REQUIRED) find_package(GEAR REQUIRED) @@ -21,7 +22,7 @@ set(G2CDArbor_srcs src/*.cpp) # Modules gaudi_add_module(G2CDArbor ${G2CDArbor_srcs} - INCLUDE_DIRS K4FWCore GaudiKernel GaudiAlgLib CLHEP DD4hep gear ${GSL_INCLUDE_DIRS} ${LCIO_INCLUDE_DIRS} - LINK_LIBRARIES K4FWCore GaudiKernel GaudiAlgLib CLHEP DD4hep ${GEAR_LIBRARIES} ${GSL_LIBRARIES} ${LCIO_LIBRARIES} + INCLUDE_DIRS K4FWCore GaudiKernel GaudiAlgLib ${CLHEP_INCLUDE_DIR} DD4hep gear ${GSL_INCLUDE_DIRS} ${LCIO_INCLUDE_DIRS} + LINK_LIBRARIES K4FWCore GaudiKernel GaudiAlgLib ${CLHEP_LIBRARIES} DD4hep ${GEAR_LIBRARIES} ${GSL_LIBRARIES} ${LCIO_LIBRARIES} EDM4HEP::edm4hep EDM4HEP::edm4hepDict ) diff --git a/Digitisers/SimpleDigi/CMakeLists.txt b/Digitisers/SimpleDigi/CMakeLists.txt index 312af33b7b84e5463efa921ec27c8fc8539fdba4..39f60902538c88005ccca1d64f9e486e9db2c0f2 100644 --- a/Digitisers/SimpleDigi/CMakeLists.txt +++ b/Digitisers/SimpleDigi/CMakeLists.txt @@ -1,5 +1,6 @@ gaudi_subdir(SimpleDigi v0r0) +find_package(CLHEP REQUIRED;CONFIG) find_package(GEAR REQUIRED) find_package(GSL REQUIRED ) find_package(LCIO REQUIRED ) @@ -18,6 +19,6 @@ set(SimpleDigi_srcs src/*.cpp) # Modules gaudi_add_module(SimpleDigi ${SimpleDigi_srcs} - INCLUDE_DIRS K4FWCore GaudiKernel GaudiAlgLib gear ${GSL_INCLUDE_DIRS} ${LCIO_INCLUDE_DIRS} - LINK_LIBRARIES K4FWCore GaudiKernel GaudiAlgLib ${GEAR_LIBRARIES} ${GSL_LIBRARIES} ${LCIO_LIBRARIES} EDM4HEP::edm4hep EDM4HEP::edm4hepDict DataHelperLib + INCLUDE_DIRS K4FWCore GaudiKernel GaudiAlgLib ${CLHEP_INCLUDE_DIR} gear ${GSL_INCLUDE_DIRS} ${LCIO_INCLUDE_DIRS} + LINK_LIBRARIES K4FWCore GaudiKernel GaudiAlgLib ${CLHEP_LIBRARIES} ${GEAR_LIBRARIES} ${GSL_LIBRARIES} ${LCIO_LIBRARIES} EDM4HEP::edm4hep EDM4HEP::edm4hepDict DataHelperLib ) diff --git a/Digitisers/SimpleDigi/src/PlanarDigiAlg.cpp b/Digitisers/SimpleDigi/src/PlanarDigiAlg.cpp index 113f2d8b30c5a77086bdc19ffedbe7323468c35e..482082a064e982a0452b93e4a945b5e2d35d13f7 100644 --- a/Digitisers/SimpleDigi/src/PlanarDigiAlg.cpp +++ b/Digitisers/SimpleDigi/src/PlanarDigiAlg.cpp @@ -91,10 +91,10 @@ StatusCode PlanarDigiAlg::initialize() return StatusCode::FAILURE; } - MarlinTrk::IMarlinTrkSystem* _trksystem = _trackSystemSvc->getTrackSystem(); + MarlinTrk::IMarlinTrkSystem* _trksystem = _trackSystemSvc->getTrackSystem(this); _trksystem->init(); - _trackSystemSvc->removeTrackSystem(); + _trackSystemSvc->removeTrackSystem(this); return GaudiAlgorithm::initialize(); } @@ -338,8 +338,8 @@ StatusCode PlanarDigiAlg::execute() float weight = 1.0; debug() <<" Set relation between " - << " sim hit " << SimTHit - << " to tracker hit " << trkHit + << " sim hit " << SimTHit.id() + << " to tracker hit " << trkHit.id() << " with a weight of " << weight << endmsg; trkHit.addToRawHits(SimTHit.getObjectID()); diff --git a/Generator/CMakeLists.txt b/Generator/CMakeLists.txt index 41f325ec59a21f336721425edae877e76a822e5b..435b84d90f96b2da6d9a9d203ab4ed9204d4fd4d 100644 --- a/Generator/CMakeLists.txt +++ b/Generator/CMakeLists.txt @@ -23,7 +23,7 @@ find_package(LCIO) find_package(podio) find_package(EDM4HEP) find_package(HepMC) -find_package(CLHEP) +find_package(CLHEP REQUIRED;CONFIG) find_package(K4FWCore REQUIRED) if(ROOT_FOUND) diff --git a/Reconstruction/Digi_Calo/CMakeLists.txt b/Reconstruction/Digi_Calo/CMakeLists.txt index 507720ca40cd718c485744ef5e95e68c0c674bd8..b61ef2955926e93f309dad0324041656fc6cb827 100644 --- a/Reconstruction/Digi_Calo/CMakeLists.txt +++ b/Reconstruction/Digi_Calo/CMakeLists.txt @@ -1,5 +1,6 @@ gaudi_subdir(Digi_Calo v0r0) +find_package(CLHEP REQUIRED;CONFIG) find_package(DD4hep COMPONENTS DDG4 REQUIRED) find_package(EDM4HEP REQUIRED ) message("EDM4HEP_INCLUDE_DIRS: ${EDM4HEP_INCLUDE_DIR}") @@ -18,8 +19,8 @@ gaudi_depends_on_subdirs( ) ## Modules gaudi_add_module(Digi_Calo ${srcs} - INCLUDE_DIRS FWCore GaudiKernel GaudiAlgLib CLHEP DD4hep - LINK_LIBRARIES FWCore GaudiKernel GaudiAlgLib CLHEP DD4hep ${DD4hep_COMPONENT_LIBRARIES} DDRec + INCLUDE_DIRS FWCore GaudiKernel GaudiAlgLib ${CLHEP_INCLUDE_DIR} DD4hep + LINK_LIBRARIES FWCore GaudiKernel GaudiAlgLib ${CLHEP_LIBRARIES} DD4hep ${DD4hep_COMPONENT_LIBRARIES} DDRec -Wl,--no-as-needed EDM4HEP::edm4hep EDM4HEP::edm4hepDict ) diff --git a/Reconstruction/PFA/Pandora/GaudiPandora/CMakeLists.txt b/Reconstruction/PFA/Pandora/GaudiPandora/CMakeLists.txt index 48415cc53781d08424181497155480c0ffd445a5..be0c48674b7465d7ff1bfd607eff6b172439933c 100644 --- a/Reconstruction/PFA/Pandora/GaudiPandora/CMakeLists.txt +++ b/Reconstruction/PFA/Pandora/GaudiPandora/CMakeLists.txt @@ -4,7 +4,7 @@ find_package(LCIO REQUIRED ) find_package(GEAR REQUIRED) message("ENV GEAR: $ENV{GEAR}") - +find_package(CLHEP REQUIRED;CONFIG) find_package(EDM4HEP REQUIRED ) include_directories(${EDM4HEP_INCLUDE_DIR}) @@ -37,8 +37,8 @@ set(dir_srcs set(dir_include include) # Modules gaudi_add_module(GaudiPandora ${dir_srcs} - INCLUDE_DIRS ${dir_include} GaudiKernel FWCore CLHEP ${LCIO_INCLUDE_DIRS} ${ROOT_INCLUDE_DIRS} gear - LINK_LIBRARIES GaudiKernel FWCore CLHEP ROOT ${LCIO_LIBRARIES} ${GEAR_LIBRARIES} DataHelperLib + INCLUDE_DIRS ${dir_include} GaudiKernel FWCore ${CLHEP_INCLUDE_DIR} ${LCIO_INCLUDE_DIRS} ${ROOT_INCLUDE_DIRS} gear + LINK_LIBRARIES GaudiKernel FWCore ${CLHEP_LIBRARIES} ROOT ${LCIO_LIBRARIES} ${GEAR_LIBRARIES} DataHelperLib -Wl,--no-as-needed EDM4HEP::edm4hep EDM4HEP::edm4hepDict -Wl,--as-needed diff --git a/Reconstruction/PFA/Pandora/MatrixPandora/CMakeLists.txt b/Reconstruction/PFA/Pandora/MatrixPandora/CMakeLists.txt index f99292cbd2294b48f7841c2524b09247f128ca75..c8a76314f563fa9d8a37bc7c65e106afa018d958 100644 --- a/Reconstruction/PFA/Pandora/MatrixPandora/CMakeLists.txt +++ b/Reconstruction/PFA/Pandora/MatrixPandora/CMakeLists.txt @@ -1,7 +1,7 @@ gaudi_subdir(MatrixPandora v0r0) find_package(DD4hep COMPONENTS DDG4 REQUIRED) -find_package(CLHEP REQUIRED) +find_package(CLHEP REQUIRED;CONFIG) find_package(LCIO REQUIRED ) find_package(GEAR REQUIRED) find_package(EDM4HEP REQUIRED ) @@ -37,8 +37,8 @@ set(dir_srcs set(dir_include include) # Modules gaudi_add_module(MatrixPandora ${dir_srcs} - INCLUDE_DIRS ${dir_include} GaudiKernel FWCore CLHEP ${LCIO_INCLUDE_DIRS} ${ROOT_INCLUDE_DIRS} gear DD4hep - LINK_LIBRARIES GaudiKernel FWCore CLHEP ROOT ${LCIO_LIBRARIES} ${GEAR_LIBRARIES} DD4hep ${DD4hep_COMPONENT_LIBRARIES} DDRec DataHelperLib + INCLUDE_DIRS ${dir_include} GaudiKernel FWCore ${CLHEP_INCLUDE_DIR} ${LCIO_INCLUDE_DIRS} ${ROOT_INCLUDE_DIRS} gear DD4hep + LINK_LIBRARIES GaudiKernel FWCore ${CLHEP_LIBRARIES} ROOT ${LCIO_LIBRARIES} ${GEAR_LIBRARIES} DD4hep ${DD4hep_COMPONENT_LIBRARIES} DDRec DataHelperLib -Wl,--no-as-needed EDM4HEP::edm4hep EDM4HEP::edm4hepDict -Wl,--as-needed diff --git a/Reconstruction/SiliconTracking/src/ForwardTrackingAlg.cpp b/Reconstruction/SiliconTracking/src/ForwardTrackingAlg.cpp index 6833e4041883b5bfc42ec21a9ac053bf4b5b8d70..7583b9f3453926b6231779a94e62adc1cf133fe1 100644 --- a/Reconstruction/SiliconTracking/src/ForwardTrackingAlg.cpp +++ b/Reconstruction/SiliconTracking/src/ForwardTrackingAlg.cpp @@ -137,7 +137,7 @@ StatusCode ForwardTrackingAlg::initialize(){ error() << "Failed to find TrackSystemSvc ..." << endmsg; return StatusCode::FAILURE; } - _trkSystem = _trackSystemSvc->getTrackSystem(); + _trkSystem = _trackSystemSvc->getTrackSystem(this); if( _trkSystem == 0 ){ error() << "Cannot initialize MarlinTrkSystem of Type: KalTest" <<endmsg; @@ -150,7 +150,7 @@ StatusCode ForwardTrackingAlg::initialize(){ _trkSystem->setOption( MarlinTrk::IMarlinTrkSystem::CFG::useSmoothing, _SmoothOn) ; //smoothing // initialise the tracking system - //_trkSystem->init() ; + _trkSystem->init() ; /**********************************************************************************************/ /* Do a few checks, if the set parameters are right */ diff --git a/Reconstruction/SiliconTracking/src/SiliconTrackingAlg.cpp b/Reconstruction/SiliconTracking/src/SiliconTrackingAlg.cpp index 9dc2964076fb4020e1a99220e9d16839450b7e11..06822bb910b8138d387f1a85e24aebf690941e8b 100644 --- a/Reconstruction/SiliconTracking/src/SiliconTrackingAlg.cpp +++ b/Reconstruction/SiliconTracking/src/SiliconTrackingAlg.cpp @@ -64,7 +64,7 @@ DECLARE_COMPONENT( SiliconTrackingAlg ) SiliconTrackingAlg::SiliconTrackingAlg(const std::string& name, ISvcLocator* svcLoc) : GaudiAlgorithm(name, svcLoc) { - + //_description = "Pattern recognition in silicon trackers"; _fastfitter = new MarlinTrk::HelixFit(); @@ -111,7 +111,7 @@ StatusCode SiliconTrackingAlg::initialize() { error() << "Failed to find TrackSystemSvc ..." << endmsg; return StatusCode::FAILURE; } - _trksystem = _trackSystemSvc->getTrackSystem(); + _trksystem = _trackSystemSvc->getTrackSystem(this); if( _trksystem == 0 ){ error() << "Cannot initialize MarlinTrkSystem of Type: KalTest" <<endmsg; @@ -448,12 +448,12 @@ int SiliconTrackingAlg::InitialiseFTD() { const float eps = 1.0e-07; // V must be the global z axis if( fabs(V.dot(Z)) > eps ) { - error() << "SiliconTrackingAlg: VXD Hit measurment vectors V is not in the global X-Y plane. \n\n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << endmsg; + error() << "SiliconTrackingAlg: FTD Hit measurment vectors V is not in the global X-Y plane. \n\n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << endmsg; exit(1); } if( fabs(U.dot(Z)) > eps ) { - error() << "SiliconTrackingAlg: VXD Hit measurment vectors U is not in the global X-Y plane. \n\n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << endmsg; + error() << "SiliconTrackingAlg: FTD Hit measurment vectors U is not in the global X-Y plane. \n\n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << endmsg; exit(1); } // SJA:FIXME Here dU and dV are almost certainly dX and dY ... should test ... @@ -491,7 +491,7 @@ int SiliconTrackingAlg::InitialiseFTD() { else { layer = 2*layer + 1; } - + //debug() << "_petalBasedFTDWithOverlaps = true layer->2*layer,2layer+1 according to petalIndex = " << petalIndex << endmsg; } if (layer >= _nlayersFTD) { fatal() << "SiliconTrackingAlg => fatal error in FTD : layer is outside allowed range : " << layer << " number of layers = " << _nlayersFTD << endmsg; @@ -509,7 +509,7 @@ int SiliconTrackingAlg::InitialiseFTD() { int iCode = iSemiSphere + 2*layer + 2*_nlayersFTD*iPhi; _sectorsFTD[iCode].push_back( hitExt ); - debug() << " FTD Pixel Hit added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " iPhi = " << iPhi << " iSemiSphere " << iSemiSphere << " iCode = " << iCode << " layer = " << layer << endmsg; + debug() << " FTD Pixel Hit added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " iPhi = " << iPhi << " iSemiSphere " << iSemiSphere << " iCode = " << iCode << " layer = " << layer << " cellID = " << hit.getCellID() << endmsg; } } // Reading out FTD SpacePoint Collection @@ -579,7 +579,6 @@ int SiliconTrackingAlg::InitialiseFTD() { unsigned int petalIndex = static_cast<unsigned int>(getModuleID(hit)); if ( _petalBasedFTDWithOverlaps == true ) { - // as we are dealing with staggered petals we will use 2*nlayers in each directions +/- z // the layers will follow the even odd numbering of the petals if ( petalIndex % 2 == 0 ) { @@ -588,7 +587,7 @@ int SiliconTrackingAlg::InitialiseFTD() { else { layer = 2*layer + 1; } - + //debug() << "_petalBasedFTDWithOverlaps = true layer->2*layer,2layer+1 according to petalIndex = " << petalIndex << endmsg; } if (layer >= _nlayersFTD) { @@ -607,7 +606,7 @@ int SiliconTrackingAlg::InitialiseFTD() { int iCode = iSemiSphere + 2*layer + 2*_nlayersFTD*iPhi; _sectorsFTD[iCode].push_back( hitExt ); - debug() << " FTD SpacePoint Hit added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " iPhi = " << iPhi << " iSemiSphere " << iSemiSphere << " iCode = " << iCode << " layer = " << layer << endmsg; + debug() << " FTD SpacePoint Hit added : @ " << pos[0] << " " << pos[1] << " " << pos[2] << " drphi " << hitExt->getResolutionRPhi() << " dz " << hitExt->getResolutionZ() << " iPhi = " << iPhi << " iSemiSphere " << iSemiSphere << " iCode = " << iCode << " layer = " << layer << " cellID = " << hit.getCellID() << endmsg; } } @@ -677,7 +676,7 @@ int SiliconTrackingAlg::InitialiseVTX() { exit(1); } TrackerHitExtended * hitExt = new TrackerHitExtended(hit); - debug() << "Saved TrackerHit pointer in TrackerHitExtended " << ielem << ": " << hitExt->getTrackerHit() << std::endl; + //debug() << "Saved TrackerHit id in TrackerHitExtended " << ielem << ": " << hitExt->getTrackerHit().id() << std::endl; // SJA:FIXME: just use planar res for now hitExt->setResolutionRPhi(hit.getCovMatrix()[2]); @@ -799,7 +798,7 @@ int SiliconTrackingAlg::InitialiseVTX() { // or a PIXEL based SIT, using 2D TrackerHitPlane like the VXD above // by fucd //else if ( ( trkhit_P = dynamic_cast<TrackerHitPlane*>( hitCollection->getElementAt( ielem ) ) ) ) { - else if( UTIL::BitSet32( trkhit.getType() )[ 31 ]){ + else if( UTIL::BitSet32( trkhit.getType() )[ 3 ]){ // first we need to check if the measurement vectors are aligned with the global coordinates //gear::Vector3D U(1.0,trkhit_P->getU()[1],trkhit_P->getU()[0],gear::Vector3D::spherical); //gear::Vector3D V(1.0,trkhit_P->getV()[1],trkhit_P->getV()[0],gear::Vector3D::spherical); @@ -907,6 +906,7 @@ StatusCode SiliconTrackingAlg::finalize(){ //delete _trksystem ; _trksystem = 0; //delete _histos ; _histos = 0; info() << "Processed " << _nEvt << " events " << endmsg; + info() << lcio::ILDCellID0::encoder_string << " " << UTIL::ILDCellID0::encoder_string << endmsg; return GaudiAlgorithm::finalize(); } @@ -2052,7 +2052,7 @@ void SiliconTrackingAlg::AttachRemainingVTXHitsSlow() { // Here we are trying to find if a hits are too close i.e. closer than _minDistToDelta edm4hep::ConstTrackerHit trkhit1 = hit->getTrackerHit(); edm4hep::ConstTrackerHit trkhit2 = hitVector[IHIT]->getTrackerHit(); - + if ( trkhit1.getCellID() == trkhit2.getCellID() ){ // i.e. they are in the same sensor float distance = 0.; @@ -2691,9 +2691,9 @@ void SiliconTrackingAlg::FinalRefit(edm4hep::TrackCollection* trk_col) { continue ; } //TrackImpl* Track = new TrackImpl ; - auto track = trk_col->create(); + //auto track = trk_col->create(); //fucd - //edm4hep::Track track;// = new edm4hep::Track; + edm4hep::Track track;// = new edm4hep::Track; // setup initial dummy covariance matrix //std::vector<float> covMatrix; //covMatrix.resize(15); @@ -2840,7 +2840,7 @@ void SiliconTrackingAlg::FinalRefit(edm4hep::TrackCollection* trk_col) { //trk_col->addElement(Track); //fucd - //trk_col->push_back(track); + trk_col->push_back(track); for(int i=0;i<track.trackStates_size();i++){ // 1 = lcio::EVENT::TrackState::AtIP edm4hep::TrackState trkStateIP = track.getTrackStates(i); diff --git a/Reconstruction/SiliconTracking/src/SpacePointBuilderAlg.cpp b/Reconstruction/SiliconTracking/src/SpacePointBuilderAlg.cpp index e907b93acee81c58fc62febdb471e62107c5b13a..956b425122cb65927647ffa885bb3e3524b7e96e 100644 --- a/Reconstruction/SiliconTracking/src/SpacePointBuilderAlg.cpp +++ b/Reconstruction/SiliconTracking/src/SpacePointBuilderAlg.cpp @@ -59,10 +59,10 @@ StatusCode SpacePointBuilderAlg::initialize() { return StatusCode::FAILURE; } - MarlinTrk::IMarlinTrkSystem* _trksystem = _trackSystemSvc->getTrackSystem(); + MarlinTrk::IMarlinTrkSystem* _trksystem = _trackSystemSvc->getTrackSystem(this); _trksystem->init(); - _trackSystemSvc->removeTrackSystem(); + _trackSystemSvc->removeTrackSystem(this); return GaudiAlgorithm::initialize(); } @@ -204,7 +204,8 @@ StatusCode SpacePointBuilderAlg::execute(){ spacePoint.setType( UTIL::set_bit( spacePoint.getType() , UTIL::ILDTrkHitTypeBit::COMPOSITE_SPACEPOINT ) ) ; - spCol->push_back( spacePoint ); + spCol->push_back( spacePoint ); + debug() << "Hit accepted id = " << spacePoint.id() << " cellID = " << spacePoint.getCellID() << endmsg; //debug() << "push_back space point's id=" << spCol->at(spCol->size()-1).id() << endmsg; createdSpacePoints++; @@ -433,12 +434,12 @@ edm4hep::TrackerHit SpacePointBuilderAlg::createSpacePoint( edm4hep::ConstTracke cov_plane(1,1) = (0.5 * du2) / cos2_alpha; cov_plane(2,2) = (0.5 * du2) / sin2_alpha; - debug() << "\t cov_plane = " << cov_plane << endmsg; - debug() << "\tstrip_angle = " << VA.angle(VB)/(M_PI/180) / 2.0 << " degrees " << endmsg; + debug() << "cov_plane = " << cov_plane << endmsg; + debug() << "strip_angle = " << VA.angle(VB)/(M_PI/180) / 2.0 << " degrees " << endmsg; CLHEP::HepSymMatrix cov_xyz= cov_plane.similarity(rot_sensor_matrix); - debug() << "\t cov_xyz = " << cov_xyz << endmsg; + debug() << "cov_xyz = " << cov_xyz << endmsg; std::array<float, 6> cov; //EVENT::FloatVec cov( 9 ) ; @@ -453,9 +454,7 @@ edm4hep::TrackerHit SpacePointBuilderAlg::createSpacePoint( edm4hep::ConstTracke } spacePoint.setCovMatrix(cov); - - debug() << "\tHit accepted id=" << spacePoint.id() << endmsg; - + return spacePoint; } /* diff --git a/Reconstruction/SiliconTracking/src/TrackSubsetAlg.cpp b/Reconstruction/SiliconTracking/src/TrackSubsetAlg.cpp index 2a241cf569d90b1a80953f53b959014d7e937bec..2efe0f52c1fae8cbb0f154b88844469b568d4540 100644 --- a/Reconstruction/SiliconTracking/src/TrackSubsetAlg.cpp +++ b/Reconstruction/SiliconTracking/src/TrackSubsetAlg.cpp @@ -63,7 +63,7 @@ StatusCode TrackSubsetAlg::initialize() { error() << "Failed to find TrackSystemSvc ..." << endmsg; return StatusCode::FAILURE; } - _trkSystem = _trackSystemSvc->getTrackSystem(); + _trkSystem = _trackSystemSvc->getTrackSystem(this); if( _trkSystem == 0 ){ error() << "Cannot initialize MarlinTrkSystem of Type: KalTest" <<endmsg; @@ -76,7 +76,7 @@ StatusCode TrackSubsetAlg::initialize() { _trkSystem->setOption( MarlinTrk::IMarlinTrkSystem::CFG::useSmoothing, _SmoothOn) ; //smoothing // initialise the tracking system - //_trkSystem->init() ; + _trkSystem->init() ; return GaudiAlgorithm::initialize(); } diff --git a/Service/TrackSystemSvc/CMakeLists.txt b/Service/TrackSystemSvc/CMakeLists.txt index c44d75470fe19657c4ac72eeb622bc8f08401793..d12bfb8758945e3b5b43083f23de9245b6af2afe 100644 --- a/Service/TrackSystemSvc/CMakeLists.txt +++ b/Service/TrackSystemSvc/CMakeLists.txt @@ -1,5 +1,6 @@ gaudi_subdir(TrackSystemSvc v0r0) +find_package(CLHEP REQUIRED;CONFIG) find_package(ROOT 6.14 REQUIRED COMPONENTS Matrix Physics) find_package(GEAR REQUIRED) find_package(LCIO REQUIRED) @@ -20,8 +21,8 @@ gaudi_install_headers(TrackSystemSvc) gaudi_add_library(TrackSystemSvcLib ${TrackSystemSvcLib_srcs} PUBLIC_HEADERS TrackSystemSvc - INCLUDE_DIRS GaudiKernel ROOT CLHEP gear ${LCIO_INCLUDE_DIRS} ${EDM4HEP_INCLUDE_DIRS} - LINK_LIBRARIES DataHelperLib KalTestLib KalDetLib GaudiKernel ROOT CLHEP ${GEAR_LIBRARIES} ${LCIO_LIBRARIES} + INCLUDE_DIRS GaudiKernel ROOT ${CLHEP_INCLUDE_DIR} gear ${LCIO_INCLUDE_DIRS} ${EDM4HEP_INCLUDE_DIRS} + LINK_LIBRARIES DataHelperLib KalTestLib KalDetLib GaudiKernel ROOT ${CLHEP_LIBRARIES} ${GEAR_LIBRARIES} ${LCIO_LIBRARIES} -Wl,--no-as-needed EDM4HEP::edm4hep EDM4HEP::edm4hepDict -Wl,--as-needed diff --git a/Service/TrackSystemSvc/TrackSystemSvc/ITrackSystemSvc.h b/Service/TrackSystemSvc/TrackSystemSvc/ITrackSystemSvc.h index 699d26bab97dae97c77f5bf5969f5bd68631e835..88a2d31dd7a85403efdfd81c790e7f93bc876b27 100644 --- a/Service/TrackSystemSvc/TrackSystemSvc/ITrackSystemSvc.h +++ b/Service/TrackSystemSvc/TrackSystemSvc/ITrackSystemSvc.h @@ -13,9 +13,9 @@ class ITrackSystemSvc: virtual public IService { virtual ~ITrackSystemSvc() = default; //Get the track manager - virtual MarlinTrk::IMarlinTrkSystem* getTrackSystem() = 0; + virtual MarlinTrk::IMarlinTrkSystem* getTrackSystem(void* address=0) = 0; - virtual void removeTrackSystem() = 0; + virtual void removeTrackSystem(void* address=0) = 0; }; #endif diff --git a/Service/TrackSystemSvc/src/TrackSystemSvc.cpp b/Service/TrackSystemSvc/src/TrackSystemSvc.cpp index bdcef20f3b120bc0aad489f693e84dcee6948965..238c12b8a34e758fe31cccd681be2f91939c76fc 100644 --- a/Service/TrackSystemSvc/src/TrackSystemSvc.cpp +++ b/Service/TrackSystemSvc/src/TrackSystemSvc.cpp @@ -8,15 +8,15 @@ DECLARE_COMPONENT(TrackSystemSvc) TrackSystemSvc::TrackSystemSvc(const std::string& name, ISvcLocator* svc) - : base_class(name, svc), - m_trackSystem(nullptr){ + : base_class(name, svc){ } TrackSystemSvc::~TrackSystemSvc(){ } -MarlinTrk::IMarlinTrkSystem* TrackSystemSvc::getTrackSystem(){ - if(!m_trackSystem){ +MarlinTrk::IMarlinTrkSystem* TrackSystemSvc::getTrackSystem(void* address){ + std::map<void*, MarlinTrk::IMarlinTrkSystem*>::iterator it=m_trackSystems.find(address); + if(it==m_trackSystems.end()){ gear::GearMgr* mgr=0; auto _gear = service<IGearSvc>("GearSvc"); if ( !_gear ) { @@ -34,45 +34,40 @@ MarlinTrk::IMarlinTrkSystem* TrackSystemSvc::getTrackSystem(){ fatal() << "Both GearSvc and GeomSvc invalid!" << endmsg; return 0; } - - m_trackSystem = new MarlinTrk::MarlinKalTest( *mgr, _geoSvc ) ; + debug() << "GearMgr=" << mgr << " GeomSvc=" << _geoSvc << endmsg; + MarlinTrk::IMarlinTrkSystem* sys = new MarlinTrk::MarlinKalTest( *mgr, _geoSvc ); + m_trackSystems[address] = sys; + debug() << "Track system created successfully for " << address << endmsg; + return sys; } - return m_trackSystem; + return it->second; } StatusCode TrackSystemSvc::initialize(){ - - auto _gear = service<IGearSvc>("GearSvc"); - if ( !_gear ) { - error() << "Failed to find GearSvc ..." << endmsg; - return StatusCode::FAILURE; + for(std::map<void*, MarlinTrk::IMarlinTrkSystem*>::iterator it=m_trackSystems.begin();it!=m_trackSystems.end();it++){ + delete it->second; } - gear::GearMgr* mgr = _gear->getGearMgr(); + m_trackSystems.clear(); + + m_trackSystems[0] = getTrackSystem(0); - auto _geoSvc = service<IGeomSvc>("GeomSvc"); - if ( !_geoSvc ) { - error() << "Failed to find GeomSvc ..." << endmsg; - return StatusCode::FAILURE; - } - m_trackSystem = new MarlinTrk::MarlinKalTest( *mgr, _geoSvc ) ; return StatusCode::SUCCESS; } -void TrackSystemSvc::removeTrackSystem(){ - if ( m_trackSystem ) { - delete m_trackSystem; - m_trackSystem = nullptr; +void TrackSystemSvc::removeTrackSystem(void* address){ + std::map<void*, MarlinTrk::IMarlinTrkSystem*>::iterator it=m_trackSystems.find(address); + if ( it!=m_trackSystems.end() ) { + delete it->second; + m_trackSystems.erase(it); } return; } StatusCode TrackSystemSvc::finalize(){ - - // if ( m_trackSystem ) { - // delete m_trackSystem; - // m_trackSystem = nullptr; - // } + for(std::map<void*, MarlinTrk::IMarlinTrkSystem*>::iterator it=m_trackSystems.begin();it!=m_trackSystems.end();it++){ + delete it->second; + } return StatusCode::SUCCESS; } diff --git a/Service/TrackSystemSvc/src/TrackSystemSvc.h b/Service/TrackSystemSvc/src/TrackSystemSvc.h index 45830a4f966682582a380093f22b4fd24a6d0854..b03dbd9b26d74f16c954206867ab976a26e39398 100644 --- a/Service/TrackSystemSvc/src/TrackSystemSvc.h +++ b/Service/TrackSystemSvc/src/TrackSystemSvc.h @@ -9,14 +9,14 @@ class TrackSystemSvc : public extends<Service, ITrackSystemSvc>{ TrackSystemSvc(const std::string& name, ISvcLocator* svc); ~TrackSystemSvc(); - MarlinTrk::IMarlinTrkSystem* getTrackSystem() override; - void removeTrackSystem() override; + MarlinTrk::IMarlinTrkSystem* getTrackSystem(void* address=0) override; + void removeTrackSystem(void* address=0) override; StatusCode initialize() override; StatusCode finalize() override; private: - MarlinTrk::IMarlinTrkSystem* m_trackSystem; + std::map<void*, MarlinTrk::IMarlinTrkSystem*> m_trackSystems; }; #endif diff --git a/Utilities/DataHelper/src/TrackExtended.cc b/Utilities/DataHelper/src/TrackExtended.cc index 2c2d1412a8d23968af301373655048ce25d25a3d..762247e6e3ca7a3b115fc8bacb69ced2fb6133c1 100644 --- a/Utilities/DataHelper/src/TrackExtended.cc +++ b/Utilities/DataHelper/src/TrackExtended.cc @@ -33,7 +33,11 @@ TrackExtended::TrackExtended( TrackerHitExtended * trackerhit) { TrackExtended::~TrackExtended() {} ConstTrack TrackExtended::getTrack() { - return _track; + if(!_track.isAvailable()){ + std::cout << "Error: track not available" << _track.isAvailable() << " id= " << _track.id() << std::endl; + throw std::runtime_error("Error: track not available"); + } + return _track; } const float * TrackExtended::getSeedPosition() { diff --git a/Utilities/KalDet/CMakeLists.txt b/Utilities/KalDet/CMakeLists.txt index d95a5cc27cffd698cea0cb37100e81e99e086e1b..f7843bd2b474f63625b0316b3f11e7587370f590 100644 --- a/Utilities/KalDet/CMakeLists.txt +++ b/Utilities/KalDet/CMakeLists.txt @@ -5,6 +5,7 @@ gaudi_subdir(KalDet v0r0) +find_package(CLHEP REQUIRED;CONFIG) find_package(LCIO) find_package(GEAR) find_package(ROOT COMPONENTS MathCore) @@ -76,5 +77,5 @@ set( KalDetLib_srcs ${LIB_SOURCES} ${COMMON_SOURCES} ) gaudi_add_library(KalDetLib ${KalDetLib_srcs} PUBLIC_HEADERS kaldet - LINK_LIBRARIES GaudiKernel ROOT CLHEP LCIO ${GEAR_LIBRARIES} KalTestLib EDM4HEP::edm4hep EDM4HEP::edm4hepDict ${DD4hep_COMPONENT_LIBRARIES} + LINK_LIBRARIES GaudiKernel ROOT ${CLHEP_LIBRARIES} LCIO ${GEAR_LIBRARIES} KalTestLib EDM4HEP::edm4hep EDM4HEP::edm4hepDict ${DD4hep_COMPONENT_LIBRARIES} ) diff --git a/Utilities/KalDet/src/ild/common/MaterialDataBase.cc b/Utilities/KalDet/src/ild/common/MaterialDataBase.cc index 5d10f39b9c5f90f2afcc8e021b8eda403bc6765e..a6b8f6e63488dd7026ad2a2f958ecf8292374966 100644 --- a/Utilities/KalDet/src/ild/common/MaterialDataBase.cc +++ b/Utilities/KalDet/src/ild/common/MaterialDataBase.cc @@ -204,36 +204,25 @@ void MaterialDataBase::createMaterials(const gear::GearMgr& gearMgr, IGeomSvc* g // VXD Support Material - - try{ - - //const gear::SimpleMaterial& vxd_sup_mat = gearMgr.getSimpleMaterial("VXDSupportMaterial"); - /* - const gear::GearParametersImpl* vxd_sup_mat = geoSvc->getDetParameters("VXDSupportMaterial"); - - A = vxd_sup_mat->getDoubleVal("A"); - Z = vxd_sup_mat->getDoubleVal("Z"); - density = vxd_sup_mat->getDoubleVal("Density"); - radlen = vxd_sup_mat->getDoubleVal("RadL"); - name = vxd_sup_mat->getStringVal("Name"); - */ - //A = vxd_sup_mat.getA(); - //Z = vxd_sup_mat.getZ(); - //density = vxd_sup_mat.getDensity() * (1000.0/ 1000000.0); // kg/m^3 -> g/cm^3 - //radlen = vxd_sup_mat.getRadLength() / 10.0 ; // mm -> cm - //name = vxd_sup_mat.getName() ; - //std::cout << "debug fucd: " << "==================" << geoSvc << std::endl; - //TMaterial &vxdsupport = *new TMaterial(name.c_str(), "", A, Z, density, radlen, 0.); - //this->addMaterial(&vxdsupport, name); + if(geoSvc){ TMaterial* vxdsupport = geoSvc->getMaterial("VXDSupportMaterial"); - //std::cout << "debug fucd: " << "==================" << std::endl; if(vxdsupport) this->addMaterial(vxdsupport, "VXDSupportMaterial"); - else std::cout << "Material VXDSupportMaterial not found" << std::endl; + else std::cout << "Material VXDSupportMaterial not found" << std::endl; } - catch( gear::UnknownParameterException& e){ - std::cout << "Error while read material from GeomSvc!" << std::endl; + else{ + try{ + const gear::SimpleMaterial& vxd_sup_mat = gearMgr.getSimpleMaterial("VXDSupportMaterial"); + A = vxd_sup_mat.getA(); + Z = vxd_sup_mat.getZ(); + density = vxd_sup_mat.getDensity() * (1000.0/ 1000000.0); // kg/m^3 -> g/cm^3 + radlen = vxd_sup_mat.getRadLength() / 10.0 ; // mm -> cm + name = vxd_sup_mat.getName() ; + TMaterial &vxdsupport = *new TMaterial(name.c_str(), "", A, Z, density, radlen, 0.); + this->addMaterial(&vxdsupport, name); + } + catch( gear::UnknownParameterException& e){ + std::cout << "Error while read material from GeomSvc!" << std::endl; + } } - - } diff --git a/Utilities/KalDet/src/ild/ftd/ILDFTDKalDetector.cc b/Utilities/KalDet/src/ild/ftd/ILDFTDKalDetector.cc index 18fed5e0d3109fb224ff48054d0c5c147983fe7a..b87a9d0f8015cbd07434c15ecc5d12e44ca97900 100644 --- a/Utilities/KalDet/src/ild/ftd/ILDFTDKalDetector.cc +++ b/Utilities/KalDet/src/ild/ftd/ILDFTDKalDetector.cc @@ -8,6 +8,7 @@ #include "DetInterface/IGeomSvc.h" #include "DD4hep/Detector.h" #include "DDRec/DetectorData.h" +#include "DD4hep/DD4hepUnits.h" #include "gear/GEAR.h" #include "gear/BField.h" @@ -29,7 +30,6 @@ ILDFTDKalDetector::ILDFTDKalDetector( 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 { - // streamlog_out(DEBUG1) << "ILDFTDKalDetector building FTD detector using GEAR " << std::endl ; MaterialDataBase::Instance().registerForService(gearMgr, geoSvc); @@ -41,15 +41,11 @@ TVKalDetector(300), _nDisks(0) // SJA:FIXME initial size, 300 looks reasonable f } this->build_staggered_design(); - SetOwner(); } - - - void ILDFTDKalDetector::build_staggered_design() { // streamlog_out(DEBUG) << "ILDFTDKalDetector::build_staggered_design " << std::endl; @@ -532,7 +528,7 @@ void ILDFTDKalDetector::setupGearGeom( IGeomSvc* geoSvc ){ } const dd4hep::Direction& field = geoSvc->lcdd()->field().magneticField(dd4hep::Position(0,0,0)); - _bZ = field.z(); + _bZ = field.z()/dd4hep::tesla; double strip_angle_deg = ftdData->angleStrip/CLHEP::degree; bool strip_angle_present = true; diff --git a/Utilities/KalDet/src/ild/set/ILDSETKalDetector.cc b/Utilities/KalDet/src/ild/set/ILDSETKalDetector.cc index e14aa74e6a922f9697756a57501528df81c6b6c3..3da10d5c34815b8cfe9df023c6078bdbd99df3b4 100644 --- a/Utilities/KalDet/src/ild/set/ILDSETKalDetector.cc +++ b/Utilities/KalDet/src/ild/set/ILDSETKalDetector.cc @@ -20,6 +20,7 @@ #include "gearimpl/Util.h" #include "CLHEP/Units/SystemOfUnits.h" +#include "DD4hep/DD4hepUnits.h" #include "TMath.h" #include "math.h" @@ -46,7 +47,7 @@ ILDSETKalDetector::ILDSETKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* ge else{ setupGearGeom( gearMgr ); } - + if (_isStripDetector) { // streamlog_out(DEBUG4) << "\t\t building SET detector as STRIP Detector." << std::endl ; } else { @@ -311,7 +312,7 @@ void ILDSETKalDetector::setupGearGeom( IGeomSvc* geoSvc ){ } const dd4hep::Direction& field = geoSvc->lcdd()->field().magneticField(dd4hep::Position(0,0,0)); - _bZ = field.z(); + _bZ = field.z()/dd4hep::tesla; std::vector<dd4hep::rec::ZPlanarData::LayerLayout>& setlayers = setData->layers; _nLayers = setlayers.size(); diff --git a/Utilities/KalDet/src/ild/sit/ILDSITKalDetector.cc b/Utilities/KalDet/src/ild/sit/ILDSITKalDetector.cc index ded33d842629e4ccfd258250f2f13d0b064f436d..5d10c55f505d6c9186baa4cb526a931e27dc733a 100644 --- a/Utilities/KalDet/src/ild/sit/ILDSITKalDetector.cc +++ b/Utilities/KalDet/src/ild/sit/ILDSITKalDetector.cc @@ -13,6 +13,7 @@ #include "DD4hep/Detector.h" #include "DDRec/DetectorData.h" #include "CLHEP/Units/SystemOfUnits.h" +#include "DD4hep/DD4hepUnits.h" #include <gear/GEAR.h> #include "gear/BField.h" @@ -30,28 +31,26 @@ ILDSITKalDetector::ILDSITKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* geoSvc ) : TVKalDetector(300) // SJA:FIXME initial size, 300 looks reasonable for ILD, though this would be better stored as a const somewhere { + // std::cout << "ILDSITKalDetector building SIT detector using GEAR " << std::endl ; - - // streamlog_out(DEBUG4) << "ILDSITKalDetector building SIT detector using GEAR " << std::endl ; - MaterialDataBase::Instance().registerForService(gearMgr, geoSvc); TMaterial & air = *MaterialDataBase::Instance().getMaterial("air"); TMaterial & silicon = *MaterialDataBase::Instance().getMaterial("silicon"); TMaterial & carbon = *MaterialDataBase::Instance().getMaterial("carbon"); + if(geoSvc){ this->setupGearGeom(geoSvc); } else{ this->setupGearGeom(gearMgr) ; } - + if (_isStripDetector) { // streamlog_out(DEBUG4) << "\t\t building SIT detector as STRIP Detector." << std::endl ; } else { // streamlog_out(DEBUG4) << "\t\t building SIT detector as PIXEL Detector." << std::endl ; } - //--The Ladder structure (realistic ladder)-- int nLadders; @@ -141,9 +140,7 @@ ILDSITKalDetector::ILDSITKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* ge } - // streamlog_out(DEBUG0) << "ILDSITKalDetector add surface with CellID = " - // << CellID - // << std::endl ; + //std::cout << "ILDSITKalDetector add surface with CellID = " << CellID << std::endl ; } @@ -187,10 +184,8 @@ ILDSITKalDetector::ILDSITKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* ge Add(new ILDParallelPlanarMeasLayer(silicon, silicon, sensitive_distance+sensitive_thickness*0.5, currPhi, _bZ, measurement_plane_sorting_policy, width, sensor_length, offset, z_centre_sensor, offset, true, CellID, "SITMeaslayer" )) ; } - - // streamlog_out(DEBUG0) << "ILDSITKalDetector add surface with CellID = " - // << CellID - // << std::endl ; + + //std::cout << "ILDSITKalDetector add surface with CellID = " << CellID << std::endl ; } @@ -312,7 +307,7 @@ void ILDSITKalDetector::setupGearGeom( IGeomSvc* geoSvc ){ } const dd4hep::Direction& field = geoSvc->lcdd()->field().magneticField(dd4hep::Position(0,0,0)); - _bZ = field.z(); + _bZ = field.z()/dd4hep::tesla; std::vector<dd4hep::rec::ZPlanarData::LayerLayout>& sitlayers = sitData->layers; _nLayers = sitlayers.size(); diff --git a/Utilities/KalDet/src/ild/support/ILDSupportKalDetector.cc b/Utilities/KalDet/src/ild/support/ILDSupportKalDetector.cc index 8bf0126bc3108e68b95efd935da8c16b3c45adca..fe92415eb745cdfcbf883be2fd2c9ac3dc878bdb 100644 --- a/Utilities/KalDet/src/ild/support/ILDSupportKalDetector.cc +++ b/Utilities/KalDet/src/ild/support/ILDSupportKalDetector.cc @@ -19,6 +19,7 @@ #include "DD4hep/Detector.h" #include "DDRec/DetectorData.h" #include "CLHEP/Units/SystemOfUnits.h" +#include "DD4hep/DD4hepUnits.h" #include "gear/GEAR.h" #include "gear/BField.h" @@ -37,11 +38,11 @@ TVKalDetector(10) const dd4hep::rec::ConicalSupportData* pBeamPipeData = geoSvc->getBeamPipeData(); const std::vector<dd4hep::rec::ConicalSupportData::Section>& sections = pBeamPipeData->sections; const dd4hep::Direction& field = geoSvc->lcdd()->field().magneticField(dd4hep::Position(0,0,0)); - bz = field.z(); + bz = field.z()/dd4hep::tesla; for(int i=0;i<sections.size();i++){ - z.push_back(sections[i].zPos); - rInner.push_back(sections[i].rInner); - rOuter.push_back(sections[i].rOuter); + z.push_back(sections[i].zPos*CLHEP::cm ); + rInner.push_back(sections[i].rInner*CLHEP::cm ); + rOuter.push_back(sections[i].rOuter*CLHEP::cm ); } } else{ @@ -54,8 +55,8 @@ TVKalDetector(10) rInner = pBeamPipe.getDoubleVals("RInner"); //inner radius of the cone rOuter = pBeamPipe.getDoubleVals("ROuter"); //outer radius of the cone } - //for(int i=0;i<sections.size();i++){ - // std::cout << z[i] << " " << rInner[i] << " " << rOuter[i] << std::endl; + //for(int i=0;i<z.size();i++){ + //std::cout << z[i] << " " << rInner[i] << " " << rOuter[i] << std::endl; //} MaterialDataBase::Instance().registerForService(gearMgr, geoSvc); TMaterial & beam = *MaterialDataBase::Instance().getMaterial("beam"); diff --git a/Utilities/KalDet/src/ild/tpc/ILDTPCKalDetector.cc b/Utilities/KalDet/src/ild/tpc/ILDTPCKalDetector.cc index 5f15894c7a838bad6f1a116296b66bc17ac7fa93..61af3eedbd6421b78c3d9189de7daea2d203160a 100644 --- a/Utilities/KalDet/src/ild/tpc/ILDTPCKalDetector.cc +++ b/Utilities/KalDet/src/ild/tpc/ILDTPCKalDetector.cc @@ -14,6 +14,7 @@ #include "DD4hep/Detector.h" #include "DDRec/DetectorData.h" #include "CLHEP/Units/SystemOfUnits.h" +#include "DD4hep/DD4hepUnits.h" #include "gear/GEAR.h" #include "gear/BField.h" @@ -52,7 +53,7 @@ TVKalDetector(250) // SJA:FIXME initial size, 250 looks reasonable for ILD, thou } const dd4hep::Direction& field = geoSvc->lcdd()->field().magneticField(dd4hep::Position(0,0,0)); - bz = field.z(); + bz = field.z()/dd4hep::tesla; nlayers = tpcData->maxRow; lhalf = tpcData->driftLength*CLHEP::cm; rstep = tpcData->padHeight*CLHEP::cm; diff --git a/Utilities/KalDet/src/ild/vxd/ILDVXDKalDetector.cc b/Utilities/KalDet/src/ild/vxd/ILDVXDKalDetector.cc index 89bc88d716efd820dec047ca13f5bd33f7e6646a..5cacd3d76205a515295102ec6c6f865bd24776a4 100644 --- a/Utilities/KalDet/src/ild/vxd/ILDVXDKalDetector.cc +++ b/Utilities/KalDet/src/ild/vxd/ILDVXDKalDetector.cc @@ -14,6 +14,7 @@ #include "DD4hep/Detector.h" #include "DDRec/DetectorData.h" #include "CLHEP/Units/SystemOfUnits.h" +#include "DD4hep/DD4hepUnits.h" #include <gear/GEAR.h> #include "gear/BField.h" @@ -48,13 +49,14 @@ ILDVXDKalDetector::ILDVXDKalDetector( const gear::GearMgr& gearMgr, IGeomSvc* ge TMaterial & aluminium = *MaterialDataBase::Instance().getMaterial("aluminium"); _vxd_Cryostat.exists = false; - + if(geoSvc){ this->setupGearGeom(geoSvc) ; } else{ this->setupGearGeom(gearMgr) ; } + //--The Ladder structure (realistic ladder)-- int nLadders; @@ -398,7 +400,7 @@ void ILDVXDKalDetector::setupGearGeom( IGeomSvc* geoSvc){ } */ const dd4hep::Direction& field = geoSvc->lcdd()->field().magneticField(dd4hep::Position(0,0,0)); - _bZ = field.z(); + _bZ = field.z()/dd4hep::tesla; const gear::ZPlanarParametersImpl* pVXDDetMain = geoSvc->getVXDParameters(); const gear::VXDLayerLayout& pVXDLayerLayout = pVXDDetMain->getVXDLayerLayout(); diff --git a/Utilities/KiTrack/CMakeLists.txt b/Utilities/KiTrack/CMakeLists.txt index 22d2f25da413adfae3d0013bd3d5136046f79ce5..7a20e00dc7f34a62d6a97bad8e0c4886d6db62a9 100644 --- a/Utilities/KiTrack/CMakeLists.txt +++ b/Utilities/KiTrack/CMakeLists.txt @@ -1,5 +1,6 @@ gaudi_subdir(KiTrack v0r0) +find_package(CLHEP REQUIRED;CONFIG) find_package(ROOT REQUIRED) #find_package(DD4hep REQUIRED) find_package(GSL REQUIRED) @@ -16,6 +17,6 @@ include_directories(src) gaudi_add_library(KiTrackLib ${KiTrackLib_srcs} PUBLIC_HEADERS KiTrack - LINK_LIBRARIES DataHelperLib TrackSystemSvcLib ROOT CLHEP GSL EDM4HEP::edm4hep LCIO + LINK_LIBRARIES DataHelperLib TrackSystemSvcLib ROOT ${CLHEP_LIBRARIES} GSL EDM4HEP::edm4hep LCIO # DD4hep )