ProtonReconstructionAlgorithm Class Reference

#include <ProtonReconstructionAlgorithm.h>

Classes
class	ChiSquareCalculator
	class for calculation of chi^2 More...
struct	RPOpticsData
	optics data associated with 1 RP More...

Public Member Functions
void	init (const LHCInterpolatedOpticalFunctionsSetCollection &opticalFunctions)
	ProtonReconstructionAlgorithm (bool fit_vtx_y, bool improved_estimate, const std::string &multiRPAlgorithm, unsigned int verbosity)
reco::ForwardProton	reconstructFromMultiRP (const CTPPSLocalTrackLiteRefVector &tracks, const LHCInfo &lhcInfo, std::ostream &os) const
	run proton reconstruction using multiple-RP strategy More...
reco::ForwardProton	reconstructFromSingleRP (const CTPPSLocalTrackLiteRef &track, const LHCInfo &lhcInfo, std::ostream &os) const
	run proton reconstruction using single-RP strategy More...
void	release ()
	~ProtonReconstructionAlgorithm ()=default

Private Types
enum	{ mrpaUndefined, mrpaChi2, mrpaNewton, mrpaAnalIter }

Static Private Member Functions
static void	doLinearFit (const std::vector< double > &vx, const std::vector< double > &vy, double &b, double &a)
static double	newtonGoalFcn (double xi, double x_N, double x_F, const RPOpticsData &i_N, const RPOpticsData &i_F)

Private Attributes
std::unique_ptr< ChiSquareCalculator >	chiSquareCalculator_
	object to calculate chi^2 More...
std::unique_ptr< ROOT::Fit::Fitter >	fitter_
	fitter object More...
bool	fitVtxY_
bool	initialized_
std::map< unsigned int, RPOpticsData >	m_rp_optics_
	map: RP id –> optics data More...
enum ProtonReconstructionAlgorithm:: { ... }	multi_rp_algorithm_
bool	useImprovedInitialEstimate_
unsigned int	verbosity_

Detailed Description

Definition at line 26 of file ProtonReconstructionAlgorithm.h.

Member Enumeration Documentation

anonymous enum

anonymous enum

private

Enumerator
mrpaUndefined
mrpaChi2
mrpaNewton
mrpaAnalIter

Definition at line 51 of file ProtonReconstructionAlgorithm.h.

{ mrpaUndefined, mrpaChi2, mrpaNewton, mrpaAnalIter } multi_rp_algorithm_;

Constructor & Destructor Documentation

ProtonReconstructionAlgorithm()

ProtonReconstructionAlgorithm::ProtonReconstructionAlgorithm	(	bool	fit_vtx_y,
		bool	improved_estimate,
		const std::string &	multiRPAlgorithm,
		unsigned int	verbosity
	)

Definition at line 22 of file ProtonReconstructionAlgorithm.cc.

References chiSquareCalculator_, Exception, fitter_, mrpaAnalIter, mrpaChi2, mrpaNewton, mrpaUndefined, and multi_rp_algorithm_.

    : verbosity_(verbosity),
      fitVtxY_(fit_vtx_y),
      useImprovedInitialEstimate_(improved_estimate),
      multi_rp_algorithm_(mrpaUndefined),
      initialized_(false),
      fitter_(new ROOT::Fit::Fitter),
      chiSquareCalculator_(new ChiSquareCalculator) {
  // needed for thread safety
  TMinuitMinimizer::UseStaticMinuit(false);

  // check and set multi-RP algorithm
  if (multiRPAlgorithm == "chi2")
    multi_rp_algorithm_ = mrpaChi2;
  else if (multiRPAlgorithm == "newton")
    multi_rp_algorithm_ = mrpaNewton;
  else if (multiRPAlgorithm == "anal-iter")
    multi_rp_algorithm_ = mrpaAnalIter;

  if (multi_rp_algorithm_ == mrpaUndefined)
    throw cms::Exception("ProtonReconstructionAlgorithm") << "Algorithm '" << multiRPAlgorithm << "' not understood.";

  // initialise fitter
  double pStart[] = {0, 0, 0, 0};
  fitter_->SetFCN(4, *chiSquareCalculator_, pStart, 0, true);
}

