NuFCExperimentFactory Class Reference

Factory to generate pseudo-experiments for the nubar FC analysis. More...

#include <NuFCExperimentFactory.h>

List of all members.

Public Member Functions

 NuFCExperimentFactory (TRandom3 *rgen, TString nd_path, TString fd_path, TString tau_path)
 Construct the factory.
virtual ~NuFCExperimentFactory ()
const std::string & GetName () const
void GenerateNewExperiment (NuMMRunNuBar *mmRunSource, const NuMMParameters &mmPars, bool reroll=true)
 Generates a new experiment.
const NuMatrixSpectrumGetPQSpectrum () const
 Get a reference to the FD PQ Spectrum.
const NuMatrixSpectrumGetNewPQSpectrum () const
 Get a new copy (pointer) to the FD PQ Spectrum.
const NuMatrixSpectrumGetNQSpectrum () const
 Get a reference to the FD NQ Spectrum (blank, but with the right PoTs).
const NuMatrixSpectrumGetNewNQSpectrum () const
 Get a new copy (pointer) to the FD NQ Spectrum (blank, but with the right PoTs).
const NuMatrixSpectrumGetPQSpectrumND () const
 Get a reference to the ND PQ Spectrum.
const NuMatrixSpectrumGetNewPQSpectrumND () const
 Get a new copy (pointer) to the ND PQ Spectrum.
const NuMatrixSpectrumGetNQSpectrumND () const
 Get a reference to the ND NQ Spectrum (this is static).
const NuMatrixSpectrumGetNewNQSpectrumND () const
 Get a new copy (pointer) to the ND NQ Spectrum (this is static).
void SetGausSyst (bool _gs)
 Set to true to have gaussian distributed syst values and false to have them be uniform.
void SetSystNorm (double size)
 Normalization systematic. Shift value is 1 +/- size. 0 = off.
void SetSystCurv (double size)
 Absolute Track Energy from curvature far systematic. Shift value is 1 +/- size. 0 = off.
void SetSystRelCurv (double size)
 Relative Track Energy from curvature far systematic. Shift value is 1 +/- size. 0 = off.
void SetSystRange (double size)
 Track Energy from range systematic. Shift value is 1 +/- size. 0 = off.
void SetSystShower (double size)
 Shower Energy overall systematic. Shift value is +/- size sigma. 0 = off.
void SetSystRelShower (double size)
 Relative shower energy (far only). Shift value is 1 +/- size. 0 = off.
bool IsDP (const NuFCEvent *nu) const
void SetDPCorr (bool _cor)
 Set to true to do the fully correlated DP systematic.
void SetOldDP (bool _old)
 Set to true to do the the old symmetric DP syst.
void SetSystDP (double size)
 Decay systematic. Shift value is 1 +/- size. 0 = off.
void SetSystNC (double size)
 NC and WS background systematic. Shift value is 1 +/- size. 0 = off.
void SetSystWS (double size)
 WS background systematic. Shift value is 1 +/- size. 0 = off.
void SetSystXSec (double size)
 All combined and NuBar cross-section systematics. Shift value is +/- size sigma. 0 = off.
void SetSystSKZP (double size)
 SKZP flux systematics. Shift value is +/- size sigma. 0 = off.
void SetSystAccept (double size)
 Acceptance systematics. Shift value is +/- size sigma. 0 = off.
void NoSystematics ()
 Turn off all systematics.
void NoND (bool n=true)
 Don't generate a new ND spectrum each time. This means no ND systematics.
void DebugMode (bool deb=true)
 Fix systematics so the same is used in every experiment.
void PrintLoopProgress (Int_t entry, Float_t nEntries) const
 Print out a progress report in 5% increments.
int Counts (Sample_t s)
 Return the number of events in Sample s in last experiment.

Private Member Functions

void ScaleRejection (double &reject, int &nevents, double rand, double shift)
void FillMCEvents (TString mc_path)
void FillFromFilter (TTree *tin, int det)
void FillFromFCTree (TTree *tin, int det)
void FillFromDST (TString mc_path, TString xmlFileName, int det)
int Int (Sample_t s) const
void GenerateND (NuMMRunNuBar *mmRunSource)
int EBin (double E) const
void Distro (Sample_t s, TString name, TString title)
TH1D * Distro (Sample_t s) const
int N (Sample_t s) const
void StoreEvent (Sample_t s, float *vals, int det)
void StoreEvent (Sample_t s, NuEvent ev, int det)
void StoreEvent (Sample_t s, NuFCEvent ev, int det)
double GetReject (Sample_t s, int i) const
void MakeSampler (Sample_t s, NuMatrixSpectrum pred)
double GetNormShift () const
double GetTrkEnergy (const NuFCEvent *ev) const
double GetShowerEnergy (const NuFCEvent *ev) const
double GetDPShift (const NuFCEvent *ev) const
double GetBGShift (const NuFCEvent *ev) const
double GetBGShift (const Sample_t s) const
double GetXSecShift (const NuFCEvent *ev) const
double GetSKZPShift (const NuFCEvent *ev) const
double GetAcceptShift (const NuFCEvent *ev) const
double GetShift (double size) const
void RollSystematics ()
void InvertShifts ()
 ClassDef (NuFCExperimentFactory, 0)

Private Attributes

bool noND
const int IntNDPQ
const int IntNDNQ
double ndPoT
TH1D * spect [12]
TH1D * hSampler [12]
int counts [12]
std::vector< NuFCEvent * > events [12]
NuMatrixSpectrum fPQHist
NuMatrixSpectrum fNQHist
NuMatrixSpectrum fPQHistND
NuMatrixSpectrum fNQHistND
Int_t numRecoBins
Int_t numTrueBins
std::vector< Double_t > vReco
std::vector< Double_t > vTrue
TRandom * fRandy
std::string fName
 Pointer to the random number generator.
bool fixedSyst
bool gausSyst
double norm_size
double norm_shift
double curv_size
double curv_shift
double relcurv_size
double relcurv_shift
double range_size
double range_shift
double shower_size
double shower_shift
double relshow_size
double relshow_shift
bool dp_corr
bool dp_old
double dp_size
double dp_shift
double nc_size
double nc_shift
double ws_size
double ws_shift
double xsec_size
double xsec_shift
double skzp_size
double skzp_shift
double accept_size
double accept_shift

Static Private Attributes

static int counter = 0
 String name of this experiment.

Detailed Description

Factory to generate pseudo-experiments for the nubar FC analysis.

Author:
Alex Himmel, Last checkin
Author
mmathis
Version:
Revision
1.29
Date:
Date
2011/07/25 18:20:35

Created on: Sat Aug 16, 2008

Definition at line 36 of file NuFCExperimentFactory.h.


Constructor & Destructor Documentation

NuFCExperimentFactory::NuFCExperimentFactory ( TRandom3 *  rgen,
TString  nd_path,
TString  fd_path,
TString  tau_path 
)

Construct the factory.

The library files can be in one of three forms. They can be trees of NuFCEvent objects (these are available in the NuBar common area and can be generated with NuDSTAna::MakeFCTree), trees of the right variables made with NuFarNear, or simple NuDSTs.

Parameters:
rgen TRandom3* random engine.
nd_path Location of the near events library .root file.
fd_path Location of the far events library .root file.
tau_path Location of the tau events library .root file.

NuDSTAna::MakeFCTree

Definition at line 49 of file NuFCExperimentFactory.cxx.

References accept_size, counter, curv_size, Distro(), dp_size, FillMCEvents(), Plot::Format(), fRandy, Sample::kAppearedNQ, Sample::kAppearedPQ, Msg::kInfo, Sample::kNCNQ, Sample::kNCPQ, Sample::kSignalNQ, Sample::kSignalPQ, Sample::kTauNQ, Sample::kTauPQ, Sample::kWrongSignNQ, Sample::kWrongSignPQ, MAXMSG, nc_size, noND, norm_size, numRecoBins, numTrueBins, range_size, NuUtilities::RecoBins(), relcurv_size, relshow_size, shower_size, skzp_size, NuUtilities::TrueBins(), vReco, vTrue, ws_size, and xsec_size.

