EventFEBStatus_module.cc

Go to the documentation of this file.

///////////////////////////////////////////////////////////////////////
// \file    EventFEBStatus_module.cc
// \brief   Analyzer module that counts hits on each FEB to determine if the FEB was live for a given event time.
// \author  lackey32@fnal.gov
////////////////////////////////////////////////////////////////////////

// C/C++ includes
#include <map>
#include <set>
#include <string>

// ROOT includes
#include "TTree.h"
#include "TH2F.h"
#include "TFile.h"

// Framework includes
#include "art/Framework/Core/EDAnalyzer.h"
#include "art/Framework/Core/ModuleMacros.h"
#include "art/Framework/Principal/Event.h"
#include "art/Framework/Principal/Handle.h"
#include "art/Framework/Principal/Run.h"
#include "art/Framework/Principal/SubRun.h"
#include "art_root_io/TFileService.h"
#include "canvas/Persistency/Common/FindManyP.h"
#include "canvas/Persistency/Common/FindOneP.h"
#include "canvas/Utilities/InputTag.h"
#include "fhiclcpp/ParameterSet.h"
#include "messagefacility/MessageLogger/MessageLogger.h"

// Novasoft DAQ includes
#include "DAQChannelMap/DAQChannelMap.h"
#include "NovaDAQConventions/DAQConventions.h"
// Novasoft includes
#include "CMap/service/CMapService.h"
#include "NovaTimingUtilities/TimingUtilities.h"
#include "RawData/RawDigit.h"
#include "RawData/RawTrigger.h"

namespace febstat {
  class EventFEBStatus : public art::EDAnalyzer {
  public:
    explicit EventFEBStatus(fhicl::ParameterSet const& pset);
    ~EventFEBStatus();
    void ReadEventFile();
    void beginJob() override;
    void analyze(const art::Event& evt) override;
    void ClearMapAndVectors();

  private:

    std::string  fRawDataLabel;  ///< Module that produced raw digits (to get trigger/event time)
    std::string  fCellHitLabel;  ///< Module that produced cell hits
    std::string  fEventFile;     ///< List of run/subrun/event/evt_times
    double       fWindow;        ///< Time window to check before/after event
    unsigned int fNhitCut;       ///< If # hits on an FEB in fWindow before or after the event time is below this, the feb is marked bad
    
    int run;                     
    int subrun;
    int spillEvent;              ///< Event number of current art event (last event with window in it is the number that goes in the TTree)
    int beamlineEvent;           ///< Event number of corresponding event in Beamline stream (provided in event-list file)
    Double_t spillEvtTime;       ///< Trigger time of current art event
    Double_t beamlineEvtTime;    ///< Trigger time of corresponding Beamline event
    std::vector<uint32_t> dcm;
    std::vector<uint32_t> feb;
    std::vector<int> nhit_pre;   ///< Number of hits in fWindow before beamline event time
    std::vector<int> nhit_post;  ///< Number of hits in fWindow after beamline event time
    std::vector<bool> febStatus; ///< 0 if nhit_pre or nhit_post is less than fNhitCut
    std::vector<bool> febShutoffBit; ///< 1 if shutoff bit set prior to event time and FEB recorded no hits after
    bool eventStatus;            ///< 0 if any FEB has a status of 0

    TTree *febshutoff;

    std::set<unsigned int> fRuns;
    std::set< std::pair<unsigned int, unsigned int> > fSubRuns;
    std::set< std::tuple<unsigned int, unsigned int, double> > fTimes;

    std::map< double, unsigned int > fEventMap;          ///< Map of beamline event time to event number

    std::map< std::pair<unsigned int, unsigned int> , 
              std::pair<unsigned int, unsigned int> > fFEBhitsNULL; ///< Null map to set hit count to zero for every FEB
    std::map< std::pair<unsigned int, unsigned int> , 
              std::pair<unsigned int, unsigned int> > fFEBhits;     ///< Map of hits counts before and after event for each dcm,feb combination