~ProtonReconstructionAlgorithm()

ProtonReconstructionAlgorithm::~ProtonReconstructionAlgorithm ( )

default

Member Function Documentation

doLinearFit()

void ProtonReconstructionAlgorithm::doLinearFit ( const std::vector< double > &  vx, const std::vector< double > &  vy, double &  b, double &  a  )

staticprivate

Definition at line 99 of file ProtonReconstructionAlgorithm.cc.

References a, b, ztail::d, and mps_fire::i.

Referenced by init().

                                                           {
  double s_1 = 0., s_x = 0., s_xx = 0., s_y = 0., s_xy = 0.;
  for (unsigned int i = 0; i < vx.size(); ++i) {
    s_1 += 1.;
    s_x += vx[i];
    s_xx += vx[i] * vx[i];
    s_y += vy[i];
    s_xy += vx[i] * vy[i];
  }

  const double d = s_xx * s_1 - s_x * s_x;
  a = (s_1 * s_xy - s_x * s_y) / d;
  b = (-s_x * s_xy + s_xx * s_y) / d;
}

init()

void ProtonReconstructionAlgorithm::init

(

const LHCInterpolatedOpticalFunctionsSetCollection &

opticalFunctions

)

Definition at line 54 of file ProtonReconstructionAlgorithm.cc.

References ProtonReconstructionAlgorithm::RPOpticsData::ch0, ProtonReconstructionAlgorithm::RPOpticsData::ch1, doLinearFit(), LHCOpticalFunctionsSet::eLx, LHCOpticalFunctionsSet::eLy, LHCOpticalFunctionsSet::evy, LHCOpticalFunctionsSet::exd, LHCOpticalFunctionsSet::eyd, LHCOpticalFunctionsSet::getFcnValues(), LHCOpticalFunctionsSet::getXiValues(), initialized_, ProtonReconstructionAlgorithm::RPOpticsData::la0, ProtonReconstructionAlgorithm::RPOpticsData::la1, m_rp_optics_, eostools::move(), profile_base_cff::opticalFunctions, ProtonReconstructionAlgorithm::RPOpticsData::optics, AlCaHLTBitMon_ParallelJobs::p, release(), ProtonReconstructionAlgorithm::RPOpticsData::s_L_x_vs_xi, ProtonReconstructionAlgorithm::RPOpticsData::s_L_y_vs_xi, ProtonReconstructionAlgorithm::RPOpticsData::s_v_y_vs_xi, ProtonReconstructionAlgorithm::RPOpticsData::s_x_d_vs_xi, ProtonReconstructionAlgorithm::RPOpticsData::s_xi_vs_x_d, ProtonReconstructionAlgorithm::RPOpticsData::s_y_d_vs_xi, LHCInterpolatedOpticalFunctionsSet::splines(), LHCInterpolatedOpticalFunctionsSet::transport(), LHCInterpolatedOpticalFunctionsSet::Kinematics::x, ProtonReconstructionAlgorithm::RPOpticsData::x0, LHCInterpolatedOpticalFunctionsSet::Kinematics::y, and ProtonReconstructionAlgorithm::RPOpticsData::y0.

