sPHENIX MAPS tracker discussion

[Yasser Corrales Morales <ycmorales AT rcf.rhic.bnl.gov>]

Hi Joe,

could you test this patch?

thank you,

Yasser.

On 8/31/20 12:40, Anthony Frawley wrote:

        I am looking at it. Give me a bit.
      
        
      
      
        Tony

---

From: sPHENIX-MAPS-l <sphenix-maps-l-bounces AT lists.bnl.gov> on behalf of Osborn, Joe <osbornjd AT ornl.gov>
Sent: Monday, August 31, 2020 1:55 PM
To: Cameron Thomas Dean <cameron.dean AT cern.ch>; sPHENIX-tracking-l <sphenix-tracking-l-bounces AT lists.bnl.gov>; sphenix-maps-l AT lists.bnl.gov <sphenix-maps-l AT lists.bnl.gov>
Subject: Re: [Sphenix-maps-l] [EXTERNAL] Re: MVTX TGeoNode Name Change?

 

Hi Cameron,

 

Thanks for the email – it looks like that doesn’t fix it. I’ve been looking through the G4_Mvtx history and I can’t seem to find anything obvious that would be causing this issue. I’m wondering if there was a calibration file or something similar that was responsible for the Mvtx geometry that got changed? I think Tony originally wrote the code that is now causing Acts to crash with the most up to date macros, so maybe he knows where “av_1” actually came from and this can help us track down if/where it got changed.

 

---------------------------   Joe Osborn, Ph.D. Postdoctoral Research Associate Oak Ridge National Laboratory osbornjd AT ornl.gov (859)-433-8738

 

 

From: Cameron Thomas Dean <cameron.dean AT cern.ch>
Date: Monday, August 31, 2020 at 1:44 PM
To: "Osborn, Joe" <osbornjd AT ornl.gov>, sPHENIX-tracking-l <sphenix-tracking-l-bounces AT lists.bnl.gov>, "sphenix-maps-l AT lists.bnl.gov" <sphenix-maps-l AT lists.bnl.gov>
Subject: [EXTERNAL] Re: [Sphenix-maps-l] MVTX TGeoNode Name Change?

 

Hi Joe,

Sorry about the issues here. To try and trace this down, would it be possible to run Acts with the following set to false?
https://github.com/sPHENIX-Collaboration/macros/blob/master/common/G4_Mvtx.C#L31

Line 31: bool MVTX_SERVICE = false;

I don’t recall any changes to objects with “av_1” in their name but if the problem disappears with this change then the issue is related to the service barrel. I’ll be able to debug further if this is the case.

 

Cheers,
Cameron

From: sPHENIX-MAPS-l <sphenix-maps-l-bounces AT lists.bnl.gov> on behalf of "Osborn, Joe" <osbornjd AT ornl.gov>
Date: Monday, August 31, 2020 at 8:53 AM
To: sPHENIX-tracking-l <sphenix-tracking-l-bounces AT lists.bnl.gov>, "sphenix-maps-l AT lists.bnl.gov" <sphenix-maps-l AT lists.bnl.gov>
Subject: [Sphenix-maps-l] MVTX TGeoNode Name Change?

 

Hi MVTX experts,

 

I updated my macros repository this morning and it seems that the Acts tracking no longer works because of a change in the MVTX geometry. Previously we were accessing the TGeoNodes in the geometry by iterating over the TGeoVolume and accessing the MVTX objects by getting the nodes that had “av_1” in their name. This now returns no valid MVTX nodes, which means I’m guessing that the node names were changed. Perhaps this is related to the recent addition of the MVTX services into the geometry?

 

Does anyone have a pointer as to what happened here? Just for clarity, the way we were accessing the MVTX nodes previously can be found here in lines 807-826.

 

Thanks for any help you may have,

 

Joe Osborn

 

          ---------------------------
           
          Joe
              Osborn, Ph.D.
          Postdoctoral
              Research Associate
          Oak Ridge
              National Laboratory
          osbornjd AT ornl.gov
          (859)-433-8738

 

_______________________________________________
sPHENIX-MAPS-l mailing list
sPHENIX-MAPS-l AT lists.bnl.gov
https://lists.bnl.gov/mailman/listinfo/sphenix-maps-l
    