00050   : IntNDPQ(0), IntNDNQ(11), ndPoT(0), 
00051     fPQHist(), fNQHist(), fPQHistND(), fNQHistND(), 
00052     fRandy(rgen),  
00053     fixedSyst(false), gausSyst(true), dp_corr(false), dp_old(false)
00054 {
00055   MAXMSG("NuFCExperimentFactory",Msg::kInfo,1) << "Creating a NuFCExperimentFactory" << endl;
00056     
00057   // Check if our random number generator is correctly setup
00058   if (!fRandy) {
00059     // If not, just use gRandom
00060     fRandy = gRandom;
00061   }
00062   noND = false;
00063   
00064   norm_size    = 0.016;
00065 
00066   curv_size    = 0.03;
00067   relcurv_size = 0.00;
00068   range_size   = 0.02;
00069   shower_size  = 1.0;
00070   relshow_size = 0.0212; // 1.85 q+ 1.05
00071   
00072   dp_size      = 0.00;
00073   //bg_size      = 0.50;
00074   nc_size      = 0.20;
00075   ws_size      = 0.30;
00076   
00077   xsec_size    = 1;
00078   skzp_size    = 0;
00079   accept_size  = 0;
00080 
00081   // Get the binning scheme set up
00082   const NuUtilities cuts;
00083   NuBinningScheme::NuBinningScheme_t binningScheme =
00084   static_cast<NuBinningScheme::NuBinningScheme_t>(4);
00085   vReco = cuts.RecoBins(binningScheme);
00086   vTrue = cuts.TrueBins(binningScheme);
00087   
00088   numRecoBins = vReco.size() - 1;
00089   numTrueBins = vTrue.size() - 1;
00090   //recoBinsArray = &(vReco[0]);
00091   //trueBinsArray = &(vTrue[0]);
00092   
00093   TString suffix = TString::Format("_%i",counter);
00094 
00095   Distro(kSignalPQ, "signalSpect"+suffix, "Signal Distro");
00096   Distro(kWrongSignPQ, "wsSpect"+suffix, "Wrong Sign Distro");
00097   Distro(kNCPQ, "ncSpect"+suffix, "NC Distro");
00098   Distro(kTauPQ, "tauSpect"+suffix, "Tau Distro");
00099   Distro(kAppearedPQ, "appearedSpect"+suffix, "Appeared Distro");
00100 
00101   Distro(kSignalNQ, "signalSpectNQ"+suffix, "Signal Distro NQ");
00102   Distro(kWrongSignNQ, "wsSpectNQ"+suffix, "Wrong Sign Distro NQ");
00103   Distro(kNCNQ, "ncSpectNQ"+suffix, "NC Distro NQ");
00104   Distro(kTauNQ, "tauSpectNQ"+suffix, "Tau Distro NQ");
00105   Distro(kAppearedNQ, "appearedSpectNQ"+suffix, "Appeared Distro NQ");
00106   
00107   counter++;
00108   
00109   // Set up the FD MC chain
00110   FillMCEvents(nd_path);
00111   FillMCEvents(fd_path);
00112   FillMCEvents(tau_path);
00113     
00114   MAXMSG("NuFCExperimentFactory",Msg::kInfo,1) << "Done creating a NuFCExperimentFactory" << endl;
00115 }

virtual NuFCExperimentFactory::~NuFCExperimentFactory (  )  [inline, virtual]

Definition at line 42 of file NuFCExperimentFactory.h.

00042 {};


Member Function Documentation

NuFCExperimentFactory::ClassDef ( NuFCExperimentFactory  ,
 
) [private]
int NuFCExperimentFactory::Counts ( Sample_t  s  )  [inline]

Return the number of events in Sample s in last experiment.

Definition at line 120 of file NuFCExperimentFactory.h.

References counts, and Int().

00120 {return counts[Int(s)];};

void NuFCExperimentFactory::DebugMode ( bool  deb = true  ) 

Fix systematics so the same is used in every experiment.

Instead of randomly drawing the systematics for each new experiment, every systematic is always taken at its passed (maximum) value.

Parameters:
deb true = debug on, false = debug off. true by default.

Definition at line 917 of file NuFCExperimentFactory.cxx.

References fixedSyst, Msg::kInfo, and MSG.

00918 {
00919   fixedSyst = deb;
00920   
00921   if (fixedSyst) {
00922     MSG("NuFCExperimentFactory",Msg::kInfo) << "Now in debugging mode.  Systematics are not random -- they are fixed at the input values." << endl;
00923   }
00924   else {
00925     MSG("NuFCExperimentFactory",Msg::kInfo) << "Now out of debugging mode.  Systematics are randomly chosen based on input values." << endl;
00926   }
00927 }

TH1D * NuFCExperimentFactory::Distro ( Sample_t  s  )  const [private]

Return the unweighted spectrum of the library events for sample s

Definition at line 797 of file NuFCExperimentFactory.cxx.

References Int(), and spect.

00798 {
00799   return spect[Int(s)];
00800 }

void NuFCExperimentFactory::Distro ( Sample_t  s,
TString  name,
TString  title 
) [inline, private]

Initialize the histograms needed for sample s. This is the nweighted spectrum of events in the library and the rejection probability sampling historgram.

Definition at line 787 of file NuFCExperimentFactory.cxx.

References hSampler, Int(), numTrueBins, spect, and vTrue.

Referenced by MakeSampler(), NuFCExperimentFactory(), and StoreEvent().

00788 {
00789   spect[Int(s)] = new TH1D(name, title, numTrueBins, &(vTrue[0]));
00790   hSampler[Int(s)] = new TH1D();  // Prevents seg faults
00791 }

int NuFCExperimentFactory::EBin ( double  E  )  const [private]

A faster way of filling a histogram. Retuns the bin number that energy E refers to in binning scheme 4. Faster than Root's bin finder by more than a factor of 2.

Definition at line 717 of file NuFCExperimentFactory.cxx.

Referenced by GenerateND().

00718 {
00719   if (E < 0) return 0;
00720   if (E < 0.5) return 1;
00721   if (E < 20.) return static_cast<int>(E*4.);
00722   if (E < 30.) return static_cast<int>(E)+60;
00723   if (E < 50.) return static_cast<int>(E)/2+75;
00724   return 100;
00725 }

void NuFCExperimentFactory::FillFromDST ( TString  mc_path,
TString  xmlFileName,
int  det 
) [private]

Fill the event library for detector det from the NuDSTs in mc_path.

Definition at line 478 of file NuFCExperimentFactory.cxx.

References NuXMLConfig::AnaVersionString(), NuEvent::charge, NuInputEvents::GetEntriesNuEvent(), NuInputEvents::GetNextNuEvent(), NuEvent::iaction, NuInputEvents::InitialiseChains(), NuInputEvents::InitialiseNuEventBranch(), NuInputEvents::InputFileName(), NuEvent::inu, NuEvent::inunoosc, Sample::kAppearedNQ, Sample::kAppearedPQ, Msg::kInfo, Sample::kNCNQ, Sample::kNCPQ, Sample::kSignalNQ, Sample::kSignalPQ, Sample::kTauNQ, Sample::kTauPQ, Sample::kWrongSignNQ, Sample::kWrongSignPQ, MSG, NuInputEvents::ResetNuEventLoopPositionToStart(), and StoreEvent().

00479 {
00480   MSG("NuFCExperimentFactory",Msg::kInfo) << "Entering FillFromDST(" << mc_path << ")" << endl;
00481   
00482   // TChain holding all the FD MC
00483   NuInputEvents* mc_input = new NuInputEvents();
00484   mc_input->InputFileName(mc_path.Data());
00485   mc_input->InitialiseChains();
00486   mc_input->InitialiseNuEventBranch();  
00487   mc_input->ResetNuEventLoopPositionToStart();
00488   
00489   NuXMLConfig xmlConfig(xmlFileName);
00490   NuCutter cutter(xmlConfig.AnaVersionString(), 0);
00491   
00492   
00493   for (Int_t i=0; i < mc_input->GetEntriesNuEvent(); ++i) {
00494     NuEvent &nu=const_cast<NuEvent&>(mc_input->GetNextNuEvent(Msg::kInfo));
00495     
00496     //this->PrintLoopProgress(i,mc_input->GetEntriesNuEvent());
00497     
00498     //if (!IsGoodStdCuts(nu)) continue;
00499     cutter.MakeCuts(nu);
00500     if (cutter.Failed()) continue;
00501     
00502     if (nu.charge == 1) {
00503       if (nu.inu == 16 || nu.inu == -16) {
00504         // Ignore taus from Nue's, NC's
00505         if (nu.inunoosc == 12 || nu.inunoosc==-12 || nu.iaction == 0) continue;
00506         StoreEvent(kTauPQ, nu, det);
00507       }
00508       else if (nu.iaction == 0) StoreEvent(kNCPQ, nu, det);
00509       else if (nu.inu == -14)   {
00510         StoreEvent(kAppearedPQ, nu, det);
00511         StoreEvent(kSignalPQ, nu, det);
00512       }
00513       else StoreEvent(kWrongSignPQ, nu, det);
00514     }
00515     else {  
00516       if (nu.inu == 16 || nu.inu == -16) {
00517         // Ignore taus from Nue's, NC's
00518         if (nu.inunoosc == 12 || nu.inunoosc==-12 || nu.iaction == 0) continue;
00519         StoreEvent(kTauNQ, nu, det);
00520       }
00521       else if (nu.iaction == 0) StoreEvent(kNCNQ, nu, det);
00522       else if (nu.inu == 14)   {
00523         StoreEvent(kAppearedNQ, nu, det);
00524         StoreEvent(kSignalNQ, nu, det);
00525       }
00526       else StoreEvent(kWrongSignNQ, nu, det);
00527     }
00528   }
00529   
00530   MSG("NuFCExperimentFactory",Msg::kInfo) << "Leaving FillFromDST" << endl;  
00531 }

void NuFCExperimentFactory::FillFromFCTree ( TTree *  tin,
int  det 
) [private]

Fill the event library for detector det from the NuFCEvent tree tin.

Definition at line 432 of file NuFCExperimentFactory.cxx.

References NuFCEvent::charge, NuFCEvent::iaction, NuFCEvent::inu, NuFCEvent::inunoosc, Sample::kAppearedNQ, Sample::kAppearedPQ, Msg::kInfo, Sample::kNCNQ, Sample::kNCPQ, Sample::kSignalNQ, Sample::kSignalPQ, Sample::kTauNQ, Sample::kTauPQ, Sample::kWrongSignNQ, Sample::kWrongSignPQ, MSG, and StoreEvent().

