g4simulation/g4main/PHG4HitContainer.cc

0001 #include "PHG4HitContainer.h"
0002
0003 #include "PHG4Hit.h"
0004 #include "PHG4Hitv1.h"
0005
0006 #include <phool/phool.h>
0007
0008 #include <TSystem.h>
0009
0010 #include <cstdlib>
0011
0012 using namespace std;
0013
0014 PHG4HitContainer::PHG4HitContainer()
0015   : id(-1)
0016   , hitmap()
0017   , layers()
0018 {
0019 }
0020
0021 PHG4HitContainer::PHG4HitContainer(const std::string &nodename)
0022   : id(PHG4HitDefs::get_volume_id(nodename))
0023   , hitmap()
0024   , layers()
0025 {
0026 }
0027
0028 void PHG4HitContainer::Reset()
0029 {
0030   while (hitmap.begin() != hitmap.end())
0031   {
0032     delete hitmap.begin()->second;
0033     hitmap.erase(hitmap.begin());
0034   }
0035   return;
0036 }
0037
0038 void PHG4HitContainer::identify(ostream &os) const
0039 {
0040   ConstIterator iter;
0041   os << "Number of hits: " << size() << endl;
0042   for (iter = hitmap.begin(); iter != hitmap.end(); ++iter)
0043   {
0044     os << "hit key 0x" << hex << iter->first << dec << endl;
0045     (iter->second)->identify();
0046   }
0047   set<unsigned int>::const_iterator siter;
0048   os << "Number of layers: " << num_layers() << endl;
0049   for (siter = layers.begin(); siter != layers.end(); ++siter)
0050   {
0051     os << "layer : " << *siter << endl;
0052   }
0053   return;
0054 }
0055
0056 PHG4HitDefs::keytype
0057 PHG4HitContainer::getmaxkey(const unsigned int detid)
0058 {
0059   ConstRange miter = getHits(detid);
0060   // first handle 2 special cases where there is no hit in the current layer
0061   // no hits in this layer and higher layers (lower layers can contain hits)
0062   if (miter.first == hitmap.end())
0063   {
0064     return 0;
0065   }
0066   // no hits in this layer - but hits in higher layers
0067   if (miter.first == miter.second)
0068   {
0069     return 0;
0070   }
0071   PHG4HitDefs::keytype detidlong = detid;
0072   PHG4HitDefs::keytype shiftval = detidlong << PHG4HitDefs::hit_idbits;
0073   ConstIterator lastlayerentry = miter.second;
0074   --lastlayerentry;
0075   PHG4HitDefs::keytype iret = lastlayerentry->first - shiftval;  // subtract layer mask
0076   return iret;
0077 }
0078
0079 PHG4HitDefs::keytype
0080 PHG4HitContainer::genkey(const unsigned int detid)
0081 {
0082   PHG4HitDefs::keytype detidlong = detid;
0083   if ((detidlong >> PHG4HitDefs::keybits) > 0)
0084   {
0085     cout << PHWHERE << " detector id too large: " << detid << endl;
0086     gSystem->Exit(1);
0087   }
0088   PHG4HitDefs::keytype shiftval = detidlong << PHG4HitDefs::hit_idbits;
0089   //  cout << "max index: " << (detminmax->second)->first << endl;
0090   // after removing hits with no energy deposition, we have holes
0091   // in our hit ranges. This construct will get us the last hit in
0092   // a layer and return it's hit id. Adding 1 will put us at the end of this layer
0093   PHG4HitDefs::keytype hitid = getmaxkey(detid);
0094   hitid++;
0095   PHG4HitDefs::keytype newkey = hitid | shiftval;
0096   if (hitmap.find(newkey) != hitmap.end())
0097   {
0098     cout << PHWHERE << " duplicate key: 0x"
0099          << hex << newkey << dec
0100          << " for detector " << detid
0101          << " hitmap.size: " << hitmap.size()
0102          << " hitid: " << hitid << " exiting now" << endl;
0103     exit(1);
0104   }
0105   return newkey;
0106 }
0107
0108 PHG4HitContainer::ConstIterator
0109 PHG4HitContainer::AddHit(PHG4Hit *newhit)
0110 {
0111   PHG4HitDefs::keytype key = newhit->get_hit_id();
0112   if (hitmap.find(key) != hitmap.end())
0113   {
0114     cout << "hit with id  0x" << hex << key << dec << " exists already" << endl;
0115     return hitmap.find(key);
0116   }
0117   PHG4HitDefs::keytype detidlong = key >> PHG4HitDefs::hit_idbits;
0118   unsigned int detid = detidlong;
0119   layers.insert(detid);
0120   return hitmap.insert(std::make_pair(key, newhit)).first;
0121 }
0122
0123 PHG4HitContainer::ConstIterator
0124 PHG4HitContainer::AddHit(const unsigned int detid, PHG4Hit *newhit)
0125 {
0126   PHG4HitDefs::keytype key = genkey(detid);
0127   layers.insert(detid);
0128   newhit->set_hit_id(key);
0129   return hitmap.insert(std::make_pair(key, newhit)).first;
0130 }
0131
0132 PHG4HitContainer::ConstRange PHG4HitContainer::getHits(const unsigned int detid) const
0133 {
0134   PHG4HitDefs::keytype detidlong = detid;
0135   if ((detidlong >> PHG4HitDefs::keybits) > 0)
0136   {
0137     cout << " detector id too large: " << detid << endl;
0138     exit(1);
0139   }
0140   PHG4HitDefs::keytype keylow = detidlong << PHG4HitDefs::hit_idbits;
0141   PHG4HitDefs::keytype keyup = ((detidlong + 1) << PHG4HitDefs::hit_idbits) - 1;
0142   ConstRange retpair;
0143   retpair.first = hitmap.lower_bound(keylow);
0144   retpair.second = hitmap.upper_bound(keyup);
0145   return retpair;
0146 }
0147
0148 PHG4HitContainer::ConstRange PHG4HitContainer::getHits() const
0149 {
0150   return std::make_pair(hitmap.begin(), hitmap.end());
0151 }
0152
0153 PHG4HitContainer::Iterator PHG4HitContainer::findOrAddHit(PHG4HitDefs::keytype key)
0154 {
0155   PHG4HitContainer::Iterator it = hitmap.find(key);
0156   if (it == hitmap.end())
0157   {
0158     hitmap[key] = new PHG4Hitv1();
0159     it = hitmap.find(key);
0160     PHG4Hit *mhit = it->second;
0161     mhit->set_hit_id(key);
0162     mhit->set_edep(0.);
0163     layers.insert(mhit->get_layer());  // add layer to our set of layers
0164   }
0165   return it;
0166 }
0167
0168 PHG4Hit *PHG4HitContainer::findHit(PHG4HitDefs::keytype key)
0169 {
0170   PHG4HitContainer::ConstIterator it = hitmap.find(key);
0171   if (it != hitmap.end())
0172   {
0173     return it->second;
0174   }
0175
0176   return nullptr;
0177 }
0178
0179 void PHG4HitContainer::RemoveZeroEDep()
0180 {
0181   //  unsigned int hitsbef = hitmap.size();
0182   Iterator itr = hitmap.begin();
0183   Iterator last = hitmap.end();
0184   for (; itr != last;)
0185   {
0186     PHG4Hit *hit = itr->second;
0187     if (hit->get_edep() == 0)
0188     {
0189       delete hit;
0190       hitmap.erase(itr++);
0191     }
0192     else
0193     {
0194       ++itr;
0195     }
0196   }
0197   //   unsigned int hitsafter = hitmap.size();
0198   //   cout << "hist before: " << hitsbef
0199   //        << ", hits after: " << hitsafter << endl;
0200   return;
0201 }