EnergyCorrelatorGeneralizedTool.cxx

Go to the documentation of this file.

/*
  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
*/
 
#include "JetSubStructureMomentTools/EnergyCorrelatorGeneralizedTool.h"
#include "JetSubStructureUtils/EnergyCorrelatorGeneralized.h"
#include "JetSubStructureUtils/EnergyCorrelator.h"
#include "AthContainers/ConstAccessor.h"
 
EnergyCorrelatorGeneralizedTool::EnergyCorrelatorGeneralizedTool(const std::string& name) : 
  JetSubStructureMomentToolsBase(name)
{
  declareProperty("Beta",       m_Beta = 1.0);
  declareProperty("BetaList",   m_rawBetaVals = {});
  declareProperty("DoM3",       m_doM3 = false);
  declareProperty("DoN3",       m_doN3 = false);
  declareProperty("DoLSeries",  m_doLSeries = false);
  declareProperty("DoDichroic", m_doDichroic = false);
}
 
StatusCode EnergyCorrelatorGeneralizedTool::initialize() {
 
  ATH_CHECK( JetSubStructureMomentToolsBase::initialize() );
 
  m_moments.emplace( 1.0, moments_t(1.0, m_prefix) );
 
  if( std::abs(m_Beta-1.0) > 1.0e-5 ) {
 
    ATH_MSG_WARNING( "The Beta property is deprecated, please use the BetaList property to provide a list of values" );
 
    m_moments.emplace( m_Beta, moments_t(m_Beta, m_prefix) );
 
  }
 
  for( float beta : m_rawBetaVals ) {
 
    float betaFix = round( beta * 10.0 ) / 10.0;
    if( std::abs(beta-betaFix) > 1.0e-5 ) ATH_MSG_DEBUG( "beta = " << beta << " has been rounded to " << betaFix );
 
    if( betaFix < 0.0 ) {
      ATH_MSG_WARNING( "beta must be positive. Skipping beta = " << beta );
      continue;
    }
 
    m_moments.emplace( betaFix, moments_t(betaFix, m_prefix) );
 
  }
 
  for( auto const& moment : m_moments ) {
    ATH_MSG_DEBUG( "Including beta = " << moment.first );
  }
 
  m_dec_ECFG_2_1_2 = std::make_unique< SG::AuxElement::Decorator<float> >(m_prefix+"ECFG_2_1_2");
  m_dec_ECFG_3_1_1 = std::make_unique< SG::AuxElement::Decorator<float> >(m_prefix+"ECFG_3_1_1");
  m_dec_ECFG_3_2_1 = std::make_unique< SG::AuxElement::Decorator<float> >(m_prefix+"ECFG_3_2_1");
  m_dec_ECFG_3_2_2 = std::make_unique< SG::AuxElement::Decorator<float> >(m_prefix+"ECFG_3_2_2");
  m_dec_ECFG_3_3_1 = std::make_unique< SG::AuxElement::Decorator<float> >(m_prefix+"ECFG_3_3_1");
  m_dec_ECFG_4_2_2 = std::make_unique< SG::AuxElement::Decorator<float> >(m_prefix+"ECFG_4_2_2");
  m_dec_ECFG_4_4_1 = std::make_unique< SG::AuxElement::Decorator<float> >(m_prefix+"ECFG_4_4_1");
 
  m_dec_ECFG_3_3_2 = std::make_unique< SG::AuxElement::Decorator<float> >(m_prefix+"ECFG_3_3_2");
 
  return StatusCode::SUCCESS;
 
}
 
