From d88420f571e27aeb2210e3c39cceb9838066258a Mon Sep 17 00:00:00 2001
From: myliu <201916234@mail.sdu.edu.cn>
Date: Wed, 21 Oct 2020 14:59:03 +0800
Subject: [PATCH] Put the initialization of phi0 into the drift chamber
 construction stage

---
 .../src/driftchamber/DriftChamber.cpp         |  6 +-
 .../DetSegmentation/GridDriftChamber.h        | 61 +++++++++++++------
 .../DetSegmentation/src/GridDriftChamber.cpp  | 28 ++++++---
 3 files changed, 63 insertions(+), 32 deletions(-)

diff --git a/Detector/DetDriftChamber/src/driftchamber/DriftChamber.cpp b/Detector/DetDriftChamber/src/driftchamber/DriftChamber.cpp
index fe219aa3..07b95b98 100644
--- a/Detector/DetDriftChamber/src/driftchamber/DriftChamber.cpp
+++ b/Detector/DetDriftChamber/src/driftchamber/DriftChamber.cpp
@@ -93,7 +93,7 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
 
     // - layer
       for(int layer_id = 0; layer_id < (inner_chamber_layer_number+outer_chamber_layer_number); layer_id++) {
-        double rmin,rmax;
+        double rmin,rmax,offset;
         std::string layer_name;
         dd4hep::Volume* current_vol_ptr = nullptr;
 
@@ -117,8 +117,10 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
         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);
+        DCHseg->setGeomParams(layer_id, layer_Phi, rmid, epsilon, offset);
         DCHseg->setWiresInLayer(layer_id, numWire);
 
         dd4hep::Tube layer_solid(rmin,rmax,chamber_length*0.5);
diff --git a/Detector/DetSegmentation/DetSegmentation/GridDriftChamber.h b/Detector/DetSegmentation/DetSegmentation/GridDriftChamber.h
index 64a219d7..fc01bf11 100644
--- a/Detector/DetSegmentation/DetSegmentation/GridDriftChamber.h
+++ b/Detector/DetSegmentation/DetSegmentation/GridDriftChamber.h
@@ -15,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 {
@@ -42,8 +57,9 @@ public:
     return std::atan2(aposition.Y, aposition.X) + M_PI ;
   }
 
-  inline void setGeomParams(int layer, double layerphi, double R, double eps ) {
-     layer_params[layer] = {layerphi,R,eps};
+  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)
@@ -52,20 +68,21 @@ public:
     updateParams(layer);
     for (int i = 0; i<numWires; ++i) {
 
-    if(layer % 2 == 0) { phi0 = 0.; }
-    else { phi0 = 0.5 * _currentLayerphi; }
+//    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;
+      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 Wstart = returnWirePosition(phi_start, 1);
+      TVector3 Wend = returnWirePosition(phi_end, -1);
 
-    TVector3 Wmid = (Wstart+Wend)*(1/2.0);
-    TVector3 Wdirection = (Wend - Wstart);
+      TVector3 Wmid = (Wstart+Wend)*(1/2.0);
+      TVector3 Wdirection = (Wend - Wstart);
 
-    m_wiresPositions[layer].push_back(std::make_pair(Wmid, Wdirection));
+      m_wiresPositions[layer].push_back(std::make_pair(Wmid, Wdirection));
       }
   }
 
@@ -82,20 +99,23 @@ public:
   void updateParams(int layer)  const{
     auto it_end = layer_params.cend();
     --it_end;
-    double layerphi = it_end->second[0];
-    double radius = it_end->second[1];
-    double eps = it_end->second[2];
+    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[0];
-     radius = map_it->second[1];
-     eps = map_it->second[2];
+     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));
@@ -105,7 +125,7 @@ public:
 protected:
   /* *** nalipour *** */
   double phi(const CellID& cID) const;
-  std::map<int,std::vector<double>> layer_params; // <layer, {layerphi, R, eps}>
+  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;
@@ -118,6 +138,7 @@ protected:
   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 3355482f..57a6d7eb 100644
--- a/Detector/DetSegmentation/src/GridDriftChamber.cpp
+++ b/Detector/DetSegmentation/src/GridDriftChamber.cpp
@@ -44,23 +44,31 @@ CellID GridDriftChamber::cellID(const Vector3D& /*localPosition*/, const Vector3
   double phi_hit = phiFromXY(globalPosition);
   double posx = globalPosition.X;
   double posy = globalPosition.Y;
-  double offsetphi;
+  double offsetphi= m_offset;
   int _lphi;
-  if(layerID % 2 == 0) {
-      offsetphi = 0.;
-     _lphi = (int) (phi_hit / _currentLayerphi);
-   }
-  else {
-    offsetphi = _currentLayerphi / 2.;
+//  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);
+      _lphi = (int) ((phi_hit - offsetphi)/ _currentLayerphi);
     }
     else {
-    _lphi = (int) ((phi_hit - offsetphi + 2 * M_PI)/ _currentLayerphi);
+      _lphi = (int) ((phi_hit - offsetphi + 2 * M_PI)/ _currentLayerphi);
     }
-  }
   int lphi = _lphi;
   _decoder->set(cID, m_phiID, lphi);
+
+//std::cout << "#######################################: " 
+//          <<  " offset : " << m_offset
+//          << " offsetphi: " << offsetphi
+//          << " layerID: " << layerID
+//          << " r: " << _currentRadius
+//          << " layerphi: " << _currentLayerphi
+//          << std::endl;
+
   return cID;
 }
 
-- 
GitLab