Referenced by FillMCEvents().

00433 {
00434   MSG("NuFCExperimentFactory",Msg::kInfo) << "Entering FillFromFCTree" << endl;
00435   
00436   NuFCEvent *nu = new NuFCEvent();
00437   tin->SetBranchAddress("NuFCEvent", &nu);
00438 
00439   for (Int_t i=0; i < tin->GetEntries(); ++i) {
00440     tin->GetEntry(i);
00441     //this->PrintLoopProgress(i,tin->GetEntries());
00442 
00443     
00444     if (nu->charge == 1) {
00445       if (nu->inu == 16 || nu->inu == -16) {
00446         // Ignore taus from Nue's, NC's
00447         if (nu->inunoosc == 12 || nu->inunoosc==-12 || nu->iaction == 0) continue;
00448         StoreEvent(kTauPQ, *nu, det);
00449       }
00450       else if (nu->iaction == 0) StoreEvent(kNCPQ, *nu, det);
00451       else if (nu->inu == -14)   {
00452         StoreEvent(kAppearedPQ, *nu, det);
00453         StoreEvent(kSignalPQ, *nu, det);
00454       }
00455       else StoreEvent(kWrongSignPQ, *nu, det);
00456     }
00457     else {
00458       if (nu->inu == 16 || nu->inu == -16) {
00459         // Ignore taus from Nue's, NC's
00460         if (nu->inunoosc == 12 || nu->inunoosc==-12 || nu->iaction == 0) continue;
00461         StoreEvent(kTauNQ, *nu, det);
00462       }
00463       else if (nu->iaction == 0) StoreEvent(kNCNQ, *nu, det);
00464       else if (nu->inu == 14)   {
00465         StoreEvent(kAppearedNQ, *nu, det);
00466         StoreEvent(kSignalNQ, *nu, det);
00467       }
00468       else StoreEvent(kWrongSignNQ, *nu, det);
00469     }
00470   }
00471   MSG("NuFCExperimentFactory",Msg::kInfo) << "Leaving FillFromFCTree" << endl;  
00472 }

void NuFCExperimentFactory::FillFromFilter ( TTree *  tin,
int  det 
) [private]

Fill the event library for detector det from the NuFarNear Filter tree tin.

Definition at line 369 of file NuFCExperimentFactory.cxx.

References Sample::kAppearedNQ, Sample::kAppearedPQ, Msg::kInfo, Sample::kNCNQ, Sample::kNCPQ, Sample::kSignalNQ, Sample::kSignalPQ, Sample::kTauNQ, Sample::kTauPQ, Sample::kWrongSignNQ, Sample::kWrongSignPQ, MSG, and StoreEvent().

Referenced by FillMCEvents().

00370 {
00371   MSG("NuFCExperimentFactory",Msg::kInfo) << "Entering FillFromFilter" << endl;
00372   
00373   Float_t energy;       tin->SetBranchAddress("energy",       &energy);
00374   Float_t trkEn;        tin->SetBranchAddress("trkEn",        &trkEn);
00375   Float_t shwEn;        tin->SetBranchAddress("shwEn",        &shwEn);
00376   Float_t neuEnMC;      tin->SetBranchAddress("neuEnMC",      &neuEnMC);  
00377   Int_t   iaction;      tin->SetBranchAddress("iaction",      &iaction);    
00378   Int_t   inu;          tin->SetBranchAddress("inu",          &inu);
00379   Float_t charge;       tin->SetBranchAddress("charge",       &charge);
00380   Float_t ppvz;         tin->SetBranchAddress("ppvz",         &ppvz); 
00381   Int_t   containedTrk; tin->SetBranchAddress("containedTrk", &containedTrk); 
00382   Float_t fluxErr;      tin->SetBranchAddress("fluxErr",      &fluxErr); 
00383   
00384   float vals[10];
00385   for (Int_t i=0; i < tin->GetEntries(); ++i) {
00386     tin->GetEntry(i);
00387     //this->PrintLoopProgress(i,tin->GetEntries());
00388    
00389     vals[0] = energy;
00390     vals[1] = trkEn;
00391     vals[2] = shwEn;
00392     vals[3] = neuEnMC;
00393     vals[4] = (float)iaction;
00394     vals[5] = (float)inu;
00395     vals[6] = charge;
00396     vals[7] = ppvz;    
00397     vals[8] = (float)containedTrk;
00398     vals[9] = fluxErr;
00399 
00400     if (charge == 1) {
00401       if (inu == 16 || inu == -16) {
00402         if (iaction == 0) continue;
00403         StoreEvent(kTauPQ, vals, det);
00404       }
00405       else if (iaction == 0) StoreEvent(kNCPQ, vals, det);
00406       else if (inu == -14)   {
00407         StoreEvent(kAppearedPQ, vals, det);
00408         StoreEvent(kSignalPQ, vals, det);
00409       }
00410       else StoreEvent(kWrongSignPQ, vals, det);
00411     }
00412     else {
00413       if (inu == 16 || inu == -16) {
00414         if (iaction == 0) continue;
00415         StoreEvent(kTauNQ, vals, det);
00416       }
00417       else if (iaction == 0) StoreEvent(kNCNQ, vals, det);
00418       else if (inu == 14)   {
00419         StoreEvent(kAppearedNQ, vals, det);
00420         StoreEvent(kSignalNQ, vals, det);
00421       }
00422       else StoreEvent(kWrongSignNQ, vals, det);
00423     }
00424   }
00425   MSG("NuFCExperimentFactory",Msg::kInfo) << "Leaving FillFromFilter" << endl;  
00426 }

void NuFCExperimentFactory::FillMCEvents ( TString  mc_path  )  [private]

Fill the event library for the files in mc_path. What type of files and what detector they come from is determined automatically.

Definition at line 322 of file NuFCExperimentFactory.cxx.

References det, FillFromFCTree(), FillFromFilter(), Msg::kError, Msg::kInfo, MSG, and ndPoT.

Referenced by NuFCExperimentFactory().

00323 {
00324   MSG("NuFCExperimentFactory",Msg::kInfo) << "Entering FillMCEvents(\"" << mc_path << "\")" << endl;
00325   TFile *fin = new TFile(mc_path,"READ");
00326 
00327   // Test the file opened, and fail otherwise
00328   if (fin->IsZombie()) {
00329     MSG("NuFCExperimentFactory",Msg::kError) << "Failed to read in events file " << mc_path << endl;
00330     throw runtime_error("Failed to read MC events file");
00331   }
00332   
00333   // Determine which Detector this file is for
00334   TH1D *hDetector = (TH1D*)gROOT->FindObject("hDetector");
00335   int det = (int)hDetector->GetMean();
00336  
00337   // Get the ND Pot
00338   if (det == 1) {
00339     TH1D *hTotalPot = (TH1D*)gROOT->FindObject("hTotalPot");
00340     ndPoT += hTotalPot->Integral();
00341   }
00342   
00343   // Attempt to load the known tree types
00344   TTree *t1 = (TTree*)fin->Get("FDFilter");
00345   TTree *t2 = (TTree*)fin->Get("TAUFilter");
00346   TTree *t3 = (TTree*)fin->Get("NDFilter");
00347   TTree *t4 = (TTree*)fin->Get("FCTree");
00348   TTree *t5 = (TTree*)fin->Get("s");
00349 
00350   // If the tree exists in the file, load it
00351   if (t1)      FillFromFilter(t1, det);
00352   else if (t2) FillFromFilter(t2, det);
00353   else if (t3) FillFromFilter(t3, det);
00354   else if (t4) FillFromFCTree(t4, det);
00355   else if (t5) {
00356     MSG("NuFCExperimentFactory",Msg::kError) << "Use FillFromDST(path, xmlfile, detector) directly -- cannot use FillMCEvents for DST sources." << endl;
00357   }
00358   else {
00359     MSG("NuFCExperimentFactory",Msg::kError) << "No recognized tree in file: " << mc_path << endl;      
00360   }
00361   
00362   MSG("NuFCExperimentFactory",Msg::kInfo) << "Leaving FillMCEvents" << endl;  
00363 }  

void NuFCExperimentFactory::GenerateND ( NuMMRunNuBar mmRunSource  )  [private]

Just a copy of the standard cuts used in NuDSTana. It is needed in order to get the library events from NuDSTs. Generate the ND spectra. mmRunSource is needed to get the ND NQ spectrum since those events are not stored.

Definition at line 611 of file NuFCExperimentFactory.cxx.

References counter, EBin(), events, fNQHistND, fPQHistND, GetAcceptShift(), GetBGShift(), GetDPShift(), NuMMRunNuBar::GetNDBarData(), NuMMRunNuBar::GetNDNuData(), GetShowerEnergy(), GetSKZPShift(), GetTrkEnergy(), GetXSecShift(), IntNDNQ, IntNDPQ, it, Msg::kInfo, MAXMSG, ndPoT, noND, numRecoBins, NuMMRunNuBar::PredictNuBars(), NuMMRunNuBar::PredictNus(), NuMatrix::ResetPOT(), NuMatrixSpectrum::ResetSpectrum(), and vReco.

Referenced by GenerateNewExperiment().

