g4simulation/g4histos/G4ScintillatorTowerTTree.cc

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