Cherry-pick MR !221 to stripping21-patches branch

Added angle calculator for decay B+ -> K*+ MuMu, accepting both subdecays of the…

See merge request !221 (merged)

(cherry picked from commit c48ae388)

32bc8956 Added angle calculator for decay B+ -> K*+ MuMu, accepting both subdecays of the… 0519ae96 AngleCalculator for Bu2KstarMuMu: use nullptr instead of NULL, output if K+ and… a2027610 Flipped sign of phi for Bminus case

Phys/DaVinciAngleCalculators/src/Bu2KstarMuMuAngleCalculator.cpp

+// $Id: Bu2KstarMuMuAngleCalculator.cpp,v 1.0 2017-12-07 09:41:00 dgerick Exp $
+// Include files 
+// from Gaudi
+
+#include "GaudiKernel/PhysicalConstants.h"
+#include "GaudiKernel/DeclareFactoryEntries.h" 
+#include "Event/Particle.h" 
+#include "Kernel/IParticleDescendants.h"
+
+// local
+#include "Bu2KstarMuMuAngleCalculator.h"
+#include "DaVinciP2VVAngles.h"
+
+//-----------------------------------------------------------------------------
+// Implementation file for class : Bu2KstarMuMuAngleCalculator
+//
+// 2007-08-02 : Thomas Blake
+// 2017-12-07 : David Gerick
+//-----------------------------------------------------------------------------
+
+// Declaration of the Tool Factory
+
+DECLARE_TOOL_FACTORY( Bu2KstarMuMuAngleCalculator )
+
+using namespace DaVinci::P2VVAngles;
+
+//=============================================================================
+// Standard constructor, initializes variables
+//=============================================================================
+
+Bu2KstarMuMuAngleCalculator::Bu2KstarMuMuAngleCalculator( const std::string& type,
+                                                          const std::string& name,
+                                                          const IInterface* parent )
+  : GaudiTool ( type, name , parent )
+{
+  declareInterface<IP2VVPartAngleCalculator>(this);    
+}
+//=============================================================================
+// Destructor
+//=============================================================================
+Bu2KstarMuMuAngleCalculator::~Bu2KstarMuMuAngleCalculator() {} 
+
+//=============================================================================
+
+
+StatusCode Bu2KstarMuMuAngleCalculator::daughters( const LHCb::Particle* mother ){
+
+  m_pMuMinus = nullptr ;
+  m_pMuPlus = nullptr ;
+  m_pK = nullptr ;
+  m_pPi = nullptr ;
+
+  LHCb::Particle::ConstVector descendants = m_descendants->descendants( mother, 2 ) ; 
+  LHCb::Particle::ConstVector::const_iterator iter;
+  
+    
+  if ( 4 > descendants.size() ) return StatusCode::FAILURE;
+ 
+  int absidPi = 0, absidK = 0;
+
+  for ( iter = descendants.begin(); iter != descendants.end(); ++iter ){
+    int pid = (*iter)->particleID().pid() ;
+    int absid = (*iter)->particleID().abspid() ;
+
+    if ( 13 == pid || 11 == pid ) m_pMuMinus = (*iter);
+    else if ( -13 == pid || -11 == pid ) m_pMuPlus = (*iter);
+    else if ( 321 == absid || 310 == absid ){
+	m_pK = (*iter);
+	absidK = absid;
+    }
+    else if ( 111 == absid || 211 == absid ){
+	m_pPi = (*iter);
+	absidPi = absid;
+    }
+  }
+
+  if ( nullptr == m_pK || nullptr == m_pMuPlus || 
+       nullptr == m_pMuMinus || nullptr == m_pPi ){
+    return StatusCode::FAILURE ;
+  }
+
+  if(321 == absidK && 111 != absidPi){
+    error() << "Cannot use K+ and pi+ as K*+ daughters for angle calculation!" << endmsg;
+    return StatusCode::FAILURE ; //allow only K+ and pi0...
+  }
+  if(310 == absidK && 211 != absidPi){
+    error() << "Cannot use KS0 and pi0 as K*+ daughters for angle calculation!" << endmsg;
+    return StatusCode::FAILURE ; //... or KS0 and pi+ combination
+  }
+
+  return StatusCode::SUCCESS ;
+}
+
+
+StatusCode Bu2KstarMuMuAngleCalculator::initialize() {
+  
+  StatusCode sc = GaudiTool::initialize() ;
+  if (!sc) return sc ;
+  
+  if (msgLevel(MSG::DEBUG)) debug() << "Initializing Angle Calculator Tool" << endmsg ;
+  
+  m_descendants = tool<IParticleDescendants>("ParticleDescendants", this);
+  
+  return sc ;
+}
+
+
+
+StatusCode Bu2KstarMuMuAngleCalculator::calculateAngles( 
+   const LHCb::Particle* particle, double& thetal, 
+   double& thetak, double& phi)
+{
+  StatusCode sc = daughters( particle ) ; 
+  
+  if ( sc.isFailure() ){
+    error() << "Could not find required particles !" << endmsg;
+    return sc; 
+  }
+  
+  bool isBplus = (particle->particleID().pid() > 0);
+
+  
+  if ( isBplus ) {
+      thetal = calculateHelicityPolarAngle( m_pK->momentum(), m_pPi->momentum(),
+                                          m_pMuPlus->momentum(), m_pMuMinus->momentum() );  
+    
+    thetak = calculateHelicityPolarAngle( m_pMuPlus->momentum(), m_pMuMinus->momentum(),
+                                          m_pK->momentum(), m_pPi->momentum() );
+  
+    phi = calculatePlaneAngle( m_pMuPlus->momentum(), m_pMuMinus->momentum(),
+                               m_pK->momentum(), m_pPi->momentum() ); 
+  } 
+  else {
+    thetal = calculateHelicityPolarAngle( m_pK->momentum(), m_pPi->momentum(), 
+                                          m_pMuMinus->momentum(), m_pMuPlus->momentum() );
+    
+    thetak = calculateHelicityPolarAngle( m_pMuMinus->momentum(), m_pMuPlus->momentum(),
+                                          m_pK->momentum(), m_pPi->momentum() );
+   
+    //flip sign for phi in Bminus case 
+    phi = - calculatePlaneAngle( m_pMuMinus->momentum(), m_pMuPlus->momentum(),
+                               m_pK->momentum(), m_pPi->momentum() );  
+  }
+  
+  return sc ;
+}
+
+
+StatusCode Bu2KstarMuMuAngleCalculator::calculateTransversityAngles( 
+   const LHCb::Particle* particle,
+   double& Theta_tr, double& Phi_tr, double& Theta_V )
+{
+
+  StatusCode sc = daughters( particle ); 
+    
+  if ( sc.isFailure() ){
+    error() << "Could not find required particles !" << endmsg;
+    return sc; 
+  }
+ 
+  bool isBplus = (particle->particleID().pid() > 0);
+  if ( isBplus ){
+    Theta_V =  calculateHelicityPolarAngle( m_pMuPlus->momentum(), m_pMuMinus->momentum(), 
+                                            m_pK->momentum(), m_pPi->momentum() );  
+  
+    Theta_tr = calculateThetaTr( m_pMuPlus->momentum(), m_pMuMinus->momentum(), 
+                                 m_pK->momentum(), m_pPi->momentum() );
+                                 
+    
+    Phi_tr = calculatePhiTr( m_pMuPlus->momentum(), m_pMuMinus->momentum(), 
+                             m_pK->momentum(), m_pPi->momentum() ) ;  
+  }
+  else {
+    
+    Theta_V =  calculateHelicityPolarAngle( m_pMuMinus->momentum(), m_pMuPlus->momentum(),
+                                            m_pK->momentum(), m_pPi->momentum() );  
+    
+    Theta_tr = calculateThetaTr( m_pMuMinus->momentum(), m_pMuPlus->momentum(),
+                                 m_pK->momentum(), m_pPi->momentum() );
+    
+    
+    Phi_tr = calculatePhiTr( m_pMuMinus->momentum(), m_pMuPlus->momentum(),
+                             m_pK->momentum(), m_pPi->momentum() ) ;  
+  }
+  
+  return StatusCode::SUCCESS ;
+}
+
+
+
+double Bu2KstarMuMuAngleCalculator::calculatePhi( const LHCb::Particle* particle )
+{
+  
+  StatusCode sc = daughters( particle ) ; 
+  
+  if ( sc.isFailure() ){
+     Exception("Could not find required particles !");
+  }
+
+  bool isBplus = (particle->particleID().pid() > 0);
+  
+  double  phi = -999;
+
+  if ( isBplus ){
+    phi = calculatePlaneAngle( m_pMuPlus->momentum(), m_pMuMinus->momentum(),
+                               m_pK->momentum(), m_pPi->momentum() );  
+  }
+  else {
+    phi = calculatePlaneAngle( m_pMuMinus->momentum(), m_pMuPlus->momentum(),
+                               m_pK->momentum(), m_pPi->momentum() );
+  }
+
+  return phi ;
+}
+
+
+// as far as here
+
+double Bu2KstarMuMuAngleCalculator::calculateThetaL( const LHCb::Particle* particle )
+{
+  
+  StatusCode sc = daughters( particle ); 
+
+  if ( sc.isFailure() ){
+     Exception("Could not find required particles !");
+  }
+  
+  bool isBplus = (particle->particleID().pid() > 0);
+  
+  double theta = -999;
+  
+  if ( isBplus ){
+    theta = calculateHelicityPolarAngle( m_pK->momentum(), m_pPi->momentum(), 
+                                         m_pMuPlus->momentum(), m_pMuMinus->momentum() );  
+  }
+  else {
+    
+    theta = calculateHelicityPolarAngle( m_pK->momentum(), m_pPi->momentum(), 
+                                         m_pMuMinus->momentum(), m_pMuPlus->momentum() );  
+  }
+  return theta ;
+}
+
+
+double Bu2KstarMuMuAngleCalculator::calculateThetaK( const LHCb::Particle* particle ) 
+{
+  StatusCode sc = daughters( particle ); 
+    
+  if ( sc.isFailure() ){
+     Exception("Could not find required particles !");
+  }
+  
+  double theta = calculateHelicityPolarAngle( m_pK->momentum(), m_pPi->momentum(), 
+                                              m_pMuPlus->momentum(), m_pMuMinus->momentum() );  
+  return theta ;
+}
+
+
+
+double Bu2KstarMuMuAngleCalculator::calculateTransThetaTr( 
+   const LHCb::Particle* particle ) 
+{
+  StatusCode sc = daughters( particle );
+  
+  if ( sc.isFailure() ){
+     Exception("Could not find required particles !");
+  }
+
+  bool isBplus = (particle->particleID().pid() > 0);
+  double Theta_tr = -999.;
+  
+  if ( isBplus ){
+    Theta_tr = calculateThetaTr( m_pMuPlus->momentum(), m_pMuMinus->momentum(),
+                                 m_pK->momentum(), m_pPi->momentum() );
+  }
+  else {
+    Theta_tr = calculateThetaTr( m_pMuMinus->momentum(), m_pMuPlus->momentum(),
+                                 m_pK->momentum(), m_pPi->momentum() );  
+  }  
+  
+  return Theta_tr;   
+}
+
+double Bu2KstarMuMuAngleCalculator::calculateTransPhiTr( 
+   const LHCb::Particle* particle ) 
+{
+  StatusCode sc = daughters( particle ) ;
+    
+  if ( sc.isFailure() ){
+    Exception("Could not find required particles !");
+  }
+  
+  bool isBplus = (particle->particleID().pid() > 0);
+  double Phi_tr = -999.;
+  
+  if ( isBplus ){
+    Phi_tr = calculatePhiTr( m_pMuPlus->momentum(), m_pMuMinus->momentum(),
+                             m_pK->momentum(), m_pPi->momentum() );
+  }
+  else {
+    Phi_tr = calculatePhiTr( m_pMuMinus->momentum(), m_pMuPlus->momentum(),
+                             m_pK->momentum(), m_pPi->momentum() );
+  }
+  
+  return Phi_tr;
+}
+
+double Bu2KstarMuMuAngleCalculator::calculateTransThetaV( const LHCb::Particle* particle ) 
+{
+  return calculateThetaK( particle ) ;
+}
+