CTPPS/src/TotemT2Segmentation.cc

0001
0002 /****************************************************************************
0003  *
0004  * This is a part of TOTEM offline software.
0005  * Author:
0006  *   Laurent Forthomme
0007  *   Arkadiusz Cwikla
0008  *
0009  ****************************************************************************/
0010
0011 #include "DQM/CTPPS/interface/TotemT2Segmentation.h"
0012 #include "FWCore/Utilities/interface/Exception.h"
0013
0014 #include "TH2D.h"
0015
0016 TotemT2Segmentation::TotemT2Segmentation(size_t nbinsx, size_t nbinsy) : nbinsx_(nbinsx), nbinsy_(nbinsy) {
0017   for (unsigned short arm = 0; arm <= CTPPSDetId::maxArm; ++arm)
0018     for (unsigned short plane = 0; plane <= TotemT2DetId::maxPlane; ++plane)
0019       for (unsigned short id = 0; id <= TotemT2DetId::maxChannel; ++id) {
0020         const TotemT2DetId detid(arm, plane, id);
0021         bins_map_[detid] = computeBins(detid);
0022       }
0023 }
0024
0025 void TotemT2Segmentation::fill(TH2D* hist, const TotemT2DetId& detid, double value) {
0026   if (bins_map_.count(detid) == 0)
0027     throw cms::Exception("TotemT2Segmentation") << "Failed to retrieve list of bins for TotemT2DetId " << detid << ".";
0028   if (static_cast<size_t>(hist->GetXaxis()->GetNbins()) != nbinsx_ ||
0029       static_cast<size_t>(hist->GetYaxis()->GetNbins()) != nbinsy_)
0030     throw cms::Exception("TotemT2Segmentation")
0031         << "Trying to fill a summary plot with invalid number of bins. "
0032         << "Should be of dimension (" << nbinsx_ << ", " << nbinsy_ << "), but has dimension ("
0033         << hist->GetXaxis()->GetNbins() << ", " << hist->GetYaxis()->GetNbins() << ").";
0034   for (const auto& bin : bins_map_.at(detid))
0035     hist->Fill(bin.first, bin.second, value);
0036 }
0037
0038 std::vector<std::pair<short, short> > TotemT2Segmentation::computeBins(const TotemT2DetId& detid) const {
0039   std::vector<std::pair<short, short> > bins;
0040   // find the histogram centre
0041   const auto ox = floor(nbinsx_ * 0.5), oy = floor(nbinsy_ * 0.5);
0042   // compute the ellipse parameters
0043   const auto ax = ceil(nbinsx_ * 0.5), by = ceil(nbinsy_ * 0.5);
0044
0045   const float max_half_angle_rad = 0.3;
0046
0047   // Without use of the Totem Geometry, follow an octant-division of channels,
0048   // with even and odd planes alternating around the circle. Starting at 9AM
0049   // and going clockwise around in increments of 1h30min=45deg, we have
0050   // Ch0 on even planes, Ch0+odd, Ch1+even, Ch1+odd, up to Ch3+odd at 7:30PM
0051
0052   const float tile_angle_rad = (180 - 45. / 2 - (detid.plane() % 2 ? 45 : 0) - (detid.channel() * 90)) * M_PI / 180.;
0053   // Geometric way of associating a DetId to a vector<ix, iy> of bins given the size (nx_, ny_) of
0054   // the TH2D('s) to be filled
0055   for (size_t ix = 0; ix < nbinsx_; ++ix)
0056     for (size_t iy = 0; iy < nbinsy_; ++iy) {
0057       const auto ell_rad_norm = std::pow((ix - ox) / ax, 2) + std::pow((iy - oy) / by, 2);
0058       if (ell_rad_norm < 1. && ell_rad_norm >= 0.1 &&
0059           fabs(std::atan2(iy - oy, ix - ox) - tile_angle_rad) < max_half_angle_rad)
0060         bins.emplace_back(ix, iy);
0061     }
0062
0063   return bins;
0064 }