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 {
0027 }
0028
0029 int G4RawTowerTTree::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("TOWERHIST");
0039   etot_hist = new TH1F("etot", "total deposited energy", 200, 0, 20);
0040   hm->registerHisto(etot_hist);
0041   PHNodeIterator iter(topNode);
0042   PHCompositeNode *dstNode = static_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "DST"));
0043   G4RootRawTowerContainer *towers = new G4RootRawTowerContainer();
0044   PHIODataNode<PHObject> *node = new PHIODataNode<PHObject>(towers, _outnodename, "PHObject");
0045   dstNode->addNode(node);
0046   evtno = 0;
0047   return 0;
0048 }
0049
0050 int G4RawTowerTTree::process_event(PHCompositeNode *topNode)
0051 {
0052   evtno++;
0053   G4RootRawTowerContainer *towers = findNode::getClass<G4RootRawTowerContainer>(topNode, _outnodename);
0054   RawTowerGeomContainer *rawtowergeom = findNode::getClass<RawTowerGeomContainer>(topNode, _towergeomnodename);
0055
0056   // RawTowerContainer *g4towers = findNode::getClass<RawTowerContainer>(topNode, _towernodename);
0057   TowerInfoContainer *g4towers = findNode::getClass<TowerInfoContainer>(topNode, _towernodename);
0058   if (!g4towers)
0059   {
0060     std::cout << "could not find " << _towernodename << std::endl;
0061     gSystem->Exit(1);
0062   }
0063
0064   double etot = 0;
0065
0066   unsigned int nchannels = g4towers->size();
0067   for (unsigned int channel = 0; channel < nchannels; channel++)
0068   {
0069     TowerInfo *intower = g4towers->get_tower_at_channel(channel);
0070     if (savetowers)
0071     {
0072       unsigned int towerkey = g4towers->encode_key(channel);
0073       int ieta = g4towers->getTowerEtaBin(towerkey);
0074       int iphi = g4towers->getTowerPhiBin(towerkey);
0075
0076       G4RootRawTower roottwr(rawtowergeom->get_etacenter(ieta), rawtowergeom->get_phicenter(iphi), intower->get_energy());
0077       towers->AddG4RootRawTower(roottwr);
0078     }
0079     etot += intower->get_energy();
0080   }
0081   etot_hist->Fill(etot);
0082   towers->set_etotal(etot);
0083   towers->set_event(evtno);
0084   return 0;
0085 }
0086
0087 int G4RawTowerTTree::End(PHCompositeNode * /*topNode*/)
0088 {
0089   hm->dumpHistos(_histofilename);
0090   delete hm;
0091   return 0;
0092 }
0093
0094 void G4RawTowerTTree::Detector(const std::string &det)
0095 {
0096   _detector = det;
0097   _outnodename = "G4RootRawTower_" + det;
0098   _towernodename = "TOWERINFO_CALIB_" + det;
0099   _towergeomnodename = "TOWERGEOM_" + det;
0100   if (_histofilename.empty())
0101   {
0102     _histofilename = "RawTowerHistos_" + det + ".root";
0103   }
0104 }