00612 {
00613   if (noND) {
00614     MAXMSG("NuFCExperimentFactory",Msg::kInfo,1) << "Not Generating the ND" << endl;
00615     fPQHistND = *(mmRunSource->GetNDBarData());
00616     fNQHistND = *(mmRunSource->GetNDNuData());
00617     return;
00618   }
00619   MAXMSG("NuFCExperimentFactory",Msg::kInfo,1) << "Entering GenerateND" << endl;
00620 
00621   if (mmRunSource->PredictNuBars()) {
00622     TString name = "ndSpectPQ";
00623     name += counter;
00624     TH1D newND(name,"ND Spect PQ", numRecoBins, &(vReco[0]));
00625     TH1D newDP(name+"_dp","ND Spect PQ, DP", numRecoBins, &(vReco[0]));      
00626     /*
00627      unsigned int max = event[IntND].size();
00628      if (max == 0) {
00629      MSG("NuFCExperimentFactor",Msg::kError) << "No ND events loaded.  Cannot make an ND spectrum."
00630      }
00631      */
00632     
00633     double E;
00634     double rw;
00635     vector<NuFCEvent*>::const_iterator it;
00636     for (it = events[IntNDPQ].begin(); it != events[IntNDPQ].end(); ++it) {
00637       E = GetTrkEnergy(*it) + GetShowerEnergy(*it);
00638       rw = (*it)->rw;
00639       rw *= GetBGShift(*it);
00640       rw *= GetXSecShift(*it);
00641       rw *= GetSKZPShift(*it);
00642       rw *= GetAcceptShift(*it);
00643       
00644 //      if (dp_corr) {
00645 //        if (IsDP(*it)) newDP.AddBinContent(EBin(E), rw);
00646 //        else           newND.AddBinContent(EBin(E), rw);
00647 //      }
00648 //      else {
00649       rw *= GetDPShift(*it);
00650       newND.AddBinContent(EBin(E), rw);
00651 //      }
00652     }
00653     
00654 //    if (dp_corr) {
00655 //      // Set the decay pipe systematic size
00656 //      int lowb = 1;
00657 //      int highb = newND.GetXaxis()->FindFixBin(13.0);
00658 //      double U = newND.Integral(lowb, highb) / ndPoT * 2.92878e+20;
00659 //      double D = newDP.Integral(lowb, highb)/ ndPoT * 2.92878e+20;
00660 //      double N = 198064.;
00661 //      dp_shift = (N - U) / D;
00662 //      if (dp_shift < 0) dp_shift = 0;
00663 //      newDP.Scale(dp_shift);
00664 //      newND.Add(&newDP);
00665 //    }
00666     
00667     fPQHistND.ResetSpectrum(newND);
00668     fPQHistND.ResetPOT(ndPoT);
00669     MAXMSG("NuFCExperimentFactory",Msg::kInfo,1) << "Leaving GenerateND NuBar" << endl;
00670   }
00671   else {
00672     fPQHistND = *(mmRunSource->GetNDBarData());
00673   }
00674 
00675 
00676   if (mmRunSource->PredictNus()) {
00677     TString name = "ndSpectNQ";
00678     name += counter;
00679     TH1D newND(name,"ND Spect NQ", numRecoBins, &(vReco[0]));
00680     
00681     double E;
00682     double rw;
00683     
00684     vector<NuFCEvent*>::const_iterator it;
00685     for (it = events[IntNDNQ].begin(); it != events[IntNDNQ].end(); ++it) {
00686       E = GetTrkEnergy(*it) + GetShowerEnergy(*it);
00687       rw = (*it)->rw;
00688       rw *= GetBGShift(*it);
00689       rw *= GetSKZPShift(*it);
00690       rw *= GetAcceptShift(*it);
00691 
00692 
00693 //      cout << "Adding event: E=" << E 
00694 //      << ", GetTrkEnergy=" << GetTrkEnergy(*it)
00695 //      << ", trkEn=" << (*it)->trkEn
00696 //      << ", GetShowerEnergy=" << GetTrkEnergy(*it)
00697 //      << ", shwEn=" << (*it)->shwEn
00698 //      << ", rw=" << rw << endl;
00699       newND.AddBinContent(EBin(E), rw);
00700     }
00701     
00702     fNQHistND.ResetSpectrum(newND);
00703     fNQHistND.ResetPOT(ndPoT);
00704     MAXMSG("NuFCExperimentFactory",Msg::kInfo,1) << "Leaving GenerateND NuMu" << endl;
00705   }
00706   else {
00707     fNQHistND = *(mmRunSource->GetNDNuData());
00708   }
00709 }

void NuFCExperimentFactory::GenerateNewExperiment ( NuMMRunNuBar mmRunSource,
const NuMMParameters mmPars,
bool  reroll = true 
)

Generates a new experiment.

The internal NQ, PQ, Near, Far spectra are replaced with newly generated histograms with a new set of random systematics. The parameters are used to create the FD prediction that is the target distribution of the rejection sampling.

Parameters:
mmRunSource A NuMMRunFC or NuMMRunTransition object, as need. It also provides the ND NQ spectrum.
mmPars A set of parameters with the desired oscillations and/or transitions.

Definition at line 128 of file NuFCExperimentFactory.cxx.

References counter, counts, events, fNQHist, fPQHist, fRandy, GenerateND(), GetAcceptShift(), GetBGShift(), GetDPShift(), GetNormShift(), NuMatrix::GetPOT(), GetReject(), GetShowerEnergy(), GetSKZPShift(), GetTrkEnergy(), GetXSecShift(), Int(), Sample::kAppearedPQ, Msg::kInfo, Sample::kWrongSignPQ, MakeSampler(), MAXMSG, MSG, NuMatrixSpectrum::Multiply(), N(), numRecoBins, NuMMRunNuBar::PredictNuBars(), NuMMRunNuBar::PredictNus(), NuMMRun::QuietModeOn(), NuMatrix::ResetPOT(), NuMatrixSpectrum::ResetSpectrum(), RollSystematics(), ScaleRejection(), NuMatrixSpectrum::Spectrum(), NuMMRun::TrueComponents(), and vReco.

Referenced by NuFCGridPointNSINu::Run(), NuFCGridPointNSINubar::Run(), and NuFCGridPoint::Run().