Referenced by CTPPSProtonProducer::produce().

                                                                                                             {
  // reset cache
  release();

  // build optics data for each object
  for (const auto &p : opticalFunctions) {
    const LHCInterpolatedOpticalFunctionsSet &ofs = p.second;

    // make record
    RPOpticsData rpod;
    rpod.optics = &p.second;
    rpod.s_x_d_vs_xi = ofs.splines()[LHCOpticalFunctionsSet::exd];
    rpod.s_L_x_vs_xi = ofs.splines()[LHCOpticalFunctionsSet::eLx];
    rpod.s_y_d_vs_xi = ofs.splines()[LHCOpticalFunctionsSet::eyd];
    rpod.s_v_y_vs_xi = ofs.splines()[LHCOpticalFunctionsSet::evy];
    rpod.s_L_y_vs_xi = ofs.splines()[LHCOpticalFunctionsSet::eLy];

    vector<double> xiValues =
        ofs.getXiValues();  // local copy made since the TSpline constructor needs non-const parameters
    vector<double> xDValues = ofs.getFcnValues()[LHCOpticalFunctionsSet::exd];
    rpod.s_xi_vs_x_d = make_shared<TSpline3>("", xDValues.data(), xiValues.data(), xiValues.size());

    // calculate auxiliary data
    LHCInterpolatedOpticalFunctionsSet::Kinematics k_in = {0., 0., 0., 0., 0.};
    LHCInterpolatedOpticalFunctionsSet::Kinematics k_out;
    rpod.optics->transport(k_in, k_out);
    rpod.x0 = k_out.x;
    rpod.y0 = k_out.y;

    doLinearFit(ofs.getXiValues(), ofs.getFcnValues()[LHCOpticalFunctionsSet::exd], rpod.ch0, rpod.ch1);
    rpod.ch0 -= rpod.x0;

    doLinearFit(ofs.getXiValues(), ofs.getFcnValues()[LHCOpticalFunctionsSet::eLx], rpod.la0, rpod.la1);

    // insert record
    const CTPPSDetId rpId(p.first);
    m_rp_optics_.emplace(rpId, std::move(rpod));
  }

  // update settings
  initialized_ = true;
}

newtonGoalFcn()

double ProtonReconstructionAlgorithm::newtonGoalFcn ( double  xi, double  x_N, double  x_F, const RPOpticsData &  i_N, const RPOpticsData &  i_F  )

staticprivate

Definition at line 160 of file ProtonReconstructionAlgorithm.cc.

References ProtonReconstructionAlgorithm::RPOpticsData::s_L_x_vs_xi, ProtonReconstructionAlgorithm::RPOpticsData::s_x_d_vs_xi, and protons_cff::xi.

Referenced by reconstructFromMultiRP().

                                                                                                          {
  return (x_N - i_N.s_x_d_vs_xi->Eval(xi)) * i_F.s_L_x_vs_xi->Eval(xi) -
         (x_F - i_F.s_x_d_vs_xi->Eval(xi)) * i_N.s_L_x_vs_xi->Eval(xi);
}

reconstructFromMultiRP()

reco::ForwardProton ProtonReconstructionAlgorithm::reconstructFromMultiRP	(	const CTPPSLocalTrackLiteRefVector &	tracks,
		const LHCInfo &	lhcInfo,
		std::ostream &	os
	)		const

run proton reconstruction using multiple-RP strategy

Definition at line 171 of file ProtonReconstructionAlgorithm.cc.

References a, funct::abs(), protons_cff::arm, b, c, ProtonReconstructionAlgorithm::RPOpticsData::ch0, ProtonReconstructionAlgorithm::RPOpticsData::ch1, hltPixelTracks_cff::chi2, chiSquareCalculator_, LHCInfo::energy(), fitter_, fitVtxY_, g, runTauDisplay::gp, mps_fire::i, createfilelist::int, dqmiolumiharvest::j, ProtonReconstructionAlgorithm::RPOpticsData::la0, ProtonReconstructionAlgorithm::RPOpticsData::la1, m_rp_optics_, HLT_2022v14_cff::maxIterations, mrpaAnalIter, mrpaChi2, mrpaNewton, multi_rp_algorithm_, newtonGoalFcn(), AlCaHLTBitMon_ParallelJobs::p, submitPVValidationJobs::params, mps_fire::result, ProtonReconstructionAlgorithm::RPOpticsData::s_L_x_vs_xi, ProtonReconstructionAlgorithm::RPOpticsData::s_L_y_vs_xi, ProtonReconstructionAlgorithm::RPOpticsData::s_v_y_vs_xi, ProtonReconstructionAlgorithm::RPOpticsData::s_x_d_vs_xi, ProtonReconstructionAlgorithm::RPOpticsData::s_y_d_vs_xi, mathSSE::sqrt(), HLT_2022v14_cff::track, tracks, useImprovedInitialEstimate_, validateGeometry_cfg::valid, verbosity_, bphysicsOniaDQM_cfi::vertex, protons_cff::xi, and y.

