packages/CaloReco/CaloTowerStatus.cc

0001 #include "CaloTowerStatus.h"
0002 #include "CaloTowerDefs.h"
0003
0004 #include <calobase/TowerInfo.h>  // for TowerInfo
0005 #include <calobase/TowerInfoContainer.h>
0006 #include <calobase/TowerInfoContainerv1.h>
0007 #include <calobase/TowerInfoContainerv2.h>
0008 #include <calobase/TowerInfov1.h>
0009 #include <calobase/TowerInfov2.h>
0010
0011 #include <cdbobjects/CDBTTree.h>  // for CDBTTree
0012
0013 #include <ffamodules/CDBInterface.h>
0014
0015 #include <ffaobjects/EventHeader.h>
0016
0017 #include <fun4all/Fun4AllReturnCodes.h>
0018 #include <fun4all/SubsysReco.h>  // for SubsysReco
0019
0020 #include <phool/PHCompositeNode.h>
0021 #include <phool/PHIODataNode.h>    // for PHIODataNode
0022 #include <phool/PHNode.h>          // for PHNode
0023 #include <phool/PHNodeIterator.h>  // for PHNodeIterator
0024 #include <phool/PHObject.h>        // for PHObject
0025 #include <phool/getClass.h>
0026 #include <phool/phool.h>
0027 #include <phool/recoConsts.h>
0028
0029 #include <TSystem.h>
0030
0031 #include <cmath>
0032 #include <cstdlib>    // for exit
0033 #include <exception>  // for exception
0034 #include <iostream>   // for operator<<, basic_ostream
0035 #include <stdexcept>  // for runtime_error
0036
0037 //____________________________________________________________________________..
0038 CaloTowerStatus::CaloTowerStatus(const std::string &name)
0039   : SubsysReco(name)
0040 {
0041   if (Verbosity() > 0)
0042   {
0043     std::cout << "CaloTowerStatus::CaloTowerStatus(const std::string &name) Calling ctor" << std::endl;
0044   }
0045 }
0046
0047 //____________________________________________________________________________..
0048 CaloTowerStatus::~CaloTowerStatus()
0049 {
0050   if (Verbosity() > 0)
0051   {
0052     std::cout << "CaloTowerStatus::~CaloTowerStatus() Calling dtor" << std::endl;
0053   }
0054   delete m_cdbttree_chi2;
0055   delete m_cdbttree_time;
0056   delete m_cdbttree_hotMap;
0057 }
0058
0059 //____________________________________________________________________________..
0060 int CaloTowerStatus::InitRun(PHCompositeNode *topNode)
0061 {
0062   PHNodeIterator nodeIter(topNode);
0063
0064   if (m_dettype == CaloTowerDefs::CEMC)
0065   {
0066     m_detector = "CEMC";
0067   }
0068   else if (m_dettype == CaloTowerDefs::HCALIN)
0069   {
0070     m_detector = "HCALIN";
0071   }
0072   else if (m_dettype == CaloTowerDefs::HCALOUT)
0073   {
0074     m_detector = "HCALOUT";
0075   }
0076   else if (m_dettype == CaloTowerDefs::ZDC)
0077   {
0078     m_detector = "ZDC";
0079   }
0080
0081   else if (m_dettype == CaloTowerDefs::SEPD)
0082   {
0083     m_detector = "SEPD";
0084   }
0085
0086   m_calibName_chi2 = m_detector + "_hotTowers_fracBadChi2";
0087   m_fieldname_chi2 = "fraction";
0088
0089   std::string calibdir = CDBInterface::instance()->getUrl(m_calibName_chi2);
0090   if (!calibdir.empty())
0091   {
0092      m_cdbttree_chi2 = new CDBTTree(calibdir);
0093     if (Verbosity() > 0)
0094     {
0095        std::cout << "CaloTowerStatus::InitRun Found " << m_calibName_chi2 << "  Doing isHot for frac bad chi2" << std::endl;
0096     }
0097   }
0098   else
0099   {
0100     if (use_directURL_chi2)
0101     {
0102       calibdir = m_directURL_chi2;
0103       std::cout << "CaloTowerStatus::InitRun: Using default hotBadChi2" << calibdir << std::endl;
0104       m_cdbttree_chi2 = new CDBTTree(calibdir);
0105     }
0106     else
0107     {
0108       m_doHotChi2 = false;
0109       if (Verbosity() > 0)
0110       {
0111         std::cout << "CaloTowerStatus::InitRun No masking file for domain " << m_calibName_chi2 << " found, not doing isHot from isBadChi2" << std::endl;
0112       }
0113     }
0114   }
0115
0116   m_calibName_time = m_detector + "_meanTime";
0117   m_fieldname_time = "time";
0118
0119   calibdir = CDBInterface::instance()->getUrl(m_calibName_time);
0120   if (!calibdir.empty())
0121   {
0122     m_cdbttree_time = new CDBTTree(calibdir);
0123     if (Verbosity() > 0)
0124     {
0125       std::cout << "CaloTowerStatus::InitRun Found " << m_calibName_time << " not Doing isBadTime" << std::endl;
0126     }
0127   }
0128   else
0129   {
0130     if (use_directURL_time)
0131     {
0132       calibdir = m_directURL_time;
0133       std::cout << "CaloTowerStatus::InitRun: Using default time  " << calibdir << std::endl;
0134       m_cdbttree_time = new CDBTTree(calibdir);
0135     }
0136     else
0137     {
0138       m_doTime = false;
0139       if (Verbosity() > 1)
0140       {
0141         std::cout << "CaloTowerStatus::InitRun no timing info, " << m_calibName_time << " not found, not doing isBadTime" << std::endl;
0142       }
0143     }
0144   }
0145
0146   m_calibName_hotMap = m_detector + "nome";
0147   if (m_dettype == CaloTowerDefs::CEMC)
0148   {
0149     m_calibName_hotMap = m_detector + "_BadTowerMap";
0150   }
0151   m_fieldname_hotMap = "status";
0152   m_fieldname_z_score = m_detector + "_sigma";
0153
0154   calibdir = CDBInterface::instance()->getUrl(m_calibName_hotMap);
0155   if (!calibdir.empty())
0156   {
0157     m_cdbttree_hotMap = new CDBTTree(calibdir);
0158     if (Verbosity() > 1)
0159     {
0160       std::cout << "CaloTowerStatus::Init " << m_detector << "  hot map found " << m_calibName_hotMap << " Ddoing isHot" << std::endl;
0161     }
0162   }
0163   else
0164   {
0165     if (m_doAbortNoHotMap)
0166     {
0167       std::cout << "CaloTowerStatus::InitRun: No hot map.. exiting" << std::endl;
0168       gSystem->Exit(1);
0169     }
0170     if (use_directURL_hotMap)
0171     {
0172       calibdir = m_directURL_hotMap;
0173       std::cout << "CaloTowerStatus::InitRun: Using default map " << calibdir << std::endl;
0174       m_cdbttree_hotMap = new CDBTTree(calibdir);
0175     }
0176     else
0177     {
0178       m_doHotMap = false;
0179       if (Verbosity() > 1)
0180       {
0181         std::cout << "CaloTowerStatus::InitRun hot map info, " << m_calibName_hotMap << " not found, not doing isHot" << std::endl;
0182       }
0183     }
0184   }
0185
0186   if (Verbosity() > 0)
0187   {
0188     std::cout << "CaloTowerStatus::Init " << m_detector << "  doing time status =" <<  std::boolalpha << m_doTime << "  doing hotBadChi2=" <<  std::boolalpha << m_doHotChi2 << " doing hot map=" << std::boolalpha << m_doHotMap << std::endl;
0189   }
0190
0191   PHNodeIterator iter(topNode);
0192
0193   // Looking for the DST node
0194   PHCompositeNode *dstNode;
0195   dstNode = dynamic_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "DST"));
0196   if (!dstNode)
0197   {
0198     std::cout << Name() << "::" << m_detector << "::" << __PRETTY_FUNCTION__
0199               << "DST Node missing, doing nothing." << std::endl;
0200     exit(1);
0201   }
0202   try
0203   {
0204     CreateNodeTree(topNode);
0205     LoadCalib();
0206   }
0207   catch (std::exception &e)
0208   {
0209     std::cout << e.what() << std::endl;
0210     return Fun4AllReturnCodes::ABORTRUN;
0211   }
0212   if (Verbosity() > 0)
0213   {
0214     topNode->print();
0215   }
0216   return Fun4AllReturnCodes::EVENT_OK;
0217 }
0218
0219 void CaloTowerStatus::LoadCalib()
0220 {
0221   unsigned int ntowers = m_raw_towers->size();
0222   m_cdbInfo_vec.resize(ntowers);
0223
0224   for (unsigned int channel = 0; channel < ntowers; channel++)
0225   {
0226     unsigned int key = m_raw_towers->encode_key(channel);
0227
0228     if (m_doHotChi2)
0229     {
0230       m_cdbInfo_vec[channel].fraction_badChi2 = m_cdbttree_chi2->GetFloatValue(key, m_fieldname_chi2);
0231     }
0232     if (m_doTime)
0233     {
0234       m_cdbInfo_vec[channel].mean_time = m_cdbttree_time->GetFloatValue(key, m_fieldname_time);
0235     }
0236     if (m_doHotMap)
0237     {
0238       m_cdbInfo_vec[channel].hotMap_val = m_cdbttree_hotMap->GetIntValue(key, m_fieldname_hotMap);
0239       m_cdbInfo_vec[channel].z_score = m_cdbttree_hotMap->GetFloatValue(key, m_fieldname_z_score);
0240     }
0241   }
0242 }
0243
0244 //____________________________________________________________________________..
0245 int CaloTowerStatus::process_event(PHCompositeNode * /*topNode*/)
0246 {
0247   unsigned int ntowers = m_raw_towers->size();
0248   float fraction_badChi2 = 0;
0249   float mean_time = 0;
0250   int hotMap_val = 0;
0251   float z_score = 0;
0252   for (unsigned int channel = 0; channel < ntowers; channel++)
0253   {
0254     // only reset what we will set
0255     m_raw_towers->get_tower_at_channel(channel)->set_isHot(false);
0256     m_raw_towers->get_tower_at_channel(channel)->set_isBadTime(false);
0257     m_raw_towers->get_tower_at_channel(channel)->set_isBadChi2(false);
0258
0259     if (m_doHotChi2)
0260     {
0261       fraction_badChi2 = m_cdbInfo_vec[channel].fraction_badChi2;
0262     }
0263     if (m_doTime)
0264     {
0265       mean_time = m_cdbInfo_vec[channel].mean_time;
0266     }
0267     if (m_doHotMap)
0268     {
0269       hotMap_val = m_cdbInfo_vec[channel].hotMap_val;
0270       z_score = m_cdbInfo_vec[channel].z_score;
0271     }
0272     float chi2 = m_raw_towers->get_tower_at_channel(channel)->get_chi2();
0273     float time = m_raw_towers->get_tower_at_channel(channel)->get_time();
0274     float adc = m_raw_towers->get_tower_at_channel(channel)->get_energy();
0275
0276     if (fraction_badChi2 > fraction_badChi2_threshold && m_doHotChi2)
0277     {
0278       m_raw_towers->get_tower_at_channel(channel)->set_isHot(true);
0279     }
0280     if (!m_raw_towers->get_tower_at_channel(channel)->get_isZS() && std::fabs(time - mean_time) > time_cut && m_doTime)
0281     {
0282       m_raw_towers->get_tower_at_channel(channel)->set_isBadTime(true);
0283     }
0284     if (( hotMap_val == 1 || // dead
0285           std::fabs(z_score) > z_score_threshold || // hot or cold
0286           (hotMap_val == 3 && z_score >= -1 * z_score_threshold_default)) // cold part 2
0287           && m_doHotMap)
0288     {
0289       m_raw_towers->get_tower_at_channel(channel)->set_isHot(true);
0290     }
0291     if (chi2 > std::min(std::max(badChi2_treshold_const, adc * adc * badChi2_treshold_quadratic),badChi2_treshold_max))
0292     {
0293       m_raw_towers->get_tower_at_channel(channel)->set_isBadChi2(true);
0294     }
0295   }
0296   return Fun4AllReturnCodes::EVENT_OK;
0297 }
0298
0299 void CaloTowerStatus::CreateNodeTree(PHCompositeNode *topNode)
0300 {
0301   std::string RawTowerNodeName = m_inputNodePrefix + m_detector;
0302   m_raw_towers = findNode::getClass<TowerInfoContainer>(topNode, RawTowerNodeName);
0303   if (!m_raw_towers)
0304   {
0305     std::cout << Name() << "::" << m_detector.c_str() << "::" << __PRETTY_FUNCTION__
0306               << " " << RawTowerNodeName << " Node missing, doing bail out!"
0307               << std::endl;
0308     throw std::runtime_error(
0309         "Failed to find " + RawTowerNodeName + " node in CaloTowerStatus::CreateNodes");
0310   }
0311
0312   return;
0313 }