MakeActsGeometry.h

Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file MakeActsGeometry.h

#ifndef MAKE_ACTS_GEOMETRY_H
#define MAKE_ACTS_GEOMETRY_H

#include <fun4all/SubsysReco.h>
#include <trackbase/TrkrDefs.h>

#include <trackbase/ActsGeometry.h>

#include <Acts/Definitions/Algebra.hpp>
#include <Acts/Utilities/BinnedArray.hpp>                      
#include <Acts/Utilities/Logger.hpp>                           
#include <Acts/EventData/MeasurementHelpers.hpp> 
#include <Acts/Geometry/TrackingGeometry.hpp>
#include <Acts/MagneticField/MagneticFieldContext.hpp>
#include <Acts/Utilities/CalibrationContext.hpp>
#include <Acts/MagneticField/MagneticFieldProvider.hpp>

#include <ActsExamples/Detector/TGeoDetectorWithOptions.hpp>

#include <boost/program_options.hpp>

#include <map>
#include <memory>            
#include <string>
#include <vector>

class PHCompositeNode;
class PHG4CylinderGeomContainer;
class PHG4TpcCylinderGeomContainer;
class TGeoManager;
class TGeoNode;
class TGeoVolume;

namespace Acts 
{
  class Surface;
}

using Surface = std::shared_ptr<const Acts::Surface>;
using TrackingGeometry = std::shared_ptr<const Acts::TrackingGeometry>;
//using TrackingGeometry = std::shared_ptr<Acts::TrackingGeometry>;
using TrackingVolumePtr = std::shared_ptr<const Acts::TrackingVolume>;

class MakeActsGeometry : public SubsysReco
{
 public:

  MakeActsGeometry(const std::string& name = "MakeActsGeometry");

  ~MakeActsGeometry() override = default;

  int Init(PHCompositeNode *topNode) override;
  int InitRun(PHCompositeNode *topNode) override;

  void loadMagField(const bool field) { m_useField = field; }
  void setMagField(const std::string &magField)
    {m_magField = magField;}
  void setMagFieldRescale(double magFieldRescale)
    {m_magFieldRescale = magFieldRescale;}


  void setMvtxDev(double array[6])
  {
    m_mvtxDevs[0] = array[0];
    m_mvtxDevs[1] = array[1];
    m_mvtxDevs[2] = array[2];
    m_mvtxDevs[3] = array[3];
    m_mvtxDevs[4] = array[4];
    m_mvtxDevs[5] = array[5];

    mvtxParam = true;
  }
  void setInttDev(double array[6])
  {
    m_inttDevs[0] = array[0];
    m_inttDevs[1] = array[1];
    m_inttDevs[2] = array[2];
    m_inttDevs[3] = array[3];
    m_inttDevs[4] = array[4];
    m_inttDevs[5] = array[5];

    inttParam = true;
  }
  void setTpcDev(double array[6])
  {
    m_tpcDevs[0] = array[0];
    m_tpcDevs[1] = array[1];
    m_tpcDevs[2] = array[2];
    m_tpcDevs[3] = array[3];
    m_tpcDevs[4] = array[4];
    m_tpcDevs[5] = array[5];

    tpcParam = true;
  }
  void setMmDev(double array[6])
  {
    m_mmDevs[0] = array[0];
    m_mmDevs[1] = array[1];
    m_mmDevs[2] = array[2];
    m_mmDevs[3] = array[3];
    m_mmDevs[4] = array[4];
    m_mmDevs[5] = array[5];

    mmParam = true;
  }


  void misalignmentFactor(uint8_t layer, const double misalignment)
  {
    auto it = m_misalignmentFactor.find(layer);
    if(it != m_misalignmentFactor.end())
      {
        it->second = misalignment;
        return;
      }
  }

  double getSurfStepPhi() {return m_surfStepPhi;}
  double getSurfStepZ() {return m_surfStepZ;}

  void set_drift_velocity(double vd){m_drift_velocity = vd;}

  void set_nSurfPhi( unsigned int value )
  { m_nSurfPhi = value; }
  void set_intt_survey(bool surv) { m_inttSurvey = surv; }

 private:
  int buildAllGeometry(PHCompositeNode *topNode);
 
  int getNodes(PHCompositeNode*);
  
  int createNodes(PHCompositeNode*);
  
  void setPlanarSurfaceDivisions();
  void editTPCGeometry(PHCompositeNode *topNode);
  void addActsTpcSurfaces(TGeoVolume *tpc_gas_vol, 
                          TGeoManager *geoManager);

  void buildActsSurfaces();