    std::map< std::pair<unsigned int, unsigned int>,
              std::pair<bool, double> > fShutoffMapNULL; ///< Null map to reset shutoff state for every FEB
    std::map< std::pair<unsigned int, unsigned int>,
              std::pair<bool, double> > fShutoffMap;     ///< Map of shutoff status and last hit time for each dcm,feb combination

    // Completeness booleans to determine if map needs to be reset
    bool completePre;           ///< Have we looked at the full fWindow before the event?
    bool completePost;          ///< Have we looked at the full fWindow after the event?

  };

  //......................................................................

  EventFEBStatus::EventFEBStatus(fhicl::ParameterSet const& pset)
    : EDAnalyzer(pset),
      fRawDataLabel(pset.get<std::string>("RawDataLabel")),
      fCellHitLabel(pset.get<std::string>("CellHitLabel")),
      fEventFile(pset.get<std::string>("EventFile")),
      fWindow(pset.get<double>("Window")),
      fNhitCut(pset.get<unsigned int>("NhitCut"))
  {
    ReadEventFile();
  }

  //......................................................................

  EventFEBStatus::~EventFEBStatus()
  {
  }

  //......................................................................
  ///
  /// Read in beamline event information from text file.
  /// File format should be (lines starting with # are ignored):
  /// #  run  subrun     event         event time
  /// 100792       1      3860  1583451706.619207 
  /// 100792       2     13403  1583452785.790757 
  /// 100792       5      9158  1583454586.687062 
  /// 100792       7     11686  1583455849.647850
  ///
  void EventFEBStatus::ReadEventFile()
  {
    fRuns.    clear();
    fSubRuns. clear();
    fTimes.  clear();
    std::ifstream ifs;
    ifs.open(fEventFile.c_str());
    if (ifs.is_open()) {
      while (ifs.good()) {
        std::string buff;
        std::getline(ifs, buff);
        if (buff[0]=='#') continue;
      
        int run, subrun, event;
        double time;
              
        sscanf(buff.c_str(), "%d %d %d %lf", &run, &subrun, &event, &time);

        fRuns.insert(run);
        fSubRuns.insert(std::make_pair(run, subrun));
        fEventMap.insert( std::pair<double, unsigned int>(time, event) );
        fTimes.insert(std::make_tuple(run, subrun, time));
      } // close while loop
      ifs.close();
    } // close input file
    else {
      std::cerr << __FILE__ << ":" << __LINE__ 
                << " Failed to open file " << fEventFile << " for read" 
                << std::endl;
    }

  }

  //......................................................................
  ///
  /// Clear vectors and map that are filled for each beamline event.
  /// Since these are filled by beamline event, but this module runs over the spill stream,
  /// they should not be cleared after each art event, but rather once we know we have the 
  /// full window before/after the time of interest.
  ///
  void EventFEBStatus::ClearMapAndVectors()
  {
    fFEBhits.clear();
    // Cleanup vectors
    dcm.clear();
    feb.clear();
    nhit_pre.clear();
    nhit_post.clear();
    febStatus.clear();
    febShutoffBit.clear();

    // Reset map to zeros and completion flags to false
    fFEBhits = fFEBhitsNULL;
    fShutoffMap = fShutoffMapNULL;
    completePre  = false;
    completePost = false;

  }

  //......................................................................

  void EventFEBStatus::beginJob()
  {

    // Initialize null map
    for (unsigned int dcmId=1; dcmId<4; ++dcmId){
      for (unsigned int febId=0; febId<64; ++febId){
        if ( (dcmId<3&&(febId==8||febId==24||febId==40||febId==56)) || (dcmId==2&&(febId==47||febId==63)) || (dcmId==3&&febId>7) )
          continue;
        fFEBhitsNULL.insert(std::pair <std::pair <unsigned int, unsigned int>, std::pair <unsigned int, unsigned int> >
                           (std::make_pair(dcmId,febId),std::make_pair(0,0)));
        // Start with assuming all FEBs are shutoff - only set to false if a hit is detected in 110 ms before
        // This takes care of FEBs that didn't have a single hit during this time and therefore did not report the shutoff bit.
        fShutoffMapNULL.insert(std::pair <std::pair <unsigned int, unsigned int>, std::pair <bool, double> >
                           (std::make_pair(dcmId,febId),std::make_pair(true,0)));
      }
    }
    ClearMapAndVectors();

    art::ServiceHandle<art::TFileService> tfs;
    febshutoff = tfs->make<TTree>("evts","");
    febshutoff->Branch("run",&run,"run/I");
    febshutoff->Branch("subrun",&subrun,"subrun/I");
    febshutoff->Branch("spillEvent",&spillEvent,"spillEvent/I");    
    febshutoff->Branch("beamlineEvent",&beamlineEvent,"beamlineEvent/I");    
    febshutoff->Branch("dcm",&dcm);
    febshutoff->Branch("feb",&feb);
    febshutoff->Branch("spillEvtTime",&spillEvtTime,"spillEvtTime/D");
    febshutoff->Branch("beamlineEvtTime",&beamlineEvtTime,"beamlineEvtTime/D");
    febshutoff->Branch("nhit_pre",&nhit_pre);
    febshutoff->Branch("nhit_post",&nhit_post);
    febshutoff->Branch("febStatus",&febStatus);
    febshutoff->Branch("eventStatus",&eventStatus,"eventStatus/O");
    
    febshutoff->Branch("febShutoffBit",&febShutoffBit);

  }

  //......................................................................

