Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-12-17 09:21:29 

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