utilities.h

/**************************************************************************
 * BASF2 (Belle Analysis Framework 2)                                     *
 * Copyright(C) 2016 - Belle II Collaboration                             *
 *                                                                        *
 * Author: The Belle II Collaboration                                     *
 * Contributors: Nils Braun                                               *
 *                                                                        *
 * This software is provided "as is" without any warranty.                *
 **************************************************************************/
#pragma once
 
#include <analysis/dataobjects/ParticleList.h>
 
#include <mdst/dataobjects/ECLCluster.h>
#include <tracking/dataobjects/RecoTrack.h>
 
#include <analysis/utility/PCmsLabTransform.h>
#include <analysis/ClusterUtility/ClusterUtils.h>
#include <TLorentzVector.h>
 
namespace Belle2 {
  namespace SoftwareTrigger {
    template<class T>
    inline TLorentzVector getFourVector(const T& item);
 
    template<>
    inline TLorentzVector getFourVector(const std::reference_wrapper<const ECLCluster>& cluster)
    {
      ClusterUtils C;
      const TLorentzVector& v = C.Get4MomentumFromCluster(&(cluster.get()), ECLCluster::EHypothesisBit::c_nPhotons);
      return TLorentzVector(v.Px(), v.Py(), v.Pz(), v.E());
    }
 
    template<>
    inline TLorentzVector getFourVector(const RecoTrack& track)
    {
      const TVector3& positionSeed = track.getPositionSeed();
      return TLorentzVector(positionSeed.X(), positionSeed.Y(), positionSeed.Z(),
                            sqrt(positionSeed.Mag2() + Const::pionMass * Const::pionMass));
    }
 
    template<class AnEntityType>
    inline const AnEntityType* getElementFromParticle(const Particle& particle);
 
    template<>
    inline const Particle* getElementFromParticle(const Particle& particle)
    {
      return &particle;
    }
 
    template<>
    inline const ECLCluster* getElementFromParticle(const Particle& particle)
    {
      return particle.getECLCluster();
    }
 
    inline double BeamEnergyCMS()
    {
      PCmsLabTransform T;
      return T.getCMSEnergy() / 2.0;
    }
 
    template<class AnEntity>
    inline double getRho(const AnEntity* entity)
    {
      if (not entity) {
        return -1;
      }
 
      const TLorentzVector& fourVector = entity->get4Vector();
      return PCmsLabTransform::labToCms(fourVector).Rho();
    }
 
    // Template specialization for ECLCluster
    template<>
    inline double getRho(const ECLCluster* entity)
    {
      if (not entity) {
        return -1;
      }
 
      ClusterUtils C;
      return PCmsLabTransform::labToCms(C.Get4MomentumFromCluster(entity, ECLCluster::EHypothesisBit::c_nPhotons)).Rho();
    }
 
    template<class AReturnType>
    static const AReturnType* getElementWithMaximumRhoBelow(const StoreObjPtr<ParticleList>& particles,
                                                            const double belowLimit)
    {
      const AReturnType* elementMaximumRho = nullptr;
      double maximumRho = -1.;
      for (const Particle& currentParticle : *particles) {
        const auto currentElement = getElementFromParticle<AReturnType>(currentParticle);
        const double& currentRho = getRho(currentElement);
        if (currentRho >= belowLimit) {
          continue;
        }
        if (currentRho > maximumRho) {
          maximumRho = currentRho;
          elementMaximumRho = currentElement;
        }
      }
      return elementMaximumRho;
    }
 
    template<class AReturnType>
    inline const AReturnType* getElementWithMaximumRho(const StoreObjPtr<ParticleList>& particles)
    {
      return getElementWithMaximumRhoBelow<AReturnType>(particles, std::nan(""));
    }
 
    extern double getRhoOfECLClusterWithMaximumRhoBelow(const StoreObjPtr<ParticleList>& pions,
                                                        const StoreObjPtr<ParticleList>& gammas,
                                                        const double belowLimit);
 
    inline double getRhoOfECLClusterWithMaximumRho(const StoreObjPtr<ParticleList>& pionshlt,
                                                   const StoreObjPtr<ParticleList>& gammahlt)
    {
      return getRhoOfECLClusterWithMaximumRhoBelow(pionshlt, gammahlt, std::nan(""));
    }
  }
}