DTMonitorClient/src/DTOccupancyClusterBuilder.cc

0001
0002 /*
0003  *  See header file for a description of this class.
0004  *
0005  *  \author G. Cerminara - INFN Torino
0006  */
0007
0008 #include "DTOccupancyClusterBuilder.h"
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010
0011 #include "TCanvas.h"
0012 #include "TH2F.h"
0013
0014 #include <algorithm>
0015 #include <sstream>
0016 #include <iostream>
0017
0018 using namespace std;
0019 using namespace edm;
0020
0021 DTOccupancyClusterBuilder::DTOccupancyClusterBuilder() : maxMean(-1.), maxRMS(-1.) {}
0022
0023 DTOccupancyClusterBuilder::~DTOccupancyClusterBuilder() {}
0024
0025 void DTOccupancyClusterBuilder::addPoint(const DTOccupancyPoint& point) {
0026   // loop over points already stored
0027   for (set<DTOccupancyPoint>::const_iterator pt = thePoints.begin(); pt != thePoints.end(); ++pt) {
0028     theDistances[(*pt).distance(point)] = make_pair(*pt, point);
0029   }
0030   thePoints.insert(point);
0031 }
0032
0033 void DTOccupancyClusterBuilder::buildClusters() {
0034   while (buildNewCluster()) {
0035     if (thePoints.size() <= 1)
0036       break;
0037   }
0038
0039   // build single point clusters with the remaining points
0040   for (set<DTOccupancyPoint>::const_iterator pt = thePoints.begin(); pt != thePoints.end(); ++pt) {
0041     DTOccupancyCluster clusterCandidate(*pt);
0042     theClusters.push_back(clusterCandidate);
0043     // store the range for building the histograms later
0044     if (clusterCandidate.maxMean() > maxMean)
0045       maxMean = clusterCandidate.maxMean();
0046     if (clusterCandidate.maxRMS() > maxRMS)
0047       maxRMS = clusterCandidate.maxRMS();
0048   }
0049   LogTrace("DTDQM|DTMonitorClient|DTOccupancyTest|DTOccupancyClusterBuilder")
0050       << " # of valid clusters: " << theClusters.size() << endl;
0051   sortClusters();
0052 }
0053
0054 void DTOccupancyClusterBuilder::drawClusters(std::string canvasName) {
0055   int nBinsX = 100;
0056   int nBinsY = 100;
0057   int colorMap[12] = {632, 600, 800, 400, 820, 416, 432, 880, 616, 860, 900, 920};
0058
0059   TCanvas* canvas = new TCanvas(canvasName.c_str(), canvasName.c_str());
0060   canvas->cd();
0061   for (vector<DTOccupancyCluster>::const_iterator cluster = theClusters.begin(); cluster != theClusters.end();
0062        ++cluster) {
0063     stringstream stream;
0064     stream << canvasName << "_" << cluster - theClusters.begin();
0065     string histoName = stream.str();
0066     TH2F* histo = (*cluster).getHisto(histoName,
0067                                       nBinsX,
0068                                       0,
0069                                       maxMean + 3 * maxMean / 100.,
0070                                       nBinsY,
0071                                       0,
0072                                       maxRMS + 3 * maxRMS / 100.,
0073                                       colorMap[cluster - theClusters.begin()]);
0074     if (cluster == theClusters.begin())
0075       histo->Draw("box");
0076     else
0077       histo->Draw("box,same");
0078   }
0079 }
0080
0081 std::pair<DTOccupancyPoint, DTOccupancyPoint> DTOccupancyClusterBuilder::getInitialPair() {
0082   return theDistances.begin()->second;
0083 }
0084
0085 void DTOccupancyClusterBuilder::computePointToPointDistances() {
0086   theDistances.clear();
0087   for (set<DTOccupancyPoint>::const_iterator pt_i = thePoints.begin(); pt_i != thePoints.end(); ++pt_i) {    // i loopo
0088     for (set<DTOccupancyPoint>::const_iterator pt_j = thePoints.begin(); pt_j != thePoints.end(); ++pt_j) {  // j loop
0089       if (*pt_i != *pt_j) {
0090         theDistances[pt_i->distance(*pt_j)] = make_pair(*pt_i, *pt_j);
0091       }
0092     }
0093   }
0094 }
0095
0096 void DTOccupancyClusterBuilder::computeDistancesToCluster(const DTOccupancyCluster& cluster) {
0097   theDistancesFromTheCluster.clear();
0098   for (set<DTOccupancyPoint>::const_iterator pt = thePoints.begin(); pt != thePoints.end(); ++pt) {
0099     theDistancesFromTheCluster[cluster.distance(*pt)] = *pt;
0100   }
0101 }
0102
0103 bool DTOccupancyClusterBuilder::buildNewCluster() {
0104   LogTrace("DTDQM|DTMonitorClient|DTOccupancyTest|DTOccupancyClusterBuilder")
0105       << "--------- New Cluster Candidate ----------------------" << endl;
0106   pair<DTOccupancyPoint, DTOccupancyPoint> initialPair = getInitialPair();
0107   LogTrace("DTDQM|DTMonitorClient|DTOccupancyTest|DTOccupancyClusterBuilder")
0108       << "   Initial Pair: " << endl
0109       << "           point1: mean " << initialPair.first.mean() << " rms " << initialPair.first.rms() << endl
0110       << "           point2: mean " << initialPair.second.mean() << " rms " << initialPair.second.rms() << endl;
0111   DTOccupancyCluster clusterCandidate(initialPair.first, initialPair.second);
0112   if (clusterCandidate.isValid()) {
0113     // remove already used pair
0114     thePoints.erase(initialPair.first);
0115     thePoints.erase(initialPair.second);
0116     if (!thePoints.empty()) {
0117       computeDistancesToCluster(clusterCandidate);
0118       while (clusterCandidate.addPoint(theDistancesFromTheCluster.begin()->second)) {
0119         thePoints.erase(theDistancesFromTheCluster.begin()->second);
0120         if (thePoints.empty())
0121           break;
0122         computeDistancesToCluster(clusterCandidate);
0123       }
0124     }
0125   } else {
0126     return false;
0127   }
0128   LogTrace("DTDQM|DTMonitorClient|DTOccupancyTest|DTOccupancyClusterBuilder")
0129       << "   # of layers: " << clusterCandidate.nPoints() << " avrg. mean: " << clusterCandidate.averageMean()
0130       << " avrg. rms: " << clusterCandidate.averageRMS() << endl;
0131   theClusters.push_back(clusterCandidate);
0132   // store the range for building the histograms later
0133   if (clusterCandidate.maxMean() > maxMean)
0134     maxMean = clusterCandidate.maxMean();
0135   if (clusterCandidate.maxRMS() > maxRMS)
0136     maxRMS = clusterCandidate.maxRMS();
0137   computePointToPointDistances();
0138   return true;
0139 }
0140
0141 void DTOccupancyClusterBuilder::sortClusters() {
0142   LogTrace("DTDQM|DTMonitorClient|DTOccupancyTest|DTOccupancyClusterBuilder") << " sorting" << endl;
0143   sort(theClusters.begin(), theClusters.end(), clusterIsLessThan);
0144   // we save the detid of the clusters which are not the best one
0145   for (vector<DTOccupancyCluster>::const_iterator cluster = ++(theClusters.begin()); cluster != theClusters.end();
0146        ++cluster) {  // loop over clusters skipping the first
0147     set<DTLayerId> clusterLayers = (*cluster).getLayerIDs();
0148     LogTrace("DTDQM|DTMonitorClient|DTOccupancyTest|DTOccupancyClusterBuilder")
0149         << "     # layers in the cluster: " << clusterLayers.size() << endl;
0150     theProblematicLayers.insert(clusterLayers.begin(), clusterLayers.end());
0151   }
0152   LogTrace("DTDQM|DTMonitorClient|DTOccupancyTest|DTOccupancyClusterBuilder")
0153       << " # of problematic layers: " << theProblematicLayers.size() << endl;
0154 }
0155
0156 DTOccupancyCluster DTOccupancyClusterBuilder::getBestCluster() const { return theClusters.front(); }
0157
0158 bool DTOccupancyClusterBuilder::isProblematic(DTLayerId layerId) const {
0159   if (theProblematicLayers.find(layerId) != theProblematicLayers.end()) {
0160     return true;
0161   }
0162   return false;
0163 }