Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-08-05 08:17:39 

0001 #include "RawTowerBuilder.h"
0002 
0003 #include <calobase/RawTower.h>  // for RawTower
0004 #include <calobase/RawTowerContainer.h>
0005 #include <calobase/RawTowerDefs.h>           // for encode_towerid
0006 #include <calobase/RawTowerGeom.h>           // for RawTowerGeom
0007 #include <calobase/RawTowerGeomContainer.h>  // for RawTowerGeomC...
0008 #include <calobase/RawTowerGeomContainer_Cylinderv1.h>
0009 #include <calobase/RawTowerGeomv1.h>
0010 #include <calobase/RawTowerv1.h>
0011 #include <calobase/TowerInfoContainer.h>
0012 #include <calobase/TowerInfoContainerv1.h>
0013 #include <calobase/TowerInfo.h>
0014 
0015 #include <g4detectors/PHG4CylinderCellGeom.h>
0016 #include <g4detectors/PHG4CylinderCellGeomContainer.h>
0017 #include <g4detectors/PHG4Cell.h>
0018 #include <g4detectors/PHG4CellContainer.h>
0019 #include <g4detectors/PHG4CellDefs.h>
0020 
0021 #include <g4main/PHG4Utils.h>
0022 
0023 #include <fun4all/Fun4AllReturnCodes.h>
0024 #include <fun4all/SubsysReco.h>  // for SubsysReco
0025 
0026 #include <phool/PHCompositeNode.h>
0027 #include <phool/PHIODataNode.h>
0028 #include <phool/PHNode.h>  // for PHNode
0029 #include <phool/PHNodeIterator.h>
0030 #include <phool/PHObject.h>  // for PHObject
0031 #include <phool/getClass.h>
0032 #include <phool/phool.h>  // for PHWHERE
0033 
0034 #include <boost/io/ios_state.hpp>
0035 
0036 #include <cmath>      // for fabs, tan, atan2
0037 #include <cstdlib>    // for exit
0038 #include <exception>  // for exception
0039 #include <iostream>
0040 #include <map>
0041 #include <stdexcept>
0042 #include <utility>  // for pair, make_pair
0043 
0044 RawTowerBuilder::RawTowerBuilder(const std::string &name)
0045   : SubsysReco(name)
0046 {
0047 }
0048 
0049 int RawTowerBuilder::InitRun(PHCompositeNode *topNode)
0050 {
0051   std::string geonodename = "CYLINDERCELLGEOM_" + m_Detector;
0052   PHG4CylinderCellGeomContainer *cellgeos = findNode::getClass<PHG4CylinderCellGeomContainer>(topNode, geonodename);
0053   if (!cellgeos)
0054   {
0055     std::cout << PHWHERE << " " << geonodename
0056               << " Node missing, doing nothing." << std::endl;
0057     throw std::runtime_error(
0058         "Failed to find " + geonodename + " node in RawTowerBuilder::CreateNodes");
0059   }
0060 
0061   // fill the number of layers in the calorimeter
0062   m_NumLayers = cellgeos->get_NLayers();
0063 
0064   PHNodeIterator iter(topNode);
0065 
0066   // Looking for the DST node
0067   PHCompositeNode *dstNode = dynamic_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "DST"));
0068   if (!dstNode)
0069   {
0070     std::cout << PHWHERE << "DST Node missing, doing nothing." << std::endl;
0071     exit(1);
0072   }
0073 
0074   try
0075   {
0076     CreateNodes(topNode);
0077   }
0078   catch (std::exception &e)
0079   {
0080     std::cout << e.what() << std::endl;
0081     //exit(1);
0082   }
0083 
0084   if (Verbosity() >= 1)
0085   {
0086     std::cout << "RawTowerBuilder::InitRun :";
0087     if (m_TowerEnergySrcEnum == kEnergyDeposition)
0088     {
0089       std::cout << "save Geant4 energy deposition as the weight of the cells"
0090                 << std::endl;
0091     }
0092     else if (m_TowerEnergySrcEnum == kLightYield)
0093     {
0094       std::cout << "save light yield as the weight of the cells" << std::endl;
0095     }
0096   }
0097   return Fun4AllReturnCodes::EVENT_OK;
0098 }
0099 
0100 int RawTowerBuilder::process_event(PHCompositeNode *topNode)
0101 {
0102   if (Verbosity())
0103   {
0104     std::cout << PHWHERE << "Process event entered" << std::endl;
0105   }
0106 
0107   //load get TowerInfoContainer node from node tree:
0108   TowerInfoContainer*  m_TowerInfoContainer = findNode::getClass<TowerInfoContainer>(topNode,m_TowerInfoNodeName);
0109   if (!m_TowerInfoContainer && m_UseTowerInfo > 0)
0110     {
0111       std::cout << PHWHERE << "TowerInfoContainer Node missing, doing nothing." << std::endl;
0112       exit(1);
0113     }
0114 
0115     // get cells
0116   std::string cellnodename = "G4CELL_" + m_Detector;
0117   PHG4CellContainer *cells = findNode::getClass<PHG4CellContainer>(topNode, cellnodename);
0118   if (!cells)
0119   {
0120     std::cout << PHWHERE << " " << cellnodename
0121               << " Node missing, doing nothing." << std::endl;
0122     return Fun4AllReturnCodes::ABORTEVENT;
0123   }
0124 
0125   // loop over all cells in an event
0126   PHG4CellContainer::ConstIterator cell_iter;
0127   PHG4CellContainer::ConstRange cell_range = cells->getCells();
0128   for (cell_iter = cell_range.first; cell_iter != cell_range.second; ++cell_iter)
0129   {
0130     PHG4Cell *cell = cell_iter->second;
0131 
0132     if (Verbosity() > 2)
0133     {
0134       std::cout << PHWHERE << " print out the cell:" << std::endl;
0135       cell->identify();
0136     }
0137 
0138 
0139     //Calculate the cell weight
0140     float cell_weight = 0;
0141     if (m_TowerEnergySrcEnum == kEnergyDeposition)
0142     {
0143       cell_weight = cell->get_edep();
0144     }
0145     else if (m_TowerEnergySrcEnum == kLightYield)
0146     {
0147       cell_weight = cell->get_light_yield();
0148     }
0149 
0150 
0151     // add the energy to the corresponding tower
0152     if (m_UseTowerInfo  != 1)
0153       {
0154     RawTower *tower = nullptr;
0155     short int firstpar;
0156     short int secondpar;
0157     if (cell->has_binning(PHG4CellDefs::sizebinning))
0158       {
0159         firstpar = PHG4CellDefs::SizeBinning::get_zbin(cell->get_cellid());
0160         secondpar = PHG4CellDefs::SizeBinning::get_phibin(cell->get_cellid());
0161       }
0162     else if (cell->has_binning(PHG4CellDefs::spacalbinning))
0163       {
0164         firstpar = PHG4CellDefs::SpacalBinning::get_etabin(cell->get_cellid());
0165         secondpar = PHG4CellDefs::SpacalBinning::get_phibin(cell->get_cellid());
0166       }
0167     else
0168       {
0169         boost::io::ios_flags_saver ifs(std::cout);
0170         std::cout << "unknown cell binning, implement 0x" << std::hex << PHG4CellDefs::get_binning(cell->get_cellid()) << std::dec << std::endl;
0171         exit(1);
0172       }
0173     tower = m_TowerContainer->getTower(firstpar, secondpar);
0174     if (!tower)
0175       {
0176         tower = new RawTowerv1();
0177         tower->set_energy(0);
0178         m_TowerContainer->AddTower(firstpar, secondpar, tower);
0179       }
0180     
0181     tower->add_ecell(cell->get_cellid(), cell_weight);
0182     
0183     PHG4Cell::ShowerEdepConstRange range = cell->get_g4showers();
0184     for (PHG4Cell::ShowerEdepConstIterator shower_iter = range.first;
0185          shower_iter != range.second;
0186          ++shower_iter)
0187       {
0188         tower->add_eshower(shower_iter->first, shower_iter->second);
0189       }
0190     
0191     tower->set_energy(tower->get_energy() + cell_weight);
0192     
0193     if (Verbosity() > 2)
0194       {
0195         m_RawTowerGeomContainer = findNode::getClass<RawTowerGeomContainer>(topNode, m_TowerGeomNodeName);
0196         tower->identify();
0197       }
0198       }
0199     if (m_UseTowerInfo > 0)
0200       {
0201     TowerInfo *towerinfo;
0202     unsigned int etabin;
0203     unsigned int phibin;
0204     if (cell->has_binning(PHG4CellDefs::spacalbinning))
0205       {
0206         etabin = PHG4CellDefs::SpacalBinning::get_etabin(cell->get_cellid());
0207         phibin = PHG4CellDefs::SpacalBinning::get_phibin(cell->get_cellid());
0208       }
0209     else
0210       {
0211         boost::io::ios_flags_saver ifs(std::cout);
0212         std::cout << "unknown cell binning, implement 0x" << std::hex << PHG4CellDefs::get_binning(cell->get_cellid()) << std::dec << std::endl;
0213         exit(1);
0214       }
0215     unsigned int towerkey = (etabin << 16U) + phibin;
0216     // unsigned int towerindex = m_TowerInfoContainer->decode_key(towerkey);
0217     towerinfo = m_TowerInfoContainer->get_tower_at_key(towerkey);
0218         if (!towerinfo)
0219         {
0220           std::cout << __PRETTY_FUNCTION__ << ": missing towerkey = " << towerkey << " in m_TowerInfoContainer!";
0221           exit(1);
0222         }
0223         else
0224         {
0225           towerinfo->set_energy(towerinfo->get_energy() + cell_weight);
0226         }
0227       }
0228   }
0229 
0230   if (m_UseTowerInfo != 1 )
0231     {
0232       double towerE = 0;
0233       if (m_ChkEnergyConservationFlag)
0234     {
0235       double cellE = cells->getTotalEdep();
0236       towerE = m_TowerContainer->getTotalEdep();
0237       if (fabs(cellE - towerE) / cellE > 1e-5)
0238         {
0239           std::cout << "towerE: " << towerE << ", cellE: " << cellE << ", delta: "
0240             << cellE - towerE << std::endl;
0241         }
0242     }
0243       if (Verbosity())
0244     {
0245       towerE = m_TowerContainer->getTotalEdep();
0246     }
0247       
0248       m_TowerContainer->compress(m_Emin);
0249       if (Verbosity())
0250     {
0251       std::cout << "Energy lost by dropping towers with less than " << m_Emin
0252             << " GeV energy, lost energy: " << towerE - m_TowerContainer->getTotalEdep()
0253             << std::endl;
0254       m_TowerContainer->identify();
0255       RawTowerContainer::ConstRange begin_end = m_TowerContainer->getTowers();
0256       RawTowerContainer::ConstIterator iter;
0257       for (iter = begin_end.first; iter != begin_end.second; ++iter)
0258         {
0259           iter->second->identify();
0260         }
0261     }
0262     }
0263   return Fun4AllReturnCodes::EVENT_OK;
0264 }
0265 
0266 void RawTowerBuilder::CreateNodes(PHCompositeNode *topNode)
0267 {
0268   PHNodeIterator iter(topNode);
0269   PHCompositeNode *runNode = dynamic_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "RUN"));
0270   if (!runNode)
0271   {
0272     std::cout << PHWHERE << "Run Node missing, doing nothing." << std::endl;
0273     throw std::runtime_error("Failed to find Run node in RawTowerBuilder::CreateNodes");
0274   }
0275 
0276   PHNodeIterator runIter(runNode);
0277   PHCompositeNode *RunDetNode = dynamic_cast<PHCompositeNode *>(runIter.findFirst("PHCompositeNode", m_Detector));
0278   if (!RunDetNode)
0279   {
0280     RunDetNode = new PHCompositeNode(m_Detector);
0281     runNode->addNode(RunDetNode);
0282   }
0283 
0284   const RawTowerDefs::CalorimeterId caloid = RawTowerDefs::convert_name_to_caloid(m_Detector);
0285 
0286   // get the cell geometry and build up the tower geometry object
0287   std::string geonodename = "CYLINDERCELLGEOM_" + m_Detector;
0288   PHG4CylinderCellGeomContainer *cellgeos = findNode::getClass<PHG4CylinderCellGeomContainer>(topNode, geonodename);
0289   if (!cellgeos)
0290   {
0291     std::cout << PHWHERE << " " << geonodename
0292               << " Node missing, doing nothing." << std::endl;
0293     throw std::runtime_error(
0294         "Failed to find " + geonodename + " node in RawTowerBuilder::CreateNodes");
0295   }
0296   m_TowerGeomNodeName = "TOWERGEOM_" + m_Detector;
0297   m_RawTowerGeomContainer = findNode::getClass<RawTowerGeomContainer>(topNode, m_TowerGeomNodeName);
0298   if (!m_RawTowerGeomContainer)
0299   {
0300     m_RawTowerGeomContainer = new RawTowerGeomContainer_Cylinderv1(caloid);
0301     PHIODataNode<PHObject> *newNode = new PHIODataNode<PHObject>(m_RawTowerGeomContainer, m_TowerGeomNodeName, "PHObject");
0302     RunDetNode->addNode(newNode);
0303   }
0304   // fill the number of layers in the calorimeter
0305   m_NumLayers = cellgeos->get_NLayers();
0306 
0307   // Create the tower nodes on the tree
0308   PHG4CylinderCellGeomContainer::ConstIterator miter;
0309   PHG4CylinderCellGeomContainer::ConstRange begin_end =
0310       cellgeos->get_begin_end();
0311   int ifirst = 1;
0312   int first_layer = -1;
0313   PHG4CylinderCellGeom *first_cellgeo = nullptr;
0314   double inner_radius = 0;
0315   double thickness = 0;
0316   for (miter = begin_end.first; miter != begin_end.second; ++miter)
0317   {
0318     PHG4CylinderCellGeom *cellgeo = miter->second;
0319 
0320     if (Verbosity())
0321     {
0322       cellgeo->identify();
0323     }
0324     thickness += cellgeo->get_thickness();
0325     if (ifirst)
0326     {
0327       first_cellgeo = miter->second;
0328       m_CellBinning = cellgeo->get_binning();
0329       m_NumPhiBins = cellgeo->get_phibins();
0330       m_PhiMin = cellgeo->get_phimin();
0331       m_PhiStep = cellgeo->get_phistep();
0332       if (m_CellBinning == PHG4CellDefs::etaphibinning || m_CellBinning == PHG4CellDefs::etaslatbinning)
0333       {
0334         m_NumEtaBins = cellgeo->get_etabins();
0335         m_EtaMin = cellgeo->get_etamin();
0336         m_EtaStep = cellgeo->get_etastep();
0337       }
0338       else if (m_CellBinning == PHG4CellDefs::sizebinning)
0339       {
0340         m_NumEtaBins = cellgeo->get_zbins();  // bin eta in the same number of z bins
0341       }
0342       else if (m_CellBinning == PHG4CellDefs::spacalbinning)
0343       {
0344         // use eta definiton for each row of towers
0345         m_NumEtaBins = cellgeo->get_etabins();
0346       }
0347       else
0348       {
0349         std::cout << "RawTowerBuilder::CreateNodes::" << Name()
0350                   << " - Fatal Error - unsupported cell binning method "
0351                   << m_CellBinning << std::endl;
0352       }
0353       inner_radius = cellgeo->get_radius();
0354       first_layer = cellgeo->get_layer();
0355       ifirst = 0;
0356     }
0357     else
0358     {
0359       if (m_CellBinning != cellgeo->get_binning())
0360       {
0361         std::cout << "inconsistent binning method from " << m_CellBinning
0362                   << " layer " << cellgeo->get_layer() << ": "
0363                   << cellgeo->get_binning() << std::endl;
0364         exit(1);
0365       }
0366       if (inner_radius > cellgeo->get_radius())
0367       {
0368         std::cout << "radius of layer " << cellgeo->get_layer() << " is "
0369                   << cellgeo->get_radius() << " which smaller than radius "
0370                   << inner_radius << " of first layer in list " << first_layer
0371                   << std::endl;
0372       }
0373       if (m_NumPhiBins != cellgeo->get_phibins())
0374       {
0375         std::cout << "mixing number of phibins, fisrt layer: " << m_NumPhiBins
0376                   << " layer " << cellgeo->get_layer() << ": "
0377                   << cellgeo->get_phibins() << std::endl;
0378         exit(1);
0379       }
0380       if (m_PhiMin != cellgeo->get_phimin())
0381       {
0382         std::cout << "mixing number of phimin, fisrt layer: " << m_PhiMin
0383                   << " layer " << cellgeo->get_layer() << ": "
0384                   << cellgeo->get_phimin() << std::endl;
0385         exit(1);
0386       }
0387       if (m_PhiStep != cellgeo->get_phistep())
0388       {
0389         std::cout << "mixing phisteps first layer: " << m_PhiStep << " layer "
0390                   << cellgeo->get_layer() << ": " << cellgeo->get_phistep()
0391                   << " diff: " << m_PhiStep - cellgeo->get_phistep() << std::endl;
0392         exit(1);
0393       }
0394       if (m_CellBinning == PHG4CellDefs::etaphibinning || m_CellBinning == PHG4CellDefs::etaslatbinning)
0395       {
0396         if (m_NumEtaBins != cellgeo->get_etabins())
0397         {
0398           std::cout << "mixing number of EtaBins , first layer: "
0399                     << m_NumEtaBins << " layer " << cellgeo->get_layer() << ": "
0400                     << cellgeo->get_etabins() << std::endl;
0401           exit(1);
0402         }
0403         if (fabs(m_EtaMin - cellgeo->get_etamin()) > 1e-9)
0404         {
0405           std::cout << "mixing etamin, fisrt layer: " << m_EtaMin << " layer "
0406                     << cellgeo->get_layer() << ": " << cellgeo->get_etamin()
0407                     << " diff: " << m_EtaMin - cellgeo->get_etamin() << std::endl;
0408           exit(1);
0409         }
0410         if (fabs(m_EtaStep - cellgeo->get_etastep()) > 1e-9)
0411         {
0412           std::cout << "mixing eta steps first layer: " << m_EtaStep
0413                     << " layer " << cellgeo->get_layer() << ": "
0414                     << cellgeo->get_etastep() << " diff: "
0415                     << m_EtaStep - cellgeo->get_etastep() << std::endl;
0416           exit(1);
0417         }
0418       }
0419 
0420       else if (m_CellBinning == PHG4CellDefs::sizebinning)
0421       {
0422         if (m_NumEtaBins != cellgeo->get_zbins())
0423         {
0424           std::cout << "mixing number of z bins , first layer: " << m_NumEtaBins
0425                     << " layer " << cellgeo->get_layer() << ": "
0426                     << cellgeo->get_zbins() << std::endl;
0427           exit(1);
0428         }
0429       }
0430     }
0431   }
0432   m_RawTowerGeomContainer->set_radius(inner_radius);
0433   m_RawTowerGeomContainer->set_thickness(thickness);
0434   m_RawTowerGeomContainer->set_phibins(m_NumPhiBins);
0435   //  m_RawTowerGeomContainer->set_phistep(m_PhiStep);
0436   //  m_RawTowerGeomContainer->set_phimin(m_PhiMin);
0437   m_RawTowerGeomContainer->set_etabins(m_NumEtaBins);
0438 
0439   if (!first_cellgeo)
0440   {
0441     std::cout << "RawTowerBuilder::CreateNodes - ERROR - can not find first layer of cells "
0442               << std::endl;
0443 
0444     exit(1);
0445   }
0446 
0447   for (int ibin = 0; ibin < first_cellgeo->get_phibins(); ibin++)
0448   {
0449     const std::pair<double, double> range = first_cellgeo->get_phibounds(ibin);
0450 
0451     m_RawTowerGeomContainer->set_phibounds(ibin, range);
0452   }
0453 
0454   if (m_CellBinning == PHG4CellDefs::etaphibinning || m_CellBinning == PHG4CellDefs::etaslatbinning || m_CellBinning == PHG4CellDefs::spacalbinning)
0455   {
0456     const double r = inner_radius;
0457 
0458     for (int ibin = 0; ibin < first_cellgeo->get_etabins(); ibin++)
0459     {
0460       const std::pair<double, double> range = first_cellgeo->get_etabounds(ibin);
0461 
0462       m_RawTowerGeomContainer->set_etabounds(ibin, range);
0463     }
0464 
0465     // setup location of all towers
0466     for (int iphi = 0; iphi < m_RawTowerGeomContainer->get_phibins(); iphi++)
0467     {
0468       for (int ieta = 0; ieta < m_RawTowerGeomContainer->get_etabins(); ieta++)
0469       {
0470         const RawTowerDefs::keytype key =
0471             RawTowerDefs::encode_towerid(caloid, ieta, iphi);
0472 
0473         const double x(r * cos(m_RawTowerGeomContainer->get_phicenter(iphi)));
0474         const double y(r * sin(m_RawTowerGeomContainer->get_phicenter(iphi)));
0475         const double z(r / tan(PHG4Utils::get_theta(m_RawTowerGeomContainer->get_etacenter(ieta))));
0476 
0477         RawTowerGeom *tg = m_RawTowerGeomContainer->get_tower_geometry(key);
0478         if (tg)
0479         {
0480           if (Verbosity() > 0)
0481           {
0482             std::cout << "RawTowerBuilder::CreateNodes - Tower geometry " << key << " already exists" << std::endl;
0483           }
0484 
0485           if (fabs(tg->get_center_x() - x) > 1e-4)
0486           {
0487             std::cout << "RawTowerBuilder::CreateNodes - Fatal Error - duplicated Tower geometry " << key << " with existing x = " << tg->get_center_x() << " and expected x = " << x
0488                       << std::endl;
0489 
0490             exit(1);
0491           }
0492           if (fabs(tg->get_center_y() - y) > 1e-4)
0493           {
0494             std::cout << "RawTowerBuilder::CreateNodes - Fatal Error - duplicated Tower geometry " << key << " with existing y = " << tg->get_center_y() << " and expected y = " << y
0495                       << std::endl;
0496             exit(1);
0497           }
0498           if (fabs(tg->get_center_z() - z) > 1e-4)
0499           {
0500             std::cout << "RawTowerBuilder::CreateNodes - Fatal Error - duplicated Tower geometry " << key << " with existing z= " << tg->get_center_z() << " and expected z = " << z
0501                       << std::endl;
0502             exit(1);
0503           }
0504         }
0505         else
0506         {
0507           if (Verbosity() > 0)
0508           {
0509             std::cout << "RawTowerBuilder::CreateNodes - building tower geometry " << key << "" << std::endl;
0510           }
0511 
0512           tg = new RawTowerGeomv1(key);
0513 
0514           tg->set_center_x(x);
0515           tg->set_center_y(y);
0516           tg->set_center_z(z);
0517           m_RawTowerGeomContainer->add_tower_geometry(tg);
0518         }
0519       }
0520     }
0521   }
0522   else if (m_CellBinning == PHG4CellDefs::sizebinning)
0523   {
0524     const double r = inner_radius;
0525 
0526     for (int ibin = 0; ibin < first_cellgeo->get_zbins(); ibin++)
0527     {
0528       const std::pair<double, double> z_range = first_cellgeo->get_zbounds(ibin);
0529       //          const double r = first_cellgeo->get_radius();
0530       const double eta1 = -log(tan(atan2(r, z_range.first) / 2.));
0531       const double eta2 = -log(tan(atan2(r, z_range.second) / 2.));
0532       m_RawTowerGeomContainer->set_etabounds(ibin, std::make_pair(eta1, eta2));
0533     }
0534 
0535     // setup location of all towers
0536     for (int iphi = 0; iphi < m_RawTowerGeomContainer->get_phibins(); iphi++)
0537     {
0538       for (int ibin = 0; ibin < first_cellgeo->get_zbins(); ibin++)
0539       {
0540         const RawTowerDefs::keytype key = RawTowerDefs::encode_towerid(caloid, ibin, iphi);
0541 
0542         const double x(r * cos(m_RawTowerGeomContainer->get_phicenter(iphi)));
0543         const double y(r * sin(m_RawTowerGeomContainer->get_phicenter(iphi)));
0544         const double z(first_cellgeo->get_zcenter(ibin));
0545 
0546         RawTowerGeom *tg = m_RawTowerGeomContainer->get_tower_geometry(key);
0547         if (tg)
0548         {
0549           if (Verbosity() > 0)
0550           {
0551             std::cout << "RawTowerBuilder::CreateNodes - Tower geometry " << key << " already exists" << std::endl;
0552           }
0553 
0554           if (fabs(tg->get_center_x() - x) > 1e-4)
0555           {
0556             std::cout << "RawTowerBuilder::CreateNodes - Fatal Error - duplicated Tower geometry " << key << " with existing x = " << tg->get_center_x() << " and expected x = " << x
0557                       << std::endl;
0558 
0559             exit(1);
0560           }
0561           if (fabs(tg->get_center_y() - y) > 1e-4)
0562           {
0563             std::cout << "RawTowerBuilder::CreateNodes - Fatal Error - duplicated Tower geometry " << key << " with existing y = " << tg->get_center_y() << " and expected y = " << y
0564                       << std::endl;
0565             exit(1);
0566           }
0567           if (fabs(tg->get_center_z() - z) > 1e-4)
0568           {
0569             std::cout << "RawTowerBuilder::CreateNodes - Fatal Error - duplicated Tower geometry " << key << " with existing z= " << tg->get_center_z() << " and expected z = " << z
0570                       << std::endl;
0571             exit(1);
0572           }
0573         }
0574         else
0575         {
0576           if (Verbosity() > 0)
0577           {
0578             std::cout << "RawTowerBuilder::CreateNodes - building tower geometry " << key << "" << std::endl;
0579           }
0580 
0581           tg = new RawTowerGeomv1(key);
0582 
0583           tg->set_center_x(x);
0584           tg->set_center_y(y);
0585           tg->set_center_z(z);
0586           m_RawTowerGeomContainer->add_tower_geometry(tg);
0587         }
0588       }
0589     }
0590   }
0591   else
0592   {
0593     std::cout << "RawTowerBuilder::CreateNodes - ERROR - unsupported cell geometry "
0594               << m_CellBinning << std::endl;
0595     exit(1);
0596   }
0597 
0598   if (Verbosity() >= 1)
0599   {
0600     m_RawTowerGeomContainer->identify();
0601   }
0602 
0603   PHCompositeNode *dstNode = dynamic_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "DST"));
0604   if (!dstNode)
0605   {
0606     std::cout << PHWHERE << "DST Node missing, doing nothing." << std::endl;
0607     throw std::runtime_error(
0608         "Failed to find DST node in RawTowerBuilder::CreateNodes");
0609   }
0610 
0611   PHNodeIterator dstiter(dstNode);
0612   PHCompositeNode *DetNode = dynamic_cast<PHCompositeNode *>(dstiter.findFirst("PHCompositeNode", m_Detector));
0613   if (!DetNode)
0614   {
0615     DetNode = new PHCompositeNode(m_Detector);
0616     dstNode->addNode(DetNode);
0617   }
0618 
0619   
0620   if (m_UseTowerInfo != 1)
0621     {
0622       // Create the tower nodes on the tree
0623       if (m_SimTowerNodePrefix.empty())
0624     {
0625       // no prefix, consistent with older convention
0626       m_TowerNodeName = "TOWER_" + m_Detector;
0627     }
0628       else
0629     {
0630       m_TowerNodeName = "TOWER_" + m_SimTowerNodePrefix + "_" + m_Detector;
0631     }
0632       m_TowerContainer = findNode::getClass<RawTowerContainer>(DetNode, m_TowerNodeName.c_str());
0633       if (m_TowerContainer == nullptr)
0634     {
0635       m_TowerContainer = new RawTowerContainer(caloid);
0636       
0637       PHIODataNode<PHObject> *towerNode = new PHIODataNode<PHObject>(m_TowerContainer, m_TowerNodeName, "PHObject");
0638       DetNode->addNode(towerNode);
0639     }
0640     }
0641 
0642 
0643 
0644   if (m_UseTowerInfo > 0 )
0645     {
0646      if (m_SimTowerNodePrefix.empty())
0647     {
0648       m_TowerInfoNodeName = "TOWERINFO_" + m_Detector;
0649     }
0650       else
0651     {
0652       m_TowerInfoNodeName = "TOWERINFO_" + m_SimTowerNodePrefix + "_" + m_Detector;
0653     }
0654      TowerInfoContainer* m_TowerInfoContainer = findNode::getClass<TowerInfoContainer>(DetNode,m_TowerInfoNodeName);
0655      if (m_TowerInfoContainer == nullptr)
0656        {
0657      TowerInfoContainer::DETECTOR detec;
0658      if (caloid == RawTowerDefs::CalorimeterId::CEMC)
0659        {
0660          detec = TowerInfoContainer::DETECTOR::EMCAL;
0661        }
0662      else if (caloid == RawTowerDefs::CalorimeterId::HCALIN || caloid == RawTowerDefs::CalorimeterId::HCALOUT)
0663        {
0664          detec = TowerInfoContainer::DETECTOR::HCAL;
0665        }
0666      else
0667        {
0668          std::cout << PHWHERE << "Detector not implemented into the TowerInfoContainer object, defaulting to HCal implementation." << std::endl;
0669          detec = TowerInfoContainer::DETECTOR::HCAL;
0670        }
0671      m_TowerInfoContainer = new TowerInfoContainerv1(detec);
0672      PHIODataNode<PHObject> *TowerInfoNode = new PHIODataNode<PHObject>(m_TowerInfoContainer, m_TowerInfoNodeName, "PHObject");
0673      DetNode->addNode(TowerInfoNode);
0674        }
0675     }
0676   
0677   return;
0678 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!