g4simulation/g4ihcal/PHG4IHCalSubsystem.cc

0001 #include "PHG4IHCalSubsystem.h"
0002
0003 #include "PHG4IHCalDetector.h"
0004 #include "PHG4IHCalDisplayAction.h"
0005 #include "PHG4IHCalSteppingAction.h"
0006
0007 #include <g4detectors/PHG4DetectorSubsystem.h>  // for PHG4DetectorSubsystem
0008 #include <g4detectors/PHG4HcalDefs.h>
0009
0010 #include <phparameter/PHParameters.h>
0011
0012 #include <g4main/PHG4DisplayAction.h>  // for PHG4DisplayAction
0013 #include <g4main/PHG4HitContainer.h>
0014 #include <g4main/PHG4SteppingAction.h>  // for PHG4SteppingAction
0015
0016 #include <phool/PHCompositeNode.h>
0017 #include <phool/PHIODataNode.h>    // for PHIODataNode
0018 #include <phool/PHNode.h>          // for PHNode
0019 #include <phool/PHNodeIterator.h>  // for PHNodeIterator
0020 #include <phool/PHObject.h>        // for PHObject
0021 #include <phool/getClass.h>
0022
0023 #include <cstdlib>
0024 #include <iostream>  // for operator<<, basic_ostream
0025 #include <limits>
0026 #include <set>
0027
0028 class PHG4Detector;
0029
0030 //_______________________________________________________________________
0031 PHG4IHCalSubsystem::PHG4IHCalSubsystem(const std::string &name, const int lyr)
0032   : PHG4DetectorSubsystem(name, lyr)
0033 {
0034   InitializeParameters();
0035 }
0036
0037 //_______________________________________________________________________
0038 PHG4IHCalSubsystem::~PHG4IHCalSubsystem()
0039 {
0040   delete m_DisplayAction;
0041 }
0042
0043 //_______________________________________________________________________
0044 int PHG4IHCalSubsystem::InitRunSubsystem(PHCompositeNode *topNode)
0045 {
0046   PHNodeIterator iter(topNode);
0047   PHCompositeNode *dstNode = dynamic_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "DST"));
0048
0049   // create display settings before detector
0050   m_DisplayAction = new PHG4IHCalDisplayAction(Name());
0051
0052   // create detector
0053   m_Detector = new PHG4IHCalDetector(this, topNode, GetParams(), Name());
0054   m_Detector->SuperDetector(SuperDetector());
0055   m_Detector->OverlapCheck(CheckOverlap());
0056   std::set<std::string> nodes;
0057   if (GetParams()->get_int_param("active"))
0058   {
0059     PHNodeIterator dstIter(dstNode);
0060     PHCompositeNode *DetNode = dstNode;
0061
0062     if (SuperDetector() != "NONE" && !SuperDetector().empty())
0063     {
0064       PHNodeIterator iter_dst(dstNode);
0065       DetNode = dynamic_cast<PHCompositeNode *>(iter_dst.findFirst("PHCompositeNode", SuperDetector()));
0066
0067       if (!DetNode)
0068       {
0069         DetNode = new PHCompositeNode(SuperDetector());
0070         dstNode->addNode(DetNode);
0071       }
0072     }
0073     std::string detector_suffix = SuperDetector();
0074     if (detector_suffix == "NONE" || detector_suffix.empty())
0075     {
0076       detector_suffix = Name();
0077     }
0078
0079     m_HitNodeName = "G4HIT_" + detector_suffix;
0080     nodes.insert(m_HitNodeName);
0081     m_AbsorberNodeName = "G4HIT_ABSORBER_" + detector_suffix;
0082     if (GetParams()->get_int_param("absorberactive"))
0083     {
0084       nodes.insert(m_AbsorberNodeName);
0085     }
0086     for (const auto &nodename : nodes)
0087     {
0088       PHG4HitContainer *g4_hits = findNode::getClass<PHG4HitContainer>(topNode, nodename);
0089       if (!g4_hits)
0090       {
0091         g4_hits = new PHG4HitContainer(nodename);
0092         DetNode->addNode(new PHIODataNode<PHObject>(g4_hits, nodename, "PHObject"));
0093       }
0094     }
0095
0096     // create stepping action
0097     m_SteppingAction = new PHG4IHCalSteppingAction(m_Detector, GetParams());
0098     m_SteppingAction->InitWithNode(topNode);
0099     m_SteppingAction->SetHitNodeName("G4HIT", m_HitNodeName);
0100     m_SteppingAction->SetHitNodeName("G4HIT_ABSORBER", m_AbsorberNodeName);
0101   }
0102   else
0103   {
0104     // if this is a black hole it does not have to be active
0105     if (GetParams()->get_int_param("blackhole"))
0106     {
0107       m_SteppingAction = new PHG4IHCalSteppingAction(m_Detector, GetParams());
0108       m_SteppingAction->InitWithNode(topNode);
0109     }
0110   }
0111   return 0;
0112 }
0113
0114 //_______________________________________________________________________
0115 int PHG4IHCalSubsystem::process_event(PHCompositeNode *topNode)
0116 {
0117   // pass top node to stepping action so that it gets
0118   // relevant nodes needed internally
0119   if (m_SteppingAction)
0120   {
0121     m_SteppingAction->SetInterfacePointers(topNode);
0122   }
0123   return 0;
0124 }
0125
0126 void PHG4IHCalSubsystem::Print(const std::string &what) const
0127 {
0128   std::cout << Name() << " Parameters: " << std::endl;
0129   GetParams()->Print();
0130   if (m_Detector)
0131   {
0132     m_Detector->Print(what);
0133   }
0134   if (m_SteppingAction)
0135   {
0136     m_SteppingAction->Print(what);
0137   }
0138
0139   return;
0140 }
0141
0142 //_______________________________________________________________________
0143 PHG4Detector *PHG4IHCalSubsystem::GetDetector() const
0144 {
0145   return m_Detector;
0146 }
0147
0148 void PHG4IHCalSubsystem::SetDefaultParameters()
0149 {
0150   set_default_double_param(PHG4HcalDefs::innerrad, 115);
0151   set_default_double_param("light_balance_inner_corr", std::numeric_limits<double>::quiet_NaN());
0152   set_default_double_param("light_balance_inner_radius", std::numeric_limits<double>::quiet_NaN());
0153   set_default_double_param("light_balance_outer_corr", std::numeric_limits<double>::quiet_NaN());
0154   set_default_double_param("light_balance_outer_radius", std::numeric_limits<double>::quiet_NaN());
0155   set_default_double_param("phistart", std::numeric_limits<double>::quiet_NaN());
0156   set_default_double_param("scinti_eta_coverage_neg", 1.1);
0157   set_default_double_param("scinti_eta_coverage_pos", 1.1);
0158   set_default_double_param(PHG4HcalDefs::outerrad, 274.010 / 2 + 3);
0159   set_default_double_param("place_x", 0.);
0160   set_default_double_param("place_y", 0.);
0161   set_default_double_param("place_z", 0.);
0162   set_default_double_param("rot_x", 0.);
0163   set_default_double_param("rot_y", 180.);
0164   set_default_double_param("rot_z", 0.);
0165   set_default_double_param("size_z", 435.000 + 10);
0166   set_default_double_param("Birk_const", 0.07943);
0167   set_default_double_param("tmin", -20.);
0168   set_default_double_param("tmax", 60.);
0169   set_default_double_param("dt", 100.);
0170
0171   set_default_int_param("light_scint_model", 1);
0172   set_default_int_param(PHG4HcalDefs::n_towers, 64);
0173   set_default_int_param(PHG4HcalDefs::scipertwr, 4);
0174   set_default_int_param(PHG4HcalDefs::n_scinti_tiles, 12);
0175   set_default_int_param("etabins", 24);
0176   set_default_int_param("saveg4hit", 1);
0177
0178   set_default_string_param("GDMPath", "HCALIN_GDML");
0179   set_default_string_param("MapFileName", "HCALIN_MEPHI_MAP");
0180   set_default_string_param("MapHistoName", "ihcalcombinedgdmlnormtbyt");
0181 }
0182
0183 void PHG4IHCalSubsystem::SetLightCorrection(const double inner_radius, const double inner_corr, const double outer_radius, const double outer_corr)
0184 {
0185   set_double_param("light_balance_inner_corr", inner_corr);
0186   set_double_param("light_balance_inner_radius", inner_radius);
0187   set_double_param("light_balance_outer_corr", outer_corr);
0188   set_double_param("light_balance_outer_radius", outer_radius);
0189   return;
0190 }