RealisticFacetFilter.cc

/**************************************************************************
 * basf2 (Belle II Analysis Software Framework)                           *
 * Author: The Belle II Collaboration                                     *
 *                                                                        *
 * See git log for contributors and copyright holders.                    *
 * This file is licensed under LGPL-3.0, see LICENSE.md.                  *
 **************************************************************************/
#include <tracking/trackFindingCDC/filters/facet/RealisticFacetFilter.h>
 
#include <tracking/trackingUtilities/eventdata/hits/CDCFacet.h>
 
#include <tracking/trackingUtilities/numerics/Quadratic.h>
 
#include <tracking/trackingUtilities/utilities/StringManipulation.h>
 
#include <framework/core/ModuleParamList.templateDetails.h>
#include <framework/geometry/VectorUtil.h>
 
#include <Math/Vector2D.h>
 
#include <cmath>
 
using namespace Belle2;
using namespace TrackFindingCDC;
using namespace TrackingUtilities;
 
RealisticFacetFilter::RealisticFacetFilter()
  : m_param_phiPullCut(11)
{
}
 
RealisticFacetFilter::RealisticFacetFilter(double phiPullCut)
  : m_param_phiPullCut(phiPullCut)
{
}
 
void RealisticFacetFilter::exposeParameters(ModuleParamList* moduleParamList,
                                            const std::string& prefix)
{
  Super::exposeParameters(moduleParamList, prefix);
  moduleParamList->addParameter(prefixed(prefix, "phiPullCut"),
                                m_param_phiPullCut,
                                "Acceptable angle pull in the angle of adjacent tangents to the "
                                "drift circles.",
                                m_param_phiPullCut);
}
 
Weight RealisticFacetFilter::operator()(const CDCFacet& facet)
{
  facet.adjustFitLine();
 
  const CDCRLWireHit& startRLWirehit = facet.getStartRLWireHit();
  const double startDriftLengthVar = startRLWirehit.getRefDriftLengthVariance();
  const double startDriftLengthStd = sqrt(startDriftLengthVar);
 
  const CDCRLWireHit& middleRLWirehit = facet.getMiddleRLWireHit();
  const double middleDriftLengthVar = middleRLWirehit.getRefDriftLengthVariance();
  const double middleDriftLengthStd = sqrt(middleDriftLengthVar);
 
  const CDCRLWireHit& endRLWirehit = facet.getEndRLWireHit();
  const double endDriftLengthVar = endRLWirehit.getRefDriftLengthVariance();
  const double endDriftLengthStd = sqrt(endDriftLengthVar);
 
  const ParameterLine2D& startToMiddleLine = facet.getStartToMiddleLine();
  const ParameterLine2D& startToEndLine = facet.getStartToEndLine();
  const ParameterLine2D& middleToEndLine = facet.getMiddleToEndLine();
 
  const ROOT::Math::XYVector& startToMiddleTangentialVector = startToMiddleLine.tangential();
  const ROOT::Math::XYVector& startToEndTangentialVector = startToEndLine.tangential();
  const ROOT::Math::XYVector& middleToEndTangentialVector = middleToEndLine.tangential();
 
  const double startToMiddleLength = startToMiddleTangentialVector.R();
  const double startToEndLength = startToEndTangentialVector.R();
  const double middleToEndLength = middleToEndTangentialVector.R();
 
  const double startCos = VectorUtil::CosPhi(startToMiddleTangentialVector, startToEndTangentialVector);
  const double middleCos = VectorUtil::CosPhi(startToMiddleTangentialVector, middleToEndTangentialVector);
  const double endCos = VectorUtil::CosPhi(startToEndTangentialVector, middleToEndTangentialVector);
 
  const double startPhi = acos(startCos);
  const double middlePhi = acos(middleCos);
  const double endPhi = acos(endCos);
 
  const double startToMiddleSigmaPhi = startDriftLengthStd / startToMiddleLength;
  const double startToEndSigmaPhi = startDriftLengthStd / startToEndLength;
 
  const double middleToStartSigmaPhi = middleDriftLengthStd / startToMiddleLength;
  const double middleToEndSigmaPhi = middleDriftLengthStd / middleToEndLength;
 
  const double endToStartSigmaPhi = endDriftLengthStd / startToEndLength;
  const double endToMiddleSigmaPhi = endDriftLengthStd / middleToEndLength;
 
  const double startPhiSigma = hypot3(startToEndSigmaPhi - startToMiddleSigmaPhi,
                                      middleToStartSigmaPhi,
                                      endToStartSigmaPhi);
 
  const double middlePhiSigma = hypot3(startToMiddleSigmaPhi,
                                       middleToStartSigmaPhi + middleToEndSigmaPhi,
                                       endToMiddleSigmaPhi);
 
  const double endPhiSigma = hypot3(startToEndSigmaPhi,
                                    middleToEndSigmaPhi,
                                    endToStartSigmaPhi - endToMiddleSigmaPhi);
 
  const double startPhiPull = startPhi / startPhiSigma;
  const double middlePhiPull = middlePhi / middlePhiSigma;
  const double endPhiPull = endPhi / endPhiSigma;
 
  if (startPhiPull < m_param_phiPullCut and
      middlePhiPull < m_param_phiPullCut and
      endPhiPull < m_param_phiPullCut) {
 
    // Introducing a mini penalty to distinguish straighter facets as better
    // This is important since otherwise the ordering of the wires takes precedence and biases
    // for counterclockwise tracks.
    double miniPenalty =
      std::fmin(0.1, (startPhiPull + middlePhiPull + endPhiPull) / m_param_phiPullCut / 1000);
 
    // Good facet contains three points of the track
    // the amount carried by this facet can the adjusted more realistically
    return 3 - miniPenalty;
  } else {
    return NAN;
  }
}