-- 
-------
Yasser Corrales Morales

Postdoctoral Research Associate
PO Box 1663
MS H846
LANL
Los Alamos, NM 87544
  

diff --git a/offline/packages/trackreco/MakeActsGeometry.cc 
b/offline/packages/trackreco/MakeActsGeometry.cc
index e0094056..37d22e1d 100644
--- a/offline/packages/trackreco/MakeActsGeometry.cc
+++ b/offline/packages/trackreco/MakeActsGeometry.cc
@@ -128,10 +128,10 @@ int MakeActsGeometry::buildAllGeometry(PHCompositeNode 
*topNode)
 
 void MakeActsGeometry::editTPCGeometry(PHCompositeNode *topNode)
 {
-  
+
   PHGeomTGeo *geomNode = PHGeomUtility::GetGeomTGeoNode(topNode, true);
   assert(geomNode);
-  
+
   /// Reset the geometry node, which we will recreate with the TPC edits
   if(geomNode->isValid())
     {
@@ -140,11 +140,11 @@ void MakeActsGeometry::editTPCGeometry(PHCompositeNode 
*topNode)
   PHGeomIOTGeo *dstGeomIO = PHGeomUtility::GetGeomIOTGeoNode(topNode, false);
   assert(dstGeomIO);
   assert(dstGeomIO->isValid());
-  
+
   TGeoManager *geoManager = dstGeomIO->ConstructTGeoManager();
   geomNode->SetGeometry(geoManager);
   assert(geoManager);
-  
+
   TGeoVolume *World_vol = geoManager->GetTopVolume();
   TGeoNode *tpc_envelope_node = nullptr;
   TGeoNode *tpc_gas_north_node = nullptr;
@@ -221,7 +221,7 @@ void MakeActsGeometry::addActsTpcSurfaces(TGeoVolume 
*tpc_gas_vol, TGeoManager *
 
   TGeoVolume *tpc_gas_measurement_vol[m_nTpcLayers];
   double tan_half_phi = tan(m_surfStepPhi / 2.0);
-  
+
   for(unsigned int ilayer = 0; ilayer < m_nTpcLayers; ++ilayer)
     {
       // make a box for this layer
@@ -229,12 +229,12 @@ void MakeActsGeometry::addActsTpcSurfaces(TGeoVolume 
*tpc_gas_vol, TGeoManager *
       sprintf(bname,"tpc_gas_measurement_%u",ilayer);
 
       // Because we use a box, not a section of a cylinder, we need this to 
prevent overlaps
-      // set the nominal r*phi dimension of the box so they just touch at 
the inner edge when placed 
+      // set the nominal r*phi dimension of the box so they just touch at 
the inner edge when placed
       double box_r_phi = 2.0 * tan_half_phi * (m_layerRadius[ilayer] - 
m_layerThickness[ilayer] / 2.0);
 
-      tpc_gas_measurement_vol[ilayer] = geoManager->MakeBox(bname, 
tpc_gas_medium, 
-                                                           
m_layerThickness[ilayer]*half_width_clearance_thick, 
-                                                           
box_r_phi*half_width_clearance_phi, 
+      tpc_gas_measurement_vol[ilayer] = geoManager->MakeBox(bname, 
tpc_gas_medium,
+                                                           
m_layerThickness[ilayer]*half_width_clearance_thick,
+                                                           
box_r_phi*half_width_clearance_phi,
                                                            
m_surfStepZ*half_width_clearance_z);
 
       tpc_gas_measurement_vol[ilayer]->SetLineColor(kBlack);
@@ -243,18 +243,18 @@ void MakeActsGeometry::addActsTpcSurfaces(TGeoVolume 
*tpc_gas_vol, TGeoManager *
 
       if(m_verbosity > 30)
        {
-         cout << m_verbosity << " Made box for layer " << ilayer << " with 
dx " << m_layerThickness[ilayer] << " dy " 
+         cout << m_verbosity << " Made box for layer " << ilayer << " with 
dx " << m_layerThickness[ilayer] << " dy "
               << box_r_phi << " ref arc " << 
m_surfStepPhi*m_layerRadius[ilayer] << " dz " << m_surfStepZ << endl;
          tpc_gas_measurement_vol[ilayer]->Print();
        }
     }
 
-  int copy = 0;              
+  int copy = 0;
   for (unsigned int iz = 0; iz < m_nSurfZ; ++iz)
     {
       // The (half) tpc gas volume is 105.5 cm long and is symmetric around 
(x,y,z) = (0,0,0) in its frame
       double z_center = -105.5/2.0 + m_surfStepZ / 2.0 + (double) iz * 
m_surfStepZ;
-      
+
       for (unsigned int imod = 0; imod < m_nTpcModulesPerLayer; ++imod)
        {
          for (unsigned int iphi = 0; iphi < m_nSurfPhi; ++iphi)
@@ -263,31 +263,31 @@ void MakeActsGeometry::addActsTpcSurfaces(TGeoVolume 
*tpc_gas_vol, TGeoManager *
              double min_phi = m_modulePhiStart + (double) imod * 
m_moduleStepPhi + (double) iphi * m_surfStepPhi;
              double phi_center = min_phi + m_surfStepPhi / 2.0;
              double phi_center_degrees = phi_center * 180.0 / M_PI;
-             
+
              for (unsigned int ilayer = 0; ilayer < m_nTpcLayers; ++ilayer)
                {
                  copy++;
-                 
+
                  // place copies of the gas volume to fill up the layer
-                 
+
                  double x_center = m_layerRadius[ilayer] * cos(phi_center);
                  double y_center = m_layerRadius[ilayer] * sin(phi_center);
-                 
+
                  char rot_name[500];
                  sprintf(rot_name,"tpc_gas_rotation_%i", copy);
                  TGeoCombiTrans *tpc_gas_measurement_location = new 
TGeoCombiTrans(x_center, y_center, z_center,
                                                                              
      new TGeoRotation(rot_name,phi_center_degrees, 0, 0));
-                 
+
                  tpc_gas_vol->AddNode(tpc_gas_measurement_vol[ilayer], copy, 
tpc_gas_measurement_location);
-                 
-                 if(m_verbosity > 30 && ilayer == 30) 
+
+                 if(m_verbosity > 30 && ilayer == 30)
                    {
                      cout << " Made copy " << copy << " iz " << iz << " imod 
" << imod << " ilayer " << ilayer << " iphi " << iphi << endl;
                      cout << "    x_center " << x_center << " y_center " << 
y_center << " z_center " << z_center << " phi_center_degrees " << 
phi_center_degrees << endl;
                    }
                }
            }
-       }      
+       }
     }
 }
 
@@ -302,7 +302,7 @@ void MakeActsGeometry::buildActsSurfaces()
 
   // Response file contains arguments necessary for geometry building
   const std::string argstr[argc]{
-    "-n1", "-l0", 
+    "-n1", "-l0",
       "--response-file=tgeo-sphenix.response",
       "--bf-values", "0", "0", "1.4"
       };
@@ -313,15 +313,15 @@ void MakeActsGeometry::buildActsSurfaces()
     // need a copy, since .c_str() returns a const char * and process 
geometry will not take a const
     arg[i] = strdup(argstr[i].c_str());
   }
-  
-  // We replicate the relevant functionality of  
+
+  // We replicate the relevant functionality of
   //acts/Examples/Run/Common/src/GeometryExampleBase::ProcessGeometry() in 
MakeActsGeometry()
   // so we get access to the results. The layer builder magically gets the 
TGeoManager
 
   makeGeometry(argc, arg, m_detector);
 }
 
-void MakeActsGeometry::makeGeometry(int argc, char *argv[], 
+void MakeActsGeometry::makeGeometry(int argc, char *argv[],
                                    FW::IBaseDetector &detector)
 {
   /// setup and parse options
@@ -373,7 +373,7 @@ void MakeActsGeometry::makeGeometry(int argc, char 
*argv[],
 
 
   /// m_tGeometry is a TrackingGeometry pointer
-  /// vol is a TrackingVolume pointer  
+  /// vol is a TrackingVolume pointer
   auto vol = m_tGeometry->highestTrackingVolume();
 
   if(m_verbosity > 10 )
@@ -390,23 +390,23 @@ void MakeActsGeometry::makeGeometry(int argc, char 
*argv[],
   /// We have several volumes to walk through with the tpc and silicon
   auto firstVolumes = volumeVector.at(0)->confinedVolumes();
   auto topVolumesVector = firstVolumes->arrayObjects();
-  
+
   if(m_verbosity > 10 )
     {
       for(long unsigned int i = 0; i<topVolumesVector.size(); i++)
        {
-         std::cout<< "TopVolume name: " 
+         std::cout<< "TopVolume name: "
                   << topVolumesVector.at(i)->volumeName() << std::endl;
        }
     }
 
   /// This actually contains the silicon volumes
   auto siliconVolumes = topVolumesVector.at(1)->confinedVolumes();
-  
+
   if(m_verbosity > 10 )
     {
       for(long unsigned int i =0; i<siliconVolumes->arrayObjects().size(); 
i++){
-       std::cout << "SiliconVolumeName: " 
+       std::cout << "SiliconVolumeName: "
                  << siliconVolumes->arrayObjects().at(i)->volumeName()
                  << std::endl;
       }
@@ -414,7 +414,7 @@ void MakeActsGeometry::makeGeometry(int argc, char 
*argv[],
 
   /// siliconVolumes is a shared_ptr<TrackingVolumeArray>
   /// Now get the individual TrackingVolumePtrs corresponding to each 
silicon volume
-  
+
   auto mvtxVolumes = siliconVolumes->arrayObjects().at(0);
   auto mvtxConfinedVolumes = mvtxVolumes->confinedVolumes();
   auto mvtxBarrel = mvtxConfinedVolumes->arrayObjects().at(1);
@@ -429,7 +429,7 @@ void MakeActsGeometry::makeGeometry(int argc, char 
*argv[],
 
   /// Same for the TPC - only one volume
   auto tpcVolume = volumeVector.at(1);
-  
+
   makeTpcMapPairs(tpcVolume);
 
   return;
@@ -455,20 +455,20 @@ void 
MakeActsGeometry::makeTpcMapPairs(TrackingVolumePtr &tpcVolume)
        {
          auto surf = surfaceVector.at(j)->getSharedPtr();
          auto vec3d = surf->center(m_geoCtxt);
-         
+
 
          /// convert to cm
-         std::vector<double> world_center = {vec3d(0) / 10.0, 
+         std::vector<double> world_center = {vec3d(0) / 10.0,
                                              vec3d(1) / 10.0,
                                              vec3d(2) / 10.0};
-       
+
          TrkrDefs::hitsetkey hitsetkey = 
getTpcHitSetKeyFromCoords(world_center);
 
          /// If there is already an entry for this hitsetkey, add the surface
          /// to its corresponding vector
          std::map<TrkrDefs::hitsetkey, std::vector<Surface>>::iterator 
mapIter;
          mapIter = m_clusterSurfaceMapTpcEdit.find(hitsetkey);
-         
+
          if(mapIter != m_clusterSurfaceMapTpcEdit.end())
            {
              mapIter->second.push_back(surf);
@@ -478,11 +478,11 @@ void 
MakeActsGeometry::makeTpcMapPairs(TrackingVolumePtr &tpcVolume)
              /// Otherwise make a new map entry
              std::vector<Surface> dumvec;
              dumvec.push_back(surf);
-             std::pair<TrkrDefs::hitsetkey, std::vector<Surface>> tmp = 
+             std::pair<TrkrDefs::hitsetkey, std::vector<Surface>> tmp =
                std::make_pair(hitsetkey, dumvec);
              m_clusterSurfaceMapTpcEdit.insert(tmp);
            }
-         
+
        }
     }
 
@@ -490,7 +490,7 @@ void MakeActsGeometry::makeTpcMapPairs(TrackingVolumePtr 
&tpcVolume)
 
 void MakeActsGeometry::makeInttMapPairs(TrackingVolumePtr &inttVolume)
 {
-  
+
   if(m_verbosity > 10)
     {
       std::cout << "intt volume name: "  << inttVolume->volumeName()
@@ -522,7 +522,7 @@ void MakeActsGeometry::makeInttMapPairs(TrackingVolumePtr 
&inttVolume)
       double ref_rad[4] = {8.987, 9.545, 10.835, 11.361};
 
       std::vector<double> world_center = {vec3d(0) / 10.0, vec3d(1) / 10.0, 
vec3d(2) / 10.0};  // convert from mm to cm
-      /// The Acts geometry builder combines layers 4 and 5 together, 
+      /// The Acts geometry builder combines layers 4 and 5 together,
       /// and layers 6 and 7 together. We need to use the radius to figure
       /// out which layer to use to get the layergeom
       double layer_rad = sqrt(pow(world_center[0], 2) + pow(world_center[1], 
2));
@@ -573,10 +573,10 @@ void 
MakeActsGeometry::makeInttMapPairs(TrackingVolumePtr &inttVolume)
 
 void MakeActsGeometry::makeMvtxMapPairs(TrackingVolumePtr &mvtxVolume)
 {
-  
+
   if(m_verbosity > 10)
     {
-      std::cout << "MVTX Barrel name to step surfaces through is " 
+      std::cout << "MVTX Barrel name to step surfaces through is "
              << mvtxVolume->volumeName() << std::endl;
     }
   // Now get the LayerArrays corresponding to each volume
@@ -654,7 +654,7 @@ Surface 
MakeActsGeometry::getTpcSurfaceFromCoords(TrkrDefs::hitsetkey hitsetkey,
 {
   std::map<TrkrDefs::hitsetkey, std::vector<Surface>>::iterator mapIter;
   mapIter = m_clusterSurfaceMapTpcEdit.find(hitsetkey);
-  
+
   if(mapIter == m_clusterSurfaceMapTpcEdit.end())
     {
       cout << PHWHERE << "Error: hitsetkey not found in clusterSurfaceMap, 
hitsetkey = " << hitsetkey << endl;
@@ -663,13 +663,13 @@ Surface 
MakeActsGeometry::getTpcSurfaceFromCoords(TrkrDefs::hitsetkey hitsetkey,
 
   double world_phi = atan2(world[1], world[0]);
   double world_z = world[2];
-  
+
   std::vector<Surface> surf_vec = mapIter->second;
   unsigned int surf_index = 999;
   for(unsigned int i=0;i<surf_vec.size(); ++i)
     {
       Surface this_surf = surf_vec[i];
-  
+
       auto vec3d = this_surf->center(m_geoCtxt);
       std::vector<double> surf_center = {vec3d(0) / 10.0, vec3d(1) / 10.0, 
vec3d(2) / 10.0};  // convert from mm to cm
       double surf_phi = atan2(surf_center[1], surf_center[0]);
@@ -677,7 +677,7 @@ Surface 
MakeActsGeometry::getTpcSurfaceFromCoords(TrkrDefs::hitsetkey hitsetkey,
       if( (world_phi > surf_phi - m_surfStepPhi / 2.0 && world_phi < 
surf_phi + m_surfStepPhi / 2.0 ) &&
          (world_z > surf_z -m_surfStepZ / 2.0 && world_z < surf_z + 
m_surfStepZ / 2.0) )
        {
-         surf_index = i;         
+         surf_index = i;
          break;
        }
     }
@@ -686,7 +686,7 @@ Surface 
MakeActsGeometry::getTpcSurfaceFromCoords(TrkrDefs::hitsetkey hitsetkey,
       cout << PHWHERE << "Error: TPC surface index not defined, skipping 
cluster!" << endl;
       return nullptr;
     }
- 
+
   return surf_vec[surf_index];
 
 }
@@ -708,21 +708,21 @@ TrkrDefs::hitsetkey 
MakeActsGeometry::getTpcHitSetKeyFromCoords(std::vector<doub
          break;
        }
     }
-  if(layer >= m_nTpcLayers) 
+  if(layer >= m_nTpcLayers)
     {
       cout << PHWHERE << "Error: undefined layer, do nothing world =  " << 
world[0] << "  " << world[1] << "  " << world[2] << " layer " << layer << 
endl;
       return Fun4AllReturnCodes::ABORTEVENT;
     }
   layer += 7;
 
-  // side -  from world z sign 
+  // side -  from world z sign
   unsigned int side;
   if(world[2] < 0)
     side = 0;
   else
     side = 1;
 
-  // readout module 
+  // readout module
   unsigned int readout_mod = 999;
   double phi_world = atan2(world[1], world[0]);
   for(unsigned int imod=0; imod<m_nTpcModulesPerLayer; ++imod)
@@ -745,7 +745,7 @@ TrkrDefs::hitsetkey 
MakeActsGeometry::getTpcHitSetKeyFromCoords(std::vector<doub
   if(m_verbosity > 3)
     if(layer == 30)
       cout << "   world = " << world[0] << "  " << world[1] << "  " << 
world[2] << " phi_world " << phi_world*180/3.14159 << " layer " << layer << " 
readout_mod " << readout_mod << " side " << side << " hitsetkey " << 
hitset_key<< endl;
-  
+
   return hitset_key;
 }
 
@@ -813,16 +813,27 @@ void MakeActsGeometry::makeTGeoNodeMap(PHCompositeNode 
*topNode)
     TGeoNode *node = dynamic_cast<TGeoNode *>(obj);
     std::string node_str = node->GetName();
 
-    std::string mvtx("av_1");
+    std::string mvtx("MVTX_Wrapper");
     std::string intt("ladder");
     std::string intt_ext("ladderext");
     std::string tpc("tpc_envelope");
 
     if (node_str.compare(0, mvtx.length(), mvtx) == 0)  // is it in the MVTX?
     {
-      if (m_verbosity > 2) 
-       cout << " node " << node->GetName() << " is in the MVTX" << endl;
-      getMvtxKeyFromNode(node);
+      std:string mvtx_stave_patt("av_1");
+      TIter iMVTX(node->GetNodes());
+      while (TObject *obj = iter_mvtx())
+      {
+        TGeoNode *mvtx_node = dynamic_cast<TGeoNode *>(obj);
+        std::string mvtx_node_str = node->GetName();
+        if ( (mvtx_node_str.compare(0, mvtx_stave_patt.length(), 
mvtx_stave_patt) == 0) && \
+             (mvtx_node_str.find("MVTXHalfStave") != std::string::npos) )
+        {
+          if (m_verbosity > 2)
+                 cout << " node " << node->GetName() << " is in the MVTX" << 
endl;
+          getMvtxKeyFromNode(mvtx_node);
+        }
+      }
     }
     else if (node_str.compare(0, intt.length(), intt) == 0)  // is it in the 
INTT?
     {
@@ -830,7 +841,7 @@ void MakeActsGeometry::makeTGeoNodeMap(PHCompositeNode 
*topNode)
       if (node_str.compare(0, intt_ext.length(), intt_ext) == 0)
         continue;
 
-      if (m_verbosity > 2) 
+      if (m_verbosity > 2)
        cout << " node " << node->GetName() << " is in the INTT" << endl;
       getInttKeyFromNode(node);
     }
@@ -1018,17 +1029,17 @@ void MakeActsGeometry::isActive(TGeoNode *gnode, int 
nmax_print)
   {
     cout << "          ******* Found INTT active volume,  node is " << 
gnode->GetName()
         << " volume name is " << gnode->GetVolume()->GetName() << endl;
-    
+
     //const TGeoMatrix* tgMatrix = gnode->GetMatrix();
     //tgMatrix->Print();
-    
+
     return;
   }
   else if (node_str.compare(0, mvtx_refactive.length(), mvtx_refactive) == 0)
   {
     cout << "          ******* Found MVTX active volume,  node is " << 
gnode->GetName()
         << " volume name is " << gnode->GetVolume()->GetName() << endl;
- 
+
      return;
   }
   else if (node_str.compare(0, tpc_refactive.length(), tpc_refactive) == 0)
@@ -1047,10 +1058,10 @@ void MakeActsGeometry::isActive(TGeoNode *gnode, int 
nmax_print)
 
   for (int i = 0; i < ndaught; ++i)
   {
-    std::cout << "     " << gnode->GetVolume()->GetName() 
-             << "  daughter " << i << " has name " 
+    std::cout << "     " << gnode->GetVolume()->GetName()
+             << "  daughter " << i << " has name "
              << gnode->GetDaughter(i)->GetVolume()->GetName() << std::endl;
-    
+
     isActive(gnode->GetDaughter(i), nmax_print);
   }
 }