g4simulation/g4histos/G4RawTowerTTree.cc

0001 #include "G4RawTowerTTree.h"
0002
0003 #include "G4RootRawTower.h"
0004 #include "G4RootRawTowerContainer.h"
0005
0006 #include <calobase/RawTowerGeomContainer.h>
0007 #include <calobase/TowerInfo.h>
0008 #include <calobase/TowerInfoContainer.h>
0009
0010 #include <fun4all/Fun4AllHistoManager.h>
0011 #include <fun4all/SubsysReco.h>  // for SubsysReco
0012
0013 #include <phool/PHCompositeNode.h>
0014 #include <phool/PHIODataNode.h>    // for PHIODataNode
0015 #include <phool/PHNodeIterator.h>  // for PHNodeIterator
0016 #include <phool/PHObject.h>        // for PHObject
0017 #include <phool/getClass.h>
0018
0019 #include <TH1.h>
0020 #include <TSystem.h>
0021
0022 #include <iostream>  // for operator<<, endl, basic_...
0023
0024 G4RawTowerTTree::G4RawTowerTTree(const std::string &name)
0025   : SubsysReco(name)
0026   , savetowers(1)
0027   , evtno(0)
0028   , hm(nullptr)
0029   , etot_hist(nullptr)
0030 {
0031 }
0032
0033 int G4RawTowerTTree::Init(PHCompositeNode *topNode)
0034 {
0035   if (_detector.empty())
0036   {
0037     std::cout << "Detector not set via Detector(<name>) method" << std::endl;
0038     std::cout << "(it is the name appended to the G4TOWER_<name> nodename)" << std::endl;
0039     std::cout << "you do not want to run like this, exiting now" << std::endl;
0040     gSystem->Exit(1);
0041   }
0042   hm = new Fun4AllHistoManager("TOWERHIST");
0043   etot_hist = new TH1F("etot", "total deposited energy", 200, 0, 20);
0044   hm->registerHisto(etot_hist);
0045   PHNodeIterator iter(topNode);
0046   PHCompositeNode *dstNode = static_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "DST"));
0047   G4RootRawTowerContainer *towers = new G4RootRawTowerContainer();
0048   PHIODataNode<PHObject> *node = new PHIODataNode<PHObject>(towers, _outnodename, "PHObject");
0049   dstNode->addNode(node);
0050   evtno = 0;
0051   return 0;
0052 }
0053
0054 int G4RawTowerTTree::process_event(PHCompositeNode *topNode)
0055 {
0056   evtno++;
0057   G4RootRawTowerContainer *towers = findNode::getClass<G4RootRawTowerContainer>(topNode, _outnodename);
0058   RawTowerGeomContainer *rawtowergeom = findNode::getClass<RawTowerGeomContainer>(topNode, _towergeomnodename);
0059
0060   // RawTowerContainer *g4towers = findNode::getClass<RawTowerContainer>(topNode, _towernodename);
0061   TowerInfoContainer *g4towers = findNode::getClass<TowerInfoContainer>(topNode, _towernodename);
0062   if (!g4towers)
0063   {
0064     std::cout << "could not find " << _towernodename << std::endl;
0065     gSystem->Exit(1);
0066   }
0067
0068   double etot = 0;
0069
0070   unsigned int nchannels = g4towers->size();
0071   for (unsigned int channel = 0; channel < nchannels; channel++)
0072   {
0073     TowerInfo *intower = g4towers->get_tower_at_channel(channel);
0074     if (savetowers)
0075     {
0076       unsigned int towerkey = g4towers->encode_key(channel);
0077       int ieta = g4towers->getTowerEtaBin(towerkey);
0078       int iphi = g4towers->getTowerPhiBin(towerkey);
0079
0080       G4RootRawTower roottwr(rawtowergeom->get_etacenter(ieta), rawtowergeom->get_phicenter(iphi), intower->get_energy());
0081       towers->AddG4RootRawTower(roottwr);
0082     }
0083     etot += intower->get_energy();
0084   }
0085   etot_hist->Fill(etot);
0086   towers->set_etotal(etot);
0087   towers->set_event(evtno);
0088   return 0;
0089 }
0090
0091 int G4RawTowerTTree::End(PHCompositeNode * /*topNode*/)
0092 {
0093   hm->dumpHistos(_histofilename);
0094   delete hm;
0095   return 0;
0096 }
0097
0098 void G4RawTowerTTree::Detector(const std::string &det)
0099 {
0100   _detector = det;
0101   _outnodename = "G4RootRawTower_" + det;
0102   _towernodename = "TOWERINFO_CALIB_" + det;
0103   _towergeomnodename = "TOWERGEOM_" + det;
0104   if (!_histofilename.size())
0105   {
0106     _histofilename = "RawTowerHistos_" + det + ".root";
0107   }
0108 }