LocalTrackFitter Class Reference

Performs straight-line fit and outlier rejection. More...

#include <LocalTrackFitter.h>

Public Member Functions
bool	fit (HitCollection &, const AlignmentGeometry &, LocalTrackFit &) const
	LocalTrackFitter ()
	dummy constructor (not to be used) More...
	LocalTrackFitter (const edm::ParameterSet &)
	normal constructor More...
virtual	~LocalTrackFitter ()

Protected Member Functions
void	fitAndRemoveOutliers (HitCollection &, const AlignmentGeometry &, LocalTrackFit &, bool &failed, bool &selectionChanged) const
void	removeInsufficientPots (HitCollection &, bool &selectionChanged) const
	removes the hits of pots with too few planes active More...

Protected Attributes
double	maxResidualToSigma
	hits with higher ratio residual/sigma will be dropped More...
unsigned int	minimumHitsPerProjectionPerRP
	minimum of hits to accept data from a RP More...
unsigned int	verbosity
	verbosity level More...

Detailed Description

Performs straight-line fit and outlier rejection.

Definition at line 20 of file LocalTrackFitter.h.

Constructor & Destructor Documentation

LocalTrackFitter() [1/2]

LocalTrackFitter::LocalTrackFitter ( )

inline

dummy constructor (not to be used)

Definition at line 23 of file LocalTrackFitter.h.

{}

LocalTrackFitter() [2/2]

LocalTrackFitter::LocalTrackFitter

(

const edm::ParameterSet &

ps

)

normal constructor

Definition at line 24 of file LocalTrackFitter.cc.

    : verbosity(ps.getUntrackedParameter<unsigned int>("verbosity", 0)),
      minimumHitsPerProjectionPerRP(ps.getParameter<unsigned int>("minimumHitsPerProjectionPerRP")),
      maxResidualToSigma(ps.getParameter<double>("maxResidualToSigma")) {}

~LocalTrackFitter()

virtual LocalTrackFitter::~LocalTrackFitter ( )

inlinevirtual

Definition at line 28 of file LocalTrackFitter.h.

{}

Member Function Documentation

fit()

bool LocalTrackFitter::fit	(	HitCollection &	selection,
		const AlignmentGeometry &	geometry,
		LocalTrackFit &	trackFit
	)		const

runs the fit and outlier-removal loop returns true in case of success

Definition at line 31 of file LocalTrackFitter.cc.

References fitAndRemoveOutliers(), removeInsufficientPots(), corrVsCorr::selection, and verbosity.

Referenced by trackingPlots.Iteration::modules(), and StraightTrackAlignment::processEvent().

                                                                                                                     {
  if (verbosity > 5)
    printf(">> LocalTrackFitter::Fit\n");

  bool selectionChanged = true;
  unsigned int loopCounter = 0;
  while (selectionChanged) {
    // fit/outlier-removal loop
    while (selectionChanged) {
      if (verbosity > 5)
        printf("* fit loop %u\n", loopCounter++);

      bool fitFailed = false;
      fitAndRemoveOutliers(selection, geometry, trackFit, fitFailed, selectionChanged);

      if (fitFailed) {
        if (verbosity > 5)
          printf("\tFIT FAILED\n");
        return false;
      }
    }

    // remove pots with too few active planes
    removeInsufficientPots(selection, selectionChanged);
  }

  return true;
}

fitAndRemoveOutliers()

void LocalTrackFitter::fitAndRemoveOutliers ( HitCollection &  selection, const AlignmentGeometry &  geometry, LocalTrackFit &  trackFit, bool &  failed, bool &  selectionChanged  ) const

protected

fits the collection of hits and removes hits with too high residual/sigma ratio

Parameters

failed	whether the fit has failed
selectionChanged	whether some hits have been removed

Definition at line 62 of file LocalTrackFitter.cc.

References A, LocalTrackFit::ax, LocalTrackFit::ay, LocalTrackFit::bx, LocalTrackFit::by, LocalTrackFit::chi_sq, ztail::d, MillePedeFileConverter_cfg::e, runEdmFileComparison::failed, mps_fire::i, dqmiolumiharvest::j, maxResidualToSigma, LocalTrackFit::ndf, printId(), dttmaxenums::R, corrVsCorr::selection, theta(), verbosity, and LocalTrackFit::z0.