Referenced by CTPPSProtonProducer::produce().

                                                                                                {
  // make sure optics is available for all tracks
  for (const auto &it : tracks) {
    auto oit = m_rp_optics_.find(it->rpId());
    if (oit == m_rp_optics_.end())
      throw cms::Exception("ProtonReconstructionAlgorithm")
          << "Optics data not available for RP " << it->rpId() << ", i.e. " << CTPPSDetId(it->rpId()) << ".";
  }

  // initial estimate of xi and th_x
  double xi_init = 0., th_x_init = 0.;

  if (useImprovedInitialEstimate_) {
    double x_N = tracks[0]->x() * 1E-1,  // conversion: mm --> cm
        x_F = tracks[1]->x() * 1E-1;

    const RPOpticsData &i_N = m_rp_optics_.find(tracks[0]->rpId())->second,
                       &i_F = m_rp_optics_.find(tracks[1]->rpId())->second;

    const double a = i_F.ch1 * i_N.la1 - i_N.ch1 * i_F.la1;
    const double b =
        i_F.ch0 * i_N.la1 - i_N.ch0 * i_F.la1 + i_F.ch1 * i_N.la0 - i_N.ch1 * i_F.la0 + x_N * i_F.la1 - x_F * i_N.la1;
    const double c = x_N * i_F.la0 - x_F * i_N.la0 + i_F.ch0 * i_N.la0 - i_N.ch0 * i_F.la0;
    const double D = b * b - 4. * a * c;
    const double sqrt_D = (D >= 0.) ? sqrt(D) : 0.;

    xi_init = (-b + sqrt_D) / 2. / a;
    th_x_init = (x_N - i_N.ch0 - i_N.ch1 * xi_init) / (i_N.la0 + i_N.la1 * xi_init);
  } else {
    double s_xi0 = 0., s_1 = 0.;
    for (const auto &track : tracks) {
      auto oit = m_rp_optics_.find(track->rpId());
      double xi = oit->second.s_xi_vs_x_d->Eval(track->x() * 1E-1 + oit->second.x0);  // conversion: mm --> cm

      s_1 += 1.;
      s_xi0 += xi;
    }

    xi_init = s_xi0 / s_1;
  }

  if (!std::isfinite(xi_init))
    xi_init = 0.;
  if (!std::isfinite(th_x_init))
    th_x_init = 0.;

  // initial estimate of th_y and vtx_y
  double y[2] = {0}, v_y[2] = {0}, L_y[2] = {0};
  unsigned int y_idx = 0;
  for (const auto &track : tracks) {
    if (y_idx >= 2)
      break;

    auto oit = m_rp_optics_.find(track->rpId());

    y[y_idx] = track->y() * 1E-1 - oit->second.s_y_d_vs_xi->Eval(xi_init);  // track y: mm --> cm
    v_y[y_idx] = oit->second.s_v_y_vs_xi->Eval(xi_init);
    L_y[y_idx] = oit->second.s_L_y_vs_xi->Eval(xi_init);

    y_idx++;
  }

  double vtx_y_init = 0.;
  double th_y_init = 0.;

  if (fitVtxY_) {
    const double det_y = v_y[0] * L_y[1] - L_y[0] * v_y[1];
    vtx_y_init = (det_y != 0.) ? (L_y[1] * y[0] - L_y[0] * y[1]) / det_y : 0.;
    th_y_init = (det_y != 0.) ? (v_y[0] * y[1] - v_y[1] * y[0]) / det_y : 0.;
  } else {
    vtx_y_init = 0.;
    th_y_init = (y[1] / L_y[1] + y[0] / L_y[0]) / 2.;
  }

  if (!std::isfinite(vtx_y_init))
    vtx_y_init = 0.;
  if (!std::isfinite(th_y_init))
    th_y_init = 0.;

  unsigned int armId = CTPPSDetId((*tracks.begin())->rpId()).arm();

  if (verbosity_)
    os << "ProtonReconstructionAlgorithm::reconstructFromMultiRP(" << armId << ")" << std::endl
       << "    initial estimate: xi_init = " << xi_init << ", th_x_init = " << th_x_init
       << ", th_y_init = " << th_y_init << ", vtx_y_init = " << vtx_y_init << "." << std::endl;

  // prepare result containers
  bool valid = false;
  double chi2 = 0.;
  double xi = 0., th_x = 0., th_y = 0., vtx_y = 0.;

  using FP = reco::ForwardProton;
  FP::CovarianceMatrix cm;

  if (multi_rp_algorithm_ == mrpaChi2) {
    // minimisation
    fitter_->Config().ParSettings(0).Set("xi", xi_init, 0.005);
    fitter_->Config().ParSettings(1).Set("th_x", th_x_init, 2E-6);
    fitter_->Config().ParSettings(2).Set("th_y", th_y_init, 2E-6);
    fitter_->Config().ParSettings(3).Set("vtx_y", vtx_y_init, 10E-6);

    if (!fitVtxY_)
      fitter_->Config().ParSettings(3).Fix();

    chiSquareCalculator_->tracks = &tracks;
    chiSquareCalculator_->m_rp_optics = &m_rp_optics_;

    fitter_->FitFCN();
    fitter_->FitFCN();  // second minimisation in case the first one had troubles

    // extract proton parameters
    const ROOT::Fit::FitResult &result = fitter_->Result();
    const double *params = result.GetParams();

    valid = result.IsValid();
    chi2 = result.Chi2();
    xi = params[0];
    th_x = params[1];
    th_y = params[2];
    vtx_y = params[3];

    map<unsigned int, signed int> index_map = {
        {(unsigned int)FP::Index::xi, 0},
        {(unsigned int)FP::Index::th_x, 1},
        {(unsigned int)FP::Index::th_y, 2},
        {(unsigned int)FP::Index::vtx_y, ((fitVtxY_) ? 3 : -1)},
        {(unsigned int)FP::Index::vtx_x, -1},
    };

    for (unsigned int i = 0; i < (unsigned int)FP::Index::num_indices; ++i) {
      signed int fit_i = index_map[i];

      for (unsigned int j = 0; j < (unsigned int)FP::Index::num_indices; ++j) {
        signed int fit_j = index_map[j];

        cm(i, j) = (fit_i >= 0 && fit_j >= 0) ? result.CovMatrix(fit_i, fit_j) : 0.;
      }
    }
  }

  else if (multi_rp_algorithm_ == mrpaNewton || multi_rp_algorithm_ == mrpaAnalIter) {
    // settings
    const unsigned int maxIterations = 100;
    const double maxXiDiff = 1E-6;
    const double xi_ep = 1E-5;

    // collect input
    double x_N = tracks[0]->x() * 1E-1,  // conversion: mm --> cm
        x_F = tracks[1]->x() * 1E-1;

    double y_N = tracks[0]->y() * 1E-1,  // conversion: mm --> cm
        y_F = tracks[1]->y() * 1E-1;

    const RPOpticsData &i_N = m_rp_optics_.find(tracks[0]->rpId())->second,
                       &i_F = m_rp_optics_.find(tracks[1]->rpId())->second;

    // horizontal reconstruction - run iterations
    valid = false;
    double xi_prev = xi_init;
    for (unsigned int it = 0; it < maxIterations; ++it) {
      if (multi_rp_algorithm_ == mrpaNewton) {
        const double g = newtonGoalFcn(xi_prev, x_N, x_F, i_N, i_F);
        const double gp = (newtonGoalFcn(xi_prev + xi_ep, x_N, x_F, i_N, i_F) - g) / xi_ep;

        xi = xi_prev - g / gp;
      }

      if (multi_rp_algorithm_ == mrpaAnalIter) {
        const double d_x_eff_N = i_N.s_x_d_vs_xi->Eval(xi_prev) / xi_prev;
        const double d_x_eff_F = i_F.s_x_d_vs_xi->Eval(xi_prev) / xi_prev;

        const double l_x_N = i_N.s_L_x_vs_xi->Eval(xi_prev);
        const double l_x_F = i_F.s_L_x_vs_xi->Eval(xi_prev);

        xi = (l_x_N * x_F - l_x_F * x_N) / (l_x_N * d_x_eff_F - l_x_F * d_x_eff_N);
      }

      if (abs(xi - xi_prev) < maxXiDiff) {
        valid = true;
        break;
      }

      xi_prev = xi;
    }

    th_x = (x_F - i_F.s_x_d_vs_xi->Eval(xi)) / i_F.s_L_x_vs_xi->Eval(xi);

    // vertical reconstruction
    const double y_eff_N = y_N - i_N.s_y_d_vs_xi->Eval(xi);
    const double y_eff_F = y_F - i_F.s_y_d_vs_xi->Eval(xi);

    const double l_y_N = i_N.s_L_y_vs_xi->Eval(xi);
    const double l_y_F = i_F.s_L_y_vs_xi->Eval(xi);

    const double v_y_N = i_N.s_v_y_vs_xi->Eval(xi);
    const double v_y_F = i_F.s_v_y_vs_xi->Eval(xi);

    const double det = l_y_N * v_y_F - l_y_F * v_y_N;

    if (fitVtxY_) {
      th_y = (y_eff_N * v_y_F - y_eff_F * v_y_N) / det;
      vtx_y = (l_y_N * y_eff_F - l_y_F * y_eff_N) / det;
    } else {
      th_y = (y_eff_N / l_y_N + y_eff_F / l_y_F) / 2.;
      vtx_y = 0.;
    }

    // covariance matrix kept empty - not to compromise the speed of these special algorithms
  }

  // print results
  if (verbosity_)
    os << "    fit valid=" << valid << std::endl
       << "    xi=" << xi << ", th_x=" << th_x << ", th_y=" << th_y << ", vtx_y=" << vtx_y << ", chiSq = " << chi2
       << std::endl;

  // save reco candidate
  const double sign_z = (armId == 0) ? +1. : -1.;  // CMS convention
  const FP::Point vertex(0., vtx_y, 0.);
  const double cos_th_sq = 1. - th_x * th_x - th_y * th_y;
  const double cos_th = (cos_th_sq > 0.) ? sqrt(cos_th_sq) : 1.;
  const double p = lhcInfo.energy() * (1. - xi);
  const FP::Vector momentum(-p * th_x,  // the signs reflect change LHC --> CMS convention
                            +p * th_y,
                            sign_z * p * cos_th);
  signed int ndf = 2. * tracks.size() - ((fitVtxY_) ? 4. : 3.);

  return reco::ForwardProton(chi2, ndf, vertex, momentum, xi, cm, FP::ReconstructionMethod::multiRP, tracks, valid);
}

