CorrSpec_MichelDecomp.cxx
Go to the documentation of this file.
1 #include <bitset>
2 #include "TH2.h"
3 
5 
7 #include "FNEX/core/VarBinning.h"
9 
10 namespace fnex{
11 
12  //----------------------------------------------------------------------------
14  : CorrectedSpectrum(pset)
15  {
16  this->Reconfigure(pset);
17  }
18 
19  //----------------------------------------------------------------------------
21  {
22  fVarName = pset.get< std::string >("VarName");
23  fIsConfigured = true;
24  fFarNormalizationType = pset.get< std::string >("FarNormalizationType");
25  fNearNormalizationType = pset.get< std::string >("NearNormalizationType");
26  fDoOscillate[novadaq::cnv::kNEARDET] = pset.get< bool >("OscillateND", false);
27  fDoOscillate[novadaq::cnv::kFARDET] = pset.get< bool >("OscillateFD", true );
28 
29  return;
30  }
31 
32  //----------------------------------------------------------------------------
34  {
35  std::string r = "CorrSpec_MichelDecomp_";
36  r += fVarName;
37  return r;
38  }
39 
40  //----------------------------------------------------------------------------
41  const TH1D & CorrSpec_MichelDecomp::DataForFit() const
42  {
44  }
45 
46  //----------------------------------------------------------------------------
47  const TH1D & CorrSpec_MichelDecomp::MCForFit() const
48  {
50  }
51 
52  //----------------------------------------------------------------------------
53  const std::vector<double> & CorrSpec_MichelDecomp::DataPointsForFit() const
54  {
55  return fDataPointsForFit;
56  }
57 
58  //----------------------------------------------------------------------------
59  const std::vector<double> & CorrSpec_MichelDecomp::WeightsForFit() const
60  {
61  return fWeightsForFit;
62  }
63 
64  //----------------------------------------------------------------------------
66  {
67 
69  throw cet::exception("CorrSpec_MichelDecomp")
70  << "SpectrumList not generated! Cannot SetPlotStyles.";
71 
72  // Both the Corrected and UnCorrected maps are expected to have
73  // the same MetaData keys present
74  int color = 0;
76  auto uncorritr = fUncorrectedHistMap.begin();
77  for(auto corritr : fCorrectedHistMap){
78  md = corritr.first;
79  uncorritr = fUncorrectedHistMap.find(md);
80 
81  if(uncorritr == fUncorrectedHistMap.end())
82  throw cet::exception("CorrSpec_MichelDecomp")
83  << "Unable to find metadata "
84  << md.ToString()
85  << " in the uncorrected histogram map";
86 
87  if(md.isMC){
88  if (md.fileType == fnex::kBeam ) color = 9;
89  else if(md.fileType == fnex::kTauSwap) color = 8;
90  else if(md.fileType == fnex::kSwap ) color = 7;
91 
92  corritr.second.SetFillColor(color);
93  corritr.second.SetFillStyle(1001);
94  corritr.second.SetOption("HIST");
95 
96  uncorritr->second.SetFillColor(color);
97  uncorritr->second.SetFillStyle(1001);
98  uncorritr->second.SetOption("HIST");
99  } // end if MC
100  else{
101  corritr.second.SetFillStyle(0);
102  corritr.second.SetMarkerColor(1);
103  corritr.second.SetMarkerStyle(2);
104  corritr.second.SetLineWidth(1);
105  corritr.second.SetOption("E1");
106 
107  uncorritr->second.SetFillStyle(0);
108  uncorritr->second.SetMarkerColor(1);
109  uncorritr->second.SetMarkerStyle(2);
110  uncorritr->second.SetLineWidth(1);
111  uncorritr->second.SetOption("E1");
112  } // end if data
113  } // end loop over spectrum metadata
114 
115 
116  // DetHistMap
117  std::set<novadaq::cnv::DetId> dets({novadaq::cnv::kNEARDET, novadaq::cnv::kFARDET});
118 
119  for(auto const& det : dets){
120 
121  // Corrected
122  fCorrectedHist_All_Data.at(det).SetMarkerColor(1);
123  fCorrectedHist_All_Data.at(det).SetMarkerStyle(2);
124  fCorrectedHist_All_Data.at(det).SetLineWidth(1);
125  fCorrectedHist_All_Data.at(det).SetOption("E1");
126 
127  fCorrectedHist_All_MC.at(det).SetFillStyle(0);
128  fCorrectedHist_All_MC.at(det).SetLineColor(2);
129  fCorrectedHist_All_MC.at(det).SetLineWidth(1);
130  fCorrectedHist_All_MC.at(det).SetOption("HIST");
131 
132  fCorrectedHist_All_Data_Over_MC.at(det).SetMarkerColor(1);
133  fCorrectedHist_All_Data_Over_MC.at(det).SetMarkerStyle(2);
134  fCorrectedHist_All_Data_Over_MC.at(det).SetLineWidth(1);
135  fCorrectedHist_All_Data_Over_MC.at(det).SetOption("E1");
136 
137  // Uncorrected
138  fUncorrectedHist_All_Data.at(det).SetMarkerColor(1);
139  fUncorrectedHist_All_Data.at(det).SetMarkerStyle(2);
140  fUncorrectedHist_All_Data.at(det).SetLineWidth(1);
141  fUncorrectedHist_All_Data.at(det).SetOption("E1");
142 
143  fUncorrectedHist_All_MC.at(det).SetFillStyle(0);
144  fUncorrectedHist_All_MC.at(det).SetLineColor(2);
145  fUncorrectedHist_All_MC.at(det).SetLineWidth(1);
146  fUncorrectedHist_All_MC.at(det).SetOption("HIST");
147 
148  fUncorrectedHist_All_Data_Over_MC.at(det).SetMarkerColor(1);
149  fUncorrectedHist_All_Data_Over_MC.at(det).SetMarkerStyle(2);
150  fUncorrectedHist_All_Data_Over_MC.at(det).SetLineWidth(1);
151  fUncorrectedHist_All_Data_Over_MC.at(det).SetOption("E1");
152 
153  }
154 
155  return;
156  }
157 
158  //----------------------------------------------------------------------------
160  {
161 
162  // Clear out existing SpectrumList if it's filled
163  fSpectrumListGenerated = false;
164  fSpectrumList.clear();
165 
166  // Keep a constant reference to the EventListMap
167  fEventListMap = &pEventListMap;
168 
169  fnex::Binning const* binning;
170 
171  if(fFitEvalType == "UBL" && (fVarName == "Nu_RecoE" ||
172  fVarName == "Lep_RecoE" ||
173  fVarName == "Had_RecoE" ))
174  binning = fnex::VarBinning::Instance()->BinInfo(fVarName + "_UBL");
175  else
177 
178  // keep track of the types of MetaData which are present, at a minimum
179  // require ND data, ND beam MC, ND flux swap MC, FD data, FD beam MC, FD flux swap MC,
180  // FD tau swap MC, set the bits in that order
181  // the bitset has all locations initialized to 0 by default.
182  std::bitset<7> requiredMD;
183 
184  // Look through pEventListMap, check for basic types required;
185  // if found, assign to class-wide MetaData definition for use as keys
187  size_t bit;
188  for(auto const& eventList : *fEventListMap){
189  md = eventList.first;
190 
191  // only use the nue selected events for this extrapolation
194 
195  fSpectrumList.emplace(md, Spectrum(eventList.second, binning, md.ToString()));
196 
197  // reset the bit
198  bit = 0;
199 
200  if(md.detector == novadaq::cnv::kFARDET) bit = 3;
201 
202  switch(md.fileType){
203  case fnex::kDataFile:
204  bit += 0;
205  break;
206  case fnex::kBeam:
207  bit += 1;
208  break;
209  case fnex::kSwap:
210  bit += 2;
211  break;
212  case fnex::kTauSwap:
213  bit += 3;
214  break;
216  break;
217  default:
218  throw cet::exception("CorrSpec_MichelDecomp")
219  << "Unknown FileType: " << md.ToString();
220  }
221 
222  LOG_DEBUG("CorrSpec_MichelDecomp")
223  << md.ToString()
224  << " " << requiredMD.to_string()
225  << " " << bit
226  << " " << md.fileType;
227 
228  requiredMD.set(bit);
229  }
230 
231  // check that the expected metadata are present.
232  if( !requiredMD.all() ){
233  //TODO: need to remove this check when we have ND flux swap events
234  requiredMD.set(2);
235  if( !requiredMD.all() )
236  throw cet::exception("CorrSpec_MichelDecomp")
237  << "Some expected metadata are missing: "
238  << requiredMD.to_string() << std::endl;
239  }
240 
241 
242  // we got here, so all the spectra are accounted for
243  fSpectrumListGenerated = true;
244 
245  // Now make the corrected histograms
246  LOG_VERBATIM("CorrSpec_MichelDecomp")
247  << "Making CorrSpec_MichelDecomp Histograms"
248  << " Spectrum list size: "
249  << fSpectrumList.size();
250 
252  this->MakeCorrectedHistograms();
253 
254  // We're done here
255  return fSpectrumListGenerated;
256  }
257 
258 
259  //----------------------------------------------------------------------------
261  {
262 
263  // Check that spectrum list is generated
265  throw cet::exception("CorrSpec_MichelDecomp")
266  << "Spectrum list not generated!";
267 
268  // Make histograms, appropriately scaled to data those factors
269  std::string histname;
270  std::string histtitle;
271  std::string xaxisname = fVarName;
272  std::string yaxisname = "Events";
273 
274  // make the HistMaps (uncombined hists)
275  // use the same binning as in the histograms for each metadata type
276  auto hist = fSpectrumList.begin()->second.Histogram(fVarName);
277  int bins = hist->GetNbinsX();
278  double low = hist->GetXaxis()->GetXmin();
279  double high = hist->GetXaxis()->GetXmax();
280 
282  for(auto & spectrum : fSpectrumList){
283  md = spectrum.first;
284  histname = md.ToString() + "_Uncorr_MichelDecomp_" + fVarName;
285  histtitle = md.ToStringLatex() + " (Uncorr::MichelDecomp) (" + fVarName + ")";
286  histtitle += ";" + xaxisname + ";" + yaxisname;
287  fUncorrectedHistMap.emplace(md, TH1D(histname.c_str(), histtitle.c_str(), bins, low, high) );
288  this->AddSpectrumHistogram(spectrum.second, fUncorrectedHistMap.at(md), fVarName);
289  }
290 
291  return true;
292  }
293 
294  //----------------------------------------------------------------------------
296  {
297 
298  // Check that spectrum list is generated
300  throw cet::exception("CorrSpec_MichelDecomp")
301  << "Spectrum list not generated!";
302 
303  // Make histograms, appropriately scaled to data those factors
304  std::string histname;
305  std::string histtitle;
306  std::string xaxisname = fVarName;
307  std::string yaxisname = "Events";
308 
309  // make the HistMaps (uncombined hists)
310  // use the same binning as in the histograms for each metadata type
311  auto hist = fSpectrumList.begin()->second.Histogram(fVarName);
312  int bins = hist->GetNbinsX();
313  double low = hist->GetXaxis()->GetXmin();
314  double high = hist->GetXaxis()->GetXmax();
315 
316  // First the uncorrected HistMap (uncombined hists)
318  for(auto & spectrum : fSpectrumList){
319  md = spectrum.first;
320  histname = md.ToString() + "_Corr_MichelDecomp_" + fVarName;
321  histtitle = md.ToStringLatex() + " (Corr::MichelDecomp) (" + fVarName + ")";
322  histtitle += ";" + xaxisname + ";" + yaxisname;
323  fCorrectedHistMap.emplace(md, TH1D(histname.c_str(), histtitle.c_str(), bins, low, high) );
324  this->AddSpectrumHistogram(spectrum.second, fCorrectedHistMap.at(md), fVarName);
325  }
326 
327  // Now name and title all of the histograms
328  // (ROOT's the worst; you'll get segfaults and other unpleasantries
329  // without unique names and titles for all TH1D instances)
330  std::set<std::pair<novadaq::cnv::DetId, std::string> > detectors({std::make_pair(novadaq::cnv::kNEARDET, "ND"),
332 
333  for(auto const& det : detectors){
334 
335  histname = det.second + "_Uncorr_MichelDecomp_Data_" + fVarName;
336  histtitle = det.second + "::Uncorr::MichelDecomp::Data::" + fVarName;
337  histtitle += ";" + xaxisname + ";" + yaxisname;
338  fUncorrectedHist_All_Data.emplace(det.first, TH1D(histname.c_str(), histtitle.c_str(), bins, low, high));
339 
340  fUncorrectedHist_All_Data.at(det.first).Sumw2();
341 
342  histname = det.second + "_Uncorr_MichelDecomp_MC_" + fVarName;
343  histtitle = det.second + "::Uncorr::MichelDecomp::MC::" + fVarName;
344  histtitle += ";" + xaxisname + ";" + yaxisname;
345  fUncorrectedHist_All_MC.emplace(det.first, TH1D(histname.c_str(), histtitle.c_str(), bins, low, high));
346  fUncorrectedHist_All_MC.at(det.first).Sumw2();
347 
348  histname = det.second + "_Uncorr_MichelDecomp_Data_Over_MC_" + fVarName;
349  histtitle = det.second + "::Uncorr::MichelDecomp::Data_Over_MC::" + fVarName;
350  histtitle += ";" + xaxisname + ";" + yaxisname;
351  fUncorrectedHist_All_Data_Over_MC.emplace(det.first, TH1D(histname.c_str(), histtitle.c_str(), bins, low, high));
352  fUncorrectedHist_All_Data_Over_MC.at(det.first).Sumw2();
353 
354  histname = det.second + "_Uncorr_MichelDecomp_Data_" + fVarName;
355  histtitle = det.second + "::Uncorr::MichelDecomp::Data::" + fVarName;
356  histtitle += ";" + xaxisname + ";" + yaxisname;
357  fCorrectedHist_All_Data.emplace(det.first, TH1D(histname.c_str(), histtitle.c_str(), bins, low, high));
358  fCorrectedHist_All_Data.at(det.first).Sumw2();
359 
360  histname = det.second + "_Uncorr_MichelDecomp_MC_" + fVarName;
361  histtitle = det.second + "::Uncorr::MichelDecomp::MC::" + fVarName;
362  histtitle += ";" + xaxisname + ";" + yaxisname;
363  fCorrectedHist_All_MC.emplace(det.first, TH1D(histname.c_str(), histtitle.c_str(), bins, low, high));
364  fCorrectedHist_All_MC.at(det.first).Sumw2();
365 
366  histname = det.second + "_Uncorr_MichelDecomp_Data_Over_MC_" + fVarName;
367  histtitle = det.second + "::Uncorr::MichelDecomp::Data_Over_MC::" + fVarName;
368  histtitle += ";" + xaxisname + ";" + yaxisname;
369  fCorrectedHist_All_Data_Over_MC.emplace(det.first, TH1D(histname.c_str(), histtitle.c_str(), bins, low, high));
370 
371  fCorrectedHist_All_Data_Over_MC.at(det.first).Sumw2();
372 
373  }
374 
375  LOG_VERBATIM("CorrSpec_MichelDecomp")
376  << "Finished making histograms. There are "
377  << "\nUncorrectedHistMap: " << fUncorrectedHistMap.size()
378  << "\nCorrectedHistMap: " << fCorrectedHistMap.size()
379  << "\nUncorrectedAllData: " << fUncorrectedHist_All_Data.size()
380  << "\nUncorrectedAllMC: " << fUncorrectedHist_All_MC.size()
381  << "\nUncorrectedAllData/MC: " << fUncorrectedHist_All_Data_Over_MC.size()
382  << "\nCorrectedAllData: " << fCorrectedHist_All_Data.size()
383  << "\nCorrectedAllMC: " << fCorrectedHist_All_MC.size()
384  << "\nCorrectedAllData/MC: " << fCorrectedHist_All_Data_Over_MC.size();
385 
386  //std::cout << " about to set plot styles" << histname << std::endl;
387  this->SetPlotStyles();
388 
389  // We're done here
390  return true;
391  }
392 
393 
394  //----------------------------------------------------------------------------
396  {
397 
398 
399  // KM June 6 2016
400  // Reset fDataPointsForFit/fWeightsForFit. These are used for UBL fits.
401  fDataPointsForFit.resize(0);
402  fWeightsForFit.resize(0);
403 
404  // Check that spectrum list is generated
406  throw cet::exception("CorrSpec_MichelDecomp")
407  << "Spectrum list not generated!";
408 
409  std::unordered_map<fnex::MetaData,
410  std::pair<double, double>,
411  fnex::MetaDataHasher> mdToEventsPOT;
412 
413  // store the number of events and POT for each detector's data
414  // the key is det * 1e4 + epoch
415  std::map<int, std::pair<double, double> > dataEventsPOT;
416 
417  // The key to fEventListMap and fSpectrumList is a fnex::MetaData,
418  // and both maps should have the same metadata present
419  fnex::MetaData md;
421  auto specItr = fSpectrumList.begin();
422  double events = 0.;
423  double pot = 0.;
424  int epochDetKey = 0;
425 
426  LOG_VERBATIM("CorrSpec_MichelDecompEventsPerPOT")
427  << "EVENTS/POT:";
428 
429  // loop over the event lists to ensure that we have
430  // a spectrum for each list
431 
432  for(auto evitr : *fEventListMap){
433  md = evitr.first;
434 
435  // we don't want to mess with numu selected events here
436  if(md.IsNuMuSelected()) continue;
437 
438  det = md.detector;
439  pot = evitr.second.ListSpillSummary().goodPOT;
440  events = 1. * fUncorrectedHistMap.at(md).GetEntries();
441 
442  mdToEventsPOT[md] = std::make_pair(events, pot);
443 
444  LOG_VERBATIM("CorrSpec_MichelDecompEventsPerPOT")
445  << std::setw(35) << md.ToString()
446  << " " << std::setw(15) << events
447  << " / " << std::setw(15) << pot
448  << " = " << std::setw(15) << events / pot;
449 
450  if(!md.isMC){
451  epochDetKey = 1e4 * det + md.epoch;
452  dataEventsPOT[epochDetKey] = std::make_pair(events, pot);
453  }
454 
455  specItr = fSpectrumList.find(md);
456  if(specItr == fSpectrumList.end() )
457  throw cet::exception("CorrSpec_MichelDecomp")
458  << "EventListMap and SpectrumList have different metadata keys: "
459  << md.ToString();
460 
461  // reset the uncorrected and corrected histograms, fill the uncorrected
462  // histograms now
463  fUncorrectedHistMap.at(md).Reset();
464  fCorrectedHistMap .at(md).Reset();
465  if(md.isMC){
466  this->AddSpectrumHistogram(specItr->second,
467  fUncorrectedHistMap.at(md),
468  fVarName);
469  }else{
470  if(det == novadaq::cnv::kFARDET){
471  // Grab FD data points as well as making histogram,
472  // and add to fDataPointsForFit. There will never be
473  // very many FD data events, this should not have
474  // an appreciable effect on speed.
475  this->AddSpectrumHistogram(specItr->second, fUncorrectedHistMap.at(md), fVarName);
476  std::vector<double> thisFDData = specItr->second.BaseVals(fVarName);
478  thisFDData.begin(),
479  thisFDData.end());
480 
481  // SJB Aug 3 2016 - also grab data weights in order to pass
482  // these to the UBL fit object so that UBL fit value
483  // contributions can be weighted correctly.
484  std::vector<double> thisFDWeights = specItr->second.BaseVals(fnex::KeyToVarName(fnex::kFake_Weight));
485  fWeightsForFit.insert(fWeightsForFit.end(),
486  thisFDWeights.begin(),
487  thisFDWeights.end());
488  }else{
489  this->AddSpectrumHistogram(specItr->second, fUncorrectedHistMap.at(md), fVarName);
490  }
491  }
492  } // end loop over metadata to fill maps
493 
494  // reset the (un)corrected histograms for data and MC in each detector
495  std::set<novadaq::cnv::DetId> detectors({novadaq::cnv::kNEARDET, novadaq::cnv::kFARDET});
496  for(auto const& det : detectors){
497  fUncorrectedHist_All_Data .at(det).Reset();
498  fUncorrectedHist_All_MC .at(det).Reset();
499  fUncorrectedHist_All_Data_Over_MC.at(det).Reset();
500  fCorrectedHist_All_Data .at(det).Reset();
501  fCorrectedHist_All_MC .at(det).Reset();
502  fCorrectedHist_All_Data_Over_MC .at(det).Reset();
503  }
504 
505  // Make scale factors
506  // loop over the event and pot maps to determine the scale factors, we
507  // couldn't do this above because we did not know where the data were
508  std::pair<double, double> eventsPOT;
510  double scale = 0.;
511 
512  for(auto itr : *fEventListMap){
513  md = itr.first;
514 
515  // We don't want to do anything with numu selected events in
516  // nue analysis
517  if(md.IsNuMuSelected()) continue;
518 
519 
520  det = md.detector;
521  scale = 1.;
522  eventsPOT = mdToEventsPOT.find(md)->second;
523  epochDetKey = 1e4 * det + md.epoch;
524 
525  // no need to normalize data
526  if(md.isMC){
527 
528  // figure out how we are normalizing for each detector
530  else if(det == novadaq::cnv::kFARDET ) norm = fFarNormalizationType;
531 
532  if(norm.compare("POT") == 0){
533  scale = dataEventsPOT.find(epochDetKey)->second.second/eventsPOT.second;
534  LOG_VERBATIM("CorrSpec_MichelDecompScale")
535  << "POT AND SCALES:"
536  << std::setw(50)
537  << md.ToString()
538  << " "
539  << std::setw(15)
540  << eventsPOT.second
541  << " / "
542  << scale;
543  }
544  else if(norm.compare("AREA") == 0){
545  scale = dataEventsPOT.find(epochDetKey)->second.first/eventsPOT.first;
546  LOG_VERBATIM("CorrSpec_MichelDecompScale")
547  << "AREA AND SCALES:\n"
548  << std::setw(50)
549  << md.ToString()
550  << " "
551  << std::setw(15)
552  << eventsPOT.first
553  << " / "
554  << scale;
555  }
556  else
557  LOG_VERBATIM("CorrSpec_MichelDecomp")
558  << std::setw(40) << md.ToString() << " NO KNOWN NEAR NORM TYPE REGISTERED, ALL SCALES SET = 1.0";
559  } // end if MC
560 
561  // scale the histograms for each type of metadata accordingly
562  fUncorrectedHistMap.at(md).Scale(scale);
563 
564  if( md.isMC ){
565  // combine the results for all the metadata for each detector
566  //std::cout << "Adding in PLACE 6 " << md.ToString()<< std::endl;
567  fUncorrectedHist_All_MC.at(det).Add( &(fUncorrectedHistMap.at(md)) );
568  }
569  else{
570  // there is only 1 metadata type for each detector for data, so if we
571  // are not looking at MC events, then md is the only option for this det
572  fUncorrectedHist_All_Data.at(det).Add( &(fUncorrectedHistMap.at(md) ) );
573  }
574  } // end loop over the event lists (ie MetaData) for uncorrected histograms
575 
576  // take the ratio of the data to MC for the uncorrected histograms
579 
582 
583 
584  //------------------------------------------------------
585  // Here's where the corrected histograms are made!!!
586  // #MichelDecompChampion #LindasWedding
587  //
588  // Possible interaction types are listed in FNEX/dataProducts/Constants.h
589  // Also listed here (accurate on 06/29/2016 @ 11:10 CDT)
590  //
591  // kUnknownInteraction = 0 (Also used to mean data)
592  // kNuMuCC,
593  // kNuMuBarCC,
594  // kNuECC,
595  // kNuEBarCC,
596  // kNuTauCC,
597  // kNuTauBarCC,
598  // kNC,
599  // kCosmicMuon
600  //------------------------------------------------------
601 
602  // Now make the corrected HistMap (uncombined hists)
603  for(auto & itr : *fEventListMap){
604  md = itr.first;
605  det = md.detector;
606 
607  // don't want to use numu selected events here
608  if(md.IsNuMuSelected()) continue;
609 
610  // start off with the corrected histogram being the same as the
611  // uncorrected one
612  fCorrectedHistMap.at(md).Add( &(fUncorrectedHistMap.at(md)) );
613 
614  if(md.isMC){
615 
616  // correct the FD MC histograms for the extrapolation from the ND
617  if(det == novadaq::cnv::kFARDET){
618 
619  // KM 09 07 2016 Pretty sure this is all we have to do for the backgrounds.
620  // Nue signal was already reweighted in the call to NearTrueEnergyRatio( )
621  if(md.interactionType == kNuMuCC ||
622  md.interactionType == kNuMuBarCC ||
623  md.interactionType == kNC ||
624  (md.interactionType == kNuECC && md.fileType == fnex::kBeam) ||
626  )
628 
629  } // end if FD
630 
631  fCorrectedHist_All_MC.at(det).Add( &(fCorrectedHistMap.at(md )) );
632 
633  } // end if MC
634  else{
635  // there is only 1 metadata type for each detector for data, so if we
636  // are not looking at MC events, then md is the only option for this det
637  fCorrectedHist_All_Data.at(det).Add( &(fCorrectedHistMap.at(md )) );
638  } // end if data
639 
640 
641  } // end loop over metadata for corrected histograms
642 
643  // make the ratios of data to MC
646 
649 
650  // We're done here
651  return true;
652  }
653 
654  //----------------------------------------------------------------------------
656  TH1D & hist,
657  std::string const& varName)
658  {
659  TH1D* specHist = (TH1D*)spectrum.Histogram(varName);
660 
661  if( specHist ) hist.Add(specHist);
662 
663  return;
664  }
665 
666 
667  //----------------------------------------------------------------------------
669  {
670 
671  // Assign reference input point
672  fRefInputPoint = pInputPoint;
673 
674  LOG_DEBUG("CorrSpec_MichelDecomp")
675  << "Applying point\n"
676  << " Getting shift collections:";
677 
678  // Now that we have reset the shifter and weighters for the input point
679  // in the MultiExperiment that called the Experiment that called this
680  // CorrectedSpectrum,
681  // we need to calculate the weights for the nue signal extrapolation - the
682  // Near E_true prediction to MC ratio
684 
685  // Currently just updates the histograms; no way to pass the
686  // information into the spectrum right now
687  LOG_DEBUG("CorrSpec_MichelDecomp")
688  << " Calling update on corrected histograms: ";
689 
690  LOG_DEBUG("CorrSpec_MichelDecomp")
691  << "CORRSPEC_MichelDecomp: Calling UpdateCorrectedHistograms : \n";
692 
694 
695  return;
696  }
697 
698 
699 } //end of namespace fnex
#define LOG_DEBUG(stream)
Definition: Messenger.h:149
std::vector< double > fDataPointsForFit
bool GenerateSpectrumList(EventListMap const &pEventListMap)
const EventListMap * fEventListMap
DetHistMap fUncorrectedHist_All_Data_Over_MC
bool IsNuMuSelected() const
Definition: Structs.cxx:256
const TH1 * Histogram(std::string const &varName)
Prepare a histogram using the fnex::Event objects in the fnex::EventList this Spectrum is using...
Definition: Spectrum.cxx:53
const std::vector< double > & WeightsForFit() const
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
std::map< fnex::MetaData, fnex::EventList > EventListMap
Definition: Event.h:186
std::pair< Spectrum *, CheatDecomp * > make_pair(SpectrumLoaderBase &loader_data, SpectrumLoaderBase &loader_mc, HistAxis *axis, Cut *cut, const SystShifts &shift, const Var &wei)
Definition: DataMCLoad.C:336
static const unsigned char kFake_Weight
Definition: VarVals.h:39
Create a list of fnex::Events to be used in fits.
void ApplyPoint(InputPoint const &pInputPoint)
static ShifterAndWeighter * Instance()
const std::string ToString() const
Double_t scale
Definition: plot.C:25
std::string const ToString() const
Definition: Structs.cxx:114
Far Detector at Ash River, MN.
SpectrumList fSpectrumList
Internal-use function that updates the weight collections for a new set of oscillation parameters...
The Spectrum is the means for interaction with Event collections.
Definition: Spectrum.h:29
T get(std::string const &key) const
Definition: ParameterSet.h:231
DetHistMap fCorrectedHist_All_Data_Over_MC
#define pot
Near Detector in the NuMI cavern.
const Binning bins
novadaq::cnv::DetId detector
Definition: Structs.h:50
static VarBinning * Instance()
Definition: VarBinning.cxx:15
fnex::FileType_t fileType
Definition: Structs.h:51
std::vector< double > fWeightsForFit
fnex::SelectionType_t selectionType
Definition: Structs.h:52
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
Float_t norm
void AddSpectrumHistogram(fnex::Spectrum &spectrum, TH1D &hist, std::string const &varName)
fnex::InteractionType_t interactionType
Definition: Structs.h:53
void CalcNearTrueEnergyRatioWeights(const fnex::EventListMap *pEventListMap)
void Reconfigure(fhicl::ParameterSet const &pset)
TRandom3 r(0)
CorrSpec_MichelDecomp(fhicl::ParameterSet const &pset)
std::string KeyToVarName(unsigned char const &key)
Definition: VarVals.h:445
const std::vector< double > & DataPointsForFit() const
#define LOG_VERBATIM(category)
Base class for &#39;binning&#39; objects. Pure virtual – you want either a FixedBinning or a VariableBinning...
Definition: Binning.h:20
std::string const ToStringLatex() const
Definition: Structs.cxx:144
fnex::Binning const * BinInfo(std::string const &varName)
Definition: VarBinning.cxx:65
void events(int which)
Definition: Cana.C:52