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