Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 10:56:16

0001 // -*- C++ -*-
0002 //
0003 // Package:     SiPixelPhase1Common
0004 // Class  :     HistogramManager
0005 //
0006 #include "DQM/SiPixelPhase1Common/interface/HistogramManager.h"
0007 
0008 #include <sstream>
0009 #include <boost/algorithm/string.hpp>
0010 
0011 // Geometry stuff
0012 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0013 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0014 
0015 // Logger
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 
0018 HistogramManager::HistogramManager(const edm::ParameterSet& iconfig, GeometryInterface& geo)
0019     : iConfig(iconfig),
0020       geometryInterface(geo),
0021       enabled(iconfig.getParameter<bool>("enabled")),
0022       perLumiHarvesting(iconfig.getParameter<bool>("perLumiHarvesting")),
0023       bookUndefined(iconfig.getParameter<bool>("bookUndefined")),
0024       top_folder_name(iconfig.getParameter<std::string>("topFolderName")),
0025       name(iconfig.getParameter<std::string>("name")),
0026       title(iconfig.getParameter<std::string>("title")),
0027       xlabel(iconfig.getParameter<std::string>("xlabel")),
0028       ylabel(iconfig.getParameter<std::string>("ylabel")),
0029       dimensions(iconfig.getParameter<int>("dimensions")),
0030       range_x_nbins(iconfig.getParameter<int>("range_nbins")),
0031       range_x_min(iconfig.getParameter<double>("range_min")),
0032       range_x_max(iconfig.getParameter<double>("range_max")),
0033       range_y_nbins(iconfig.getParameter<int>("range_y_nbins")),
0034       range_y_min(iconfig.getParameter<double>("range_y_min")),
0035       range_y_max(iconfig.getParameter<double>("range_y_max")),
0036       statsOverflows(iconfig.getParameter<bool>("statsOverflows")) {
0037   auto spec_configs = iconfig.getParameter<edm::VParameterSet>("specs");
0038   for (const auto& spec : spec_configs) {
0039     // this would fit better in SummationSpecification(...), but it has to
0040     // happen here.
0041     auto conf = spec.getParameter<edm::ParameterSet>("conf");
0042     if (!conf.getParameter<bool>("enabled"))
0043       continue;
0044     addSpec(SummationSpecification(spec, geometryInterface));
0045   }
0046 }
0047 
0048 void HistogramManager::addSpec(SummationSpecification spec) {
0049   specs.push_back(spec);
0050   tables.push_back(Table());
0051   counters.push_back(Table());
0052   significantvalues.push_back(GeometryInterface::Values());
0053   fastpath.push_back(nullptr);
0054 }
0055 
0056 // This is the hottest function in the HistogramManager. Make sure it does not
0057 // allocate memory or other expensive things.
0058 // Currently the GeometryInterface::extract (some virtual calls) and the map
0059 // lookups should be the most expensive things, but they should only happen if
0060 // the module changed from the last call; an optimization that fails when row/
0061 // col are used.
0062 // fillInternal (called from here) does more lookups on the geometry, if EXTEND
0063 // is used; we do not attempt the optimization there since in most cases row/
0064 // col are involved.
0065 void HistogramManager::fill(double x, double y, DetId sourceModule, const edm::Event* sourceEvent, int col, int row) {
0066   if (!enabled)
0067     return;
0068 
0069   // We only need to check the module here, since the fastpath is only used to
0070   // determine which plot is filled (not which bin inside) and fillInternal
0071   // re-extracts whatever it needs for the axis.
0072   // Since do not support booking based on ROC or time, the plot can only depend
0073   // on the module.
0074   // The sourceEvent check is not really effective (pointer is always the same)
0075   // but needed for initialisation.
0076   // NOTE that this might change if we want to support per-ROC booking
0077   // NOTE that we could even cache the bin to fill if iq and spec are identical,
0078   // also across HistogramManagers.
0079   bool cached = false;
0080   if (sourceModule == this->iq.sourceModule && sourceEvent == this->iq.sourceEvent) {
0081     cached = true;
0082   }
0083   iq = GeometryInterface::InterestingQuantities{sourceEvent, sourceModule, int16_t(col), int16_t(row)};
0084   for (unsigned int i = 0; i < specs.size(); i++) {
0085     auto& s = specs[i];
0086     auto& t = s.steps[0].type == SummationStep::COUNT ? counters[i] : tables[i];
0087     if (!cached) {  // recompute ME to fill (aka fastpath)
0088       significantvalues[i].clear();
0089       geometryInterface.extractColumns(s.steps[0].columns, iq, significantvalues[i]);
0090       auto histo = t.find(significantvalues[i]);
0091       if (histo == t.end()) {
0092         if (bookUndefined) {
0093           edm::LogError("HistogramManager") << "Missing Histogram!\n"
0094                                             << "name " << tables[i].begin()->second.th1->GetName() << "\n";
0095           assert(!"Histogram not booked! Probably inconsistent geometry description.");
0096         }
0097         fastpath[i] = nullptr;
0098       } else {
0099         fastpath[i] = &(histo->second);
0100       }
0101     }
0102     // A fastpath of nullptr means there is no ME for this iq, which can be
0103     // a valid cached result.
0104     if (fastpath[i]) {
0105       if (s.steps[0].type == SummationStep::COUNT) {
0106         fastpath[i]->count++;
0107       } else {
0108         fillInternal(x, y, this->dimensions, iq, s.steps.begin() + 1, s.steps.end(), *(fastpath[i]));
0109       }
0110     }
0111   }
0112 }
0113 void HistogramManager::fill(double x, DetId sourceModule, const edm::Event* sourceEvent, int col, int row) {
0114   assert(this->dimensions == 1);
0115   fill(x, 0.0, sourceModule, sourceEvent, col, row);
0116 }
0117 void HistogramManager::fill(DetId sourceModule, const edm::Event* sourceEvent, int col, int row) {
0118   assert(this->dimensions == 0);
0119   fill(0.0, 0.0, sourceModule, sourceEvent, col, row);
0120 }
0121 
0122 void HistogramManager::fillInternal(double x,
0123                                     double y,
0124                                     int n_parameters,
0125                                     GeometryInterface::InterestingQuantities const& iq,
0126                                     std::vector<SummationStep>::iterator first,
0127                                     std::vector<SummationStep>::iterator last,
0128                                     AbstractHistogram& dest) {
0129   double fx = 0, fy = 0, fz = 0;
0130   int tot_parameters = n_parameters;
0131   for (auto it = first; it != last; ++it) {
0132     if (it->stage != SummationStep::STAGE1)
0133       break;
0134     // The specification builder precomputes where x and y go, this loop will
0135     // always do 3 iterations to set x, y, z. The builder does not know how
0136     // many parameters we have, so we have to check that and count the total.
0137     switch (it->type) {
0138       case SummationStep::USE_X:
0139         if (it->arg[0] == '1' && n_parameters >= 1)
0140           fx = x;
0141         if (it->arg[0] == '2' && n_parameters >= 2)
0142           fx = y;
0143         break;
0144       case SummationStep::USE_Y:
0145         if (it->arg[0] == '1' && n_parameters >= 1)
0146           fy = x;
0147         if (it->arg[0] == '2' && n_parameters >= 2)
0148           fy = y;
0149         break;
0150       case SummationStep::USE_Z:
0151         if (it->arg[0] == '1' && n_parameters >= 1)
0152           fz = x;
0153         if (it->arg[0] == '2' && n_parameters >= 2)
0154           fz = y;
0155         break;
0156       case SummationStep::EXTEND_X:
0157         fx = geometryInterface.extract(it->columns[0], iq).second;
0158         tot_parameters++;
0159         break;
0160       case SummationStep::EXTEND_Y:
0161         fy = geometryInterface.extract(it->columns[0], iq).second;
0162         tot_parameters++;
0163         break;
0164       case SummationStep::PROFILE:
0165         break;  // profile does not make a difference here, only in booking
0166       default:
0167         assert(!"illegal step in STAGE1!");
0168     }
0169   }
0170 
0171   switch (tot_parameters) {
0172     case 1:
0173       dest.me->Fill(fx);
0174       break;
0175     case 2:
0176       dest.me->Fill(fx, fy);
0177       break;
0178     case 3:
0179       dest.me->Fill(fx, fy, fz);
0180       break;
0181     default:
0182       edm::LogError("HistogramManager") << "got " << tot_parameters << " dimensions\n"
0183                                         << "name " << dest.th1->GetName() << "\n";
0184       assert(!"More than 3 dimensions should never occur.");
0185   }
0186 }
0187 
0188 // This is only used for ndigis-like counting. It could be more optimized, but
0189 // is probably fine for a per-event thing.
0190 void HistogramManager::executePerEventHarvesting(const edm::Event* sourceEvent) {
0191   if (!enabled)
0192     return;
0193   for (unsigned int i = 0; i < specs.size(); i++) {
0194     auto& s = specs[i];
0195     auto& t = tables[i];
0196     auto& c = counters[i];
0197     if (s.steps[0].type != SummationStep::COUNT)
0198       continue;  // no counting, done
0199     assert((s.steps.size() >= 2 && s.steps[1].type == SummationStep::GROUPBY) ||
0200            !"Incomplete spec (but this cannot be caught in Python)");
0201     for (auto& e : c) {
0202       // the iq on the counter can only be a _sample_, since many modules
0203       // could be grouped on one counter. Even worse, the sourceEvent ptr
0204       // could be dangling if the counter was not touched in this event, so
0205       // we replace it. row/col are most likely useless as well.
0206       auto iq = e.second.iq_sample;
0207       iq.sourceEvent = sourceEvent;
0208 
0209       significantvalues[i].clear();
0210       geometryInterface.extractColumns(s.steps[1].columns, iq, significantvalues[i]);
0211       auto histo = t.find(significantvalues[i]);
0212       if (histo == t.end()) {
0213         if (!bookUndefined)
0214           continue;
0215         edm::LogError("HistogramManager") << "Histogram Missing!\n"
0216                                           << "name " << t.begin()->second.th1->GetName() << "\n"
0217                                           << "ctr "
0218                                           << " detid " << iq.sourceModule << "\n";
0219         assert(!"Histogram not booked! (per-event) Probably inconsistent geometry description.");
0220       }
0221       fillInternal(e.second.count, 0, 1, iq, s.steps.begin() + 2, s.steps.end(), histo->second);
0222       e.second.count = 0;
0223     }
0224   }
0225 }
0226 
0227 std::pair<std::string, std::string> HistogramManager::makePathName(SummationSpecification const& s,
0228                                                                    GeometryInterface::Values const& significantvalues,
0229                                                                    SummationStep const* upto) {
0230   std::ostringstream dir("");
0231   std::string suffix = "";
0232 
0233   // we omit the last value here, to get all disks next to each other etc.
0234   if (!significantvalues.empty()) {
0235     for (auto it = significantvalues.begin(); it != (significantvalues.end() - 1); ++it) {
0236       auto name = geometryInterface.formatValue(it->first, it->second);
0237       if (name.empty())
0238         continue;
0239       dir << name << "/";
0240     }
0241     auto e = significantvalues[significantvalues.size() - 1];
0242     suffix = "_" + geometryInterface.formatValue(e.first, e.second);
0243   }
0244 
0245   // PERF: this is actually independent of significantvalues and iq
0246   std::string name = this->name;
0247   for (SummationStep const& step : s.steps) {
0248     if (&step == upto)
0249       break;
0250     switch (step.type) {
0251       case SummationStep::COUNT:
0252         name = "num_" + name;
0253         break;
0254       case SummationStep::EXTEND_X:
0255       case SummationStep::EXTEND_Y: {
0256         if (step.stage != SummationStep::STAGE1)
0257           break;
0258         GeometryInterface::Column col0 = step.columns[0];
0259         std::string colname = geometryInterface.pretty(col0);
0260         name = name + "_per_" + colname;
0261         break;
0262       }
0263       case SummationStep::REDUCE: {
0264         auto red = step.arg;
0265         boost::algorithm::to_lower(red);
0266         name = red + "_" + name;
0267         break;
0268       }
0269       default:
0270         // Maybe PROFILE is worth showing.
0271         break;  // not visible in name
0272     }
0273   }
0274   return std::make_pair(top_folder_name + "/" + dir.str(), name + suffix);
0275 }
0276 
0277 void HistogramManager::book(DQMStore::IBooker& iBooker, edm::EventSetup const& iSetup) {
0278   if (!geometryInterface.loaded()) {
0279     geometryInterface.load(iSetup);
0280   }
0281   if (!enabled)
0282     return;
0283 
0284   struct MEInfo {
0285     int dimensions = 1;
0286     double range_x_min = 1e12;
0287     double range_x_max = -1e12;
0288     double range_y_min = 1e12;
0289     double range_y_max = -1e12;
0290     double range_z_min = 1e12;  // z range carried around but unused
0291     double range_z_max = -1e12;
0292     int range_x_nbins = 0;
0293     int range_y_nbins = 0;
0294     int range_z_nbins = 0;
0295     GeometryInterface::Value binwidth_x = 0;  // override nbins for geom-things
0296     GeometryInterface::Value binwidth_y = 0;
0297     std::string title, xlabel, ylabel, zlabel;
0298     bool do_profile = false;
0299     bool statsOverflows = true;
0300     ;
0301   };
0302   std::map<GeometryInterface::Values, MEInfo> toBeBooked;
0303 
0304   for (unsigned int i = 0; i < specs.size(); i++) {
0305     auto& s = specs[i];
0306     auto& t = tables[i];
0307     auto& c = counters[i];
0308     toBeBooked.clear();
0309     bool bookCounters = false;
0310 
0311     auto firststep = s.steps.begin();
0312     int n_parameters = this->dimensions;
0313 
0314     if (firststep->type != SummationStep::GROUPBY) {
0315       ++firststep;
0316       n_parameters = 1;
0317       bookCounters = true;
0318     }
0319 
0320     auto laststep = std::find_if(
0321         s.steps.begin(), s.steps.end(), [](SummationStep const& step) { return step.stage == SummationStep::STAGE2; });
0322 
0323     GeometryInterface::Values significantvalues;
0324 
0325     for (auto iq : geometryInterface.allModules()) {
0326       if (bookCounters) {
0327         // add an entry for the counting step if present
0328         geometryInterface.extractColumns(s.steps[0].columns, iq, significantvalues);
0329         c[significantvalues].iq_sample = iq;
0330       }
0331 
0332       geometryInterface.extractColumns(firststep->columns, iq, significantvalues);
0333       if (!bookUndefined) {
0334         // skip if any column is UNDEFINED
0335         bool ok = true;
0336         for (auto e : significantvalues)
0337           if (e.second == GeometryInterface::UNDEFINED)
0338             ok = false;
0339         if (!ok)
0340           continue;
0341       }
0342 
0343       auto histo = toBeBooked.find(significantvalues);
0344       if (histo == toBeBooked.end()) {
0345         // create new histo
0346         MEInfo& mei = toBeBooked[significantvalues];
0347         mei.title = this->title;
0348         mei.statsOverflows = this->statsOverflows;
0349         if (bookCounters)
0350           mei.title =
0351               "Number of " + mei.title + " per Event and " + geometryInterface.pretty(*(s.steps[0].columns.end() - 1));
0352         std::string xlabel = bookCounters ? "#" + this->xlabel : this->xlabel;
0353 
0354         // refer to fillInternal() for the actual execution
0355         // compute labels, title, type, user-set ranges here
0356         int tot_parameters = n_parameters;
0357 
0358 #define SET_AXIS(to, from)                                                          \
0359   mei.to##label = from##label;                                                      \
0360   mei.range_##to##_min = ((it->nbins == -1) ? this->range_##from##_min : it->xmin); \
0361   mei.range_##to##_max = ((it->nbins == -1) ? this->range_##from##_max : it->xmax); \
0362   mei.range_##to##_nbins = ((it->nbins == -1) ? this->range_##from##_nbins : it->nbins)
0363 
0364         for (auto it = firststep + 1; it != laststep; ++it) {
0365           switch (it->type) {
0366             case SummationStep::USE_X:
0367               if (it->arg[0] == '1' && n_parameters >= 1) {
0368                 SET_AXIS(x, x);
0369               }  // TODO: make use of current nbins, xmin, xmax if set
0370               if (it->arg[0] == '2' && n_parameters >= 2) {
0371                 SET_AXIS(x, y);
0372               }
0373               break;
0374             case SummationStep::USE_Y:
0375               if (it->arg[0] == '1' && n_parameters >= 1) {
0376                 SET_AXIS(y, x);
0377               }
0378               if (it->arg[0] == '2' && n_parameters >= 2) {
0379                 SET_AXIS(y, y);
0380               }
0381               break;
0382             case SummationStep::USE_Z:
0383               if (it->arg[0] == '1' && n_parameters >= 1) {
0384                 SET_AXIS(z, x);
0385               }
0386               if (it->arg[0] == '2' && n_parameters >= 2) {
0387                 SET_AXIS(z, y);
0388               }
0389               break;
0390             case SummationStep::EXTEND_X: {
0391               assert(mei.range_x_nbins == 0);
0392               auto col = it->columns[0];
0393               mei.xlabel = geometryInterface.pretty(col);
0394               mei.title = mei.title + " by " + mei.xlabel;
0395               if (geometryInterface.minValue(col) != GeometryInterface::UNDEFINED)
0396                 mei.range_x_min = geometryInterface.minValue(col);
0397               if (geometryInterface.maxValue(col) != GeometryInterface::UNDEFINED)
0398                 mei.range_x_max = geometryInterface.maxValue(col);
0399               mei.binwidth_x = geometryInterface.binWidth(col);
0400               tot_parameters++;
0401             } break;
0402             case SummationStep::EXTEND_Y: {
0403               auto col = it->columns[0];
0404               mei.ylabel = geometryInterface.pretty(col);
0405               mei.title = mei.title + " by " + mei.ylabel;
0406               if (geometryInterface.minValue(col) != GeometryInterface::UNDEFINED)
0407                 mei.range_y_min = geometryInterface.minValue(col);
0408               if (geometryInterface.maxValue(col) != GeometryInterface::UNDEFINED)
0409                 mei.range_y_max = geometryInterface.maxValue(col);
0410               mei.binwidth_y = geometryInterface.binWidth(col);
0411               tot_parameters++;
0412             } break;
0413             case SummationStep::PROFILE:
0414               mei.do_profile = true;
0415               break;
0416             default:
0417               assert(!"illegal step in STAGE1! (booking)");
0418           }
0419         }
0420         mei.dimensions = tot_parameters;
0421         if (mei.do_profile)
0422           mei.title = "Profile of " + mei.title;
0423         if (!mei.zlabel.empty())
0424           mei.title = mei.title + " (Z: " + mei.zlabel + ")";
0425       }
0426       // only update range
0427       MEInfo& mei = toBeBooked[significantvalues];
0428       double val;
0429 
0430       for (auto it = firststep + 1; it != laststep; ++it) {
0431         switch (it->type) {
0432           case SummationStep::EXTEND_X:
0433             val = geometryInterface.extract(it->columns[0], iq).second;
0434             if (val == GeometryInterface::UNDEFINED)
0435               break;
0436             mei.range_x_min = std::min(mei.range_x_min, val);
0437             mei.range_x_max = std::max(mei.range_x_max, val);
0438             break;
0439           case SummationStep::EXTEND_Y:
0440             val = geometryInterface.extract(it->columns[0], iq).second;
0441             if (val == GeometryInterface::UNDEFINED)
0442               break;
0443             mei.range_y_min = std::min(mei.range_y_min, val);
0444             mei.range_y_max = std::max(mei.range_y_max, val);
0445             break;
0446           default:  // ignore the rest, code above will catch bugs
0447             break;
0448         }
0449       }
0450     }
0451 
0452     // Now do the actual booking.
0453     for (auto& e : toBeBooked) {
0454       AbstractHistogram& h = t[e.first];
0455       MEInfo& mei = e.second;
0456       auto name = makePathName(s, e.first, &(*laststep));
0457       iBooker.setCurrentFolder(name.first);
0458 
0459       // determine nbins for geometry derived quantities
0460       // due to how we counted above, we need to include lower and upper bound
0461       // For Coord-values, which are not precisely aligned with the bins, we force
0462       // alignment.
0463       if (mei.binwidth_x != 0) {
0464         double range = (mei.range_x_max - mei.range_x_min) / mei.binwidth_x;
0465         if ((range - int(range)) == 0.0) {
0466           mei.range_x_min -= mei.binwidth_x / 2;
0467           mei.range_x_max += mei.binwidth_x / 2;
0468         } else {
0469           mei.range_x_min = std::floor(mei.range_x_min / mei.binwidth_x) * mei.binwidth_x;
0470           mei.range_x_max = std::ceil(mei.range_x_max / mei.binwidth_x) * mei.binwidth_x;
0471         }
0472         mei.range_x_nbins = int((mei.range_x_max - mei.range_x_min) / mei.binwidth_x);
0473       }
0474       if (mei.binwidth_y != 0) {
0475         double range = (mei.range_y_max - mei.range_y_min) / mei.binwidth_y;
0476         if ((range - int(range)) == 0.0) {
0477           mei.range_y_min -= mei.binwidth_y / 2;
0478           mei.range_y_max += mei.binwidth_y / 2;
0479         } else {
0480           mei.range_y_min = std::floor(mei.range_y_min / mei.binwidth_y) * mei.binwidth_y;
0481           mei.range_y_max = std::ceil(mei.range_y_max / mei.binwidth_y) * mei.binwidth_y;
0482         }
0483         mei.range_y_nbins = int((mei.range_y_max - mei.range_y_min) / mei.binwidth_y);
0484       }
0485 
0486       if (mei.dimensions == 1) {
0487         h.me = iBooker.book1D(
0488             name.second, (mei.title + ";" + mei.xlabel).c_str(), mei.range_x_nbins, mei.range_x_min, mei.range_x_max);
0489       } else if (mei.dimensions == 2 && !mei.do_profile) {
0490         h.me = iBooker.book2D(name.second,
0491                               (mei.title + ";" + mei.xlabel + ";" + mei.ylabel).c_str(),
0492                               mei.range_x_nbins,
0493                               mei.range_x_min,
0494                               mei.range_x_max,
0495                               mei.range_y_nbins,
0496                               mei.range_y_min,
0497                               mei.range_y_max);
0498       } else if (mei.dimensions == 2 && mei.do_profile) {
0499         h.me = iBooker.bookProfile(name.second,
0500                                    (mei.title + ";" + mei.xlabel + ";" + mei.ylabel).c_str(),
0501                                    mei.range_x_nbins,
0502                                    mei.range_x_min,
0503                                    mei.range_x_max,
0504                                    0.0,
0505                                    0.0,
0506                                    "");
0507       } else if (mei.dimensions == 3 && mei.do_profile) {
0508         h.me = iBooker.bookProfile2D(name.second,
0509                                      (mei.title + ";" + mei.xlabel + ";" + mei.ylabel).c_str(),
0510                                      mei.range_x_nbins,
0511                                      mei.range_x_min,
0512                                      mei.range_x_max,
0513                                      mei.range_y_nbins,
0514                                      mei.range_y_min,
0515                                      mei.range_y_max,
0516                                      0.0,
0517                                      0.0);  // Z range is ignored if min==max
0518       } else {
0519         edm::LogError("HistogramManager") << "Illegal Histogram!\n"
0520                                           << "name " << name.second << "\n"
0521                                           << "dim " << mei.dimensions << " profile " << mei.do_profile << "\n";
0522         assert(!"Illegal Histogram kind.");
0523       }
0524       h.th1 = h.me->getTH1();
0525       h.me->setStatOverflows(mei.statsOverflows);
0526     }
0527   }
0528 }
0529 
0530 void HistogramManager::executePerLumiHarvesting(DQMStore::IBooker& iBooker,
0531                                                 DQMStore::IGetter& iGetter,
0532                                                 edm::LuminosityBlock const& lumiBlock,
0533                                                 edm::EventSetup const& iSetup) {
0534   if (!enabled)
0535     return;
0536   // this should also give us the GeometryInterface for offline, though it is a
0537   // bit dirty and might explode.
0538   if (!geometryInterface.loaded()) {
0539     geometryInterface.load(iSetup);
0540   }
0541   if (perLumiHarvesting) {
0542     this->lumisection = &lumiBlock;  // "custom" steps can use this
0543     executeHarvesting(iBooker, iGetter);
0544     this->lumisection = nullptr;
0545   }
0546 }
0547 
0548 void HistogramManager::loadFromDQMStore(SummationSpecification& s, Table& t, DQMStore::IGetter& iGetter) {
0549   t.clear();
0550   GeometryInterface::Values significantvalues;
0551   auto firststep = s.steps.begin();
0552   if (firststep->type != SummationStep::GROUPBY)
0553     ++firststep;
0554   auto laststep = std::find_if(
0555       s.steps.begin(), s.steps.end(), [](SummationStep const& step) { return step.stage == SummationStep::STAGE2; });
0556 
0557   for (auto iq : geometryInterface.allModules()) {
0558     geometryInterface.extractColumns(firststep->columns, iq, significantvalues);
0559 
0560     auto histo = t.find(significantvalues);
0561     if (histo == t.end()) {
0562       auto name = makePathName(s, significantvalues, &(*laststep));
0563       std::string path = name.first + name.second;
0564       MonitorElement* me = iGetter.get(path);
0565       if (!me) {
0566         if (bookUndefined)
0567           edm::LogError("HistogramManager") << "ME " << path << " not found\n";
0568         // else this will happen quite often
0569       } else {
0570         AbstractHistogram& histo = t[significantvalues];
0571         histo.me = me;
0572         histo.th1 = histo.me->getTH1();
0573         histo.iq_sample = iq;
0574       }
0575     }
0576   }
0577 }
0578 
0579 void HistogramManager::executeGroupBy(SummationStep const& step,
0580                                       Table& t,
0581                                       DQMStore::IBooker& iBooker,
0582                                       SummationSpecification const& s) {
0583   // Simple regrouping, sum histos if they end up in the same place.
0584   Table out;
0585   GeometryInterface::Values significantvalues;
0586   for (auto& e : t) {
0587     TH1* th1 = e.second.th1;
0588     geometryInterface.extractColumns(step.columns, e.second.iq_sample, significantvalues);
0589     AbstractHistogram& new_histo = out[significantvalues];
0590     if (!new_histo.me) {
0591       auto name = makePathName(s, significantvalues, &step);
0592       iBooker.setCurrentFolder(name.first);
0593       if (dynamic_cast<TH1F*>(th1))
0594         new_histo.me = iBooker.book1D(name.second, (TH1F*)th1);
0595       else if (dynamic_cast<TH2F*>(th1))
0596         new_histo.me = iBooker.book2D(name.second, (TH2F*)th1);
0597       else if (dynamic_cast<TProfile*>(th1))
0598         new_histo.me = iBooker.bookProfile(name.second, (TProfile*)th1);
0599       else if (dynamic_cast<TProfile2D*>(th1))
0600         new_histo.me = iBooker.bookProfile2D(name.second, (TProfile2D*)th1);
0601       else
0602         assert(!"No idea how to book this.");
0603       new_histo.th1 = new_histo.me->getTH1();
0604       new_histo.iq_sample = e.second.iq_sample;
0605     } else {
0606       new_histo.th1->Add(th1);
0607     }
0608     new_histo.me->setStatOverflows(e.second.me->getStatOverflows());
0609   }
0610   t.swap(out);
0611 }
0612 
0613 void HistogramManager::executeExtend(SummationStep const& step,
0614                                      Table& t,
0615                                      std::string const& reduce_type,
0616                                      DQMStore::IBooker& iBooker,
0617                                      SummationSpecification const& s) {
0618   // For the moment only X.
0619   // first pass determines the range.
0620   std::map<GeometryInterface::Values, int> nbins;
0621   // separators collects meta info for the render plugin about the boundaries.
0622   // for each EXTEND, this is added to the axis label. In total this is not
0623   // fully correct since we have to assume the the substructure of each sub-
0624   // histogram is the same, which is e.g. not true for layers. It still works
0625   // since layers only contain leaves (ladders).
0626   std::map<GeometryInterface::Values, std::string> separators;
0627 
0628   GeometryInterface::Values significantvalues;
0629 
0630   for (auto& e : t) {
0631     geometryInterface.extractColumns(step.columns, e.second.iq_sample, significantvalues);
0632     int& n = nbins[significantvalues];
0633     assert(e.second.th1 || !"invalid histogram");
0634     // if we reduce, every histogram only needs one bin
0635     int bins = !reduce_type.empty() ? 1 : e.second.th1->GetXaxis()->GetNbins();
0636     if (bins > 1)
0637       separators[significantvalues] += std::to_string(n) + ",";
0638     n += bins;
0639   }
0640   for (auto& e : separators)
0641     e.second = "(" + e.second + ")/";
0642 
0643   Table out;
0644   for (auto& e : t) {
0645     geometryInterface.extractColumns(step.columns, e.second.iq_sample, significantvalues);
0646     TH1* th1 = e.second.th1;
0647     assert(th1);
0648 
0649     AbstractHistogram& new_histo = out[significantvalues];
0650     if (!new_histo.me) {
0651       // we put the name of the actual, last column of a input histo there.
0652       std::string colname = geometryInterface.pretty((e.first.end() - 1)->first);
0653 
0654       auto separator = separators[significantvalues];
0655 
0656       auto name = makePathName(s, significantvalues, &step);
0657       auto title =
0658           std::string("") + th1->GetTitle() + " per " + colname + ";" + colname + separator +
0659           (!reduce_type.empty() ? th1->GetYaxis()->GetTitle() : th1->GetXaxis()->GetTitle()) + ";" +
0660           (!reduce_type.empty() ? reduce_type + " of " + th1->GetXaxis()->GetTitle() : th1->GetYaxis()->GetTitle());
0661       iBooker.setCurrentFolder(name.first);
0662 
0663       if (th1->GetDimension() == 1) {
0664         new_histo.me =
0665             iBooker.book1D(name.second, title, nbins[significantvalues], 0.5, nbins[significantvalues] + 0.5);
0666       } else {
0667         assert(!"2D extend not implemented in harvesting.");
0668       }
0669       new_histo.th1 = new_histo.me->getTH1();
0670       new_histo.iq_sample = e.second.iq_sample;
0671       new_histo.count = 1;  // used as a fill pointer. Assumes histograms are
0672                             // ordered correctly (map should provide that)
0673     }
0674 
0675     // now add data.
0676     if (new_histo.th1->GetDimension() == 1) {
0677       if (reduce_type.empty()) {  // no reduction, concatenate
0678         for (int i = 1; i <= th1->GetXaxis()->GetNbins(); i++) {
0679           new_histo.th1->SetBinContent(new_histo.count, th1->GetBinContent(i));
0680           new_histo.th1->SetBinError(new_histo.count, th1->GetBinError(i));
0681           new_histo.count += 1;
0682         }
0683       } else if (reduce_type == "MEAN") {
0684         new_histo.th1->SetBinContent(new_histo.count, th1->GetMean());
0685         new_histo.th1->SetBinError(new_histo.count, th1->GetMeanError());
0686         new_histo.count += 1;
0687       } else {
0688         assert(!"Reduction type not supported");
0689       }
0690       new_histo.me->setStatOverflows(e.second.me->getStatOverflows());
0691     } else {
0692       assert(!"2D extend not implemented in harvesting.");
0693     }
0694   }
0695   t.swap(out);
0696 }
0697 
0698 void HistogramManager::executeHarvesting(DQMStore::IBooker& iBooker, DQMStore::IGetter& iGetter) {
0699   if (!enabled)
0700     return;
0701   // Debug output
0702   for (auto& s : specs) {
0703     edm::LogInfo log("HistogramManager");
0704     log << "Specs for " << name << " ";
0705     s.dump(log, geometryInterface);
0706   }
0707 
0708   for (unsigned int i = 0; i < specs.size(); i++) {
0709     auto& s = specs[i];
0710     auto& t = tables[i];
0711     loadFromDQMStore(s, t, iGetter);
0712 
0713     std::string reduce_type = "";
0714 
0715     // now execute step2.
0716     for (SummationStep const& step : s.steps) {
0717       if (step.stage == SummationStep::STAGE2) {
0718         switch (step.type) {
0719           case SummationStep::SAVE:
0720             // no explicit implementation atm.
0721             break;
0722           case SummationStep::GROUPBY:
0723             executeGroupBy(step, t, iBooker, s);
0724             break;
0725           case SummationStep::REDUCE:
0726             // reduction is done in the following EXTEND
0727             reduce_type = step.arg;
0728             break;
0729           case SummationStep::EXTEND_X:
0730             executeExtend(step, t, reduce_type, iBooker, s);
0731             reduce_type = "";
0732             break;
0733           case SummationStep::EXTEND_Y:
0734             assert(!"EXTEND_Y currently not supported in harvesting.");
0735             break;
0736           default:
0737             assert(!"Operation not supported in harvesting.");
0738         }
0739       }
0740     }
0741   }
0742 }