g4simulation/g4ohcal/PHG4OHCalSubsystem.cc

0001 #include "PHG4OHCalSubsystem.h"
0002
0003 #include "PHG4OHCalDetector.h"
0004 #include "PHG4OHCalDisplayAction.h"
0005 #include "PHG4OHCalSteppingAction.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>   // for getenv
0024 #include <iostream>  // for operator<<, basic_ostream
0025 #include <limits>
0026 #include <set>  // for set
0027
0028 class PHG4Detector;
0029
0030 //_______________________________________________________________________
0031 PHG4OHCalSubsystem::PHG4OHCalSubsystem(const std::string &name, const int lyr)
0032   : PHG4DetectorSubsystem(name, lyr)
0033 {
0034   InitializeParameters();
0035 }
0036
0037 //_______________________________________________________________________
0038 PHG4OHCalSubsystem::~PHG4OHCalSubsystem()
0039 {
0040   delete m_DisplayAction;
0041 }
0042
0043 //_______________________________________________________________________
0044 int PHG4OHCalSubsystem::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 PHG4OHCalDisplayAction(Name());
0051
0052   if (get_string_param("IronFieldMapPath") == "DefaultParameters-InvadPath")
0053   {
0054     std::cout << __PRETTY_FUNCTION__ << ": invalid string parameter IronFieldMapPath, where we expect a 3D field map" << std::endl;
0055     exit(1);
0056   }
0057
0058   // create detector
0059   m_Detector = new PHG4OHCalDetector(this, topNode, GetParams(), Name());
0060   m_Detector->SuperDetector(SuperDetector());
0061   m_Detector->OverlapCheck(CheckOverlap());
0062   std::set<std::string> nodes;
0063   if (GetParams()->get_int_param("active"))
0064   {
0065     PHNodeIterator dstIter(dstNode);
0066     PHCompositeNode *DetNode = dstNode;
0067     if (SuperDetector() != "NONE" && !SuperDetector().empty())
0068     {
0069       PHNodeIterator iter_dst(dstNode);
0070       DetNode = dynamic_cast<PHCompositeNode *>(iter_dst.findFirst("PHCompositeNode", SuperDetector()));
0071
0072       if (!DetNode)
0073       {
0074         DetNode = new PHCompositeNode(SuperDetector());
0075         dstNode->addNode(DetNode);
0076       }
0077     }
0078     std::string detector_suffix = SuperDetector();
0079     if (detector_suffix == "NONE" || detector_suffix.empty())
0080     {
0081       detector_suffix = Name();
0082     }
0083
0084     m_HitNodeName = "G4HIT_" + detector_suffix;
0085     nodes.insert(m_HitNodeName);
0086     m_AbsorberNodeName = "G4HIT_ABSORBER_" + detector_suffix;
0087     if (GetParams()->get_int_param("absorberactive"))
0088     {
0089       nodes.insert(m_AbsorberNodeName);
0090     }
0091     for (const auto &nodename : nodes)
0092     {
0093       PHG4HitContainer *g4_hits = findNode::getClass<PHG4HitContainer>(topNode, nodename);
0094       if (!g4_hits)
0095       {
0096         g4_hits = new PHG4HitContainer(nodename);
0097         DetNode->addNode(new PHIODataNode<PHObject>(g4_hits, nodename, "PHObject"));
0098       }
0099     }
0100     // create stepping action
0101     m_SteppingAction = new PHG4OHCalSteppingAction(m_Detector, GetParams());
0102     m_SteppingAction->InitWithNode(topNode);
0103     m_SteppingAction->SetHitNodeName("G4HIT", m_HitNodeName);
0104     m_SteppingAction->SetHitNodeName("G4HIT_ABSORBER", m_AbsorberNodeName);
0105   }
0106   else
0107   {
0108     if (GetParams()->get_int_param("blackhole"))
0109     {
0110       m_SteppingAction = new PHG4OHCalSteppingAction(m_Detector, GetParams());
0111       m_SteppingAction->InitWithNode(topNode);
0112     }
0113   }
0114
0115   return 0;
0116 }
0117
0118 //_______________________________________________________________________
0119 int PHG4OHCalSubsystem::process_event(PHCompositeNode *topNode)
0120 {
0121   // pass top node to stepping action so that it gets
0122   // relevant nodes needed internally
0123   if (m_SteppingAction)
0124   {
0125     m_SteppingAction->SetInterfacePointers(topNode);
0126   }
0127   return 0;
0128 }
0129
0130 void PHG4OHCalSubsystem::Print(const std::string &what) const
0131 {
0132   std::cout << "Outer Hcal Parameters: " << std::endl;
0133   GetParams()->Print();
0134   if (m_Detector)
0135   {
0136     m_Detector->Print(what);
0137   }
0138   return;
0139 }
0140
0141 //_______________________________________________________________________
0142 PHG4Detector *PHG4OHCalSubsystem::GetDetector() const
0143 {
0144   return m_Detector;
0145 }
0146
0147 void PHG4OHCalSubsystem::SetLightCorrection(const double inner_radius, const double inner_corr, const double outer_radius, const double outer_corr)
0148 {
0149   set_double_param("light_balance_inner_corr", inner_corr);
0150   set_double_param("light_balance_inner_radius", inner_radius);
0151   set_double_param("light_balance_outer_corr", outer_corr);
0152   set_double_param("light_balance_outer_radius", outer_radius);
0153   return;
0154 }
0155
0156 void PHG4OHCalSubsystem::SetDefaultParameters()
0157 {
0158   set_default_double_param("inner_radius", 182.423 - 5);
0159   set_default_double_param("light_balance_inner_corr", std::numeric_limits<double>::quiet_NaN());
0160   set_default_double_param("light_balance_inner_radius", std::numeric_limits<double>::quiet_NaN());
0161   set_default_double_param("light_balance_outer_corr", std::numeric_limits<double>::quiet_NaN());
0162   set_default_double_param("light_balance_outer_radius", std::numeric_limits<double>::quiet_NaN());
0163   set_default_double_param("phistart", std::numeric_limits<double>::quiet_NaN());
0164   set_default_double_param("scinti_eta_coverage_neg", 1.1);
0165   set_default_double_param("scinti_eta_coverage_pos", 1.1);
0166   set_default_double_param("outer_radius", 269.317 + 5);
0167   set_default_double_param("place_x", 0.);
0168   set_default_double_param("place_y", 0.);
0169   set_default_double_param("place_z", 0.);
0170   set_default_double_param("rot_x", 0.);
0171   set_default_double_param("rot_y", 180.);
0172   set_default_double_param("rot_z", 0.);
0173   set_default_double_param("size_z", 639.240 + 10);
0174   set_default_double_param("Birk_const", 0.07943);
0175   set_default_double_param("tmin", -20.);
0176   set_default_double_param("tmax", 60.);
0177   set_default_double_param("dt", 100.);
0178   set_default_int_param("field_check", 0);
0179   set_default_int_param("light_scint_model", 1);
0180
0181   set_default_int_param("n_towers", 64);
0182   set_default_int_param(PHG4HcalDefs::scipertwr, 5);
0183   set_default_int_param("n_scinti_tiles", 12);
0184   set_default_int_param("etabins", 24);
0185   set_default_int_param("saveg4hit", 1);
0186
0187   set_default_string_param("GDMPath", "HCALOUT_GDML");           // use CDB
0188   set_default_string_param("MapFileName", "HCALOUT_MEPHI_MAP");  // use CDB
0189   set_default_string_param("MapHistoName", "ohcal_mephi_map_towerid_");
0190   set_default_string_param("IronFieldMapPath", "HCALOUT_STEEL_MAP");  // use CDB
0191   set_default_double_param("IronFieldMapScale", 1.);
0192 }