00129 {
00130   MAXMSG("NuFCExperimentFactory",Msg::kInfo,1) << "Entering GenerateNewExperiment" << endl;
00131   
00132   TString name = "fdExp";
00133   name += (counter++);
00134 
00135   NuMatrixSpectrum fdTemp;
00136   mmRunSource->QuietModeOn();
00137 
00138   if (reroll) RollSystematics();
00139   GenerateND(mmRunSource);
00140   
00141   static int allDraws = 0;
00142   static int gthn1s = 0;
00143   static double maxProb = 0;
00144   
00145   if (mmRunSource->PredictNuBars()) {
00146     TH1D *fdExp = new TH1D(name,"Single FC Experiment", numRecoBins, &(vReco[0]));
00147     for (int i = 1; i <= 5; i++) {
00148       Sample_t s = (Sample_t)i;
00149       
00150       fdTemp = mmRunSource->TrueComponents(mmPars, s);
00151       
00152       // Apply Normalization to the Prediction
00153       fdTemp.Multiply(GetNormShift());
00154       fdTemp.Multiply(GetBGShift(s));
00155 
00156       // Get the average events by integrating the PDF
00157       Double_t average_events = fdTemp.Spectrum()->Integral();
00158       counts[Int(s)] = 0;
00159       if (average_events == 0) continue;
00160       
00161       // Get the number of events to generate
00162       Int_t nevents = fRandy->Poisson(average_events);
00163 
00164       // Create sampling distribution
00165       MakeSampler(s, fdTemp);
00166       
00167       double numdraws = 0;
00168 
00169       // Now start sampling
00170       int ev;
00171       double rand, reject, E;
00172       bool seeking;
00173       
00174       while (nevents > 0) {
00175         // Get an event -- takes care of rejection sampling
00176         seeking = true;
00177         while (seeking) {
00178           ev = fRandy->Integer(N(s));
00179           ++numdraws;
00180           
00181           rand = fRandy->Uniform();
00182           reject = GetReject(s, ev);
00183           
00184           if (s != kAppearedPQ) {
00185             ScaleRejection(reject, nevents, rand, GetDPShift(events[Int(s)][ev]));
00186           }
00187           if (s != kWrongSignPQ) {
00188             ScaleRejection(reject, nevents, rand, GetXSecShift(events[Int(s)][ev]));
00189           }
00190           ScaleRejection(reject, nevents, rand, GetSKZPShift(events[Int(s)][ev]));
00191           // MJM July 20, 2011 - add in acceptance systematic
00192           ScaleRejection(reject, nevents, rand, GetAcceptShift(events[Int(s)][ev]));
00193 
00194           ++allDraws;
00195           if (reject > maxProb) maxProb = reject;
00196           if (reject > 1.) ++gthn1s;
00197           seeking = rand > reject;
00198         }
00199 
00200         // Get the new reco energy -- apply the systematic
00201         E = GetTrkEnergy(events[Int(s)][ev]) + GetShowerEnergy(events[Int(s)][ev]);
00202 
00203         // Fill our final reco spectrum
00204         fdExp->Fill(E);
00205         ++counts[Int(s)];
00206         
00207         --nevents;
00208       }
00209     }
00210     fPQHist.ResetSpectrum(*fdExp);
00211     delete fdExp;
00212   }
00213   fPQHist.ResetPOT(fdTemp.GetPOT());
00214   
00215   
00216   
00217   if (mmRunSource->PredictNus()) {
00218     TH1D *fdExp = new TH1D(name,"Single FC Experiment", numRecoBins, &(vReco[0]));
00219     for (int i = -1; i >= -5; i--) {
00220       Sample_t s = (Sample_t)i;
00221       
00222       fdTemp = mmRunSource->TrueComponents(mmPars, s);
00223       
00224       // Apply Normalization to the Prediction
00225       fdTemp.Multiply(GetNormShift());
00226       fdTemp.Multiply(GetBGShift(s));
00227       
00228       // Get the average events by integrating the PDF
00229       Double_t average_events = fdTemp.Spectrum()->Integral();
00230       counts[Int(s)] = 0;
00231       if (average_events == 0) continue;
00232       
00233       // Get the number of events to generate
00234       Int_t nevents = fRandy->Poisson(average_events);
00235       
00236       // Create sampling distribution
00237       MakeSampler(s, fdTemp);
00238       
00239       double numdraws = 0;
00240       
00241       // Now start sampling
00242       int ev;
00243       double rand, reject, E;
00244       bool seeking;
00245       
00246       while (nevents > 0) {
00247         // Get an event -- takes care of rejection sampling
00248         seeking = true;
00249         while (seeking) {
00250           ev = fRandy->Integer(N(s));
00251           ++numdraws;
00252           
00253           rand = fRandy->Uniform();
00254           reject = GetReject(s, ev);
00255           
00256 //          if (s != kAppearedNQ) {
00257 //            ScaleRejection(reject, nevents, rand, GetDPShift(events[Int(s)][ev]));
00258 //          }
00259 //          if (s != kWrongSignNQ) {
00260 //            ScaleRejection(reject, nevents, rand, GetXSecShift(events[Int(s)][ev]));
00261 //          }
00262           ScaleRejection(reject, nevents, rand, GetSKZPShift(events[Int(s)][ev]));
00263           
00264           ++allDraws;
00265           if (reject > maxProb) maxProb = reject;
00266           if (reject > 1.) ++gthn1s;
00267           seeking = rand > reject;
00268         }
00269         
00270         // Get the new reco energy -- apply the systematic
00271         E = GetTrkEnergy(events[Int(s)][ev]) + GetShowerEnergy(events[Int(s)][ev]);
00272         
00273         // Fill our final reco spectrum
00274         fdExp->Fill(E);
00275         ++counts[Int(s)];
00276         
00277         --nevents;
00278       }
00279     }
00280     fNQHist.ResetSpectrum(*fdExp);
00281     delete fdExp;
00282   }
00283   fNQHist.ResetPOT(fdTemp.GetPOT());
00284   
00285   float frac = gthn1s;
00286   frac /= (float)allDraws;
00287   if (frac > 0.001) {
00288     MSG("NuFCExperimentFactory",Msg::kInfo) << "So far, max wt=" << maxProb
00289     << ", and " << gthn1s << "/" << allDraws << " had prob > 1" << endl;
00290   }
00291   
00292   MAXMSG("NuFCExperimentFactory",Msg::kInfo,1) << "Leaving GenerateNewExperiment" << endl;
00293 }

double NuFCExperimentFactory::GetAcceptShift ( const NuFCEvent nu  )  const [inline, private]

Return the acceptance systematic shift for this event. It is based on the acceptErr stored in nu.

Definition at line 1169 of file NuFCExperimentFactory.cxx.

References accept_shift, and NuFCEvent::acceptErr.

Referenced by GenerateND(), and GenerateNewExperiment().

01170 {
01171   return 1. + accept_shift*(nu->acceptErr - 1.);
01172 }

double NuFCExperimentFactory::GetBGShift ( const Sample_t  s  )  const [inline, private]

Return the wrong sign/nc background systematic shift for the total number of background events if s is a background sample. Return 1 (no effect) otherwise.

Definition at line 1123 of file NuFCExperimentFactory.cxx.

References Sample::kNCNQ, Sample::kNCPQ, Sample::kWrongSignNQ, Sample::kWrongSignPQ, nc_shift, and ws_shift.

01124 {
01125   //if (s == kWrongSignPQ || s == kNCPQ) return bg_shift;
01126   if (s == kNCPQ || s == kNCNQ) return nc_shift;
01127   if (s == kWrongSignPQ || s == kWrongSignNQ) return ws_shift;
01128   return 1;
01129 }

double NuFCExperimentFactory::GetBGShift ( const NuFCEvent nu  )  const [inline, private]

Return the wrong sign/nc background systematic weight shift for this event. Return 1 (no effect) if nu is not a background event.

Definition at line 1109 of file NuFCExperimentFactory.cxx.

References NuFCEvent::charge, NuFCEvent::iaction, NuFCEvent::inu, nc_shift, and ws_shift.

Referenced by GenerateND(), and GenerateNewExperiment().

01110 {
01111   //if (nu->iaction == 0 || nu->inu == 14) return bg_shift;
01112   if (nu->iaction == 0) return nc_shift;
01113   if (nu->inu == 14 && nu->charge == 1) return ws_shift;
01114   if (nu->inu == -14 && nu->charge == -1) return ws_shift;
01115   return 1;
01116 }

double NuFCExperimentFactory::GetDPShift ( const NuFCEvent nu  )  const [inline, private]

Return the decay pipe systematic weight shift. Return 1 (no effect) if nu is not a decay pipe event.

Definition at line 1099 of file NuFCExperimentFactory.cxx.

References dp_shift, and IsDP().

Referenced by GenerateND(), and GenerateNewExperiment().

01100 {
01101   if (IsDP(nu)) return dp_shift;
01102   return 1;  
01103 }

const std::string& NuFCExperimentFactory::GetName ( void   )  const [inline]

Definition at line 43 of file NuFCExperimentFactory.h.

References fName.

00043 { return fName; }

const NuMatrixSpectrum* NuFCExperimentFactory::GetNewNQSpectrum (  )  const [inline]

Get a new copy (pointer) to the FD NQ Spectrum (blank, but with the right PoTs).

Definition at line 55 of file NuFCExperimentFactory.h.

References fNQHist.

00055 { return new NuMatrixSpectrum(fNQHist); }

const NuMatrixSpectrum* NuFCExperimentFactory::GetNewNQSpectrumND (  )  const [inline]

Get a new copy (pointer) to the ND NQ Spectrum (this is static).

Definition at line 65 of file NuFCExperimentFactory.h.

References fNQHistND.

00065 { return new NuMatrixSpectrum(fNQHistND); }  

const NuMatrixSpectrum* NuFCExperimentFactory::GetNewPQSpectrum (  )  const [inline]

Get a new copy (pointer) to the FD PQ Spectrum.

Definition at line 50 of file NuFCExperimentFactory.h.

References fPQHist.

00050 { return new NuMatrixSpectrum(fPQHist); }

const NuMatrixSpectrum* NuFCExperimentFactory::GetNewPQSpectrumND (  )  const [inline]

Get a new copy (pointer) to the ND PQ Spectrum.

Definition at line 60 of file NuFCExperimentFactory.h.

References fPQHistND.

00060 { return new NuMatrixSpectrum(fPQHistND); }

double NuFCExperimentFactory::GetNormShift (  )  const [inline, private]

Return the normalization systematic shift. Use this to scale the total number of events being drawn.

Definition at line 1038 of file NuFCExperimentFactory.cxx.

References norm_shift.

Referenced by GenerateNewExperiment().

01039 {
01040   return norm_shift;
01041 }

const NuMatrixSpectrum& NuFCExperimentFactory::GetNQSpectrum (  )  const [inline]

Get a reference to the FD NQ Spectrum (blank, but with the right PoTs).

Definition at line 53 of file NuFCExperimentFactory.h.

References fNQHist.

Referenced by NuFCGridPointNSINu::Run(), NuFCGridPointNSINubar::Run(), and NuFCGridPoint::Run().

00053 { return fNQHist; }

const NuMatrixSpectrum& NuFCExperimentFactory::GetNQSpectrumND (  )  const [inline]

Get a reference to the ND NQ Spectrum (this is static).

Definition at line 63 of file NuFCExperimentFactory.h.

References fNQHistND.

Referenced by NuFCGridPointNSINu::Run(), NuFCGridPointNSINubar::Run(), and NuFCGridPoint::Run().

00063 { return fNQHistND; }

const NuMatrixSpectrum& NuFCExperimentFactory::GetPQSpectrum (  )  const [inline]

Get a reference to the FD PQ Spectrum.

Definition at line 48 of file NuFCExperimentFactory.h.

References fPQHist.

Referenced by NuFCGridPointNSINu::Run(), NuFCGridPointNSINubar::Run(), and NuFCGridPoint::Run().

00048 { return fPQHist; }

const NuMatrixSpectrum& NuFCExperimentFactory::GetPQSpectrumND (  )  const [inline]

Get a reference to the ND PQ Spectrum.

Definition at line 58 of file NuFCExperimentFactory.h.

References fPQHistND.

Referenced by NuFCGridPointNSINu::Run(), NuFCGridPointNSINubar::Run(), and NuFCGridPoint::Run().

00058 { return fPQHistND; }

