SurveyAnalysis/src/SurveyPxbDicer.cc

0001 #include "Alignment/SurveyAnalysis/interface/SurveyPxbImage.h"
0002 #include "Alignment/SurveyAnalysis/interface/SurveyPxbDicer.h"
0003
0004 //#include <stdexcept>
0005 #include <vector>
0006 #include <cmath>
0007 #include <string>
0008 #include <sstream>
0009 #include <functional>
0010 #include <algorithm>
0011 #include <fstream>
0012 #include "DataFormats/GeometryVector/interface/Point3DBase.h"
0013 #include "Math/SMatrix.h"
0014 #include "Math/SVector.h"
0015 #include "CLHEP/Random/RandGauss.h"
0016 #include "CLHEP/Random/Random.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018
0019 #include <iostream>
0020
0021 SurveyPxbDicer::SurveyPxbDicer(const std::vector<edm::ParameterSet> &pars, unsigned int seed) {
0022   CLHEP::HepRandom::setTheSeed(seed);
0023   mean_a0 = getParByName("a0", "mean", pars);
0024   sigma_a0 = getParByName("a0", "sigma", pars);
0025   mean_a1 = getParByName("a1", "mean", pars);
0026   sigma_a1 = getParByName("a1", "sigma", pars);
0027   mean_scale = getParByName("scale", "mean", pars);
0028   sigma_scale = getParByName("scale", "sigma", pars);
0029   mean_phi = getParByName("phi", "mean", pars);
0030   sigma_phi = getParByName("phi", "sigma", pars);
0031   mean_x = getParByName("x", "mean", pars);
0032   sigma_x = getParByName("x", "sigma", pars);
0033   mean_y = getParByName("y", "mean", pars);
0034   sigma_y = getParByName("y", "sigma", pars);
0035 }
0036
0037 std::string SurveyPxbDicer::doDice(const fidpoint_t &fidpointvec, const idPair_t &id, const bool rotate) {
0038   // Dice the local parameters
0039   const value_t a0 = ranGauss(mean_a0, sigma_a0);
0040   const value_t a1 = ranGauss(mean_a1, sigma_a1);
0041   const value_t scale = ranGauss(mean_scale, sigma_scale);
0042   const value_t phi = ranGauss(mean_phi, sigma_phi);
0043   const value_t a2 = scale * cos(phi);
0044   const value_t a3 = scale * sin(phi);
0045   const coord_t p0 = transform(fidpointvec[0], a0, a1, a2, a3);
0046   const coord_t p1 = transform(fidpointvec[1], a0, a1, a2, a3);
0047   const coord_t p2 = transform(fidpointvec[2], a0, a1, a2, a3);
0048   const coord_t p3 = transform(fidpointvec[3], a0, a1, a2, a3);
0049   const value_t sign = rotate ? -1 : 1;
0050   std::ostringstream oss;
0051   oss << id.first << " " << sign * ranGauss(p0.x(), sigma_x) << " " << -sign * ranGauss(p0.y(), sigma_y) << " "
0052       << sign * ranGauss(p1.x(), sigma_x) << " " << -sign * ranGauss(p1.y(), sigma_y) << " " << id.second << " "
0053       << sign * ranGauss(p2.x(), sigma_x) << " " << -sign * ranGauss(p2.y(), sigma_y) << " "
0054       << sign * ranGauss(p3.x(), sigma_x) << " " << -sign * ranGauss(p3.y(), sigma_y) << " " << sigma_x << " "
0055       << sigma_y << " " << rotate << " # MC-truth:"
0056       << " a0-a3: " << a0 << " " << a1 << " " << a2 << " " << a3 << " S: " << scale << " phi: " << phi
0057       << " x0: " << fidpointvec[0].x() << " " << fidpointvec[0].y() << " x1: " << fidpointvec[1].x() << " "
0058       << fidpointvec[1].y() << " x2: " << fidpointvec[2].x() << " " << fidpointvec[2].y()
0059       << " x3: " << fidpointvec[3].x() << " " << fidpointvec[3].y() << std::endl;
0060   return oss.str();
0061 }
0062
0063 void SurveyPxbDicer::doDice(const fidpoint_t &fidpointvec,
0064                             const idPair_t &id,
0065                             std::ofstream &outfile,
0066                             const bool rotate) {
0067   outfile << doDice(fidpointvec, id, rotate);
0068 }
0069
0070 SurveyPxbDicer::coord_t SurveyPxbDicer::transform(
0071     const coord_t &x, const value_t &a0, const value_t &a1, const value_t &a2, const value_t &a3) {
0072   return coord_t(a0 + a2 * x.x() + a3 * x.y(), a1 - a3 * x.x() + a2 * x.y());
0073 }
0074
0075 SurveyPxbDicer::value_t SurveyPxbDicer::getParByName(const std::string &name,
0076                                                      const std::string &par,
0077                                                      const std::vector<edm::ParameterSet> &pars) {
0078   std::vector<edm::ParameterSet>::const_iterator it;
0079   it = std::find_if(pars.begin(), pars.end(), [&name](auto const &c) { return findParByName()(name, c); });
0080   if (it == pars.end()) {
0081     throw std::runtime_error("Parameter not found in SurveyPxbDicer::getParByName");
0082   }
0083   return (*it).getParameter<value_t>(par);
0084 }