reconstructFromSingleRP()

reco::ForwardProton ProtonReconstructionAlgorithm::reconstructFromSingleRP	(	const CTPPSLocalTrackLiteRef &	track,
		const LHCInfo &	lhcInfo,
		std::ostream &	os
	)		const

run proton reconstruction using single-RP strategy

Definition at line 405 of file ProtonReconstructionAlgorithm.cc.

References funct::abs(), protons_cff::arm, CTPPSDetId::arm(), LHCInfo::energy(), m_rp_optics_, AlCaHLTBitMon_ParallelJobs::p, funct::pow(), edm::RefVector< C, T, F >::push_back(), mathSSE::sqrt(), HLT_2022v14_cff::track, verbosity_, bphysicsOniaDQM_cfi::vertex, and protons_cff::xi.

Referenced by CTPPSProtonProducer::produce().

                                                                                                 {
  CTPPSDetId rpId(track->rpId());

  if (verbosity_)
    os << "reconstructFromSingleRP(" << rpId.arm() * 100 + rpId.station() * 10 + rpId.rp() << ")" << std::endl;

  // make sure optics is available for the track
  auto oit = m_rp_optics_.find(track->rpId());
  if (oit == m_rp_optics_.end())
    throw cms::Exception("ProtonReconstructionAlgorithm")
        << "Optics data not available for RP " << track->rpId() << ", i.e. " << rpId << ".";

  // rough estimate of xi and th_y from each track
  const double x_full = track->x() * 1E-1 + oit->second.x0;  // conversion mm --> cm
  const double xi = oit->second.s_xi_vs_x_d->Eval(x_full);
  const double L_y = oit->second.s_L_y_vs_xi->Eval(xi);
  const double th_y = track->y() * 1E-1 / L_y;  // conversion mm --> cm

  const double ep_x = 1E-6;
  const double dxi_dx = (oit->second.s_xi_vs_x_d->Eval(x_full + ep_x) - xi) / ep_x;
  const double xi_unc = abs(dxi_dx) * track->xUnc() * 1E-1;  // conversion mm --> cm

  const double ep_xi = 1E-4;
  const double dL_y_dxi = (oit->second.s_L_y_vs_xi->Eval(xi + ep_xi) - L_y) / ep_xi;
  const double th_y_unc_sq = th_y * th_y * (pow(track->yUnc() / track->y(), 2.) + pow(dL_y_dxi * xi_unc / L_y, 2.));

  if (verbosity_)
    os << "    xi = " << xi << " +- " << xi_unc << ", th_y = " << th_y << " +- " << sqrt(th_y_unc_sq) << "."
       << std::endl;

  using FP = reco::ForwardProton;

  // save proton candidate
  const double sign_z = (CTPPSDetId(track->rpId()).arm() == 0) ? +1. : -1.;  // CMS convention
  const FP::Point vertex(0., 0., 0.);
  const double cos_th_sq = 1. - th_y * th_y;
  const double cos_th = (cos_th_sq > 0.) ? sqrt(cos_th_sq) : 1.;
  const double p = lhcInfo.energy() * (1. - xi);
  const FP::Vector momentum(0., p * th_y, sign_z * p * cos_th);

  FP::CovarianceMatrix cm;
  cm((int)FP::Index::xi, (int)FP::Index::xi) = xi_unc * xi_unc;
  cm((int)FP::Index::th_y, (int)FP::Index::th_y) = th_y_unc_sq;

  CTPPSLocalTrackLiteRefVector trk;
  trk.push_back(track);

  return reco::ForwardProton(0., 0, vertex, momentum, xi, cm, FP::ReconstructionMethod::singleRP, trk, true);
}