double NuFCExperimentFactory::GetReject ( Sample_t  s,
int  i 
) const [private]

Get the rejection probability of event number i of sample s. Essentially, determine which energy bin this event is in and return the ratio in that bin.

MakeSampler for sample s needs to have been run for the results to be sensible.

MakeSampler

Definition at line 878 of file NuFCExperimentFactory.cxx.

References events, hSampler, Int(), Sample::kNCNQ, and Sample::kNCPQ.

Referenced by GenerateNewExperiment().

00879 {
00880   double Esel, reject;
00881   if (s == kNCPQ || s == kNCNQ) Esel = events[Int(s)][i]->energy;
00882   else                          Esel = events[Int(s)][i]->neuEnMC;
00883   
00884   int bin = hSampler[Int(s)]->GetXaxis()->FindFixBin(Esel);
00885   reject = hSampler[Int(s)]->GetBinContent(bin);
00886   return reject;
00887 }

double NuFCExperimentFactory::GetShift ( double  size  )  const [inline, private]

Get a systematic shift of the size passed. This can be random or or fixed depending on the debugging mode used.

Definition at line 967 of file NuFCExperimentFactory.cxx.

References fixedSyst, fRandy, and gausSyst.

Referenced by RollSystematics().

00968 { 
00969   if (fixedSyst)     return size;
00970   else if (gausSyst) return fRandy->Gaus(0, size);
00971   else               return fRandy->Uniform(-size, size); 
00972 }

double NuFCExperimentFactory::GetShowerEnergy ( const NuFCEvent nu  )  const [inline, private]

Return the systematically shifted shower energy. Accounts for both overall and relative shower energy scales.

Definition at line 1062 of file NuFCExperimentFactory.cxx.

References NuFCEvent::detector, relshow_shift, shower_shift, and NuFCEvent::shwEn.

Referenced by GenerateND(), and GenerateNewExperiment().

01063 {
01064   // Energy-dependent systematic as described in doc-7173
01065   double offset = 6.6;
01066   double scale = 3.5;
01067   double eLife = 1.44; //GeV
01068   double enShift = 1.0 + shower_shift*0.01*
01069                    (offset + scale*TMath::Exp(-1.0*nu->shwEn/eLife));
01070 
01071   double E = nu->shwEn * enShift;
01072   if (nu->detector == 2) E *= relshow_shift;
01073 //  cout << "GetShowerEnergy:" 
01074 //  << " shwEn=" << nu->shwEn 
01075 //  << ", shower_shift=" << shower_shift
01076 //  << ", shower_size=" << shower_size
01077 //  << ", enShift=" << enShift 
01078 //  << ", relshow_shift=" << relshow_shift
01079 //  << ", returning E=" << E << endl;
01080   return E;
01081 }

double NuFCExperimentFactory::GetSKZPShift ( const NuFCEvent nu  )  const [inline, private]

Return the skzp systematic shift for this event. It is based on the fluxErr stored in nu.

Definition at line 1159 of file NuFCExperimentFactory.cxx.

References NuFCEvent::fluxErr, and skzp_shift.

Referenced by GenerateND(), and GenerateNewExperiment().

01160 {
01161   return 1. + skzp_shift*(nu->fluxErr - 1.);
01162 }

double NuFCExperimentFactory::GetTrkEnergy ( const NuFCEvent nu  )  const [inline, private]

Return the systematically shifted track energy. Accounts for both range and curvature as required for this event.

Definition at line 1047 of file NuFCExperimentFactory.cxx.

References NuFCEvent::containedTrk, curv_shift, NuFCEvent::detector, range_shift, relcurv_shift, and NuFCEvent::trkEn.

Referenced by GenerateND(), and GenerateNewExperiment().

01048 {
01049   double E = nu->trkEn;
01050   if (nu->containedTrk == 1) E *= range_shift;
01051   else {
01052     E *= curv_shift;
01053     if (nu->detector == 2) E *= relcurv_shift;
01054   }
01055   return E;
01056 }

double NuFCExperimentFactory::GetXSecShift ( const NuFCEvent ev  )  const [inline, private]

Return the cross-section systematic weight shift for CC events. This includes both combined cross-section errors and the NuBar-Nu systematic ratio errors.

Definition at line 1135 of file NuFCExperimentFactory.cxx.

References MuELoss::a, NuFCEvent::charge, NuFCEvent::iaction, NuFCEvent::neuEnMC, and xsec_shift.

Referenced by GenerateND(), and GenerateNewExperiment().

01136 {
01137   if (ev->iaction != 1) return 1;
01138   if (ev->charge == -1) return 1; // X Sec is specific to NuBars
01139   
01140   static double xp = 25.;
01141   static double a = 0.895;
01142   static double xc = 4.0;
01143   static double b = 7.59e-3;
01144   static double c = -8.05e-4;
01145   double E = ev->neuEnMC;
01146   double shift = -0.0617;
01147   
01148   if (E < 25) shift += a - 1 + 2.*(1.-a)*E/xp + (a - 1.)*E*E/(xp*xp);
01149   if (E < xc) shift += (c*(xc-E)*(xc-E)*(xc-E)+b*(xc-E)*(xc-E));
01150   
01151   shift = 1 + xsec_shift * shift;
01152   return shift;
01153 }

int NuFCExperimentFactory::Int ( Sample_t  s  )  const [private]

Return the integer value corresponding to sample s. Used to access the appropriate array elements.

Definition at line 753 of file NuFCExperimentFactory.cxx.

References Sample::kAppearedNQ, Sample::kAppearedPQ, Msg::kError, Sample::kNCNQ, Sample::kNCPQ, Sample::kSignalNQ, Sample::kSignalPQ, Sample::kTauNQ, Sample::kTauPQ, Sample::kWrongSignNQ, Sample::kWrongSignPQ, and MSG.

Referenced by Counts(), Distro(), GenerateNewExperiment(), GetReject(), MakeSampler(), N(), and StoreEvent().

00754 {
00755   switch (s) {
00756       case kSignalPQ:
00757         return 1;
00758       case kWrongSignPQ:
00759         return 2;
00760       case kNCPQ:
00761         return 3;
00762       case kTauPQ:
00763         return 4;
00764       case kAppearedPQ:
00765         return 5;
00766       case kSignalNQ:
00767         return 6;
00768       case kWrongSignNQ:
00769         return 7;
00770       case kNCNQ:
00771         return 8;
00772       case kTauNQ:
00773         return 9;
00774       case kAppearedNQ:
00775         return 10;
00776       default:
00777         MSG("NuFCExperimentFactor",Msg::kError) << s << " is not a known Sample_t.  Returning -1" << endl;
00778         return -1;
00779   }
00780 }

void NuFCExperimentFactory::InvertShifts (  )  [private]

Set every shift (the randomly drawn ones) to its negative. No longer in use.

Definition at line 1013 of file NuFCExperimentFactory.cxx.

References accept_shift, curv_shift, dp_shift, nc_shift, norm_shift, range_shift, relcurv_shift, relshow_shift, shower_shift, skzp_shift, ws_shift, and xsec_shift.

01014 {
01015   // Get a single copy of these systematics
01016   norm_shift = 2 - norm_shift;
01017   
01018   curv_shift = 2 - curv_shift;
01019   relcurv_shift = 2 - curv_shift;
01020   range_shift = 2 - range_shift;
01021   shower_shift = 2 - shower_shift;
01022   relshow_shift = 2 - relshow_shift;
01023   
01024   dp_shift = 2 - dp_shift;
01025   //bg_shift = 2 - bg_shift;
01026   nc_shift = 2 - nc_shift;
01027   ws_shift = 2 - ws_shift;
01028   
01029   xsec_shift *= -1;
01030   skzp_shift *= -1;
01031   accept_shift *= -1;
01032 }

bool NuFCExperimentFactory::IsDP ( const NuFCEvent nu  )  const

Return true if nu is a decay pipe event.

Definition at line 1087 of file NuFCExperimentFactory.cxx.

References NuFCEvent::ppvz, and NuFCEvent::ptype.

Referenced by GetDPShift().

01088 {
01089   if (nu->ppvz < 4500) return false;
01090   if (abs(nu->ptype) == 13) return false;
01091   return true;
01092   
01093 }

void NuFCExperimentFactory::MakeSampler ( Sample_t  s,
NuMatrixSpectrum  pred 
) [private]

Create a sampling distribution (the ratio of target distribution to the distribution we currently have) for rejection sampling. This ratio will be normalized so its maximum value is 1.

Parameters:
s The sample were creating the sampler for.
pred The target distribution.

Definition at line 898 of file NuFCExperimentFactory.cxx.

References Distro(), hSampler, Int(), and NuMatrixSpectrum::Spectrum().

Referenced by GenerateNewExperiment().

00899 {
00900   if (hSampler[Int(s)]) delete hSampler[Int(s)];
00901   hSampler[Int(s)] = new TH1D(*(pred.Spectrum()));
00902   hSampler[Int(s)]->Divide(Distro(s));
00903   hSampler[Int(s)]->Scale(0.87/hSampler[Int(s)]->GetMaximum());
00904   // Scale to 0.87 (empirical #) instead of 1.0 to avoid 
00905   // clipping problems from having syst. shifted prob > 1
00906 }

int NuFCExperimentFactory::N ( Sample_t  s  )  const [inline, private]