  void makeGeometry(int argc, char* argv[], 
                    ActsExamples::TGeoDetectorWithOptions& detector);
  std::pair<std::shared_ptr<const Acts::TrackingGeometry>,
          std::vector<std::shared_ptr<ActsExamples::IContextDecorator>>>
    build(const boost::program_options::variables_map& vm,
                        ActsExamples::TGeoDetectorWithOptions& detector);

  void readTGeoLayerBuilderConfigsFile(const std::string& path,
                                       ActsExamples::TGeoDetector::Config& config);
 
  void setMaterialResponseFile(std::string& responseFile,
                               std::string& materialFile);

  void getInttKeyFromNode(TGeoNode *gnode);
  void getMvtxKeyFromNode(TGeoNode *gnode);
  void getTpcKeyFromNode(TGeoNode *gnode);
  
  void makeMvtxMapPairs(TrackingVolumePtr &mvtxVolume);
  void makeInttMapPairs(TrackingVolumePtr &inttVolume);
  void makeTpcMapPairs(TrackingVolumePtr &tpcVolume);

  void makeMmMapPairs(TrackingVolumePtr &tpcVolume);
  
  TrkrDefs::hitsetkey getMvtxHitSetKeyFromCoords(unsigned int layer, 
                                                 std::vector<double> &world);
  TrkrDefs::hitsetkey getInttHitSetKeyFromCoords(unsigned int layer,
                                                 std::vector<double> &world);
  TrkrDefs::hitsetkey getTpcHitSetKeyFromCoords(std::vector<double> &world);

//   /// Helper diagnostic function for identifying active layers in subdetectors
//   void isActive(TGeoNode *gnode, int nmax_print);

//   /// Makes map of TrkrHitSetKey<-->TGeoNode
//   void makeTGeoNodeMap(PHCompositeNode *topNode);
  
  void unpackVolumes();

  PHG4CylinderGeomContainer* m_geomContainerMvtx = nullptr;
  PHG4CylinderGeomContainer* m_geomContainerIntt = nullptr;
  PHG4CylinderGeomContainer* m_geomContainerMicromegas = nullptr;
  PHG4TpcCylinderGeomContainer* m_geomContainerTpc = nullptr;
  TGeoManager* m_geoManager = nullptr;

  bool m_inttSurvey = false;

  bool m_useField = true;
  std::map<uint8_t, double> m_misalignmentFactor;

  std::map<TrkrDefs::hitsetkey, TGeoNode*> m_clusterNodeMap;
  std::map<TrkrDefs::hitsetkey, Surface> m_clusterSurfaceMapSilicon;
  std::map<unsigned int, std::vector<Surface>> m_clusterSurfaceMapTpcEdit;  // uses layer as key
  std::map<TrkrDefs::hitsetkey, Surface> m_clusterSurfaceMapMmEdit;
  
  static constexpr unsigned int m_nTpcLayers = 48;
  static constexpr unsigned int m_nTpcModulesPerLayer = 12;
  static constexpr unsigned int m_nTpcSides = 2;

  double m_minSurfZ = 0.;
  double m_maxSurfZ = 105.42;
  unsigned int m_nSurfZ = 1;
  unsigned int m_nSurfPhi = 12;
  double m_surfStepPhi = 0;
  double m_surfStepZ = 0;
  double m_moduleStepPhi = 0;
  double m_modulePhiStart = 0;

  int nprint_tpc = 0;

  const static int m_nTpcSectors = 3;
  
  double m_layerRadius[m_nTpcLayers] = {0};
  double m_layerThickness[m_nTpcLayers] = {0};

  // Spaces to prevent boxes from touching when placed
  const double half_width_clearance_thick = 0.4999;
  const double half_width_clearance_phi = 0.4999;
  const double half_width_clearance_z = 0.5;

  ActsExamples::TGeoDetectorWithOptions m_detector;

  TrackingGeometry m_tGeometry;
  std::shared_ptr<Acts::MagneticFieldProvider> m_magneticField;
  Acts::GeometryContext  m_geoCtxt;  

  ActsGeometry *m_actsGeometry = nullptr;

  int m_verbosity = 0;

  double m_drift_velocity = 8.0e-03;  // cm/ns, override from macro

  std::string m_magField ="1.4" ;
  double m_magFieldRescale = -1.;

  double m_mvtxDevs[6] = {0};
  double m_inttDevs[6] = {0};
  double m_tpcDevs[6] = {0};
  double m_mmDevs[6] = {0};

  bool mvtxParam = false;
  bool inttParam = false;
  bool tpcParam  = false;
  bool mmParam   = false;

};

#endif