MultiDetectorMultiplicity/src/MultiDetectorMultiplicity.cc

0001 #include "MultiDetectorMultiplicity.h"
0002
0003 // Centrality
0004 #include <centrality/CentralityInfov1.h>
0005
0006 // Calo includes
0007 #include <calobase/TowerInfov1.h>
0008 #include <calobase/TowerInfoContainerv1.h>
0009
0010 #include <phool/getClass.h>
0011 #include <phool/phool.h>
0012
0013 // Vertex includes
0014 #include <globalvertex/GlobalVertex.h>
0015 #include <globalvertex/GlobalVertexMap.h>
0016
0017 // Tracking includes
0018 #include <trackbase_historic/SvtxTrack.h>
0019 #include <trackbase_historic/SvtxTrackMap.h>
0020 #include <trackbase_historic/TrackSeedContainer.h>
0021 #include <trackbase/TrkrHitSet.h>
0022 #include <trackbase/TrkrHitSetContainer.h>
0023
0024 /// Fun4All includes
0025 #include <fun4all/Fun4AllReturnCodes.h>
0026
0027 /// ROOT includes
0028 #include <TFile.h>
0029 #include <TH2F.h>
0030 #include <TH1I.h>
0031 #include <TTree.h>
0032
0033 /// C++ includes
0034 #include <string>
0035
0036 /**
0037  * MultiDetectorMultiplicity is a class developed to cross-check multiplicity and energy correlation in different subsystems
0038  * Author: Antonio Silva (antonio.sphenix@gmail.com)
0039  */
0040
0041 /**
0042  * Constructor of module
0043  */
0044 MultiDetectorMultiplicity::MultiDetectorMultiplicity(const std::string &name, const std::string &filename)
0045   : SubsysReco(name)
0046   , _outfilename(filename)
0047 {
0048   /// Initialize variables and trees so we don't accidentally access
0049   /// memory that was never allocated
0050 }
0051
0052 /**
0053  * Destructor of module
0054  */
0055
0056 MultiDetectorMultiplicity::~MultiDetectorMultiplicity()
0057 {
0058
0059 }
0060
0061 /**
0062  * Initialize the module and prepare looping over events
0063  */
0064 int MultiDetectorMultiplicity::Init(PHCompositeNode *topNode)
0065 {
0066   if (Verbosity() > 5)
0067   {
0068     std::cout << "Beginning Init in MultiDetectorMultiplicity" << std::endl;
0069   }
0070
0071   initializeObjects();
0072
0073   return 0;
0074 }
0075
0076 /**
0077  * Main workhorse function where each event is looped over and
0078  * data from each event is collected from the node tree for analysis
0079  */
0080 int MultiDetectorMultiplicity::process_event(PHCompositeNode *topNode)
0081 {
0082   if (Verbosity() > 5)
0083   {
0084     std::cout << "Beginning process_event in MultiDetectorMultiplicity" << std::endl;
0085   }
0086
0087   float Ntracks = getFilteredNtracks(topNode);
0088   float MVTXclusters = getMVTXnclusters(topNode);
0089   float EMCalTotalEnergy = getEMCalTotalEnergy(topNode);
0090   float IHCalTotalEnergy = getIHCalTotalEnergy(topNode);
0091   float OHCalTotalEnergy = getOHCalTotalEnergy(topNode);
0092   float EPDTotalEnergy = getsEPDTotalEnergy(topNode);
0093
0094   bool acceptVertex = isVertexAccepted(topNode);
0095
0096   CentralityInfo* cent_node = findNode::getClass<CentralityInfo>(topNode, "CentralityInfo");
0097   if (!cent_node)
0098   {
0099     std::cout << "MultiDetectorMultiplicity::process_event - Error can not find centrality node " << std::endl;
0100     _SubsystemObjectsFail->Fill(OBJECT::CENTRALITY);
0101     return Fun4AllReturnCodes::ABORTEVENT;
0102   }
0103
0104   if(!acceptVertex) return Fun4AllReturnCodes::ABORTEVENT;
0105
0106   _mbd_centrality = cent_node->get_centile(CentralityInfo::PROP::mbd_NS);
0107   _epd_centrality = cent_node->get_centile(CentralityInfo::PROP::epd_NS);
0108
0109   _tracks_EMCALenergy->Fill(Ntracks, EMCalTotalEnergy);
0110   _tracks_IHCALenergy->Fill(Ntracks, IHCalTotalEnergy);
0111   _tracks_OHCALenergy->Fill(Ntracks, OHCalTotalEnergy);
0112   _tracks_AllCALenergy->Fill(Ntracks, EMCalTotalEnergy+IHCalTotalEnergy+OHCalTotalEnergy);
0113   _tracks_MVTXclusters->Fill(Ntracks, MVTXclusters);
0114
0115   _MVTXclusters_EMCALenergy->Fill(MVTXclusters, EMCalTotalEnergy);
0116   _MVTXclusters_IHCALenergy->Fill(MVTXclusters, IHCalTotalEnergy);
0117   _MVTXclusters_OHCALenergy->Fill(MVTXclusters, OHCalTotalEnergy);
0118   _MVTXclusters_AllCALenergy->Fill(MVTXclusters, EMCalTotalEnergy+IHCalTotalEnergy+OHCalTotalEnergy);
0119
0120   _EMCALenergy_IHCALenergy->Fill(EMCalTotalEnergy, IHCalTotalEnergy);
0121   _EMCALenergy_OHCALenergy->Fill(EMCalTotalEnergy, OHCalTotalEnergy);
0122   _IHCALenergy_OHCALenergy->Fill(IHCalTotalEnergy, OHCalTotalEnergy);
0123
0124   float centBinMBD = 100.-_mbd_centrality;
0125
0126   _MBDcentrality_tracks->Fill(centBinMBD, Ntracks);
0127   _MBDcentrality_MVTXclusters->Fill(centBinMBD, MVTXclusters);
0128   _MBDcentrality_EMCALenergy->Fill(centBinMBD, EMCalTotalEnergy);
0129   _MBDcentrality_IHCALenergy->Fill(centBinMBD, IHCalTotalEnergy);
0130   _MBDcentrality_OHCALenergy->Fill(centBinMBD, OHCalTotalEnergy);
0131   _MBDcentrality_AllCALenergy->Fill(centBinMBD, EMCalTotalEnergy+IHCalTotalEnergy+OHCalTotalEnergy);
0132   _MBDcentrality_EPDenergy->Fill(centBinMBD, EPDTotalEnergy);
0133
0134   float centBinEPD = 100.-_epd_centrality;
0135
0136   _EPDcentrality_tracks->Fill(centBinEPD, Ntracks);
0137   _EPDcentrality_MVTXclusters->Fill(centBinEPD, MVTXclusters);
0138   _EPDcentrality_EMCALenergy->Fill(centBinEPD, EMCalTotalEnergy);
0139   _EPDcentrality_IHCALenergy->Fill(centBinEPD, IHCalTotalEnergy);
0140   _EPDcentrality_OHCALenergy->Fill(centBinEPD, OHCalTotalEnergy);
0141   _EPDcentrality_AllCALenergy->Fill(centBinEPD, EMCalTotalEnergy+IHCalTotalEnergy+OHCalTotalEnergy);
0142   _EPDcentrality_EPDenergy->Fill(centBinEPD, EPDTotalEnergy);
0143
0144   _EPDcentrality_MBDcentrality->Fill(centBinEPD, centBinMBD);
0145
0146   _EPDenergy_tracks->Fill(EPDTotalEnergy, Ntracks);
0147   _EPDenergy_MVTXclusters->Fill(EPDTotalEnergy, MVTXclusters);
0148   _EPDenergy_EMCALenergy->Fill(EPDTotalEnergy, EMCalTotalEnergy);
0149   _EPDenergy_IHCALenergy->Fill(EPDTotalEnergy, IHCalTotalEnergy);
0150   _EPDenergy_OHCALenergy->Fill(EPDTotalEnergy, OHCalTotalEnergy);
0151   _EPDenergy_AllCALenergy->Fill(EPDTotalEnergy, EMCalTotalEnergy+IHCalTotalEnergy+OHCalTotalEnergy);
0152
0153   return Fun4AllReturnCodes::EVENT_OK;
0154 }
0155
0156 /**
0157  * End the module and finish any data collection. Clean up any remaining
0158  * loose ends
0159  */
0160 int MultiDetectorMultiplicity::End(PHCompositeNode *topNode)
0161 {
0162   if (Verbosity() > 1)
0163   {
0164     std::cout << "Ending MultiDetectorMultiplicity analysis package" << std::endl;
0165   }
0166
0167   /// Change to the outfile
0168   _outfile->cd();
0169
0170   /// Write and close the outfile
0171   _tracks_EMCALenergy->Write();
0172   _tracks_IHCALenergy->Write();
0173   _tracks_OHCALenergy->Write();
0174   _tracks_AllCALenergy->Write();
0175   _tracks_MVTXclusters->Write();
0176
0177   _MVTXclusters_EMCALenergy->Write();
0178   _MVTXclusters_IHCALenergy->Write();
0179   _MVTXclusters_OHCALenergy->Write();
0180   _MVTXclusters_AllCALenergy->Write();
0181
0182   _EMCALenergy_IHCALenergy->Write();
0183   _EMCALenergy_OHCALenergy->Write();
0184   _IHCALenergy_OHCALenergy->Write();
0185
0186   _MBDcentrality_tracks->Write();
0187   _MBDcentrality_MVTXclusters->Write();
0188   _MBDcentrality_EMCALenergy->Write();
0189   _MBDcentrality_IHCALenergy->Write();
0190   _MBDcentrality_OHCALenergy->Write();
0191   _MBDcentrality_AllCALenergy->Write();
0192   _MBDcentrality_EPDenergy->Write();
0193
0194   _EPDcentrality_tracks->Write();
0195   _EPDcentrality_MVTXclusters->Write();
0196   _EPDcentrality_EMCALenergy->Write();
0197   _EPDcentrality_IHCALenergy->Write();
0198   _EPDcentrality_OHCALenergy->Write();
0199   _EPDcentrality_AllCALenergy->Write();
0200   _EPDcentrality_EPDenergy->Write();
0201
0202   _EPDcentrality_MBDcentrality->Write();
0203
0204   _EPDenergy_tracks->Write();
0205   _EPDenergy_MVTXclusters->Write();
0206   _EPDenergy_EMCALenergy->Write();
0207   _EPDenergy_IHCALenergy->Write();
0208   _EPDenergy_OHCALenergy->Write();
0209   _EPDenergy_AllCALenergy->Write();
0210
0211   _SubsystemObjectsFail->Write();
0212
0213   _outfile->Close();
0214
0215   /// Clean up pointers and associated histos/trees in TFile
0216   delete _outfile;
0217
0218   if (Verbosity() > 1)
0219   {
0220     std::cout << "Finished MultiDetectorMultiplicity analysis package" << std::endl;
0221   }
0222
0223   return 0;
0224 }
0225
0226 void MultiDetectorMultiplicity::initializeObjects()
0227 {
0228   _outfile = new TFile(_outfilename.c_str(), "RECREATE");
0229
0230   _tracks_EMCALenergy = new TH2F("tracks_EMCALenergy", ";Track multiplicity;EMCal Energy (GeV)",_nbins,0.,_track_bin_max,_nbins,0.,_emcal_bin_max);
0231   _tracks_IHCALenergy = new TH2F("tracks_IHCALenergy", ";Track multiplicity;IHCal Energy (GeV)",_nbins,0.,_track_bin_max,_nbins,0.,_ihcal_bin_max);
0232   _tracks_OHCALenergy = new TH2F("tracks_OHCALenergy", ";Track multiplicity;OHCal Energy (GeV)",_nbins,0.,_track_bin_max,_nbins,0.,_ohcal_bin_max);
0233   _tracks_AllCALenergy = new TH2F("tracks_AllCALenergy", ";Track multiplicity;EMCal+IHCal+OHCal Energy (GeV)",_nbins,0.,_track_bin_max,_nbins,0.,_allcal_bin_max);
0234   _tracks_MVTXclusters = new TH2F("tracks_MVTXclusters", ";Track multiplicity;MVTX cluster multiplicity",_nbins,0.,_track_bin_max,_nbins,0.,_mvtx_bin_max);
0235
0236   _MVTXclusters_EMCALenergy = new TH2F("MVTXclusters_EMCALenergy", ";MVTX cluster multiplicity;EMCal Energy (GeV)",_nbins,0.,_mvtx_bin_max,_nbins,0.,_emcal_bin_max);
0237   _MVTXclusters_IHCALenergy = new TH2F("MVTXclusters_IHCALenergy", ";MVTX cluster multiplicity;IHCal Energy (GeV)",_nbins,0.,_mvtx_bin_max,_nbins,0.,_ihcal_bin_max);
0238   _MVTXclusters_OHCALenergy = new TH2F("MVTXclusters_OHCALenergy", ";MVTX cluster multiplicity;OHCal Energy (GeV)",_nbins,0.,_mvtx_bin_max,_nbins,0.,_ohcal_bin_max);
0239   _MVTXclusters_AllCALenergy = new TH2F("MVTXclusters_AllCALenergy", ";MVTX cluster multiplicity;EMCal+IHCal+OHCal Energy (GeV)",_nbins,0.,_mvtx_bin_max,_nbins,0.,_allcal_bin_max);
0240
0241   _EMCALenergy_IHCALenergy = new TH2F("EMCALenergy_IHCALenergy", ";EMCal Energy (GeV);IHCal Energy (GeV)",_nbins,0.,_emcal_bin_max,_nbins,0.,_ihcal_bin_max);
0242   _EMCALenergy_OHCALenergy = new TH2F("EMCALenergy_OHCALenergy", ";EMCal Energy (GeV);OHCal Energy (GeV)",_nbins,0.,_emcal_bin_max,_nbins,0.,_ohcal_bin_max);
0243   _IHCALenergy_OHCALenergy = new TH2F("IHCALenergy_OHCALenergy", ";IHCal Energy (GeV);OHCal Energy (GeV)",_nbins,0.,_ihcal_bin_max,_nbins,0.,_ohcal_bin_max);
0244
0245   _MBDcentrality_tracks = new TH2F("MBDcentrality_tracks", ";MBD Centrality;Track multiplicity",10,0.,100.,_nbins,0.,_track_bin_max);
0246   setCentralityHistoLabel(_MBDcentrality_tracks);
0247
0248   _MBDcentrality_MVTXclusters = new TH2F("MBDcentrality_MVTXclusters", ";MBD Centrality;MVTX cluster multiplicity",10,0.,100.,_nbins,0.,_mvtx_bin_max);
0249   setCentralityHistoLabel(_MBDcentrality_MVTXclusters);
0250
0251   _MBDcentrality_EMCALenergy = new TH2F("MBDcentrality_EMCALenergy", ";MBD Centrality;EMCal Energy (GeV)",10,0.,100.,_nbins,0.,_emcal_bin_max);
0252   setCentralityHistoLabel(_MBDcentrality_EMCALenergy);
0253
0254   _MBDcentrality_IHCALenergy = new TH2F("MBDcentrality_IHCALenergy", ";MBD Centrality;IHCal Energy (GeV)",10,0.,100.,_nbins,0.,_ihcal_bin_max);
0255   setCentralityHistoLabel(_MBDcentrality_IHCALenergy);
0256
0257   _MBDcentrality_OHCALenergy = new TH2F("MBDcentrality_OHCALenergy", ";MBD Centrality;OHCal Energy (GeV)",10,0.,100.,_nbins,0.,_ohcal_bin_max);
0258   setCentralityHistoLabel(_MBDcentrality_OHCALenergy);
0259
0260   _MBDcentrality_AllCALenergy = new TH2F("MBDcentrality_AllCALenergy", ";MBD Centrality;EMCal+IHCal+OHCal Energy (GeV)",10,0.,100.,_nbins,0.,_allcal_bin_max);
0261   setCentralityHistoLabel(_MBDcentrality_AllCALenergy);
0262
0263   _MBDcentrality_EPDenergy = new TH2F("MBDcentrality_EPDenergy", ";MBD Centrality;EPD Energy (GeV)",10,0.,100.,_nbins,0.,_epd_bin_max);
0264   setCentralityHistoLabel(_MBDcentrality_EPDenergy);
0265
0266   _EPDcentrality_tracks = new TH2F("EPDcentrality_tracks", ";EPD Centrality;Track multiplicity",10,0.,100.,_nbins,0.,_track_bin_max);
0267   setCentralityHistoLabel(_EPDcentrality_tracks);
0268
0269   _EPDcentrality_MVTXclusters = new TH2F("EPDcentrality_MVTXclusters", ";EPD Centrality;MVTX cluster multiplicity",10,0.,100.,_nbins,0.,_mvtx_bin_max);
0270   setCentralityHistoLabel(_EPDcentrality_MVTXclusters);
0271
0272   _EPDcentrality_EMCALenergy = new TH2F("EPDcentrality_EMCALenergy", ";EPD Centrality;EMCal Energy (GeV)",10,0.,100.,_nbins,0.,_emcal_bin_max);
0273   setCentralityHistoLabel(_EPDcentrality_EMCALenergy);
0274
0275   _EPDcentrality_IHCALenergy = new TH2F("EPDcentrality_IHCALenergy", ";EPD Centrality;IHCal Energy (GeV)",10,0.,100.,_nbins,0.,_ihcal_bin_max);
0276   setCentralityHistoLabel(_EPDcentrality_IHCALenergy);
0277
0278   _EPDcentrality_OHCALenergy = new TH2F("EPDcentrality_OHCALenergy", ";EPD Centrality;OHCal Energy (GeV)",10,0.,100.,_nbins,0.,_ohcal_bin_max);
0279   setCentralityHistoLabel(_EPDcentrality_OHCALenergy);
0280
0281   _EPDcentrality_AllCALenergy = new TH2F("EPDcentrality_AllCALenergy", ";EPD Centrality;EMCal+IHCal+OHCal Energy (GeV)",10,0.,100.,_nbins,0.,_allcal_bin_max);
0282   setCentralityHistoLabel(_EPDcentrality_AllCALenergy);
0283
0284   _EPDcentrality_EPDenergy = new TH2F("EPDcentrality_EPDenergy", ";EPD Centrality;EPD Energy (GeV)",10,0.,100.,_nbins,0.,_epd_bin_max);
0285   setCentralityHistoLabel(_EPDcentrality_EPDenergy);
0286
0287   _EPDcentrality_MBDcentrality = new TH2F("EPDcentrality_MBDcentrality", ";EPD Centrality;MBD Centrality",10,0.,100.,10,0.,100);
0288   setCentralityHistoLabel(_EPDcentrality_MBDcentrality, true);
0289
0290   _EPDenergy_tracks = new TH2F("EPDenergy_tracks", ";EPD Energy (GeV);Track multiplicity",_nbins,0.,_epd_bin_max,_nbins,0.,_track_bin_max);
0291   _EPDenergy_MVTXclusters = new TH2F("EPDenergy_MVTXclusters", ";EPD Energy (GeV);MVTX cluster multiplicity",_nbins,0.,_epd_bin_max,_nbins,0.,_mvtx_bin_max);
0292   _EPDenergy_EMCALenergy = new TH2F("EPDenergy_EMCALenergy", ";EPD Energy (GeV);EMCal Energy (GeV)",_nbins,0.,_epd_bin_max,_nbins,0.,_emcal_bin_max);
0293   _EPDenergy_IHCALenergy = new TH2F("EPDenergy_IHCALenergy", ";EPD Energy (GeV);IHCal Energy (GeV)",_nbins,0.,_epd_bin_max,_nbins,0.,_ihcal_bin_max);
0294   _EPDenergy_OHCALenergy = new TH2F("EPDenergy_OHCALenergy", ";EPD Energy (GeV);OHCal Energy (GeV)",_nbins,0.,_epd_bin_max,_nbins,0.,_ohcal_bin_max);
0295   _EPDenergy_AllCALenergy = new TH2F("EPDenergy_AllCALenergy", ";EPD Energy (GeV);EMCal+IHCal+OHCal Energy (GeV)",_nbins,0.,_epd_bin_max,_nbins,0.,_allcal_bin_max);
0296
0297   // This histo keeps track of the number of times certain objects were not found.
0298   // Ideally this histogram should be EMPTY!!!
0299   // For future updates including other objects, be sure to not reject the event before all objects are checked
0300   _SubsystemObjectsFail = new TH1I("SubsystemObjectsFail", "Times objects were NOT found;;Entries",8,0.,8.);
0301
0302   _SubsystemObjectsFail->GetXaxis()->SetBinLabel(1,"Vertex");
0303   _SubsystemObjectsFail->GetXaxis()->SetBinLabel(2,"Centrality");
0304   _SubsystemObjectsFail->GetXaxis()->SetBinLabel(3,"Tracks");
0305   _SubsystemObjectsFail->GetXaxis()->SetBinLabel(4,"MVTX clusters");
0306   _SubsystemObjectsFail->GetXaxis()->SetBinLabel(5,"EMCal towers");
0307   _SubsystemObjectsFail->GetXaxis()->SetBinLabel(6,"IHCal towers");
0308   _SubsystemObjectsFail->GetXaxis()->SetBinLabel(7,"OHCal towers");
0309   _SubsystemObjectsFail->GetXaxis()->SetBinLabel(8,"sEPD towers");
0310
0311 }
0312
0313 float MultiDetectorMultiplicity::getFilteredNtracks(PHCompositeNode *topNode)
0314 {
0315   SvtxTrackMap *trackmap = findNode::getClass<SvtxTrackMap>(topNode, "SvtxTrackMap");
0316
0317   if (!trackmap)
0318   {
0319     std::cout << PHWHERE
0320          << "SvtxTrackMap node is missing, can't collect tracks"
0321          << std::endl;
0322     _SubsystemObjectsFail->Fill(OBJECT::TRACKS);
0323     return -1;
0324   }
0325
0326   float Ntracks = 0;
0327
0328   for(SvtxTrackMap::Iter iter = trackmap->begin(); iter != trackmap->end(); ++iter)
0329     {
0330       SvtxTrack* track = iter->second;
0331     float quality = track->get_quality();
0332
0333     if (quality > 10) continue;
0334
0335     // DCA units and values have to be cross-checked
0336     //float DCAxy = track->get_dca3d_xy(); //DCAxy in micrometer
0337     //if(std::abs(DCAxy) > 20) continue;
0338     //float DCAz = track->get_dca3d_z() * 10000; //DCAz in cm?, so multipli by 10^4
0339     //if(std::abs(DCAz) > 10) continue;
0340
0341     auto silicon_seed = track->get_silicon_seed();
0342     auto tpc_seed = track->get_tpc_seed();
0343     int nsiliconhits = 0;
0344
0345     int nTPChits = 0;
0346
0347     // Getting number of TPC hits as done in SvtxEvaluator.h
0348
0349     if(tpc_seed)
0350     {
0351       for (TrackSeed::ConstClusterKeyIter iter = tpc_seed->begin_cluster_keys(); iter != tpc_seed->end_cluster_keys(); ++iter)
0352       {
0353         TrkrDefs::cluskey cluster_key = *iter;
0354         unsigned int layer = TrkrDefs::getLayer(cluster_key);
0355         if (_nlayers_tpc > 0 && layer >= (_nlayers_maps + _nlayers_intt) && layer < (_nlayers_maps + _nlayers_intt + _nlayers_tpc))
0356         {
0357           nTPChits++;
0358         }
0359       }
0360     }
0361
0362     if (silicon_seed)
0363     {
0364       for (TrackSeed::ConstClusterKeyIter iter = silicon_seed->begin_cluster_keys(); iter != silicon_seed->end_cluster_keys(); ++iter)
0365       {
0366         TrkrDefs::cluskey cluster_key = *iter;
0367         unsigned int layer = TrkrDefs::getLayer(cluster_key);
0368         if (_nlayers_tpc > 0 && layer >= (_nlayers_maps + _nlayers_intt) && layer < (_nlayers_maps + _nlayers_intt + _nlayers_tpc))
0369         {
0370           nTPChits++;
0371         }
0372       }
0373       nsiliconhits = silicon_seed->size_cluster_keys(); // Get number of hits on silicon detectors (INTT+MVTX)
0374     }
0375
0376     if(nTPChits < 20) continue; // Require at least 20 TPC hits
0377
0378     if (nsiliconhits < 2) continue; // Require at least 2 hits on silicon detectors
0379
0380     Ntracks += 1.;
0381   }
0382
0383   return Ntracks;
0384 }
0385
0386 float MultiDetectorMultiplicity::getMVTXnclusters(PHCompositeNode *topNode)
0387 {
0388   TrkrClusterContainer *MVTXclusterMap = findNode::getClass<TrkrClusterContainer>(topNode,"TRKR_CLUSTER");
0389
0390   if(!MVTXclusterMap)
0391   {
0392     std::cout << PHWHERE
0393          << "TRKR_CLUSTER node is missing, can't collect hits"
0394          << std::endl;
0395     _SubsystemObjectsFail->Fill(OBJECT::MVTXCLUSTERS);
0396     return -1;
0397   }
0398
0399   float MVTXclusters = 0.;
0400
0401   for(const auto& hitsetkey : MVTXclusterMap->getHitSetKeys(TrkrDefs::TrkrId::mvtxId))
0402   {
0403     auto range = MVTXclusterMap->getClusters(hitsetkey);
0404     for( auto clusIter = range.first; clusIter != range.second; ++clusIter )
0405     {
0406       const auto cluster = clusIter->second;
0407
0408       if(!cluster) continue;
0409
0410       MVTXclusters += 1.;
0411
0412     }
0413   }
0414
0415   return MVTXclusters;
0416
0417 }
0418
0419 float MultiDetectorMultiplicity::getEMCalTotalEnergy(PHCompositeNode *topNode)
0420 {
0421   TowerInfoContainer *towerinfosEM = findNode::getClass<TowerInfoContainer>(topNode, "TOWERINFO_CALIB_CEMC");
0422
0423   if(!towerinfosEM)
0424   {
0425     std::cout << PHWHERE
0426          << "TOWERINFO_CALIB_CEMC node is missing, can't collect EMCal towers"
0427          << std::endl;
0428     _SubsystemObjectsFail->Fill(OBJECT::EMCALTOWERS);
0429     return -1.;
0430   }
0431
0432   float EMCalTotalEnergy = 0.;
0433
0434   TowerInfo *towerInfo = nullptr;
0435   for(unsigned int i = 0; i < towerinfosEM->size(); i++)
0436   {
0437     towerInfo = towerinfosEM->get_tower_at_channel(i);
0438     EMCalTotalEnergy += towerInfo->get_energy();
0439   }
0440
0441   return EMCalTotalEnergy;
0442 }
0443
0444 float MultiDetectorMultiplicity::getIHCalTotalEnergy(PHCompositeNode *topNode)
0445 {
0446   TowerInfoContainer *towerinfosIH = findNode::getClass<TowerInfoContainer>(topNode, "TOWERINFO_CALIB_HCALIN");
0447
0448   if(!towerinfosIH)
0449   {
0450     std::cout << PHWHERE
0451          << "TOWERINFO_CALIB_HCALIN node is missing, can't collect IHCal towers"
0452          << std::endl;
0453     _SubsystemObjectsFail->Fill(OBJECT::IHCALTOWERS);
0454     return -1.;
0455   }
0456
0457   float IHCalTotalEnergy = 0.;
0458
0459   TowerInfo *towerInfo = nullptr;
0460   for(unsigned int i = 0; i < towerinfosIH->size(); i++)
0461   {
0462     towerInfo = towerinfosIH->get_tower_at_channel(i);
0463     IHCalTotalEnergy += towerInfo->get_energy();
0464   }
0465
0466   return IHCalTotalEnergy;
0467 }
0468
0469 float MultiDetectorMultiplicity::getOHCalTotalEnergy(PHCompositeNode *topNode)
0470 {
0471   TowerInfoContainer *towerinfosOH = findNode::getClass<TowerInfoContainer>(topNode, "TOWERINFO_CALIB_HCALOUT");
0472
0473   if(!towerinfosOH)
0474   {
0475     std::cout << PHWHERE
0476          << "TOWERINFO_CALIB_HCALOUT node is missing, can't collect OHCal towers"
0477          << std::endl;
0478     _SubsystemObjectsFail->Fill(OBJECT::OHCALTOWERS);
0479     return -1.;
0480   }
0481
0482   float OHCalTotalEnergy = 0.;
0483
0484   TowerInfo *towerInfo = nullptr;
0485   for(unsigned int i = 0; i < towerinfosOH->size(); i++)
0486   {
0487     towerInfo = towerinfosOH->get_tower_at_channel(i);
0488     OHCalTotalEnergy += towerInfo->get_energy();
0489   }
0490
0491   return OHCalTotalEnergy;
0492 }
0493
0494 float MultiDetectorMultiplicity::getsEPDTotalEnergy(PHCompositeNode *topNode)
0495 {
0496   TowerInfoContainer *towerinfosEPD = findNode::getClass<TowerInfoContainer>(topNode, "TOWERINFO_CALIB_EPD");
0497
0498   if(!towerinfosEPD)
0499   {
0500     std::cout << PHWHERE
0501          << "TOWERINFO_CALIB_EPD node is missing, can't collect sEPD towers"
0502          << std::endl;
0503     _SubsystemObjectsFail->Fill(OBJECT::EPDTOWERS);
0504     return -1.;
0505   }
0506
0507   float EPDTotalEnergy = 0.;
0508
0509   TowerInfo *towerInfo = nullptr;
0510   for(unsigned int i = 0; i < towerinfosEPD->size(); i++)
0511   {
0512     towerInfo = towerinfosEPD->get_tower_at_channel(i);
0513     EPDTotalEnergy += towerInfo->get_energy();
0514   }
0515
0516   return EPDTotalEnergy;
0517 }
0518
0519 void MultiDetectorMultiplicity::setCentralityHistoLabel(TH2 *histo, bool doYaxisLabels)
0520 {
0521   for(int i = 0; i < 10; i++)
0522   {
0523     histo->GetXaxis()->SetBinLabel(i+1,Form("%d-%d%%",90-(i*10),100-(i*10)));
0524     if(doYaxisLabels) histo->GetYaxis()->SetBinLabel(i+1,Form("%d-%d%%",90-(i*10),100-(i*10)));
0525   }
0526 }
0527
0528 bool MultiDetectorMultiplicity::isVertexAccepted(PHCompositeNode *topNode)
0529 {
0530   GlobalVertexMap *vertexmap = findNode::getClass<GlobalVertexMap>(topNode, "GlobalVertexMap");
0531   if (!vertexmap)
0532   {
0533     std::cout << "GlobalVertexMap node is missing. Event aborted!" << std::endl;
0534     _SubsystemObjectsFail->Fill(OBJECT::VERTEX);
0535     return false;
0536   }
0537
0538   if (vertexmap->empty())
0539   {
0540     std::cout << "GlobalVertexMap node is empty. Event aborted!" << std::endl;
0541     _SubsystemObjectsFail->Fill(OBJECT::VERTEX);
0542     return false;
0543   }
0544
0545   GlobalVertex *vtx = vertexmap->begin()->second;
0546   if (vtx == nullptr)
0547   {
0548     _SubsystemObjectsFail->Fill(OBJECT::VERTEX);
0549     return false;
0550   }
0551
0552   if(std::abs(vtx->get_z()) > _zVertexAcceptance) return false;
0553
0554   return true;
0555 }