Definition at line 152 of file NuFCExperimentFactory.h.

References events, and Int().

Referenced by GenerateNewExperiment().

00152 {return events[Int(s)].size();};

void NuFCExperimentFactory::NoND ( bool  n = true  )  [inline]

Don't generate a new ND spectrum each time. This means no ND systematics.

Definition at line 114 of file NuFCExperimentFactory.h.

References n, and noND.

00114 {noND = n;};

void NuFCExperimentFactory::NoSystematics (  ) 

Turn off all systematics.

Turns off all systematics and resets all the random shifts to 0. To test just one systematic, turn off all of them and then set the one you want to the desired value with SetSyst*.

Definition at line 937 of file NuFCExperimentFactory.cxx.

References accept_size, curv_size, dp_size, Msg::kInfo, MAXMSG, nc_size, norm_size, range_size, relcurv_size, relshow_size, RollSystematics(), shower_size, skzp_size, ws_size, and xsec_size.

Referenced by NuFCGridPoint::NuFCGridPoint(), NuFCGridPointNSINu::NuFCGridPointNSINu(), and NuFCGridPointNSINubar::NuFCGridPointNSINubar().

00938 {
00939   MAXMSG("NuFCExperimentFactory",Msg::kInfo,1) << "Turning off systematics" << endl;
00940   norm_size = 0;
00941   
00942   curv_size = 0;
00943   relcurv_size = 0;
00944   range_size = 0;
00945   shower_size = 0;
00946   relshow_size = 0;
00947   
00948   dp_size = 0;
00949   //bg_size = 0;
00950   nc_size = 0;
00951   ws_size = 0;
00952   
00953   xsec_size = 0;
00954   skzp_size = 0;
00955   accept_size = 0;
00956   
00957   // Zero out the shifts
00958   RollSystematics();
00959   MAXMSG("NuFCExperimentFactory",Msg::kInfo,1) << "Done turning off systematics" << endl;
00960 }

void NuFCExperimentFactory::PrintLoopProgress ( Int_t  entry,
Float_t  nEntries 
) const

Print out a progress report in 5% increments.

Parameters:
entry The current count.
nEntries The total count.

Definition at line 735 of file NuFCExperimentFactory.cxx.

References Msg::kInfo, and MSG.

00736 {
00737   Float_t fract=ceil(nEntries/20.);  
00738   if (ceil(((Float_t)entry)/fract)==((Float_t)entry)/fract){
00739     MSG("NuFCExperimentFactory",Msg::kInfo) 
00740     <<"Fraction of loop complete: "<<entry 
00741     <<"/"<<nEntries<<"  ("
00742     <<(Int_t)(100.*entry/nEntries)<<"%)"<<endl;
00743   }
00744 }

void NuFCExperimentFactory::RollSystematics (  )  [private]

Generate a new set of random systematic shifts.

Definition at line 978 of file NuFCExperimentFactory.cxx.

References accept_shift, accept_size, curv_shift, curv_size, dp_corr, dp_old, dp_shift, dp_size, GetShift(), nc_shift, nc_size, norm_shift, norm_size, range_shift, range_size, relcurv_shift, relcurv_size, relshow_shift, relshow_size, shower_shift, shower_size, skzp_shift, skzp_size, ws_shift, ws_size, xsec_shift, and xsec_size.

Referenced by GenerateNewExperiment(), and NoSystematics().

00979 {
00980   // Get a single copy of these systematics
00981   norm_shift = 1 + GetShift(norm_size);
00982   
00983   relcurv_shift = 1 + GetShift(relcurv_size);
00984 
00985   double tmp = GetShift(curv_size);
00986   double ratio = 0;
00987   if (curv_size != 0) ratio = range_size/curv_size;
00988   curv_shift = 1 + tmp;
00989   range_shift = 1 + tmp*ratio;
00990   
00991   shower_shift = GetShift(shower_size);
00992   relshow_shift = 1 + GetShift(relshow_size);
00993 
00994   if (!dp_corr) {
00995     if (dp_size == 0) dp_shift = 1;
00996     else if (dp_old)  dp_shift = 1.0 + GetShift(dp_size);
00997     else              dp_shift = 0.8 + GetShift(dp_size);
00998     if (dp_shift < 0) dp_shift = 0;
00999   }
01000   //bg_shift = 1 + GetShift(bg_size);
01001   nc_shift = 1 + GetShift(nc_size);
01002   ws_shift = 1 + GetShift(ws_size);
01003   
01004   xsec_shift = GetShift(xsec_size);
01005   skzp_shift = GetShift(skzp_size);
01006   accept_shift = GetShift(accept_size);
01007 }

void NuFCExperimentFactory::ScaleRejection ( double &  reject,
int &  nevents,
double  rand,
double  shift 
) [inline, private]

This is the rejection-sampling equivalent to reweighting an events. A new, shifted rejection probability is calculated -- this will be used to determine whether to keep or reject the event. But, if this new rejection probability will keep an event that would have been rejected, increase the total number we draw to account for this change in normalization. If we now reject an event we would have kept, decrease the total number of events.

Parameters:
reject The original rejection probability. This may be modified to the new rejection probability.
nevents The number of events remaining to be drawn for this experiment. This may be modified.
rand The random number drawn for rejection purposes.
shift The shift in event "weight."

Definition at line 308 of file NuFCExperimentFactory.cxx.

Referenced by GenerateNewExperiment().

00309 { 
00310   double reject_new = reject*shift;
00311   if (rand > reject && rand < reject_new) nevents++;
00312   if (rand < reject && rand > reject_new) nevents--;
00313   
00314   reject = reject_new;
00315 }

void NuFCExperimentFactory::SetDPCorr ( bool  _cor  )  [inline]

Set to true to do the fully correlated DP systematic.

Definition at line 92 of file NuFCExperimentFactory.h.

References dp_corr.

00092 {dp_corr = _cor;};

void NuFCExperimentFactory::SetGausSyst ( bool  _gs  )  [inline]

Set to true to have gaussian distributed syst values and false to have them be uniform.

Definition at line 68 of file NuFCExperimentFactory.h.

References gausSyst.

00068 {gausSyst = _gs;}

void NuFCExperimentFactory::SetOldDP ( bool  _old  )  [inline]

Set to true to do the the old symmetric DP syst.

Definition at line 94 of file NuFCExperimentFactory.h.

References dp_old.

00094 {dp_old = _old;};

void NuFCExperimentFactory::SetSystAccept ( double  size  )  [inline]

Acceptance systematics. Shift value is +/- size sigma. 0 = off.

Definition at line 109 of file NuFCExperimentFactory.h.

References accept_size.

00109 {accept_size = size;};

void NuFCExperimentFactory::SetSystCurv ( double  size  )  [inline]

Absolute Track Energy from curvature far systematic. Shift value is 1 +/- size. 0 = off.

Definition at line 79 of file NuFCExperimentFactory.h.

References curv_size.

00079 {curv_size = size;};

void NuFCExperimentFactory::SetSystDP ( double  size  )  [inline]

Decay systematic. Shift value is 1 +/- size. 0 = off.

Definition at line 96 of file NuFCExperimentFactory.h.

References dp_size.

00096 {dp_size = size;};

void NuFCExperimentFactory::SetSystNC ( double  size  )  [inline]

NC and WS background systematic. Shift value is 1 +/- size. 0 = off.

NC background systematic. Shift value is 1 +/- size. 0 = off.

Definition at line 100 of file NuFCExperimentFactory.h.

References nc_size.

00100 {nc_size = size;};

void NuFCExperimentFactory::SetSystNorm ( double  size  )  [inline]

Normalization systematic. Shift value is 1 +/- size. 0 = off.

SetSyst* Functions -- set the size of these systematic shifts The actual shifted value is varied from experiment to experiment Shift values are defined as 1 +/- size Set to 0 to turn off that systematic

Definition at line 76 of file NuFCExperimentFactory.h.

References norm_size.

00076 {norm_size = size;};

void NuFCExperimentFactory::SetSystRange ( double  size  )  [inline]

Track Energy from range systematic. Shift value is 1 +/- size. 0 = off.

Definition at line 83 of file NuFCExperimentFactory.h.

References range_size.

00083 {range_size = size;};

void NuFCExperimentFactory::SetSystRelCurv ( double  size  )  [inline]

Relative Track Energy from curvature far systematic. Shift value is 1 +/- size. 0 = off.

Definition at line 81 of file NuFCExperimentFactory.h.

References relcurv_size.

00081 {relcurv_size = size;};

void NuFCExperimentFactory::SetSystRelShower ( double  size  )  [inline]

Relative shower energy (far only). Shift value is 1 +/- size. 0 = off.

Definition at line 87 of file NuFCExperimentFactory.h.

References relshow_size.

00087 {relshow_size = size;};

void NuFCExperimentFactory::SetSystShower ( double  size  )  [inline]

Shower Energy overall systematic. Shift value is +/- size sigma. 0 = off.

Definition at line 85 of file NuFCExperimentFactory.h.

References shower_size.

00085 {shower_size = size;};

void NuFCExperimentFactory::SetSystSKZP ( double  size  )  [inline]

SKZP flux systematics. Shift value is +/- size sigma. 0 = off.

Definition at line 107 of file NuFCExperimentFactory.h.

References skzp_size.

00107 {skzp_size = size;};

