packages/CaloBase/RawTowerGeomContainer_Cylinderv1.cc

0001 #include "RawTowerGeomContainer_Cylinderv1.h"
0002
0003 #include <cassert>
0004 #include <cmath>
0005 #include <cstdlib>
0006 #include <iostream>
0007 #include <memory>
0008
0009 RawTowerGeomContainer_Cylinderv1::RawTowerGeomContainer_Cylinderv1(RawTowerDefs::CalorimeterId caloid)
0010   : RawTowerGeomContainerv1(caloid)
0011 {
0012   return;
0013 }
0014
0015 void RawTowerGeomContainer_Cylinderv1::Reset()
0016 {
0017   eta_bound_map.clear();
0018
0019   phi_bound_map.clear();
0020
0021   RawTowerGeomContainerv1::Reset();
0022 }
0023
0024 void RawTowerGeomContainer_Cylinderv1::set_etabins(const int i)
0025 {
0026   assert(i > 0);
0027   bound_t invalid_bound(std::numeric_limits<double>::signaling_NaN(), std::numeric_limits<double>::signaling_NaN());
0028   eta_bound_map.resize(i, invalid_bound);
0029 }
0030
0031 void RawTowerGeomContainer_Cylinderv1::set_phibins(const int i)
0032 {
0033   assert(i > 0);
0034   bound_t invalid_bound(std::numeric_limits<double>::signaling_NaN(), std::numeric_limits<double>::signaling_NaN());
0035   phi_bound_map.resize(i, invalid_bound);
0036 }
0037
0038 void RawTowerGeomContainer_Cylinderv1::identify(std::ostream& os) const
0039 {
0040   os << "RawTowerGeomContainer_Cylinderv1: radius: " << radius
0041      << ", thickness: " << thickness << ", etabins: " << get_etabins()
0042      << ", phibins: " << get_phibins();
0043
0044   int i = 0;
0045   for (const auto& iter : eta_bound_map)
0046   {
0047     os << "eta_bin[" << i << "](" << iter.first << ", " << iter.second
0048        << ")  ";
0049     i++;
0050   }
0051   os << std::endl;
0052   i = 0;
0053   for (const auto& iter : phi_bound_map)
0054   {
0055     os << "phi_bin[" << i << "](" << iter.first << ", " << iter.second
0056        << ")  ";
0057     i++;
0058   }
0059   os << std::endl;
0060   return;
0061 }
0062
0063 std::pair<double, double>
0064 RawTowerGeomContainer_Cylinderv1::get_etabounds(const int ibin) const
0065 {
0066   if (ibin < 0 || ibin > get_etabins())
0067   {
0068     identify();
0069     std::cout
0070         << "RawTowerGeomContainer_Cylinderv1::get_etabounds - Asking for invalid bin in eta: "
0071         << ibin << std::endl;
0072     exit(1);
0073   }
0074   return eta_bound_map[ibin];
0075 }
0076
0077 std::pair<double, double>
0078 RawTowerGeomContainer_Cylinderv1::get_phibounds(const int ibin) const
0079 {
0080   if (ibin < 0 || ibin > get_phibins())
0081   {
0082     identify();
0083     std::cout
0084         << "RawTowerGeomContainer_Cylinderv1::get_phibounds - Asking for invalid bin in phi: "
0085         << ibin << std::endl;
0086     exit(1);
0087   }
0088   return phi_bound_map[ibin];
0089 }
0090
0091 int RawTowerGeomContainer_Cylinderv1::get_etabin(const double eta) const
0092 {
0093   int ibin = -1;
0094   int i = 0;
0095
0096   // switch to search for the closest bin
0097   // since in a realistic calorimeter, there could be gaps
0098   double min_deta = 10;
0099
0100   for (const auto& iter : eta_bound_map)
0101   {
0102     const double mean_eta = 0.5 * (iter.first + iter.second);
0103
0104     if (eta >= iter.first && eta < iter.second)
0105     {
0106       // found the bin that the hit belong
0107       min_deta = 0;
0108       ibin = i;
0109       break;
0110     }
0111     else
0112     {
0113       const double deta = fabs(mean_eta - eta);
0114       if (deta < min_deta)
0115       {
0116         min_deta = deta;
0117         ibin = i;
0118       }  // keep searching
0119     }
0120
0121     i++;
0122   }
0123
0124   if (ibin < 0)
0125   {
0126     std::cout
0127         << "RawTowerGeomContainer_Cylinderv1::get_etabin - ERROR - Asking for invalid bin in eta "
0128         << eta << std::endl;
0129     exit(1);
0130   }
0131
0132   return ibin;
0133 }
0134
0135 int RawTowerGeomContainer_Cylinderv1::get_phibin(const double phi) const
0136 {
0137   int ibin = -1;
0138   int i = 0;
0139
0140   // switch to search for the closest bin
0141   // since in a realistic calorimeter, there could be gaps
0142   double min_dphi = 10;
0143
0144   for (const auto& iter : phi_bound_map)
0145   {
0146     const double mean_phi = 0.5 * (iter.first + iter.second);
0147
0148     const double phi_fold = phi - round((phi - mean_phi) / 2. / M_PI) * 2 * M_PI;
0149
0150     if (phi_fold >= iter.first && phi_fold < iter.second)
0151     {
0152       // found the bin that the hit belong
0153       min_dphi = 0;
0154       ibin = i;
0155       break;
0156     }
0157     else
0158     {
0159       const double dphi = fabs(mean_phi - phi_fold);
0160       if (dphi < min_dphi)
0161       {
0162         min_dphi = dphi;
0163         ibin = i;
0164       }  // keep searching
0165     }
0166
0167     i++;
0168   }
0169
0170   if (ibin < 0)
0171   {
0172     std::cout
0173         << "RawTowerGeomContainer_Cylinderv1::get_phibin - ERROR - Asking for invalid bin in phi "
0174         << phi << std::endl;
0175     exit(1);
0176   }
0177
0178   return ibin;
0179 }
0180
0181 double
0182 RawTowerGeomContainer_Cylinderv1::get_etacenter(const int ibin) const
0183 {
0184   if (ibin < 0 || ibin >= get_etabins())
0185   {
0186     std::cout
0187         << "RawTowerGeomContainer_Cylinderv1::get_etacenter - Asking for invalid bin in eta: "
0188         << ibin << std::endl;
0189     std::cout << "minbin: 0, maxbin " << get_etabins() << std::endl;
0190     exit(1);
0191   }
0192   return (eta_bound_map[ibin].first + eta_bound_map[ibin].second) / 2.;
0193 }
0194
0195 void RawTowerGeomContainer_Cylinderv1::set_etabounds(const int ibin,
0196                                                      const std::pair<double, double>& bounds)
0197 {
0198   if (ibin < 0 || ibin >= get_etabins())
0199   {
0200     std::cout
0201         << "RawTowerGeomContainer_Cylinderv1::set_bounds - Asking for invalid bin in eta: "
0202         << ibin << std::endl;
0203     std::cout << "minbin: 0, maxbin " << get_etabins() << std::endl;
0204     exit(1);
0205   }
0206
0207   std::pair<double, double> b_reg(bounds);
0208   if (b_reg.first > b_reg.second)
0209   {
0210     b_reg.second = bounds.first;
0211     b_reg.first = bounds.second;
0212   }
0213
0214   eta_bound_map[ibin] = b_reg;
0215 }
0216
0217 double
0218 RawTowerGeomContainer_Cylinderv1::get_phicenter(const int ibin) const
0219 {
0220   if (ibin < 0 || ibin >= get_phibins())
0221   {
0222     std::cout
0223         << "RawTowerGeomContainer_Cylinderv1::get_phicenter - Asking for invalid bin in phi: "
0224         << ibin << std::endl;
0225     std::cout << "minbin: 0, maxbin " << get_phibins() << std::endl;
0226     exit(1);
0227   }
0228   return (phi_bound_map[ibin].first + phi_bound_map[ibin].second) / 2.;
0229 }
0230
0231 void RawTowerGeomContainer_Cylinderv1::set_phibounds(const int ibin,
0232                                                      const std::pair<double, double>& bounds)
0233 {
0234   if (ibin < 0 || ibin >= get_phibins())
0235   {
0236     std::cout
0237         << "RawTowerGeomContainer_Cylinderv1::set_bounds - Asking for invalid bin in phi: "
0238         << ibin << std::endl;
0239     std::cout << "minbin: 0, maxbin " << get_phibins() << std::endl;
0240     exit(1);
0241   }
0242
0243   std::pair<double, double> b_reg(bounds);
0244   if (b_reg.first > b_reg.second)
0245   {
0246     b_reg.second = bounds.first;
0247     b_reg.first = bounds.second;
0248   }
0249
0250   phi_bound_map[ibin] = b_reg;
0251 }