Referenced by fit().

                                                                          {
  if (verbosity > 5)
    printf(" - LocalTrackFitter::FitAndRemoveOutliers\n");

  if (selection.empty()) {
    failed = true;
    return;
  }

  // check if input size is sufficient
  if (selection.size() < 4) {
    failed = true;
    return;
  }

  // build matrices and vectors
  TMatrixD A(selection.size(), 4);
  TMatrixD Vi(selection.size(), selection.size());
  TVectorD measVec(selection.size());
  unsigned int j = 0;
  for (auto it = selection.begin(); it != selection.end(); ++it, ++j) {
    const unsigned int &detId = it->id;

    const DetGeometry &d = geometry.get(detId);
    const auto &dirData = d.getDirectionData(it->dirIdx);

    A(j, 0) = it->z * dirData.dx;
    A(j, 1) = dirData.dx;
    A(j, 2) = it->z * dirData.dy;
    A(j, 3) = dirData.dy;

    measVec(j) = it->position + dirData.s - (it->z - d.z) * dirData.dz;  // in mm

    Vi(j, j) = 1. / it->sigma / it->sigma;
  }

  // evaluate local track parameter estimates (h stands for hat)
  TMatrixD AT(4, selection.size());
  AT.Transpose(A);
  TMatrixD ATViA(4, 4);
  ATViA = AT * Vi * A;
  TMatrixD ATViAI(ATViA);

  try {
    ATViAI = ATViA.Invert();
  } catch (cms::Exception &e) {
    failed = true;
    return;
  }

  TVectorD theta(4);
  theta = ATViAI * AT * Vi * measVec;

  // residuals
  TVectorD R(measVec);
  R -= A * theta;

  // save results to trackFit
  trackFit.ax = theta(0);
  trackFit.bx = theta(1);
  trackFit.ay = theta(2);
  trackFit.by = theta(3);
  trackFit.z0 = geometry.z0;
  trackFit.ndf = selection.size() - 4;
  trackFit.chi_sq = 0;
  for (int i = 0; i < R.GetNrows(); i++)
    trackFit.chi_sq += R(i) * R(i) * Vi(i, i);

  if (verbosity > 5) {
    printf("    ax = %.3f mrad, bx = %.4f mm, ay = %.3f mrad, by = %.4f mm, z0 = %.3f mm\n",
           trackFit.ax * 1E3,
           trackFit.bx,
           trackFit.ay * 1E3,
           trackFit.by,
           trackFit.z0);
    printf("    ndof = %i, chi^2/ndof/si^2 = %.3f\n", trackFit.ndf, trackFit.chi_sq / trackFit.ndf);
  }

  // check residuals
  selectionChanged = false;
  TVectorD interpolation(A * theta);
  j = 0;
  for (auto it = selection.begin(); it != selection.end(); ++j) {
    if (verbosity > 5) {
      printf("        %2u, ", j);
      printId(it->id);
      printf(", dirIdx=%u: interpol = %+8.1f um, residual = %+6.1f um, residual / sigma = %+6.2f\n",
             it->dirIdx,
             interpolation[j] * 1E3,
             R[j] * 1E3,
             R[j] / it->sigma);
    }

    double resToSigma = R[j] / it->sigma;
    if (fabs(resToSigma) > maxResidualToSigma) {
      selection.erase(it);
      selectionChanged = true;
      if (verbosity > 5)
        printf("            Removed\n");
    } else
      ++it;
  }
}

removeInsufficientPots()

void LocalTrackFitter::removeInsufficientPots ( HitCollection &  selection, bool &  selectionChanged  ) const

protected

removes the hits of pots with too few planes active

Definition at line 172 of file LocalTrackFitter.cc.

References minimumHitsPerProjectionPerRP, TotemRPDetId::plane(), CTPPSDetId::rpId(), CTPPSDetId::sdTrackingPixel, CTPPSDetId::sdTrackingStrip, corrVsCorr::selection, DetId::subdetId(), and verbosity.

Referenced by fit().

                                                                                                    {
  if (verbosity > 5)
    printf(" - RemoveInsufficientPots\n");

  selectionChanged = false;

  // map: RP id -> (active planes in projection 1, active planes in projection 2)
  map<unsigned int, pair<set<unsigned int>, set<unsigned int> > > planeMap;
  for (auto it = selection.begin(); it != selection.end(); ++it) {
    CTPPSDetId senId(it->id);
    const unsigned int rpId = senId.rpId();

    if (senId.subdetId() == CTPPSDetId::sdTrackingStrip) {
      const unsigned int plane = TotemRPDetId(it->id).plane();
      if ((plane % 2) == 0)
        planeMap[rpId].first.insert(senId);
      else
        planeMap[rpId].second.insert(senId);
    }

    if (senId.subdetId() == CTPPSDetId::sdTrackingPixel) {
      planeMap[rpId].first.insert(senId);
      planeMap[rpId].second.insert(senId);
    }
  }

  // remove RPs with insufficient information
  selectionChanged = false;

  for (auto it = planeMap.begin(); it != planeMap.end(); ++it) {
    if (it->second.first.size() < minimumHitsPerProjectionPerRP ||
        it->second.second.size() < minimumHitsPerProjectionPerRP) {
      if (verbosity > 5)
        printf("\tRP %u: projection1 = %lu, projection 2 = %lu\n",
               it->first,
               it->second.first.size(),
               it->second.second.size());

      // remove all hits from that RP
      for (auto dit = selection.begin(); dit != selection.end();) {
        if (it->first == CTPPSDetId(dit->id).rpId()) {
          if (verbosity > 5)
            printf("\t\tremoving %u\n", dit->id);
          selection.erase(dit);
          selectionChanged = true;
        } else {
          dit++;
        }
      }
    }
  }
}

Member Data Documentation

maxResidualToSigma

double LocalTrackFitter::maxResidualToSigma

protected

hits with higher ratio residual/sigma will be dropped

Definition at line 42 of file LocalTrackFitter.h.

Referenced by fitAndRemoveOutliers().

minimumHitsPerProjectionPerRP

unsigned int LocalTrackFitter::minimumHitsPerProjectionPerRP

protected

minimum of hits to accept data from a RP

Definition at line 39 of file LocalTrackFitter.h.

Referenced by removeInsufficientPots().

verbosity

unsigned int LocalTrackFitter::verbosity

protected

verbosity level

Definition at line 36 of file LocalTrackFitter.h.

Referenced by fit(), fitAndRemoveOutliers(), and removeInsufficientPots().