void NuFCExperimentFactory::SetSystWS ( double  size  )  [inline]

WS background systematic. Shift value is 1 +/- size. 0 = off.

Definition at line 102 of file NuFCExperimentFactory.h.

References ws_size.

00102 {ws_size = size;};

void NuFCExperimentFactory::SetSystXSec ( double  size  )  [inline]

All combined and NuBar cross-section systematics. Shift value is +/- size sigma. 0 = off.

Definition at line 105 of file NuFCExperimentFactory.h.

References xsec_size.

00105 {xsec_size = size;};

void NuFCExperimentFactory::StoreEvent ( Sample_t  s,
NuFCEvent  ev,
int  det 
) [private]

Store the NuFCEvent ev in the sample s cache. Set its detector type to det.

Definition at line 841 of file NuFCExperimentFactory.cxx.

References NuFCEvent::charge, NuFCEvent::detector, Distro(), NuFCEvent::energy, events, Int(), IntNDNQ, IntNDPQ, Sample::kAppearedNQ, Sample::kAppearedPQ, Sample::kNCNQ, Sample::kNCPQ, Msg::kWarning, MSG, and NuFCEvent::neuEnMC.

00842 {
00843   NuFCEvent * loc = new NuFCEvent(ev);
00844   loc->detector = det;
00845   
00846   if (loc->detector == 2) { // Far Det and Taus
00847     events[Int(s)].push_back(loc);
00848     
00849     double E;
00850     if (s == kNCPQ || s == kNCNQ) E = loc->energy;
00851     else                        E = loc->neuEnMC;
00852     
00853     Distro(s)->Fill(E);
00854   }
00855   else if (loc->detector == 1) {
00856     // Prevent double counting the appeared events.
00857     if (s == kAppearedPQ || s == kAppearedNQ) return;
00858       
00859     if (loc->charge == 1) events[IntNDPQ].push_back(loc);
00860     else                  events[IntNDNQ].push_back(loc);
00861   }
00862   else {
00863     MSG("NuFCExperimentFactory",Msg::kWarning) << "Detector " << det << " not recognized.  Event not stored." << endl;
00864   }
00865 }

void NuFCExperimentFactory::StoreEvent ( Sample_t  s,
NuEvent  ev,
int  det 
) [private]

Store the NuEvent ev in the sample s cache. Set its detector type to det. It will be converted to a NuFCEvent first.

Definition at line 830 of file NuFCExperimentFactory.cxx.

References StoreEvent().

00831 {
00832   NuFCEvent loc(ev);
00833   StoreEvent(s, loc, det);
00834 }

void NuFCExperimentFactory::StoreEvent ( Sample_t  s,
float *  vals,
int  det 
) [private]

Store the event described by the values in vals in the sample s cache. Set its detector typ to det. The expected order of fields is: energy, trkEn, shwEn, neuEnMC, iaction, inu, charge, ppvz, containedTrk, fluxErr.

Definition at line 808 of file NuFCExperimentFactory.cxx.

References NuFCEvent::charge, NuFCEvent::containedTrk, NuFCEvent::energy, NuFCEvent::fluxErr, NuFCEvent::iaction, NuFCEvent::inu, NuFCEvent::neuEnMC, NuFCEvent::ppvz, NuFCEvent::shwEn, and NuFCEvent::trkEn.

Referenced by FillFromDST(), FillFromFCTree(), FillFromFilter(), and StoreEvent().

00809 {
00810   NuFCEvent loc;
00811   loc.energy       = vals[0];
00812   loc.trkEn        = vals[1];
00813   loc.shwEn        = vals[2];
00814   loc.neuEnMC      = vals[3];
00815   loc.iaction      = (Int_t)vals[4];
00816   loc.inu          = (Int_t)vals[5];
00817   loc.charge       = (Int_t)vals[6];
00818   loc.ppvz         = vals[7];
00819   loc.containedTrk = (Int_t)vals[8];
00820   loc.fluxErr      = vals[9];    
00821   
00822   StoreEvent(s, loc, det);
00823 }


Member Data Documentation

Definition at line 224 of file NuFCExperimentFactory.h.

Referenced by GetAcceptShift(), InvertShifts(), and RollSystematics().

int NuFCExperimentFactory::counter = 0 [static, private]

String name of this experiment.

Definition at line 196 of file NuFCExperimentFactory.h.

Referenced by GenerateND(), GenerateNewExperiment(), and NuFCExperimentFactory().

int NuFCExperimentFactory::counts[12] [private]

Definition at line 174 of file NuFCExperimentFactory.h.

Referenced by Counts(), and GenerateNewExperiment().

Definition at line 209 of file NuFCExperimentFactory.h.

Referenced by GetTrkEnergy(), InvertShifts(), and RollSystematics().

Definition at line 215 of file NuFCExperimentFactory.h.

Referenced by RollSystematics(), and SetDPCorr().

Definition at line 216 of file NuFCExperimentFactory.h.

Referenced by RollSystematics(), and SetOldDP().

Definition at line 217 of file NuFCExperimentFactory.h.

Referenced by GetDPShift(), InvertShifts(), and RollSystematics().

std::vector<NuFCEvent*> NuFCExperimentFactory::events[12] [private]

Definition at line 177 of file NuFCExperimentFactory.h.

Referenced by GenerateND(), GenerateNewExperiment(), GetReject(), N(), and StoreEvent().

Definition at line 204 of file NuFCExperimentFactory.h.

Referenced by DebugMode(), and GetShift().

std::string NuFCExperimentFactory::fName [private]

Pointer to the random number generator.

Definition at line 194 of file NuFCExperimentFactory.h.

Referenced by GetName().

Definition at line 181 of file NuFCExperimentFactory.h.

Referenced by GenerateNewExperiment(), GetNewNQSpectrum(), and GetNQSpectrum().

Definition at line 183 of file NuFCExperimentFactory.h.

Referenced by GenerateND(), GetNewNQSpectrumND(), and GetNQSpectrumND().

Definition at line 180 of file NuFCExperimentFactory.h.

Referenced by GenerateNewExperiment(), GetNewPQSpectrum(), and GetPQSpectrum().

Definition at line 182 of file NuFCExperimentFactory.h.

Referenced by GenerateND(), GetNewPQSpectrumND(), and GetPQSpectrumND().

TRandom* NuFCExperimentFactory::fRandy [private]

Definition at line 193 of file NuFCExperimentFactory.h.

Referenced by GenerateNewExperiment(), GetShift(), and NuFCExperimentFactory().

Definition at line 205 of file NuFCExperimentFactory.h.

Referenced by GetShift(), and SetGausSyst().

TH1D* NuFCExperimentFactory::hSampler[12] [private]

Definition at line 173 of file NuFCExperimentFactory.h.

Referenced by Distro(), GetReject(), and MakeSampler().

const int NuFCExperimentFactory::IntNDNQ [private]

Definition at line 142 of file NuFCExperimentFactory.h.

Referenced by GenerateND(), and StoreEvent().

const int NuFCExperimentFactory::IntNDPQ [private]

Definition at line 142 of file NuFCExperimentFactory.h.

Referenced by GenerateND(), and StoreEvent().

Definition at line 219 of file NuFCExperimentFactory.h.

Referenced by GetBGShift(), InvertShifts(), and RollSystematics().

double NuFCExperimentFactory::ndPoT [private]

Definition at line 143 of file NuFCExperimentFactory.h.

Referenced by FillMCEvents(), and GenerateND().

Definition at line 141 of file NuFCExperimentFactory.h.

Referenced by GenerateND(), NoND(), and NuFCExperimentFactory().

Definition at line 207 of file NuFCExperimentFactory.h.

Referenced by GetNormShift(), InvertShifts(), and RollSystematics().

Definition at line 188 of file NuFCExperimentFactory.h.

Referenced by Distro(), and NuFCExperimentFactory().

Definition at line 211 of file NuFCExperimentFactory.h.

Referenced by GetTrkEnergy(), InvertShifts(), and RollSystematics().

Definition at line 210 of file NuFCExperimentFactory.h.

Referenced by GetTrkEnergy(), InvertShifts(), and RollSystematics().

Definition at line 213 of file NuFCExperimentFactory.h.

Referenced by GetShowerEnergy(), InvertShifts(), and RollSystematics().

Definition at line 212 of file NuFCExperimentFactory.h.

Referenced by GetShowerEnergy(), InvertShifts(), and RollSystematics().

Definition at line 223 of file NuFCExperimentFactory.h.

Referenced by GetSKZPShift(), InvertShifts(), and RollSystematics().

TH1D* NuFCExperimentFactory::spect[12] [private]

Definition at line 172 of file NuFCExperimentFactory.h.

Referenced by Distro().

std::vector<Double_t> NuFCExperimentFactory::vReco [private]
std::vector<Double_t> NuFCExperimentFactory::vTrue [private]

Definition at line 190 of file NuFCExperimentFactory.h.

Referenced by Distro(), and NuFCExperimentFactory().

Definition at line 220 of file NuFCExperimentFactory.h.

Referenced by GetBGShift(), InvertShifts(), and RollSystematics().

Definition at line 222 of file NuFCExperimentFactory.h.

Referenced by GetXSecShift(), InvertShifts(), and RollSystematics().


The documentation for this class was generated from the following files:

Generated on 2 Nov 2017 for loon by  doxygen 1.6.1