g4simulation/g4tpc/PHG4TpcSubsystem.cc

0001 #include "PHG4TpcSubsystem.h"
0002 #include "PHG4TpcDetector.h"
0003 #include "PHG4TpcDisplayAction.h"
0004 #include "PHG4TpcSteppingAction.h"
0005
0006 #include <g4detectors/PHG4DetectorSubsystem.h>  // for PHG4DetectorSubsystem
0007
0008 #include <phparameter/PHParameters.h>
0009
0010 #include <g4main/PHG4DisplayAction.h>  // for PHG4DisplayAction
0011 #include <g4main/PHG4HitContainer.h>
0012 #include <g4main/PHG4SteppingAction.h>  // for PHG4SteppingAction
0013
0014 #include <phool/PHCompositeNode.h>
0015 #include <phool/PHIODataNode.h>    // for PHIODataNode
0016 #include <phool/PHNode.h>          // for PHNode
0017 #include <phool/PHNodeIterator.h>  // for PHNodeIterator
0018 #include <phool/PHObject.h>        // for PHObject
0019 #include <phool/getClass.h>
0020
0021 #include <iostream>  // for operator<<, basic_ost...
0022 #include <set>
0023
0024 class PHG4Detector;
0025
0026 //_______________________________________________________________________
0027 PHG4TpcSubsystem::PHG4TpcSubsystem(const std::string &name, const int lyr)
0028   : PHG4DetectorSubsystem(name, lyr)
0029 {
0030   InitializeParameters();
0031 }
0032
0033 //_______________________________________________________________________
0034 PHG4TpcSubsystem::~PHG4TpcSubsystem()
0035 {
0036   delete m_DisplayAction;
0037 }
0038
0039 //_______________________________________________________________________
0040 int PHG4TpcSubsystem::InitRunSubsystem(PHCompositeNode *topNode)
0041 {
0042   PHNodeIterator iter(topNode);
0043   PHCompositeNode *dstNode = dynamic_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "DST"));
0044
0045   // create display settings before detector (detector adds its volumes to it)
0046   m_DisplayAction = new PHG4TpcDisplayAction(Name());
0047   // create detector
0048   m_Detector = new PHG4TpcDetector(this, topNode, GetParams(), Name());
0049   m_Detector->SuperDetector(SuperDetector());
0050   m_Detector->OverlapCheck(CheckOverlap());
0051   std::set<std::string> nodes;
0052   if (GetParams()->get_int_param("active"))
0053   {
0054     PHNodeIterator dstIter(dstNode);
0055     PHCompositeNode *DetNode = dstNode;
0056     if (SuperDetector() != "NONE" && !SuperDetector().empty())
0057     {
0058       PHNodeIterator iter_dst(dstNode);
0059       DetNode = dynamic_cast<PHCompositeNode *>(iter_dst.findFirst("PHCompositeNode", SuperDetector()));
0060       if (!DetNode)
0061       {
0062         DetNode = new PHCompositeNode(SuperDetector());
0063         dstNode->addNode(DetNode);
0064       }
0065     }
0066     std::string detector_suffix = SuperDetector();
0067     if (detector_suffix == "NONE" || detector_suffix.empty())
0068     {
0069       detector_suffix = Name();
0070     }
0071     m_HitNodeName = "G4HIT_" + detector_suffix;
0072     nodes.insert(m_HitNodeName);
0073     m_AbsorberNodeName = "G4HIT_ABSORBER_" + detector_suffix;
0074     if (GetParams()->get_int_param("absorberactive"))
0075     {
0076       nodes.insert(m_AbsorberNodeName);
0077     }
0078     for (const auto &nodename : nodes)
0079     {
0080       PHG4HitContainer *g4_hits = findNode::getClass<PHG4HitContainer>(topNode, nodename);
0081       if (!g4_hits)
0082       {
0083         g4_hits = new PHG4HitContainer(nodename);
0084         DetNode->addNode(new PHIODataNode<PHObject>(g4_hits, nodename, "PHObject"));
0085       }
0086     }
0087
0088     // create stepping action
0089     m_SteppingAction = new PHG4TpcSteppingAction(m_Detector, GetParams());
0090     m_SteppingAction->SetHitNodeName("G4HIT", m_HitNodeName);
0091     m_SteppingAction->SetHitNodeName("G4HIT_ABSORBER", m_AbsorberNodeName);
0092   }
0093   else
0094   {
0095     // if this is a black hole it does not have to be active
0096     if (GetParams()->get_int_param("blackhole"))
0097     {
0098       m_SteppingAction = new PHG4TpcSteppingAction(m_Detector, GetParams());
0099     }
0100   }
0101   return 0;
0102 }
0103
0104 //_______________________________________________________________________
0105 int PHG4TpcSubsystem::process_event(PHCompositeNode *topNode)
0106 {
0107   // pass top node to stepping action so that it gets
0108   // relevant nodes needed internally
0109   if (m_SteppingAction)
0110   {
0111     m_SteppingAction->SetInterfacePointers(topNode);
0112   }
0113   return 0;
0114 }
0115
0116 void PHG4TpcSubsystem::Print(const std::string &what) const
0117 {
0118   std::cout << Name() << " Parameters: " << std::endl;
0119   GetParams()->Print();
0120   if (m_Detector)
0121   {
0122     m_Detector->Print(what);
0123   }
0124   if (m_SteppingAction)
0125   {
0126     m_SteppingAction->Print(what);
0127   }
0128
0129   return;
0130 }
0131
0132 //_______________________________________________________________________
0133 PHG4Detector *PHG4TpcSubsystem::GetDetector() const
0134 {
0135   return m_Detector;
0136 }
0137
0138 void PHG4TpcSubsystem::SetDefaultParameters()
0139 {
0140   set_default_double_param("gas_inner_radius", 21.6);
0141   set_default_double_param("gas_outer_radius", 76.4);
0142   set_default_double_param("place_x", 0.);
0143   set_default_double_param("place_y", 0.);
0144   set_default_double_param("place_z", 0.);
0145   set_default_double_param("rot_x", 0.);
0146   set_default_double_param("rot_y", 0.);
0147   set_default_double_param("rot_z", 0.);
0148   set_default_double_param("tpc_length", 211.);
0149
0150   set_default_double_param("steplimits", 1);  // 1cm by default
0151
0152   // material budget:
0153   // Cu (all layers): 0.5 oz cu per square foot, 1oz == 0.0347mm --> 0.5 oz ==  0.00347cm/2.
0154   // Kapton insulation 18 layers of * 5mil = 18*0.0127=0.2286
0155   // 250 um FR4 (Substrate for Cu layers)
0156   // HoneyComb (nomex) 1/2 inch=0.5*2.54 cm
0157   set_default_string_param("cage_layer_1_material", "G4_Cu");
0158   set_default_double_param("cage_layer_1_thickness", 0.00347 / 2.);
0159
0160   set_default_string_param("cage_layer_2_material", "FR4");
0161   set_default_double_param("cage_layer_2_thickness", 0.025);
0162
0163   set_default_string_param("cage_layer_3_material", "NOMEX");
0164   set_default_double_param("cage_layer_3_thickness", 0.5 * 2.54);
0165
0166   set_default_string_param("cage_layer_4_material", "G4_Cu");
0167   set_default_double_param("cage_layer_4_thickness", 0.00347 / 2.);
0168
0169   set_default_string_param("cage_layer_5_material", "FR4");
0170   set_default_double_param("cage_layer_5_thickness", 0.025);
0171
0172   set_default_string_param("cage_layer_6_material", "G4_KAPTON");
0173   set_default_double_param("cage_layer_6_thickness", 0.2286);
0174
0175   set_default_string_param("cage_layer_7_material", "G4_Cu");
0176   set_default_double_param("cage_layer_7_thickness", 0.00347 / 2.);
0177
0178   set_default_string_param("cage_layer_8_material", "G4_KAPTON");
0179   set_default_double_param("cage_layer_8_thickness", 0.05);  // 50 um
0180
0181   set_default_string_param("cage_layer_9_material", "G4_Cu");
0182   set_default_double_param("cage_layer_9_thickness", 0.00347 / 2.);
0183
0184   // Thomas K Hemmick <Thomas.Hemmick@stonybrook.edu>
0185   //  The total thickness along Zed would be 5.6 millimeters (+/- 2.8 mm around Zed=0).
0186   //  The outer surfaces would have 0.005 inches (125 um) FR4 coated with a negligible thickness of Al. (revised to Au as below)
0187   //  The interior would be some stiffener of either honeycomb or rohacell.  The range of radiation lengths for this material are:
0188   //  Large cell honeycomb:  1450 cm  (0.028 g/cm^3 density)
0189   //  rohacell:  760 cm (0.052 g/cm^3 density)
0190   //  Close cell honeycomb:  635 cm (0.064 g/cm^3 density)
0191   //  I think a calculation just for the rohacell would be more than sufficient.
0192   set_default_string_param("window_core_material", "ROHACELL_FOAM_51");
0193   set_default_double_param("window_thickness", 0.56);  // overall thickness
0194   // I just checked with PC manufacturers and we can get 8.9 micron thick copper in reasonably large sheets.
0195   //  At normal incidence, 8.9 microns is 0.06% of a radiation length.
0196   set_default_string_param("window_surface1_material", "G4_Cu");
0197   set_default_double_param("window_surface1_thickness", 8.9e-4);  // 8.9  um outter shell thickness be default
0198   // The FR4 should be either 5 or 10 mils thick.  10 mils is 254 microns and 5 mils is 0.127 microns.  I think either of these is mechanically fine...
0199   set_default_string_param("window_surface2_material", "FR4");
0200   set_default_double_param("window_surface2_thickness", 0.0127);  // 127  um 2nd shell thickness be default
0201
0202   // for geonode initialization
0203   set_default_double_param("drift_velocity", 0.008);
0204
0205   set_default_int_param("ntpc_layers_inner", 16);
0206   set_default_int_param("ntpc_layers_mid", 16);
0207   set_default_int_param("ntpc_layers_outer", 16);
0208   set_default_int_param("tpc_minlayer_inner", 7);
0209
0210   set_default_double_param("tpc_minradius_inner", 31.105);  // 30.0);  // cm
0211   set_default_double_param("tpc_minradius_mid", 41.153);    // 40.0);
0212   set_default_double_param("tpc_minradius_outer", 58.367);  // 60.0);
0213
0214   set_default_double_param("tpc_maxradius_inner", 40.249);  // 40.0);  // cm
0215   set_default_double_param("tpc_maxradius_mid", 57.475);    // 60.0);
0216   set_default_double_param("tpc_maxradius_outer", 75.911);  // 77.0);  // from Tom
0217
0218   set_default_double_param("maxdriftlength", 105.5);       // cm
0219   set_default_double_param("extended_readout_time", 0.0);  // ns
0220   set_default_double_param("tpc_adc_clock", 53.326184);         // ns, for 18.8 MHz clock
0221
0222   set_default_double_param("tpc_sector_phi_inner", 0.5024);  // 2 * M_PI / 12 );//sector size in phi for R1 sector
0223   set_default_double_param("tpc_sector_phi_mid", 0.5087);    // 2 * M_PI / 12 );//sector size in phi for R2 sector
0224   set_default_double_param("tpc_sector_phi_outer", 0.5097);  // 2 * M_PI / 12 );//sector size in phi for R3 sector
0225
0226   set_default_int_param("ntpc_phibins_inner", 1128);  // 94 * 12
0227   set_default_int_param("ntpc_phibins_mid", 1536);    // 128 * 12
0228   set_default_int_param("ntpc_phibins_outer", 2304);  // 192 * 12
0229
0230   set_default_double_param("TPC_gas_temperature", 15.0); //in celcius
0231   set_default_double_param("TPC_gas_pressure", 1.0); //in atmospheres
0232   set_default_double_param("Ne_frac", 0.00);
0233   set_default_double_param("Ar_frac", 0.75);
0234   set_default_double_param("CF4_frac", 0.20);
0235   set_default_double_param("N2_frac", 0.00);
0236   set_default_double_param("isobutane_frac", 0.05);
0237 }