MuonID.cxx

Go to the documentation of this file.

////////////////////////////////////////////////////////////////////////
// \file    MuonID.cxx
// \brief   BDT for muon / mubar identification
// \author  Connor Johnson - Connor.Johnson@colostate.edu
////////////////////////////////////////////////////////////////////////
#include "MuonID/MuonID.h"

#include "canvas/Persistency/Common/FindOneP.h"

namespace muonid
{
  //......................................
  MuonID::MuonID()
    : rb::PID()
  {

  }

  //......................................
  MuonID::MuonID(int pdg, double val)
    : rb::PID(pdg,val)
  {

  }

  //.....................................
  unsigned int HighestPIDTrack(
                const std::vector< art::Ptr<rb::Track> >& sliceTracks,
                const std::string& muonidModuleLabel, const art::Event& e)
  {
    art::InputTag tag(muonidModuleLabel);
    return muonid::HighestPIDTrack(sliceTracks, tag, e);
  }

  //.....................................
  unsigned int HighestPIDTrack(std::vector< art::Ptr<rb::Track> > const& sliceTracks,
                               art::InputTag                      const& muonidTag,
                               art::Event                         const& e)
  {
    int    bestTrack        = -1;
    double highestPID       = -2.0;
    double bestTrackLength  = 0.0;
    double bestStartZ       = std::numeric_limits<double>::max();
    double bestADC          = std::numeric_limits<double>::max();

    art::FindOneP<muonid::MuonID> tracksToMuonID(sliceTracks, e, muonidTag);

    for(size_t c = 0; c < sliceTracks.size(); ++c){

      geo::View_t view = sliceTracks[c]->View();

      //Only use 3D tracks
      if (view != geo::kXorY){
        continue;
      } //End of loop over 2d tracks

      if(tracksToMuonID.isValid()){
        art::Ptr<muonid::MuonID> trackMuonID = tracksToMuonID.at(c);

        // If the PID is lower, we can bail, nothing to update
        if(trackMuonID->Value() < highestPID)
          continue;

        if(trackMuonID->Value() == highestPID)
        { // If the new guy is shorter, we can bail
          if (sliceTracks[c]->TotalLength() < bestTrackLength) continue;

            // If equal, need to fall to the next level, low start Z
          if (sliceTracks[c]->TotalLength() == bestTrackLength){

            // If start Z is higher, we can bail
            if (sliceTracks[c]->Start().Z() > bestStartZ) continue;

            // If equal, need to fall to the next level, total ADC
            // Use lower ADC
            if (sliceTracks[c]->Start().Z() == bestStartZ){
              
              if(sliceTracks[c]->TotalADC() > bestADC) continue;
              
              if(sliceTracks[c]->TotalADC() == bestADC)
              throw cet::exception("MuonID Matched Track Exception")<<std::endl
              << "Can't break a tie between two MuonID/tracks.  "
              << "This should have never happened, but it did.  " << std::endl;
            }
          }
        }

        highestPID = trackMuonID->Value();
        bestTrackLength = sliceTracks[c]->TotalLength();
        bestStartZ = sliceTracks[c]->Start().Z();
        bestADC    = sliceTracks[c]->TotalADC();
        bestTrack = c;

      }//End of valid MuonID object
      
    }//End of for loop over tracks

    if (bestTrack < 0)
      return 999;

    return (unsigned int)bestTrack;
  }
} // namespace muonid