  void EventFEBStatus::analyze(const art::Event& evt)
  
{
  bool inlist = false;
    // Is this run/subrun pair in the list?
    inlist = fSubRuns.count(std::make_pair(evt.run(), evt.subRun()));

    // Skip the event if the run, subrun pair is not in the list
    if(!inlist){
      return;
    }

    run = evt.run();
    subrun = evt.subRun();
    spillEvent = evt.event();
    
    // Load services
    art::ServiceHandle<cmap::CMap> cmap;
    
    // Get trigger time
    art::Handle<std::vector<rawdata::RawTrigger>> trigs;
    evt.getByLabel(fRawDataLabel, trigs);
    const rawdata::RawTrigger trig = trigs->at(0);

    struct timespec unixTime;
    novadaq::timeutils::convertNovaTimeToUnixTime(trig.fTriggerTimingMarker_TimeStart, unixTime);
    unsigned long int unixTimeSec = unixTime.tv_sec;
    unsigned long int unixTimeNanoSec = unixTime.tv_nsec;
    spillEvtTime = unixTimeSec + (1e-9)*unixTimeNanoSec;

    // Check event list for first time after current trigger time (first time potentially contained within event)
    std::set< std::tuple<unsigned int, unsigned int, double> >::iterator setIt
      = fTimes.lower_bound(std::make_tuple(run,subrun,spillEvtTime));
    beamlineEvtTime = std::get<2>(*setIt);
    //bound functiion will only find *next* time in list. Check previous time to see if diff is closer.
    setIt--;
    double beamlineEvtTime_prev = std::get<2>(*setIt);
    double diff = beamlineEvtTime - spillEvtTime;
    double diff2 = beamlineEvtTime_prev - spillEvtTime;
    if (abs(diff2)<abs(diff)){
      beamlineEvtTime = beamlineEvtTime_prev;
      diff = diff2;
    }
    // Skip event if it does not contain any of window of interest
    //if ( (diff-fWindow<0 || diff-fWindow>0.2) && (diff+fWindow<0 || diff+fWindow>0.2) )
    // Need to keep longer before for shutfoff bit stuff
    if ( (diff-0.11<0 || diff-0.11>0.2) && (diff+fWindow<0 || diff+fWindow>0.2) )
      return;

    // Lookup beamline event number in event map using beamline event time
    beamlineEvent = fEventMap.find(beamlineEvtTime)->second;

    // Check if you have the full time window of interest yet
    if (!completePre && spillEvtTime < (beamlineEvtTime-fWindow) ){
      completePre = true;
    }
    if (!completePost && (spillEvtTime + 0.2) > (beamlineEvtTime+fWindow) ) {
      completePost = true;
    }
    
    // Get CellHits
    art::Handle<std::vector<rb::CellHit>> cellhits;
    evt.getByLabel(fCellHitLabel,cellhits);
    double hitTimePrev = 0;
    for (unsigned int hitIdx=0; hitIdx<cellhits->size(); ++hitIdx){
      art::Ptr<rb::CellHit> ihit(cellhits,hitIdx);
      uint32_t dchan = ihit->DaqChannel();
      unsigned int dcmId = cmap->Map()->getDCM(dchan);
      unsigned int febId = cmap->Map()->getFEB(dchan);
      double hitTime = spillEvtTime + ihit->TNS()*1e-9;
      

      // Check for FEB shutoff in 110 ms before event
      if ( (beamlineEvtTime - hitTime) < 0.110 && (beamlineEvtTime - hitTime) > 0 ) {
        uint8_t febStatusNano = ihit->fFEBStatus;
        hitTimePrev = fShutoffMap.at(std::make_pair(dcmId,febId)).second;
        if ( febStatusNano == 32) {
          //std::cout<<std::fixed<<"Shutoff detected. dcm: "<<dcmId<<" feb: "<<febId<<" hitTime: "<<hitTime<<std::endl;
          fShutoffMap.at(std::make_pair(dcmId,febId)).first = true;
          fShutoffMap.at(std::make_pair(dcmId,febId)).second = hitTime;
        }
        else if (febStatusNano == 48 && hitTime > hitTimePrev ) {
          //if (fShutoffMap.at(std::make_pair(dcmId,febId)).first)
            //std::cout<<std::fixed<<"Hits after shutoff. dcm: "<<dcmId<<" feb: "<<febId<<" hitTime: "<<hitTime<<std::endl;
          fShutoffMap.at(std::make_pair(dcmId,febId)).first = false;
          fShutoffMap.at(std::make_pair(dcmId,febId)).second = hitTime;
        }
      }

      // count hits for fWindow before
      if ( (beamlineEvtTime - hitTime) < fWindow && (beamlineEvtTime - hitTime) > 0 ) {
        fFEBhits.at(std::make_pair(dcmId,febId)).first++;
      }
      // count hits for fWindow after
      else if ( (hitTime - beamlineEvtTime) < fWindow && (hitTime - beamlineEvtTime) > 0 ) {
        fFEBhits.at(std::make_pair(dcmId,febId)).second++;
      }

    } // end loop over cellhits


    // Fill tree if we have information from the full window and reset map to zeroes
    if (completePre && completePost) {

      // Iterate over map and fill nhit stuff, including boolean of good/bad feb
      for (auto &mapIt: fFEBhits) {
        dcm.push_back((mapIt.first).first);
        feb.push_back((mapIt.first).second);
        nhit_pre.push_back((mapIt.second).first);
        nhit_post.push_back((mapIt.second).second);
        if ( (mapIt.second).first > fNhitCut && (mapIt.second).second > fNhitCut )
          febStatus.push_back(1);
        else
          febStatus.push_back(0);
      } // end loop over map

      // Iterate over shutoff bit map and fill
      for (auto &mapIt: fShutoffMap) {
        febShutoffBit.push_back((mapIt.second).first);
      }

      if (std::count(febStatus.begin(), febStatus.end(), true) == 126)
        eventStatus = 1;
      else
        eventStatus = 0;

      febshutoff->Fill();
      ClearMapAndVectors();

    } // end filling tree

  } // end analyze

  //......................................................................
  DEFINE_ART_MODULE(febstat::EventFEBStatus)
} // end namespace febstat