/Users/jmonk/Physics/ForIA/src/AnalysisTools/ClusterSelection.cxx

Go to the documentation of this file.

#include "ForIA/AnalysisTools/ClusterSelection.hh"
#include "ForIA/Units.hh"
#include "ForIA/Utils.hh"
#include <cmath>
#include <stdexcept>

namespace ForIA {

  ClusterSelection::ClusterSelection(): 
    FourVectorSelection(), 
    m_haveSetetMin(false),
    m_haveSeteMin(false),
    m_haveSetetaMin(false),
    m_haveSetetaMax(false),
    m_haveSetEmFrac(false),
    m_haveSetEmFracSpike(false),
    m_smearE(false),
    m_scaleMCEType("NOMINAL"),
    m_useEMScale(true)
  {
    m_eScale.SetRandomSeed(4321);
  }

  const ClusterVector &ClusterSelection::acceptedClusters() const{
    if(m_freshEvent) fillVectors();
    return m_acceptedClusters;
  }
  const ClusterVector &ClusterSelection::rejectedClusters() const{
    if(m_freshEvent) fillVectors();
    return m_rejectedClusters;
  }

  bool ClusterSelection::setSelection(DefinedSelection selection){
    if(!m_canSetCuts) return false;

    m_haveSetscale = false;
    m_haveSetetMin = false;
    m_haveSeteMin  = false;
    m_haveSetetaMin= false;
    m_haveSetetaMax= false;
    m_haveSetEmFrac= false;
    m_haveSetEmFracSpike= false;

    switch (selection) {
    case LOOSE:
      setSelection(LOOSE);
      setEMin(500 * MeV);
      break;
      
    case MEDIUM:
      setSelection(MEDIUM);
      setEMin(1.0 * GeV);
      break;

    case TIGHT:
      setSelection(TIGHT);
      setEMin(5.0 * GeV);
      break;

    }

    return true;
  }

  void ClusterSelection::fillVectors() const{
    
    m_freshEvent = false;
    
    m_acceptedClusters.clear();
    m_passedFourVectors.clear();
    m_rejectedClusters.clear();
    m_rejectedFourVectors.clear();

    if(m_event == 0) return;

    ClusterVector input = m_event->clusters();

    for(ClusterVector::const_iterator cluster = input.begin();
        cluster != input.end(); ++cluster){
      

      (*cluster)->setUseEMScale(m_useEMScale);


      if (m_smearE) {
        double SF = 1.0;
        if ((*cluster)->E() > 0.0) m_eScale.getSmearingCorrection((*cluster)->eta(),(*cluster)->E()/1000. ,0,true,"2010");
        (*cluster)->scaleEnergy(SF);
      }


      if ((*cluster)->isMC() && m_scaleMCEType != "NONE" ) {
        double SF = UtilsScaleMC((*cluster)->eta(),m_scaleMCEType);
          //scaleFactor((*cluster)->eta());
                (*cluster)->scaleEnergy(SF);
      }
      bool passedCuts = true;

      if (m_haveSetetMin && passedCuts){
        passedCuts = ((*cluster)->ET() >= m_etMin);
      
      }

      if(m_haveSeteMin   && passedCuts){          
        //      std::cout << "Energy cut value " << m_eMin << std::endl;
        if(m_haveSetscale)
          passedCuts = ((*cluster)->E(m_scale) > m_eMin); //if broken, switch m_scale to Cluster::EM
        else
          passedCuts = ((*cluster)->E() > m_eMin);
      }
      if(m_haveSetetaMin && passedCuts){ passedCuts = (fabs((*cluster)->eta()) >= m_etaMin);}
      //      else continue;
      if(m_haveSetetaMax && passedCuts){ passedCuts = (fabs((*cluster)->eta()) <= m_etaMax);}
      //      else continue;
      if (m_haveSetEmFrac && passedCuts){
        float emFraction = (*cluster)->E(Cluster::EM) / ((*cluster)->E(Cluster::EM) + (*cluster)->E(Cluster::HAD));
        passedCuts = (emFraction > m_emFrac);
       }

      if (m_haveSetEmFracSpike && passedCuts){
        float emFraction = (*cluster)->E(Cluster::EM) / ((*cluster)->E(Cluster::EM) + (*cluster)->E(Cluster::HAD));
        passedCuts = (emFraction > m_emFracSpike) || (fabs((*cluster)->eta()) < 1.3) || (fabs((*cluster)->eta()) > 1.32);

      }


      if(passedCuts){
        m_acceptedClusters.push_back( *cluster );
        m_passedFourVectors.push_back( *cluster ); 
      } else {
        m_rejectedClusters.push_back( *cluster );
        m_rejectedFourVectors.push_back( *cluster );
      }
    }

    return;
  }

  bool ClusterSelection::setScale(Cluster::Energy scale){
    if(!m_canSetCuts) return false;
    m_haveSetscale = true;
    m_scale = scale;
    return true;
  }

  bool ClusterSelection::setETMin(double etMin){
    if(!m_canSetCuts) return false;
    m_haveSetetMin = true;
    m_etMin = etMin;
    return true;
  }

  bool ClusterSelection::setEMin(double eMin){
    if(!m_canSetCuts){ std::cout << "Cuts not allowed!" << std::endl; return false; }
    m_haveSeteMin = true;
    m_eMin = eMin;
    //std::cout << "This happened? " << eMin << std::endl;
    return true;
  }

  bool ClusterSelection::setEtaMin(double etaMin){
    if(!m_canSetCuts) return false;
    m_haveSetetaMin = true;
    m_etaMin = etaMin;
    return true;
  }

  bool ClusterSelection::setEtaMax(double etaMax){
    if(!m_canSetCuts) return false;
    m_haveSetetaMax = true;
    m_etaMax = etaMax;
    return true;
  }

  bool ClusterSelection::setEmFrac(double emFrac){
    if(!m_canSetCuts) return false;
    m_haveSetEmFrac = true;
    m_emFrac = emFrac;
    return true;
  }

  bool ClusterSelection::setEmFracSpike(double emFracSpike){
    if(!m_canSetCuts) return false;
    m_haveSetEmFracSpike = true;
    m_emFracSpike = emFracSpike;
    return true;
  }

void ClusterSelection::setSmearE(bool smear){
   m_smearE = smear;
   return;
  }

 void ClusterSelection::setUseEMScale(bool useEMScale){
   m_useEMScale = useEMScale;
   return;
  }

  void ClusterSelection::setScaleMCClusterEType(std::string type){
    m_scaleMCEType = type;
    return;
  }



}//end