release()

void ProtonReconstructionAlgorithm::release

(

)

Definition at line 119 of file ProtonReconstructionAlgorithm.cc.

References initialized_, and m_rp_optics_.

Referenced by init(), and CTPPSProtonProducer::produce().

                                            {
  initialized_ = false;

  m_rp_optics_.clear();
}

Member Data Documentation

chiSquareCalculator_

std::unique_ptr<ChiSquareCalculator> ProtonReconstructionAlgorithm::chiSquareCalculator_

private

object to calculate chi^2

Definition at line 84 of file ProtonReconstructionAlgorithm.h.

Referenced by ProtonReconstructionAlgorithm(), and reconstructFromMultiRP().

fitter_

std::unique_ptr<ROOT::Fit::Fitter> ProtonReconstructionAlgorithm::fitter_

private

fitter object

Definition at line 81 of file ProtonReconstructionAlgorithm.h.

Referenced by ProtonReconstructionAlgorithm(), and reconstructFromMultiRP().

fitVtxY_

bool ProtonReconstructionAlgorithm::fitVtxY_

private

Definition at line 49 of file ProtonReconstructionAlgorithm.h.

Referenced by reconstructFromMultiRP().

initialized_

bool ProtonReconstructionAlgorithm::initialized_

private

Definition at line 52 of file ProtonReconstructionAlgorithm.h.

Referenced by init(), and release().

m_rp_optics_

std::map<unsigned int, RPOpticsData> ProtonReconstructionAlgorithm::m_rp_optics_

private

map: RP id –> optics data

Definition at line 67 of file ProtonReconstructionAlgorithm.h.

Referenced by init(), reconstructFromMultiRP(), reconstructFromSingleRP(), and release().

multi_rp_algorithm_

enum { ... } ProtonReconstructionAlgorithm::multi_rp_algorithm_

Referenced by ProtonReconstructionAlgorithm(), and reconstructFromMultiRP().

useImprovedInitialEstimate_

bool ProtonReconstructionAlgorithm::useImprovedInitialEstimate_

private

Definition at line 50 of file ProtonReconstructionAlgorithm.h.

Referenced by reconstructFromMultiRP().

verbosity_

unsigned int ProtonReconstructionAlgorithm::verbosity_

private

Definition at line 48 of file ProtonReconstructionAlgorithm.h.

Referenced by reconstructFromMultiRP(), and reconstructFromSingleRP().