Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:20:23

0001 // use #include "" only for your local include and put
0002 // those in the first line(s) before any #include <>
0003 // otherwise you are asking for weird behavior
0004 // (more info - check the difference in include path search when using "" versus <>)
0005 
0006 #include "CemcMon.h"
0007 
0008 #include <onlmon/OnlMon.h>  // for OnlMon
0009 #include <onlmon/OnlMonServer.h>
0010 #include <onlmon/pseudoRunningMean.h>
0011 #include <onlmon/triggerEnum.h>
0012 
0013 #include <calobase/TowerInfoDefs.h>
0014 #include <caloreco/CaloWaveformFitting.h>
0015 
0016 #include <cdbobjects/CDBTTree.h>
0017 
0018 #include <Event/Event.h>
0019 #include <Event/EventTypes.h>
0020 #include <Event/eventReceiverClient.h>
0021 #include <Event/msg_profile.h>
0022 
0023 #include <TH1.h>
0024 #include <TH2.h>
0025 #include <TProfile.h>
0026 #include <TProfile2D.h>
0027 
0028 #include <cmath>
0029 #include <cstdio>  // for printf
0030 #include <fstream>
0031 #include <iostream>
0032 #include <sstream>
0033 #include <string>  // for allocator, string, char_traits
0034 
0035 enum
0036 {
0037   TRGMESSAGE = 1,
0038   FILLMESSAGE = 2
0039 };
0040 
0041 const int depth = 50000;
0042 // const int historyLength = 100;
0043 // const float hit_threshold = 100;
0044 const float hit_threshold = 100;
0045 const float waveform_hit_threshold = 100;
0046 const float chi2_check_threshold = 3000;
0047 
0048 using namespace std;
0049 
0050 CemcMon::CemcMon(const std::string &name)
0051   : OnlMon(name)
0052   , eventCounter(0)
0053 {
0054   // leave ctor fairly empty, its hard to debug if code crashes already
0055   // during a new CemcMon()
0056 
0057   return;
0058 }
0059 
0060 CemcMon::~CemcMon()
0061 {
0062   for (auto iter : rm_vector_twr)
0063   {
0064     delete iter;
0065   }
0066   for (auto iter : rm_vector_twrhits)
0067   {
0068     delete iter;
0069   }
0070   for (auto iter : rm_vector_twrhits_alltrig)
0071   {
0072     delete iter;
0073   }
0074 
0075   delete WaveformProcessingFast;
0076   delete WaveformProcessingTemp;
0077   delete erc;
0078   return;
0079 }
0080 
0081 int CemcMon::Init()
0082 {
0083   // read our calibrations from CemcMonData.dat
0084   const char *cemccalib = getenv("CEMCCALIB");
0085   if (!cemccalib)
0086   {
0087     std::cout << "CEMCCALIB environment variable not set" << std::endl;
0088     exit(1);
0089   }
0090   std::string fullfile = std::string(cemccalib) + "/" + "CemcMonData.dat";
0091   std::ifstream calib(fullfile);
0092   calib.close();
0093   // use printf for stuff which should go the screen but not into the message
0094   // system (all couts are redirected)
0095   printf("CemcMon::Init()\n");
0096   std::string cdbfile = std::string(cemccalib) + "/cemc_adcskipmask_40430.root";
0097   cdbttree = new CDBTTree(cdbfile.c_str());
0098   // Histograms definitions
0099   // Trigger histograms
0100   for (int itrig = 0; itrig < 64; itrig++)
0101   {
0102     h2_cemc_hits_trig[itrig] = new TH2F(Form("h2_cemc_hits_trig_bit_%d", itrig), "", 96, 0, 96, 256, 0, 256);
0103   }
0104   p2_zsFrac_etaphi = new TProfile2D("p2_zsFrac_etaphi", "", 96, 0, 96, 256, 0, 256);
0105   p2_zsFrac_etaphi_all = new TProfile2D("p2_zsFrac_etaphi_all", "", 96, 0, 96, 256, 0, 256);
0106   h1_cemc_trig = new TH1F("h1_cemc_trig", "", 64, -0.5, 63.5);
0107   h1_packet_event = new TH1F("h1_packet_event", "", 8, packetlow - 0.5, packethigh + 0.5);
0108   h2_caloPack_gl1_clock_diff = new TH2F("h2_caloPack_gl1_clock_diff", "", 8, packetlow - 0.5, packethigh + 0.5, 65536, 0, 65536);
0109   h_evtRec = new TProfile("h_evtRec", "", 1, 0, 1);
0110 
0111   // tower hit information
0112   h2_cemc_hits = new TH2F("h2_cemc_hits", "", 96, 0, 96, 256, 0, 256);
0113   h2_cemc_rm = new TH2F("h2_cemc_rm", "", 96, 0, 96, 256, 0, 256);
0114   h2_cemc_rmhits = new TH2F("h2_cemc_rmhits", "", 96, 0, 96, 256, 0, 256);
0115   h2_cemc_rmhits_alltrig = new TH2F("h2_cemc_rmhits_alltrig", "", 96, 0, 96, 256, 0, 256);
0116   h2_cemc_mean = new TH2F("h2_cemc_mean", "", 96, 0, 96, 256, 0, 256);
0117   h1_cemc_adc = new TH1F("h1_cemc_adc", "", 1000, 0.5, pow(2, 14));
0118   // event counter
0119   h1_event = new TH1F("h1_event", "", 1, 0, 1);
0120 
0121   // waveform processing
0122   // h2_waveform_twrAvg = new TH2F("h2_waveform_twrAvg", "", 16, 0.5, 16.5, 10000,0,pow(2,14));
0123   h2_waveform_twrAvg = new TH2F("h2_waveform_twrAvg", "", 12, -0.5, 11.5, 1000, 0, 15000);
0124   h1_waveform_time = new TH1F("h1_waveform_time", "", 12, -0.5, 11.5);
0125   h1_waveform_pedestal = new TH1F("h1_waveform_pedestal", "", 5000, 1, 5000);
0126 
0127   // waveform processing, template vs. fast interpolation
0128   h1_cemc_fitting_sigDiff = new TH1F("h1_fitting_sigDiff", "", 50, 0, 2);
0129   h1_cemc_fitting_pedDiff = new TH1F("h1_fitting_pedDiff", "", 50, 0, 2);
0130   h1_cemc_fitting_timeDiff = new TH1F("h1_fitting_timeDiff", "", 50, -10, 10);
0131 
0132   // packet information
0133   h1_packet_number = new TH1F("h1_packet_number", "", 128, 6000.5, 6128.5);
0134   h1_packet_length = new TH1F("h1_packet_length", "", 128, 6000.5, 6128.5);
0135   h1_packet_chans = new TH1F("h1_packet_chans", "", 128, 6000.5, 6128.5);
0136 
0137   p2_bad_chi2 = new TProfile2D("p2_bad_chi2", "", 96, 0, 96, 256, 0, 256);
0138   //s option to track rms
0139   p2_pre_post = new TProfile2D("p2_pre_post", "", 96, 0, 96, 256, 0, 256, "S");
0140 
0141   // make the per-packet running mean objects
0142   // 32 packets and 48 channels for hcal detectors
0143   for (int i = 0; i < Ntower; i++)
0144   {
0145     rm_vector_twr.push_back(new pseudoRunningMean(1, depth));
0146     rm_vector_twrhits.push_back(new pseudoRunningMean(1, depth));
0147     rm_vector_twrhits_alltrig.push_back(new pseudoRunningMean(1, depth));
0148   }
0149   
0150 
0151 
0152   OnlMonServer *se = OnlMonServer::instance();
0153   // register histograms with server otherwise client won't get them
0154 
0155   // Trigger histograms
0156   for (int itrig = 0; itrig < 64; itrig++)
0157   {
0158     se->registerHisto(this, h2_cemc_hits_trig[itrig]);
0159   }
0160   se->registerHisto(this, p2_zsFrac_etaphi);
0161   se->registerHisto(this, p2_zsFrac_etaphi_all);
0162   se->registerHisto(this, p2_bad_chi2);
0163   se->registerHisto(this, h1_cemc_trig);
0164   se->registerHisto(this, h1_packet_event);
0165   se->registerHisto(this, h2_caloPack_gl1_clock_diff);
0166   se->registerHisto(this, h_evtRec);
0167 
0168   se->registerHisto(this, h2_cemc_hits);
0169   se->registerHisto(this, h2_cemc_rm);
0170   se->registerHisto(this, h2_cemc_rmhits);
0171   se->registerHisto(this, h2_cemc_rmhits_alltrig);
0172   se->registerHisto(this, p2_pre_post);
0173   se->registerHisto(this, h2_cemc_mean);
0174   se->registerHisto(this, h1_event);
0175 
0176   se->registerHisto(this, h2_waveform_twrAvg);
0177   se->registerHisto(this, h1_waveform_time);
0178   se->registerHisto(this, h1_waveform_pedestal);
0179   se->registerHisto(this, h1_cemc_fitting_sigDiff);
0180   se->registerHisto(this, h1_cemc_fitting_pedDiff);
0181   se->registerHisto(this, h1_cemc_fitting_timeDiff);
0182   se->registerHisto(this, h1_packet_number);
0183   se->registerHisto(this, h1_packet_length);
0184   se->registerHisto(this, h1_packet_chans);
0185   se->registerHisto(this, h1_cemc_adc);
0186 
0187   // Commented until potential replacement with TProfile3D
0188   // h2_waveform=new TProfile**[nPhiIndex];
0189   // for(int iphi=0; iphi<nPhiIndex; iphi++){
0190   //   h2_waveform[iphi]=new TProfile*[nEtaIndex];
0191   //   for(int ieta=0; ieta<nEtaIndex; ieta++){
0192   //     h2_waveform[iphi][ieta]=new TProfile(Form("h2_waveform_phi%d_eta%d",iphi,ieta),Form("Profiled raw waveform for #phi %d and #eta %d",iphi,ieta),12, -0.5, 11.5, "s");
0193   //     h2_waveform[iphi][ieta]->GetXaxis()->SetTitle("sample #");
0194   //     h2_waveform[iphi][ieta]->GetYaxis()->SetTitle("ADC counts");
0195   //     h2_waveform[iphi][ieta]->SetStats(false);
0196   //     se->registerHisto(this, (TH1*)h2_waveform[iphi][ieta]);
0197   //   }
0198   // }
0199 
0200   // initialize waveform extraction tool
0201   WaveformProcessingFast = new CaloWaveformFitting();
0202 
0203   WaveformProcessingTemp = new CaloWaveformFitting();
0204 
0205   std::string cemctemplate;
0206   if (getenv("CEMCCALIB"))
0207   {
0208     cemctemplate = getenv("CEMCCALIB");
0209   }
0210   else
0211   {
0212     cemctemplate = ".";
0213   }
0214   cemctemplate += std::string("/testbeam_cemc_template.root");
0215   WaveformProcessingTemp->initialize_processing(cemctemplate);
0216 
0217   if (anaGL1)
0218   {
0219     erc = new eventReceiverClient("gl1daq");
0220   }
0221 
0222   return 0;
0223 }
0224 
0225 int CemcMon::BeginRun(const int /* runno */)
0226 {
0227   // if you need to read calibrations on a run by run basis
0228   // this is the place to do it
0229 
0230   // reset the running means
0231   std::vector<runningMean *>::iterator rm_it;
0232   for (rm_it = rm_vector_twr.begin(); rm_it != rm_vector_twr.end(); ++rm_it)
0233   {
0234     (*rm_it)->Reset();
0235   }
0236   for (rm_it = rm_vector_twrhits.begin(); rm_it != rm_vector_twrhits.end(); ++rm_it)
0237   {
0238     (*rm_it)->Reset();
0239   }
0240   for (rm_it = rm_vector_twrhits_alltrig.begin(); rm_it != rm_vector_twrhits_alltrig.end(); ++rm_it)
0241   {
0242     (*rm_it)->Reset();
0243   }
0244   if (anaGL1)
0245   {
0246     OnlMonServer *se = OnlMonServer::instance();
0247     se->UseGl1();
0248   }
0249   return 0;
0250 }
0251 
0252 // simple wavefrom analysis for possibe issues with the wavforProcessor
0253 std::vector<float> CemcMon::getSignal(Packet *p, const int channel)
0254 {
0255   double baseline = 0;
0256   for (int s = 0; s < 3; s++)
0257   {
0258     baseline += p->iValue(s, channel);
0259   }
0260   baseline /= 3.;
0261 
0262   double signal = 0;
0263   float x = 0;
0264   for (int s = 3; s < p->iValue(0, "SAMPLES"); s++)
0265   {
0266     x++;
0267     signal += p->iValue(s, channel) - baseline;
0268   }
0269 
0270   signal /= x;
0271 
0272   std::vector<float> result;
0273   result.push_back(signal);
0274   result.push_back(2);
0275   result.push_back(1);
0276   return result;
0277 }
0278 
0279 std::vector<float> CemcMon::anaWaveformFast(Packet *p, const int channel)
0280 {
0281   std::vector<float> waveform;
0282 
0283   // int nSamples = p->iValue(0, "SAMPLES");
0284   if (p->iValue(channel, "SUPPRESSED"))
0285   {
0286     waveform.push_back(p->iValue(channel, "PRE"));
0287     waveform.push_back(p->iValue(channel, "POST"));
0288   }
0289   else
0290   {
0291     int nSamples = p->iValue(0, "SAMPLES");
0292     waveform.reserve(nSamples);
0293     for (int s = 0; s < nSamples; s++)
0294     {
0295       waveform.push_back(p->iValue(s, channel));
0296     }
0297   }
0298 
0299   std::vector<std::vector<float>> multiple_wfs;
0300   multiple_wfs.push_back(waveform);
0301 
0302   std::vector<std::vector<float>> fitresults_cemc;
0303   fitresults_cemc = WaveformProcessingFast->calo_processing_fast(multiple_wfs);
0304 
0305   std::vector<float> result;
0306   result = fitresults_cemc.at(0);
0307   return result;
0308 }
0309 
0310 std::vector<float> CemcMon::anaWaveformTemp(Packet *p, const int channel)
0311 {
0312   std::vector<float> waveform;
0313 
0314   if (p->iValue(channel, "SUPPRESSED"))
0315   {
0316     waveform.push_back(p->iValue(channel, "PRE"));
0317     waveform.push_back(p->iValue(channel, "POST"));
0318   }
0319   else
0320   {
0321     waveform.reserve(p->iValue(0, "SAMPLES"));
0322     for (int s = 0; s < p->iValue(0, "SAMPLES"); s++)
0323     {
0324       waveform.push_back(p->iValue(s, channel));
0325     }
0326   }
0327   std::vector<std::vector<float>> multiple_wfs;
0328   multiple_wfs.push_back(waveform);
0329 
0330   std::vector<std::vector<float>> fitresults_cemc;
0331   fitresults_cemc = WaveformProcessingTemp->process_waveform(multiple_wfs);
0332 
0333   std::vector<float> result;
0334   result = fitresults_cemc.at(0);
0335   return result;
0336 }
0337 
0338 int CemcMon::process_event(Event *e /* evt */)
0339 {
0340   float sectorAvg[Nsector] = {0};
0341   unsigned int towerNumber = 0;
0342   bool fillhist = true;
0343   std::vector<bool> trig_bools;
0344   trig_bools.resize(64);
0345   long long int gl1_clock = 0;
0346   bool have_gl1 = false;
0347   if (anaGL1)
0348   {
0349     int evtnr = e->getEvtSequence();
0350     Event *gl1Event = erc->getEvent(evtnr);
0351     if (gl1Event)
0352     {
0353       OnlMonServer *se = OnlMonServer::instance();
0354       se->IncrementGl1FoundCounter();
0355       have_gl1 = true;
0356       Packet *p = gl1Event->getPacket(14001);
0357       h_evtRec->Fill(0.0, 1.0);
0358       if (p)
0359       {
0360         gl1_clock = p->lValue(0, "BCO");
0361         uint64_t triggervec = p->lValue(0, "ScaledVector");
0362         for (int i = 0; i < 64; i++)
0363         {
0364           bool trig_decision = ((triggervec & 0x1U) == 0x1U);
0365           trig_bools[i] = trig_decision;
0366           if (trig_decision)
0367           {
0368             h1_cemc_trig->Fill(i);
0369           }
0370           triggervec = (triggervec >> 1U);
0371         }
0372         delete p;
0373       }
0374       delete gl1Event;
0375     }
0376     else
0377     {
0378       std::cout << "GL1 event is null" << std::endl;
0379       h_evtRec->Fill(0.0, 0.0);
0380     }
0381     
0382     //this is for only process event with the MBD>=2 trigger
0383     //use MBD>=2 for AuAu
0384     if(usembdtrig){
0385       if(trig_bools.at(TriggerEnum::BitCodes::MBD_NS2_ZVRTX10) == 0 && trig_bools.at(TriggerEnum::BitCodes::MBD_NS2_ZVRTX30) == 0 && trig_bools.at(TriggerEnum::BitCodes::MBD_NS2_ZVRTX150) == 0){
0386         fillhist = false;
0387       }
0388     }
0389     
0390   }
0391 
0392   // loop over packets which contain a single sector
0393   eventCounter++;
0394   int one = 1;
0395   int zero = 0;
0396   for (int packet = packetlow; packet <= packethigh; packet++)
0397   {
0398     Packet *p = e->getPacket(packet);
0399 
0400     if (p)
0401     {
0402       unsigned int adc_skip_mask = cdbttree->GetIntValue(packet, "adcskipmask");
0403       h1_packet_number->Fill(packet);
0404       h1_packet_length->SetBinContent(packet - 6000, h1_packet_length->GetBinContent(packet - 6000) + p->getLength());
0405 
0406       h1_packet_event->SetBinContent(packet - 6000, p->lValue(0, "CLOCK"));
0407 
0408       if (have_gl1)
0409       {
0410         long long int p_clock = p->lValue(0, "CLOCK");
0411         long long int diff = (p_clock - gl1_clock) % 65536;
0412         h2_caloPack_gl1_clock_diff->Fill(packet, diff);
0413       }
0414       int nChannels = p->iValue(0, "CHANNELS");
0415       int skiped_channel = 0;
0416       if (nChannels > m_nChannels)
0417       {
0418         return -1;  // packet is corrupted, reports too many channels
0419       }
0420       //print packet and nCHannels
0421       for (int c = 0; c < nChannels; c++)
0422       {
0423         h1_packet_chans->Fill(packet);
0424         // skip masked ADC board 
0425         if(c % 64 ==0){
0426           unsigned int adcboard = (unsigned int) c / 64;
0427           if ((adc_skip_mask >> adcboard) & 0x1U){
0428             
0429              towerNumber += 64;
0430              skiped_channel += 64;
0431           }
0432         }
0433 
0434         
0435         towerNumber++;
0436         // channel mapping
0437         unsigned int key = TowerInfoDefs::encode_emcal(towerNumber - 1);
0438         unsigned int phi_bin = TowerInfoDefs::getCaloTowerPhiBin(key);
0439         unsigned int eta_bin = TowerInfoDefs::getCaloTowerEtaBin(key);
0440 
0441         // Commented until potential replacement by TProfile3D
0442         // for (int s = 0; s < p->iValue(0, "SAMPLES"); s++)
0443         //   {
0444         //     h2_waveform[phi_bin][eta_bin]->Fill(s,p->iValue(s,c));//for the moment only for good packet and with signal (potentially also bad packet later, not sure for zero suppressed)
0445         //   }
0446 
0447         ////Uninstrumented area
0448         // if ((packet==6019)||(packet==6073)){
0449         //   if(c>63&&c<128) continue;
0450         // }
0451         // if (packet==6030){
0452         //   if(c>127)continue;
0453         // }
0454 
0455         std::vector<float> resultFast = anaWaveformFast(p, c);  // fast waveform fitting
0456         float signalFast = resultFast.at(0);
0457         float timeFast = resultFast.at(1);
0458         float pedestalFast = resultFast.at(2);
0459         int bin = h2_cemc_mean->FindBin(eta_bin + 0.5, phi_bin + 0.5);
0460         float prepost = p->iValue(c, "PRE") - p->iValue(c, "POST");
0461         //________________________________for this part we only want to deal with the MBD>=1 trigger
0462         if (fillhist)
0463         {
0464           if (p->iValue(c, "SUPPRESSED"))
0465           {
0466             p2_zsFrac_etaphi->Fill(eta_bin, phi_bin, 0);
0467           }
0468           else
0469           {
0470             p2_zsFrac_etaphi->Fill(eta_bin, phi_bin, 1);
0471           }
0472 
0473           h1_waveform_pedestal->Fill(pedestalFast);
0474           
0475           rm_vector_twr[towerNumber - 1]->Add(&signalFast);
0476 
0477           if (signalFast > hit_threshold)
0478           {
0479             rm_vector_twrhits[towerNumber - 1]->Add(&one);
0480             h2_cemc_hits->SetBinContent(bin, h2_cemc_hits->GetBinContent(bin) + 1);
0481           }
0482           else
0483           {
0484             rm_vector_twrhits[towerNumber - 1]->Add(&zero);
0485           }
0486           h2_cemc_mean->SetBinContent(bin, h2_cemc_mean->GetBinContent(bin) + signalFast);
0487           h2_cemc_rm->SetBinContent(bin, rm_vector_twr[towerNumber - 1]->getMean(0));
0488           h2_cemc_rmhits->SetBinContent(bin, rm_vector_twrhits[towerNumber - 1]->getMean(0));
0489           h1_cemc_adc->Fill(signalFast);
0490         }
0491         //_______________________________________________________end of MBD trigger requirement
0492 
0493          if (p->iValue(c, "SUPPRESSED"))
0494           {
0495             p2_zsFrac_etaphi_all->Fill(eta_bin, phi_bin, 0);
0496           }
0497           else
0498           {
0499             p2_zsFrac_etaphi_all->Fill(eta_bin, phi_bin, 1);
0500           }
0501 
0502 
0503 
0504         if (signalFast > waveform_hit_threshold)
0505         {
0506           // check if hot tower and skip it
0507           if (!isHottower(packet, c))
0508           {
0509             for (int s = 0; s < p->iValue(0, "SAMPLES"); s++)
0510             {
0511               h2_waveform_twrAvg->Fill(s, p->iValue(s, c) - pedestalFast);
0512             }
0513             h1_waveform_time->Fill(timeFast);
0514           }
0515         }
0516         if (signalFast > hit_threshold)
0517         {
0518           // h2_cemc_hits->SetBinContent(bin, h2_cemc_hits->GetBinContent(bin) + signalFast);
0519           if (have_gl1)
0520           {
0521             for (int itrig = 0; itrig < 64; itrig++)
0522             {
0523               if (trig_bools[itrig])
0524               {
0525                 h2_cemc_hits_trig[itrig]->Fill(eta_bin + 0.5, phi_bin + 0.5);
0526               }
0527             }
0528           }
0529           rm_vector_twrhits_alltrig[towerNumber - 1]->Add(&one);
0530         }
0531         else
0532         {
0533           rm_vector_twrhits_alltrig[towerNumber - 1]->Add(&zero);
0534         }
0535         if(prepost > 0)
0536         {
0537           p2_pre_post->Fill(eta_bin, phi_bin, prepost);
0538           p2_pre_post->Fill(eta_bin, phi_bin, -prepost);
0539         }
0540         
0541         h2_cemc_rmhits_alltrig->SetBinContent(bin, rm_vector_twrhits_alltrig[towerNumber - 1]->getMean(0));
0542 
0543         if (signalFast > chi2_check_threshold)
0544         {
0545           std::vector<float> resultTemp = anaWaveformTemp(p, c);  // template waveform fitting
0546           float chi2 = resultTemp.at(3);
0547           float signalTemp = resultTemp.at(0);
0548           if(chi2 > signalTemp*signalTemp / 50. )
0549           {
0550             p2_bad_chi2->Fill(eta_bin, phi_bin, 1);
0551           }
0552           else
0553           {
0554             p2_bad_chi2->Fill(eta_bin, phi_bin, 0);
0555           }
0556         }
0557         
0558 
0559 
0560         /*
0561         if (!((eventCounter - 2) % 10000))
0562         {
0563           std::vector<float> resultTemp = anaWaveformTemp(p, c);  // template waveform fitting
0564           float signalTemp = resultTemp.at(0);
0565           float timeTemp = resultTemp.at(1);
0566           float pedestalTemp = resultTemp.at(2);
0567           h1_cemc_fitting_sigDiff->Fill(signalFast / signalTemp);
0568           h1_cemc_fitting_pedDiff->Fill(pedestalFast / pedestalTemp);
0569           h1_cemc_fitting_timeDiff->Fill(timeFast - timeTemp - 6);
0570         }
0571         */
0572 
0573       }  // channel loop
0574       if ((nChannels + skiped_channel) < m_nChannels)
0575       {
0576         
0577         // still need to correctly set bad channels to zero.
0578         for (int channel = 0; channel < m_nChannels - (nChannels + skiped_channel); channel++)
0579         {
0580           towerNumber++;
0581 
0582           unsigned int key = TowerInfoDefs::encode_emcal(towerNumber - 1);
0583           unsigned int phi_bin = TowerInfoDefs::getCaloTowerPhiBin(key);
0584           unsigned int eta_bin = TowerInfoDefs::getCaloTowerEtaBin(key);
0585 
0586           int sectorNumber = phi_bin / 8 + 1;
0587 
0588           int bin = h2_cemc_mean->FindBin(eta_bin + 0.5, phi_bin + 0.5);
0589 
0590           sectorAvg[sectorNumber - 1] += 0.;
0591 
0592           float signalFast = 0.0;
0593 
0594           rm_vector_twr[towerNumber - 1]->Add(&signalFast);
0595 
0596           h2_cemc_rm->SetBinContent(bin, rm_vector_twr[towerNumber - 1]->getMean(0));
0597           h2_cemc_rmhits->SetBinContent(bin, rm_vector_twrhits[towerNumber - 1]->getMean(0));
0598 
0599           h2_cemc_mean->SetBinContent(bin, h2_cemc_mean->GetBinContent(bin));
0600         }
0601       }
0602       delete p;
0603     }     // if packet good
0604     else  // packet is corrupted, treat all channels as zero suppressed
0605     {
0606       for (int channel = 0; channel < m_nChannels; channel++)
0607       {
0608         towerNumber++;
0609         unsigned int key = TowerInfoDefs::encode_emcal(towerNumber - 1);
0610         unsigned int phi_bin = TowerInfoDefs::getCaloTowerPhiBin(key);
0611         unsigned int eta_bin = TowerInfoDefs::getCaloTowerEtaBin(key);
0612 
0613         int sectorNumber = phi_bin / 8 + 1;
0614 
0615         int bin = h2_cemc_mean->FindBin(eta_bin + 0.5, phi_bin + 0.5);
0616 
0617         sectorAvg[sectorNumber - 1] += 0;
0618 
0619         float signalFast = 0;
0620 
0621         rm_vector_twr[towerNumber - 1]->Add(&signalFast);
0622 
0623         h2_cemc_rm->SetBinContent(bin, rm_vector_twr[towerNumber - 1]->getMean(0));
0624         h2_cemc_rmhits->SetBinContent(bin, rm_vector_twrhits[towerNumber - 1]->getMean(0));
0625 
0626         h2_cemc_mean->SetBinContent(bin, h2_cemc_mean->GetBinContent(bin));
0627       }
0628     }  // zero filling bad packets
0629   }    // packet loop
0630 
0631   h1_event->Fill(0);
0632 
0633   eventCounter++;
0634   return 0;
0635 }
0636 
0637 int CemcMon::Reset()
0638 {
0639   // reset our internal counters
0640   eventCounter = 0;
0641   idummy = 0;
0642   return 0;
0643 }