File indexing completed on 2024-04-06 12:09:14
0001
0002 #include "Alignment/OfflineValidation/interface/FitWithRooFit.h"
0003 #include "DQMOffline/Alignment/interface/DiMuonMassBiasClient.h"
0004 #include "DataFormats/Histograms/interface/DQMToken.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 #include "FWCore/ServiceRegistry/interface/Service.h"
0007 #include "PhysicsTools/TagAndProbe/interface/RooCMSShape.h"
0008
0009
0010 #include "TCanvas.h"
0011 #include "RooAddPdf.h"
0012 #include "RooDataHist.h"
0013 #include "RooExponential.h"
0014 #include "RooGaussian.h"
0015 #include "RooPlot.h"
0016 #include "RooRealVar.h"
0017 #include "RooVoigtian.h"
0018 #include "RooCBShape.h"
0019
0020
0021 DiMuonMassBiasClient::DiMuonMassBiasClient(edm::ParameterSet const& iConfig)
0022 : TopFolder_(iConfig.getParameter<std::string>("FolderName")),
0023 useTH1s_(iConfig.getParameter<bool>("useTH1s")),
0024 useBWtimesCB_(iConfig.getParameter<bool>("useBWtimesCB")),
0025 fitBackground_(iConfig.getParameter<bool>("fitBackground")),
0026 useRooCBShape_(iConfig.getParameter<bool>("useRooCBShape")),
0027 useRooCMSShape_(iConfig.getParameter<bool>("useRooCMSShape")),
0028 debugMode_(iConfig.getParameter<bool>("debugMode")),
0029 MEtoHarvest_(iConfig.getParameter<std::vector<std::string>>("MEtoHarvest"))
0030
0031 {
0032 edm::LogInfo("DiMuonMassBiasClient") << "DiMuonMassBiasClient::Constructing DiMuonMassBiasClient ";
0033
0034
0035 edm::ParameterSet fit_par = iConfig.getParameter<edm::ParameterSet>("fit_par");
0036 diMuonMassBias::fillArrayF(meanConfig_, fit_par, "mean_par");
0037 diMuonMassBias::fillArrayF(widthConfig_, fit_par, "width_par");
0038 diMuonMassBias::fillArrayF(sigmaConfig_, fit_par, "sigma_par");
0039
0040 if (debugMode_) {
0041 edm::LogPrint("DiMuonMassBiasClient")
0042 << "mean: " << meanConfig_[0] << " (" << meanConfig_[1] << "," << meanConfig_[2] << ") " << std::endl;
0043 edm::LogPrint("DiMuonMassBiasClient")
0044 << "width: " << widthConfig_[0] << " (" << widthConfig_[1] << "," << widthConfig_[2] << ")" << std::endl;
0045 edm::LogPrint("DiMuonMassBiasClient")
0046 << "sigma: " << sigmaConfig_[0] << " (" << sigmaConfig_[1] << "," << sigmaConfig_[2] << ")" << std::endl;
0047 }
0048 }
0049
0050
0051 DiMuonMassBiasClient::~DiMuonMassBiasClient()
0052
0053 {
0054 edm::LogInfo("DiMuonMassBiasClient") << "DiMuonMassBiasClient::Deleting DiMuonMassBiasClient ";
0055 }
0056
0057
0058 void DiMuonMassBiasClient::beginJob(void)
0059
0060 {
0061 edm::LogInfo("DiMuonMassBiasClient") << "DiMuonMassBiasClient::beginJob done";
0062 }
0063
0064
0065 void DiMuonMassBiasClient::beginRun(edm::Run const& run, edm::EventSetup const& eSetup)
0066
0067 {
0068 edm::LogInfo("DiMuonMassBiasClient") << "DiMuonMassBiasClient:: Begining of Run";
0069 }
0070
0071
0072 void DiMuonMassBiasClient::bookMEs(DQMStore::IBooker& iBooker)
0073
0074 {
0075 iBooker.setCurrentFolder(TopFolder_ + "/DiMuonMassBiasMonitor/MassBias/Profiles");
0076 for (const auto& [key, ME] : harvestTargets_) {
0077 if (ME == nullptr) {
0078 edm::LogError("DiMuonMassBiasClient") << "could not find MonitorElement for key: " << key << std::endl;
0079 continue;
0080 }
0081
0082 const auto& title = ME->getTitle();
0083 const auto& xtitle = ME->getAxisTitle(1);
0084 const auto& ytitle = ME->getAxisTitle(2);
0085
0086 const auto& nxbins = ME->getNbinsX();
0087 const auto& xmin = ME->getAxisMin(1);
0088 const auto& xmax = ME->getAxisMax(1);
0089
0090 MonitorElement* meanToBook =
0091 iBooker.book1D(("Mean" + key), (title + ";#LT M_{#mu^{-}#mu^{+}} #GT [GeV];" + ytitle), nxbins, xmin, xmax);
0092 meanHistos_.insert({key, meanToBook});
0093
0094 MonitorElement* sigmaToBook =
0095 iBooker.book1D(("Sigma" + key), (title + ";" + xtitle + ";" + "#sigma of " + ytitle), nxbins, xmin, xmax);
0096 widthHistos_.insert({key, sigmaToBook});
0097 }
0098 }
0099
0100
0101 void DiMuonMassBiasClient::getMEsToHarvest(DQMStore::IGetter& iGetter)
0102
0103 {
0104 std::string inFolder = TopFolder_ + "/DiMuonMassBiasMonitor/MassBias/";
0105
0106
0107 for (const auto& hname : MEtoHarvest_) {
0108 MonitorElement* toHarvest = iGetter.get(inFolder + hname);
0109
0110 if (toHarvest == nullptr) {
0111 edm::LogError("DiMuonMassBiasClient") << "could not find input MonitorElement: " << inFolder + hname << std::endl;
0112 continue;
0113 }
0114
0115 harvestTargets_.insert({hname, toHarvest});
0116 }
0117 }
0118
0119
0120 void DiMuonMassBiasClient::fitAndFillProfile(std::pair<std::string, MonitorElement*> toHarvest,
0121 DQMStore::IBooker& iBooker)
0122
0123 {
0124 const auto& key = toHarvest.first;
0125 const auto& ME = toHarvest.second;
0126
0127 if (debugMode_)
0128 edm::LogPrint("DiMuonMassBiasClient") << "dealing with key: " << key << std::endl;
0129
0130 if (ME == nullptr) {
0131 edm::LogError("DiMuonMassBiasClient") << "could not find MonitorElement for key: " << key << std::endl;
0132 return;
0133 }
0134
0135 const auto& title = ME->getTitle();
0136 const auto& xtitle = ME->getAxisTitle(1);
0137 const auto& ytitle = ME->getAxisTitle(2);
0138
0139 const auto& nxbins = ME->getNbinsX();
0140 const auto& xmin = ME->getAxisMin(1);
0141 const auto& xmax = ME->getAxisMax(1);
0142
0143 TProfile* p_mean = new TProfile(("Mean" + key).c_str(),
0144 (title + ";" + xtitle + ";#LT M_{#mu^{-}#mu^{+}} #GT [GeV]").c_str(),
0145 nxbins,
0146 xmin,
0147 xmax,
0148 "g");
0149
0150 TProfile* p_width = new TProfile(
0151 ("Sigma" + key).c_str(), (title + ";" + xtitle + ";#sigma of " + ytitle).c_str(), nxbins, xmin, xmax, "g");
0152
0153 p_mean->Sumw2();
0154 p_width->Sumw2();
0155
0156 TH2F* bareHisto = ME->getTH2F();
0157 for (int bin = 1; bin <= nxbins; bin++) {
0158 const auto& xaxis = bareHisto->GetXaxis();
0159 const auto& low_edge = xaxis->GetBinLowEdge(bin);
0160 const auto& high_edge = xaxis->GetBinUpEdge(bin);
0161
0162 if (debugMode_)
0163 edm::LogPrint("DiMuonMassBiasClient") << "dealing with bin: " << bin << " range: (" << std::setprecision(2)
0164 << low_edge << "," << std::setprecision(2) << high_edge << ")";
0165
0166 TH1D* Proj = bareHisto->ProjectionY(Form("%s_proj_%i", key.c_str(), bin), bin, bin);
0167 Proj->SetTitle(Form("%s #in (%.2f,%.2f), bin: %i", Proj->GetTitle(), low_edge, high_edge, bin));
0168
0169 diMuonMassBias::fitOutputs results = useBWtimesCB_ ? fitBWTimesCB(Proj) : fitLineShape(Proj, fitBackground_);
0170
0171 if (results.isInvalid()) {
0172 edm::LogWarning("DiMuonMassBiasClient") << "the current bin has invalid data" << std::endl;
0173 continue;
0174 }
0175
0176
0177 const Measurement1D& bias = results.getBias();
0178
0179
0180
0181
0182
0183
0184
0185
0186 p_mean->SetBinContent(bin, bias.value() / (bias.error() * bias.error()));
0187 p_mean->SetBinEntries(bin, 1. / (bias.error() * bias.error()));
0188
0189 if (debugMode_)
0190 LogDebug("DiMuonBassBiasClient") << " Bin: " << bin << " value: " << bias.value() << " from profile ( "
0191 << p_mean->GetBinContent(bin) << ") - error: " << bias.error()
0192 << " from profile ( " << p_mean->GetBinError(bin) << " )";
0193
0194
0195 const Measurement1D& width = results.getWidth();
0196
0197
0198 p_width->SetBinContent(bin, width.value() / (width.error() * width.error()));
0199 p_width->SetBinEntries(bin, 1. / (width.error() * width.error()));
0200 }
0201
0202
0203 iBooker.setCurrentFolder(TopFolder_ + "/DiMuonMassBiasMonitor/MassBias/Profiles");
0204 MonitorElement* meanToBook = iBooker.bookProfile(p_mean->GetName(), p_mean);
0205 meanProfiles_.insert({key, meanToBook});
0206
0207 MonitorElement* sigmaToBook = iBooker.bookProfile(p_width->GetName(), p_width);
0208 widthProfiles_.insert({key, sigmaToBook});
0209
0210 delete p_mean;
0211 delete p_width;
0212 }
0213
0214
0215 void DiMuonMassBiasClient::fitAndFillHisto(std::pair<std::string, MonitorElement*> toHarvest,
0216 DQMStore::IBooker& iBooker)
0217
0218 {
0219 const auto& key = toHarvest.first;
0220 const auto& ME = toHarvest.second;
0221
0222 if (debugMode_)
0223 edm::LogPrint("DiMuonMassBiasClient") << "dealing with key: " << key << std::endl;
0224
0225 if (ME == nullptr) {
0226 edm::LogError("DiMuonMassBiasClient") << "could not find MonitorElement for key: " << key << std::endl;
0227 return;
0228 }
0229
0230 TH2F* bareHisto = ME->getTH2F();
0231 for (int bin = 1; bin <= ME->getNbinsX(); bin++) {
0232 const auto& xaxis = bareHisto->GetXaxis();
0233 const auto& low_edge = xaxis->GetBinLowEdge(bin);
0234 const auto& high_edge = xaxis->GetBinUpEdge(bin);
0235
0236 if (debugMode_)
0237 edm::LogPrint("DiMuonMassBiasClient") << "dealing with bin: " << bin << " range: (" << std::setprecision(2)
0238 << low_edge << "," << std::setprecision(2) << high_edge << ")";
0239 TH1D* Proj = bareHisto->ProjectionY(Form("%s_proj_%i", key.c_str(), bin), bin, bin);
0240 Proj->SetTitle(Form("%s #in (%.2f,%.2f), bin: %i", Proj->GetTitle(), low_edge, high_edge, bin));
0241
0242 diMuonMassBias::fitOutputs results = useBWtimesCB_ ? fitBWTimesCB(Proj) : fitLineShape(Proj, fitBackground_);
0243
0244 if (results.isInvalid()) {
0245 edm::LogWarning("DiMuonMassBiasClient") << "the current bin has invalid data" << std::endl;
0246 continue;
0247 }
0248
0249
0250 const Measurement1D& bias = results.getBias();
0251 meanHistos_[key]->setBinContent(bin, bias.value());
0252 meanHistos_[key]->setBinError(bin, bias.error());
0253
0254
0255 const Measurement1D& width = results.getWidth();
0256 widthHistos_[key]->setBinContent(bin, width.value());
0257 widthHistos_[key]->setBinError(bin, width.error());
0258 }
0259 }
0260
0261
0262 void DiMuonMassBiasClient::dqmEndJob(DQMStore::IBooker& ibooker, DQMStore::IGetter& igetter)
0263
0264 {
0265 edm::LogInfo("DiMuonMassBiasClient") << "DiMuonMassBiasClient::endLuminosityBlock";
0266
0267 getMEsToHarvest(igetter);
0268
0269
0270 if (useTH1s_) {
0271 bookMEs(ibooker);
0272 }
0273
0274 for (const auto& element : harvestTargets_) {
0275 if (!useTH1s_) {
0276
0277 this->fitAndFillProfile(element, ibooker);
0278 } else {
0279
0280 this->fitAndFillHisto(element, ibooker);
0281 }
0282 }
0283 }
0284
0285
0286 diMuonMassBias::fitOutputs DiMuonMassBiasClient::fitLineShape(TH1* hist, const bool& fitBackground) const
0287
0288 {
0289 if (hist->GetEntries() < diMuonMassBias::minimumHits) {
0290 edm::LogWarning("DiMuonMassBiasClient") << " Input histogram:" << hist->GetName() << " has not enough entries ("
0291 << hist->GetEntries() << ") for a meaningful Voigtian fit!\n"
0292 << "Skipping!";
0293
0294 return diMuonMassBias::fitOutputs(Measurement1D(0., 0.), Measurement1D(0., 0.));
0295 }
0296
0297 TCanvas* c1 = new TCanvas();
0298 if (debugMode_) {
0299 c1->Clear();
0300 c1->SetLeftMargin(0.15);
0301 c1->SetRightMargin(0.10);
0302 }
0303
0304
0305 RooMsgService::instance().setGlobalKillBelow(RooFit::FATAL);
0306
0307 Double_t xmin = hist->GetXaxis()->GetXmin();
0308 Double_t xmax = hist->GetXaxis()->GetXmax();
0309
0310 if (debugMode_) {
0311 edm::LogPrint("DiMuonMassBiasClient") << "fitting range: (" << xmin << "-" << xmax << ")" << std::endl;
0312 }
0313
0314 RooRealVar InvMass("InvMass", "di-muon mass M(#mu^{+}#mu^{-}) [GeV]", xmin, xmax);
0315 std::unique_ptr<RooPlot> frame{InvMass.frame()};
0316 RooDataHist datahist("datahist", "datahist", InvMass, RooFit::Import(*hist));
0317 datahist.plotOn(frame.get());
0318
0319
0320 RooRealVar mean("#mu", "mean", meanConfig_[0], meanConfig_[1], meanConfig_[2]);
0321 RooRealVar width("width", "width", widthConfig_[0], widthConfig_[1], widthConfig_[2]);
0322 RooRealVar sigma("#sigma", "sigma", sigmaConfig_[0], sigmaConfig_[1], sigmaConfig_[2]);
0323 RooVoigtian voigt("voigt", "voigt", InvMass, mean, width, sigma);
0324
0325
0326 RooRealVar peakCB("peakCB", "peakCB", meanConfig_[0], meanConfig_[1], meanConfig_[2]);
0327 RooRealVar sigmaCB("#sigma", "sigma", sigmaConfig_[0], sigmaConfig_[1], sigmaConfig_[2]);
0328 RooRealVar alphaCB("#alpha", "alpha", 1., 0., 10.);
0329 RooRealVar nCB("n", "n", 1., 0., 100.);
0330 RooCBShape crystalball("crystalball", "crystalball", InvMass, peakCB, sigmaCB, alphaCB, nCB);
0331
0332
0333 RooRealVar lambda("#lambda", "slope", 0., -50., 50.);
0334 RooExponential expo("expo", "expo", InvMass, lambda);
0335
0336
0337 RooRealVar exp_alpha("#alpha", "alpha", 40.0, 20.0, 160.0);
0338 RooRealVar exp_beta("#beta", "beta", 0.05, 0.0, 2.0);
0339 RooRealVar exp_gamma("#gamma", "gamma", 0.02, 0.0, 0.1);
0340 RooRealVar exp_peak("peak", "peak", meanConfig_[0]);
0341 RooCMSShape exp_pdf("exp_pdf", "bkg shape", InvMass, exp_alpha, exp_beta, exp_gamma, exp_peak);
0342
0343
0344 RooRealVar b("N_{b}", "Number of background events", 0, hist->GetEntries() / 10.);
0345 RooRealVar s("N_{s}", "Number of signal events", 0, hist->GetEntries());
0346
0347 if (fitBackground) {
0348 RooArgList listPdf;
0349 if (useRooCBShape_) {
0350 if (useRooCMSShape_) {
0351
0352 listPdf.add(crystalball);
0353 listPdf.add(exp_pdf);
0354 } else {
0355
0356 listPdf.add(crystalball);
0357 listPdf.add(expo);
0358 }
0359 } else {
0360 if (useRooCMSShape_) {
0361
0362 listPdf.add(voigt);
0363 listPdf.add(exp_pdf);
0364 } else {
0365
0366 listPdf.add(voigt);
0367 listPdf.add(expo);
0368 }
0369 }
0370
0371 RooAddPdf fullModel("fullModel", "Signal + Background Model", listPdf, RooArgList(s, b));
0372 fullModel.fitTo(datahist, RooFit::PrintLevel(-1));
0373 fullModel.plotOn(frame.get(), RooFit::LineColor(kRed));
0374 if (useRooCMSShape_) {
0375 fullModel.plotOn(frame.get(), RooFit::Components(exp_pdf), RooFit::LineStyle(kDashed));
0376 } else {
0377 fullModel.plotOn(frame.get(), RooFit::Components(expo), RooFit::LineStyle(kDashed));
0378 }
0379 fullModel.paramOn(frame.get(), RooFit::Layout(0.65, 0.90, 0.90));
0380 } else {
0381 if (useRooCBShape_) {
0382
0383 crystalball.fitTo(datahist, RooFit::PrintLevel(-1));
0384 crystalball.plotOn(frame.get(), RooFit::LineColor(kRed));
0385 crystalball.paramOn(frame.get(),
0386 RooFit::Layout(0.65, 0.90, 0.90));
0387 } else {
0388
0389 voigt.fitTo(datahist, RooFit::PrintLevel(-1));
0390 voigt.plotOn(frame.get(), RooFit::LineColor(kRed));
0391 voigt.paramOn(frame.get(), RooFit::Layout(0.65, 0.90, 0.90));
0392 }
0393 }
0394
0395
0396 datahist.plotOn(frame.get());
0397 frame->GetYaxis()->SetTitle("n. of events");
0398 TString histName = hist->GetName();
0399 frame->SetName("frame" + histName);
0400 frame->SetTitle(hist->GetTitle());
0401 frame->Draw();
0402
0403 if (debugMode_) {
0404 c1->Print("fit_debug" + histName + ".pdf");
0405 }
0406 delete c1;
0407
0408 float mass_mean = useRooCBShape_ ? peakCB.getVal() : mean.getVal();
0409 float mass_sigma = useRooCBShape_ ? sigmaCB.getVal() : sigma.getVal();
0410
0411 float mass_mean_err = useRooCBShape_ ? peakCB.getError() : mean.getError();
0412 float mass_sigma_err = useRooCBShape_ ? sigmaCB.getError() : sigma.getError();
0413
0414 Measurement1D resultM(mass_mean, mass_mean_err);
0415 Measurement1D resultW(mass_sigma, mass_sigma_err);
0416
0417 return diMuonMassBias::fitOutputs(resultM, resultW);
0418 }
0419
0420
0421 diMuonMassBias::fitOutputs DiMuonMassBiasClient::fitBWTimesCB(TH1* hist) const
0422
0423 {
0424 if (hist->GetEntries() < diMuonMassBias::minimumHits) {
0425 edm::LogWarning("DiMuonMassBiasClient") << " Input histogram:" << hist->GetName() << " has not enough entries ("
0426 << hist->GetEntries() << ") for a meaningful Voigtian fit!\n"
0427 << "Skipping!";
0428
0429 return diMuonMassBias::fitOutputs(Measurement1D(0., 0.), Measurement1D(0., 0.));
0430 }
0431
0432 TCanvas* c1 = new TCanvas();
0433 if (debugMode_) {
0434 c1->Clear();
0435 c1->SetLeftMargin(0.15);
0436 c1->SetRightMargin(0.10);
0437 }
0438
0439
0440 RooMsgService::instance().setGlobalKillBelow(RooFit::FATAL);
0441
0442 Double_t xMean = 91.1876;
0443 Double_t xMin = hist->GetXaxis()->GetXmin();
0444 Double_t xMax = hist->GetXaxis()->GetXmax();
0445
0446 if (debugMode_) {
0447 edm::LogPrint("DiMuonMassBiasClient") << "fitting range: (" << xMin << "-" << xMax << ")" << std::endl;
0448 }
0449
0450 double sigma(2.);
0451 double sigmaMin(0.1);
0452 double sigmaMax(10.);
0453
0454 double sigma2(0.1);
0455 double sigma2Min(0.);
0456 double sigma2Max(10.);
0457
0458 std::unique_ptr<FitWithRooFit> fitter = std::make_unique<FitWithRooFit>();
0459
0460 bool useChi2(false);
0461
0462 fitter->useChi2_ = useChi2;
0463 fitter->initMean(xMean, xMin, xMax);
0464 fitter->initSigma(sigma, sigmaMin, sigmaMax);
0465 fitter->initSigma2(sigma2, sigma2Min, sigma2Max);
0466 fitter->initAlpha(1.5, 0.05, 10.);
0467 fitter->initN(1, 0.01, 100.);
0468 fitter->initFGCB(0.4, 0., 1.);
0469 fitter->initGamma(2.4952, 0., 10.);
0470 fitter->gamma()->setConstant(kTRUE);
0471 fitter->initMean2(0., -20., 20.);
0472 fitter->mean2()->setConstant(kTRUE);
0473 fitter->initSigma(1.2, 0., 5.);
0474 fitter->initAlpha(1.5, 0.05, 10.);
0475 fitter->initN(1, 0.01, 100.);
0476 fitter->initExpCoeffA0(-1., -10., 10.);
0477 fitter->initExpCoeffA1(0., -10., 10.);
0478 fitter->initExpCoeffA2(0., -2., 2.);
0479 fitter->initFsig(0.9, 0., 1.);
0480 fitter->initA0(0., -10., 10.);
0481 fitter->initA1(0., -10., 10.);
0482 fitter->initA2(0., -10., 10.);
0483 fitter->initA3(0., -10., 10.);
0484 fitter->initA4(0., -10., 10.);
0485 fitter->initA5(0., -10., 10.);
0486 fitter->initA6(0., -10., 10.);
0487
0488 fitter->fit(hist, "breitWignerTimesCB", "exponential", xMin, xMax, false);
0489
0490 TString histName = hist->GetName();
0491
0492 if (debugMode_) {
0493 c1->Print("fit_debug" + histName + ".pdf");
0494 }
0495 delete c1;
0496
0497 Measurement1D resultM(fitter->mean()->getVal(), fitter->mean()->getError());
0498 Measurement1D resultW(fitter->sigma()->getVal(), fitter->sigma()->getError());
0499
0500 return diMuonMassBias::fitOutputs(resultM, resultW);
0501 }
0502
0503
0504 void DiMuonMassBiasClient::fillDescriptions(edm::ConfigurationDescriptions& descriptions)
0505
0506 {
0507 edm::ParameterSetDescription desc;
0508 desc.add<std::string>("FolderName", "DiMuonMassBiasMonitor");
0509 desc.add<bool>("useTH1s", false);
0510 desc.add<bool>("useBWtimesCB", false);
0511 desc.add<bool>("fitBackground", false);
0512 desc.add<bool>("useRooCMSShape", false);
0513 desc.add<bool>("useRooCBShape", false);
0514 desc.add<bool>("debugMode", false);
0515
0516 edm::ParameterSetDescription fit_par;
0517 fit_par.add<std::vector<double>>("mean_par",
0518 {std::numeric_limits<float>::max(),
0519 std::numeric_limits<float>::max(),
0520 std::numeric_limits<float>::max()});
0521
0522 fit_par.add<std::vector<double>>("width_par",
0523 {std::numeric_limits<float>::max(),
0524 std::numeric_limits<float>::max(),
0525 std::numeric_limits<float>::max()});
0526
0527 fit_par.add<std::vector<double>>("sigma_par",
0528 {std::numeric_limits<float>::max(),
0529 std::numeric_limits<float>::max(),
0530 std::numeric_limits<float>::max()});
0531
0532 desc.add<edm::ParameterSetDescription>("fit_par", fit_par);
0533
0534 desc.add<std::vector<std::string>>("MEtoHarvest",
0535 {"DiMuMassVsMuMuPhi",
0536 "DiMuMassVsMuMuEta",
0537 "DiMuMassVsMuPlusPhi",
0538 "DiMuMassVsMuPlusEta",
0539 "DiMuMassVsMuMinusPhi",
0540 "DiMuMassVsMuMinusEta",
0541 "DiMuMassVsMuMuDeltaEta",
0542 "DiMuMassVsCosThetaCS"});
0543 descriptions.addWithDefaultLabel(desc);
0544 }
0545
0546 DEFINE_FWK_MODULE(DiMuonMassBiasClient);