int EnergyCorrelatorGeneralizedTool::modifyJet(xAOD::Jet &injet) const {
 
  fastjet::PseudoJet jet;
  fastjet::PseudoJet jet_ungroomed;
 
  bool calculate = SetupDecoration(jet,injet);
 
  bool calculate_ungroomed = false;
 
  if( m_doDichroic ) {
 
    static const SG::ConstAccessor<ElementLink<xAOD::JetContainer> > ParentAcc ("Parent");
    ElementLink<xAOD::JetContainer> parentLink = ParentAcc (injet);
 
    if( !parentLink.isValid() ) {
      ATH_MSG_ERROR( "Parent element link is not valid. Aborting" );
      return 1;
    }
 
    const xAOD::Jet* parentJet = *(parentLink);
    calculate_ungroomed = SetupDecoration(jet_ungroomed,*parentJet);
 
  }
 
  for( auto const& moment : m_moments ) {
 
    float beta = moment.first;
 
 
    float ECFG_2_1_value = -999.0;
    float ECFG_3_2_value = -999.0;
    
    float ECFG_2_1_ungroomed_value = -999.0;
    float ECFG_3_1_ungroomed_value = -999.0;
    float ECFG_3_2_ungroomed_value = -999.0;
 
    float ECFG_3_1_value = -999.0;
    float ECFG_4_1_value = -999.0;
    float ECFG_4_2_value = -999.0;
 
    if( calculate ) {
 
      JetSubStructureUtils::EnergyCorrelatorGeneralized ECFG_2_1(1, 2, beta, JetSubStructureUtils::EnergyCorrelator::pt_R);
      JetSubStructureUtils::EnergyCorrelatorGeneralized ECFG_3_1(1, 3, beta, JetSubStructureUtils::EnergyCorrelator::pt_R);
      JetSubStructureUtils::EnergyCorrelatorGeneralized ECFG_3_2(2, 3, beta, JetSubStructureUtils::EnergyCorrelator::pt_R);
 
      ECFG_2_1_value = ECFG_2_1.result(jet);
      ECFG_3_1_value = ECFG_3_1.result(jet);
      ECFG_3_2_value = ECFG_3_2.result(jet);
 
      if( calculate_ungroomed ) {
        ECFG_2_1_ungroomed_value = ECFG_2_1.result(jet_ungroomed);
        ECFG_3_1_ungroomed_value = ECFG_3_1.result(jet_ungroomed);
        ECFG_3_2_ungroomed_value = ECFG_3_2.result(jet_ungroomed);
      }
 
      if( m_doM3 ) {
        JetSubStructureUtils::EnergyCorrelatorGeneralized ECFG_4_1(1, 4, beta, JetSubStructureUtils::EnergyCorrelator::pt_R);
        ECFG_4_1_value = ECFG_4_1.result(jet);
      }
 
      if( m_doN3 ) {
        JetSubStructureUtils::EnergyCorrelatorGeneralized ECFG_4_2(2, 4, beta, JetSubStructureUtils::EnergyCorrelator::pt_R);
        ECFG_4_2_value = ECFG_4_2.result(jet);
      }
 
    }
 
    (*moment.second.dec_ECFG_2_1)(injet) = ECFG_2_1_value;
    (*moment.second.dec_ECFG_3_1)(injet) = ECFG_3_1_value;
    (*moment.second.dec_ECFG_3_2)(injet) = ECFG_3_2_value;
    (*moment.second.dec_ECFG_4_1)(injet) = ECFG_4_1_value;
    (*moment.second.dec_ECFG_4_2)(injet) = ECFG_4_2_value;
 
    (*moment.second.dec_ECFG_2_1_ungroomed)(injet) = ECFG_2_1_ungroomed_value;
    (*moment.second.dec_ECFG_3_1_ungroomed)(injet) = ECFG_3_1_ungroomed_value;
    (*moment.second.dec_ECFG_3_2_ungroomed)(injet) = ECFG_3_2_ungroomed_value;
 
  }
 
  float ECFG_2_1_2_value = -999;
  float ECFG_3_1_1_value = -999;
  float ECFG_3_2_1_value = -999;
  float ECFG_3_2_2_value = -999;
  float ECFG_3_3_1_value = -999;
  float ECFG_4_2_2_value = -999;
  float ECFG_4_4_1_value = -999;
 
  float ECFG_3_3_2_value = -999;
 
 
  if( calculate && m_doLSeries ) {
 
    if( m_moments.count(2.0) ) {
      ECFG_2_1_2_value = (*m_moments.at(2.0).dec_ECFG_2_1)(injet);
    }
    else {
      JetSubStructureUtils::EnergyCorrelatorGeneralized ECFG_2_1_2(1, 2, 2, JetSubStructureUtils::EnergyCorrelator::pt_R);
      ECFG_2_1_2_value = ECFG_2_1_2.result(jet);
    }
 
    if( m_doN3 ) {
      ECFG_3_1_1_value = (*m_moments.at(1.0).dec_ECFG_3_1)(injet);
    }
    else {
      JetSubStructureUtils::EnergyCorrelatorGeneralized ECFG_3_1_1(1, 3, 1, JetSubStructureUtils::EnergyCorrelator::pt_R);
      ECFG_3_1_1_value = ECFG_3_1_1.result(jet);
    }
 
    ECFG_3_2_1_value = (*m_moments.at(1.0).dec_ECFG_3_2)(injet);
 
    if( m_moments.count(2.0) ) {
      ECFG_3_2_2_value = (*m_moments.at(2.0).dec_ECFG_3_2)(injet);
    }
    else {
      JetSubStructureUtils::EnergyCorrelatorGeneralized ECFG_3_2_2(2, 3, 2, JetSubStructureUtils::EnergyCorrelator::pt_R);
      ECFG_3_2_2_value = ECFG_3_2_2.result(jet);
    }
 
    JetSubStructureUtils::EnergyCorrelatorGeneralized ECFG_3_3_1(3, 3, 1, JetSubStructureUtils::EnergyCorrelator::pt_R);
    ECFG_3_3_1_value = ECFG_3_3_1.result(jet);
 
    if( m_doN3 && m_moments.count(2.0) ) {
      ECFG_4_2_2_value = (*m_moments.at(2.0).dec_ECFG_4_2)(injet);
    }
    else {
      JetSubStructureUtils::EnergyCorrelatorGeneralized ECFG_4_2_2(2, 4, 2, JetSubStructureUtils::EnergyCorrelator::pt_R);
      ECFG_4_2_2_value = ECFG_4_2_2.result(jet);
    }
 
    JetSubStructureUtils::EnergyCorrelatorGeneralized ECFG_4_4_1(4, 4, 1, JetSubStructureUtils::EnergyCorrelator::pt_R);
    ECFG_4_4_1_value = ECFG_4_4_1.result(jet);
 
    JetSubStructureUtils::EnergyCorrelatorGeneralized ECFG_3_3_2(3, 3, 2, JetSubStructureUtils::EnergyCorrelator::pt_R);
    ECFG_3_3_2_value = ECFG_3_3_2.result(jet);
 
  }
 
  (*m_dec_ECFG_2_1_2)(injet) = ECFG_2_1_2_value;
  (*m_dec_ECFG_3_1_1)(injet) = ECFG_3_1_1_value;
  (*m_dec_ECFG_3_2_1)(injet) = ECFG_3_2_1_value;
  (*m_dec_ECFG_3_2_2)(injet) = ECFG_3_2_2_value;
  (*m_dec_ECFG_3_3_1)(injet) = ECFG_3_3_1_value;
  (*m_dec_ECFG_4_2_2)(injet) = ECFG_4_2_2_value;
  (*m_dec_ECFG_4_4_1)(injet) = ECFG_4_4_1_value;
 
  (*m_dec_ECFG_3_3_2)(injet) = ECFG_3_3_2_value;
 
  return 0;
 
}