SdHorizontalEventSelector.cc

Go to the documentation of this file.

#include "SdHorizontalEventSelector.h"

#include <fwk/CentralConfig.h>

#include <det/Detector.h>

#include <sdet/SDetector.h>
#include <sdet/Station.h>

#include <evt/Event.h>
#include <evt/ShowerRecData.h>
#include <evt/ShowerSRecData.h>

#include <sevt/SEvent.h>
#include <sevt/Station.h>

#include <utl/Vector.h>
#include <utl/ErrorLogger.h>
#include <utl/Branch.h>

using namespace fwk;
using namespace evt;
using namespace sevt;
using namespace utl;
using namespace std;
using namespace SdHorizontalEventSelectorHOG;


inline
double
Distance(const Point& p, const sdet::Station& s)
{
  return (p - s.GetPosition()).GetMag();
}


inline
bool
IsActive(const sevt::Station& s)
{
  return s.IsCandidate() ||
         s.IsSilent() ||
         (s.GetRejectionStatus() & (StationConstants::eOutOfTime |
                                    StationConstants::eRandom |
                                    StationConstants::eLonely |
                                    StationConstants::eTOTd));
}


inline
std::pair<int, std::vector<int>> 
SdHorizontalEventSelector::CalculateT5HASTrigger(const Event& event)
  const
{
  const SEvent& sEvent = event.GetSEvent();

  const ShowerRecData& shRec = event.GetRecShower();
  if (!shRec.HasSRecShower())
    return make_pair(0U, vector<int>());

  const ShowerSRecData& shSRec = shRec.GetSRecShower();
  const utl::Point& core = shSRec.GetCorePosition();

  SEvent::ConstStationIterator closestStation = sEvent.StationsEnd();
  {
    double closestDistance = 10000*meter;
    for (SEvent::ConstStationIterator sIt = sEvent.StationsBegin();
         sIt != sEvent.StationsEnd(); ++sIt) {
      if (!IsActive(*sIt))
        continue;

      const sdet::Station& ds = fSDetector->GetStation(*sIt);
      const double distanceToCore = Distance(core, ds);
      if (distanceToCore < closestDistance) {
        closestDistance = distanceToCore;
        closestStation = sIt;
      }
    }

    if (closestStation == sEvent.StationsEnd()) {
      // this should never happen
      ERROR("No station closest to core found");
      return make_pair(0U, vector<int>());
    }
  }

  pair<unsigned int, vector<int>> ret;
  ret.first = closestStation->GetId();
  const utl::Point& closestStationPos =
    fSDetector->GetStation(*closestStation).GetPosition();
  
  auto& functioningNeighbors = ret.second;
  for (SEvent::ConstStationIterator sIt = sEvent.StationsBegin();
       sIt != sEvent.StationsEnd(); ++sIt) {
    if (!IsActive(*sIt))
      continue;

    const sdet::Station& ds = fSDetector->GetStation(*sIt);
    const double dist = Distance(closestStationPos, ds);
    if (dist > fFirstCrownDistanceRange[0] && dist < fFirstCrownDistanceRange[1])
      functioningNeighbors.push_back(sIt->GetId());
  }

  if (int(functioningNeighbors.size()) >= fMinNumberOfActiveStations)
    return ret;

  return make_pair(0U, vector<int>());;
}


VModule::ResultFlag
SdHorizontalEventSelector::Init()
{
  Branch topB = CentralConfig::GetInstance()->GetTopBranch("SdHorizontalEventSelector");

  topB.GetChild("EnableT5HasTrigger").GetData(fEnableT5HasTrigger);

  topB.GetChild("RejectNonT5HasEvents").GetData(fRejectNonT5HasEvents);

  topB.GetChild("MinNumberOfActiveStations").GetData(fMinNumberOfActiveStations);

  topB.GetChild("FirstCrownDistanceRange").GetData(fFirstCrownDistanceRange);

  fSDetector = &det::Detector::GetInstance().GetSDetector();

  return eSuccess;
}


VModule::ResultFlag
SdHorizontalEventSelector::Run(evt::Event& event)
{
  if (!event.HasSEvent())
    return eContinueLoop;

  ShowerRecData& shRec = event.GetRecShower();
  if (!shRec.HasSRecShower())
    return eContinueLoop;

  if (fEnableT5HasTrigger) {

    const auto t5Stations = CalculateT5HASTrigger(event);
    const unsigned int closestStation = t5Stations.first;
    const int nActive = t5Stations.second.size();
    const bool t5Has = (closestStation && nActive >= fMinNumberOfActiveStations);
    
    bool closestIsUUB = false;
    if(t5Has)
      closestIsUUB = fSDetector->GetStation(closestStation).IsUUB();

    // Fill vector denoting which T5 neighbor stations were UUBs
    vector<int> t5NeighborsUUB;
    for (const auto sid : t5Stations.second) {
      const bool isUUB = fSDetector->GetStation(sid).IsUUB();
      if (isUUB)
        t5NeighborsUUB.push_back(sid);
    }

    ostringstream info;
    info << "T5Has = " << t5Has;
    INFO(info);

    if (fRejectNonT5HasEvents && !t5Has)
      return eContinueLoop;

    ShowerSRecData& shSRec = shRec.GetSRecShower();
    shSRec.SetT5Trigger(shSRec.GetT5Trigger() | t5Has*ShowerSRecData::eT5_Has);
    shSRec.SetT5ClosestStation(closestStation);
    shSRec.SetT5ClosestStationIsUUB(closestIsUUB);
    shSRec.SetT5PosteriorActiveNeighbors(t5Stations.second);
    shSRec.SetT5PosteriorActiveUUBNeighbors(t5NeighborsUUB);

  }

  return eSuccess;
}


VModule::ResultFlag
SdHorizontalEventSelector::Finish()
{
  return eSuccess;
}

// Configure (x)emacs for this file ...
// Local Variables:
// mode: c++
// End: