NuMMRunTransSME Class Reference

#include <NuMMRunTransSME.h>

Inheritance diagram for NuMMRunTransSME:
NuMMRunNuBar NuMMRun

List of all members.

Public Types

enum  ESyst { kShwScaleTrans, kTrkScaleTrans, kDecayPipe, kNumSyst }

Public Member Functions

 NuMMRunTransSME (NuMMHelperCPT *helper, NuMatrixSpectrum *ndNuData, NuMatrixSpectrum *ndBarData, double pot)
 NuMMRunTransSME (NuMMHelperCPT *helper, NuMatrixSpectrum *ndNuData, NuMatrixSpectrum *ndBarData, NuMatrixSpectrum *fdNuData, NuMatrixSpectrum *fdBarData)
 ~NuMMRunTransSME ()
void SetSystShifts (TFile *systFilePQ, TFile *systFileNQ)
virtual void FitSysts (Bool_t fitSysts)
Double_t ComparePredWithData (const NuMMParameters &pars)
pair< NuMatrixSpectrum,
NuMatrixSpectrum
MakeFDPred (const NuMMParameters &pars)
 Return a pair (NQ,PQ) of NuMatrixSpectra of the far detector prediction with oscillation parameters described in pars. Pure virtual.
std::vector< TH1D > WriteFDPredHistos (const NuMMParameters &pars) const
NuMatrixSpectrum TrueComponents (const NuMMParameters &pars, Sample_t s) const
 Return the individual true energy prediction of the far detector sub sample s with oscillation parameters described in pars. Default implementation just drops an error.
NuMatrixSpectrum TrueComponents (const NuMMParameters &pars, Sample_t s)
virtual NuMatrixSpectrumMakeFDPredNuMu (const NuMMParameters &pars)
 Get just the NuMu component of the prediction with oscillations parameters pars.
virtual NuMatrixSpectrumMakeFDPredNuBar (const NuMMParameters &pars)
 Get just the NuMuBar component of the prediction with oscillations parameters pars.
void OscillateNuMu (NuMatrixSpectrum &numu, const NuMMParameters &pars) const
void InverseOscillateNuTau (NuMatrixSpectrum &nutau, const NuMMParameters &pars) const
 Helper. Inverse oscillates a true energy nutau spectrum.
void OscillateNuBar (NuMatrixSpectrum &numu, const NuMMParameters &pars) const
 Helper. Oscillates, decays or decoheres a true energy nubar spectrum.
void InverseOscillateTauBar (NuMatrixSpectrum &nutau, const NuMMParameters &pars) const
 Helper. Inverse oscillates a true energy taubar spectrum.

Public Attributes

TH1D * fSystPlusPQ [kNumSyst]
TH1D * fSystPlusNQ [kNumSyst]
TH1D * fSystMinusPQ [kNumSyst]
TH1D * fSystMinusNQ [kNumSyst]
TH1D * fhOnePQ [kNumSyst]
TH1D * fhOneNQ [kNumSyst]
TH1D * fhOne [kNumSyst]
bool fFitSysts

Private Member Functions

void CacheExtrapolation (const NuMMParameters &pars)
void ResetCache ()

Private Attributes

TFile * FD_MC
TFile * FD_NuTau
TFile * FD_data
TH1D * FD_dataCC
NuMatrixSpectrumPreCalcnuPrediction
NuMatrixSpectrumPreCalcbarPrediction
NuMatrixSpectrumPreCalcnuPredictionReco
NuMatrixSpectrumPreCalcbarPredictionReco
NuMatrixSpectrumPreCalcnuAppeared
NuMatrixSpectrumPreCalcbarAppeared
NuMatrixSpectrumPreCalcnuAppearedReco
NuMatrixSpectrumPreCalcbarAppearedReco
NuMatrixSpectrumPreCalcpotentialNuTaus
NuMatrixSpectrumPreCalcpotentialTauBars
NuMatrixSpectrumPreCalcpotentialNuTausReco
NuMatrixSpectrumPreCalcpotentialTauBarsReco
NuMatrixSpectrumPreCalcpotentialNuTaus1
NuMatrixSpectrumPreCalcpotentialNuTaus2
NuMatrixSpectrumPreCalcpotentialTauBars1
NuMatrixSpectrumPreCalcpotentialTauBars2
NuMatrixSpectrumPreCalcpotentialNuTausReco1
NuMatrixSpectrumPreCalcpotentialNuTausReco2
NuMatrixSpectrumPreCalcpotentialTauBarsReco1
NuMatrixSpectrumPreCalcpotentialTauBarsReco2
NuMatrixSpectrumPreCalcnuNCBackground
NuMatrixSpectrumPreCalcbarNCBackground
NuMatrixSpectrumPreCalcwrongSignNuMus
NuMatrixSpectrumPreCalcwrongSignNuBars
NuMatrixSpectrumPreCalcwrongSignNuMusReco
NuMatrixSpectrumPreCalcwrongSignNuBarsReco
NuMatrixSpectrumPreCalcwrongSignNuMus_NuMu
NuMatrixSpectrumPreCalcwrongSignNuBars_NuMu
NuMatrixSpectrumPreCalcwrongSignNuMusReco_NuMu
NuMatrixSpectrumPreCalcwrongSignNuBarsReco_NuMu
NuMatrixSpectrumPreCalcwrongSignNuMus_NuBar
NuMatrixSpectrumPreCalcwrongSignNuBars_NuBar
NuMatrixSpectrumPreCalcwrongSignNuMusReco_NuBar
NuMatrixSpectrumPreCalcwrongSignNuBarsReco_NuBar

Detailed Description

Definition at line 25 of file NuMMRunTransSME.h.


Member Enumeration Documentation

Enumerator:
kShwScaleTrans 
kTrkScaleTrans 
kDecayPipe 
kNumSyst 

Definition at line 39 of file NuMMRunTransSME.h.

00039              {
00040           kShwScaleTrans, // Energy dependent Transitions systematics
00041           kTrkScaleTrans,
00042 //          kNuBarXSecDIS,
00043           kDecayPipe,
00044           kNumSyst    // Number of energy depended systematics
00045           };


Constructor & Destructor Documentation

NuMMRunTransSME::NuMMRunTransSME ( NuMMHelperCPT helper,
NuMatrixSpectrum ndNuData,
NuMatrixSpectrum ndBarData,
double  pot 
)

Definition at line 23 of file NuMMRunTransSME.cxx.

00027   : NuMMRunNuBar()
00028 {
00029   fPredictNeutrinos = true;
00030   fPredictAntiNeutrinos = true;
00031 
00032   TString name = "hError";
00033   name += (counter++);
00034   TH1D hBlank(name, "Something Went Wrong", 10, 0, 30);
00035   
00036   fndNuData = ndNuData;
00037   fndBarData = ndBarData;
00038   fHelper = helper;
00039   
00040   ffdNuData = new NuMatrixSpectrum(hBlank, pot);
00041   ffdBarData = new NuMatrixSpectrum(hBlank, pot);
00042   
00043   
00044   // Blank out all the Precache variables
00045   PreCalcnuPrediction   = PreCalcbarPrediction      = 0;
00046   PreCalcnuAppeared     = PreCalcbarAppeared        = 0;
00047   PreCalcpotentialNuTaus= PreCalcpotentialTauBars   = 0;
00048   PreCalcpotentialNuTaus1= PreCalcpotentialTauBars1   = 0;
00049   PreCalcpotentialNuTaus2= PreCalcpotentialTauBars2   = 0;
00050   PreCalcnuNCBackground = PreCalcbarNCBackground    = 0;
00051   PreCalcwrongSignNuMus = PreCalcwrongSignNuBars    = 0;  
00052   PreCalcwrongSignNuMus_NuMu = PreCalcwrongSignNuBars_NuMu    = 0;
00053   PreCalcwrongSignNuMus_NuBar = PreCalcwrongSignNuBars_NuBar    = 0; 
00054 
00055  
00056   fFitSysts=false;
00057 }

NuMMRunTransSME::NuMMRunTransSME ( NuMMHelperCPT helper,
NuMatrixSpectrum ndNuData,
NuMatrixSpectrum ndBarData,
NuMatrixSpectrum fdNuData,
NuMatrixSpectrum fdBarData 
)

Definition at line 60 of file NuMMRunTransSME.cxx.

References NuMMRunNuBar::ffdBarData, NuMMRunNuBar::ffdNuData, fFitSysts, NuMMRunNuBar::fHelper, NuMMRunNuBar::fndBarData, NuMMRunNuBar::fndNuData, NuMMRunNuBar::fPredictAntiNeutrinos, NuMMRunNuBar::fPredictNeutrinos, PreCalcbarAppeared, PreCalcbarNCBackground, PreCalcbarPrediction, PreCalcnuAppeared, PreCalcnuNCBackground, PreCalcnuPrediction, PreCalcpotentialNuTaus, PreCalcpotentialNuTaus1, PreCalcpotentialNuTaus2, PreCalcpotentialTauBars, PreCalcpotentialTauBars1, PreCalcpotentialTauBars2, PreCalcwrongSignNuBars, PreCalcwrongSignNuBars_NuBar, PreCalcwrongSignNuBars_NuMu, PreCalcwrongSignNuMus, PreCalcwrongSignNuMus_NuBar, and PreCalcwrongSignNuMus_NuMu.

00065   : NuMMRunNuBar()
00066 { 
00067   fPredictNeutrinos = true;
00068   fPredictAntiNeutrinos = true;
00069 
00070   fndNuData = ndNuData;
00071   fndBarData = ndBarData;
00072   ffdNuData = fdNuData;
00073   ffdBarData = fdBarData;
00074   fHelper = helper;
00075    
00076   
00077   // Blank out all the Precache variables
00078   PreCalcnuPrediction   = PreCalcbarPrediction      = 0;
00079   PreCalcnuAppeared     = PreCalcbarAppeared        = 0;
00080   PreCalcpotentialNuTaus = PreCalcpotentialTauBars = 0;
00081   PreCalcpotentialNuTaus1 = PreCalcpotentialNuTaus2 = 0;
00082   PreCalcpotentialTauBars1   = PreCalcpotentialTauBars2 = 0;
00083   PreCalcnuNCBackground = PreCalcbarNCBackground    = 0;
00084   PreCalcwrongSignNuMus = PreCalcwrongSignNuBars    = 0;  
00085   PreCalcwrongSignNuMus_NuMu = PreCalcwrongSignNuBars_NuMu    = 0;
00086   PreCalcwrongSignNuMus_NuBar = PreCalcwrongSignNuBars_NuBar    = 0;
00087 
00088   fFitSysts=false;
00089 }

NuMMRunTransSME::~NuMMRunTransSME (  ) 

Definition at line 93 of file NuMMRunTransSME.cxx.

References ResetCache().

00094 {
00095   ResetCache();
00096 }


Member Function Documentation

void NuMMRunTransSME::CacheExtrapolation ( const NuMMParameters pars  )  [private]

Pre-calculate much of the extrapolation. Most of the matrix method extrapolation does not depend on the oscillation parameters. This includes lots of complicated matrix multiplications and histogram divisions. This function seperates out this code, and is called from the constructor (as everything it needs should exist with a validly called constructor)

Definition at line 146 of file NuMMRunTransSME.cxx.

References NuMMHelperCPT::BarBeamMatrix(), NuMMHelperCPT::BarBeamMatrixNuMuCCXSec(), NuMMHelperCPT::BarBeamMatrixTauCCXSec(), NuMMHelperCPT::BarXSecGraph(), NuMatrixSpectrum::Divide(), NuMMHelperCPT::DoTaus(), NuMatrixSpectrum::ExtrapolateNDToFD(), NuMMRunNuBar::fCached, NuMMHelperCPT::FDBarEfficiency(), NuMMHelperCPT::FDBarNCContamination(), NuMMHelperCPT::FDBarPurity(), NuMMHelperCPT::FDBarRecoVsTrue(), NuMMHelperCPT::FDNuEfficiency(), NuMMHelperCPT::FDNuNCContamination(), NuMMHelperCPT::FDNuPurity(), NuMMHelperCPT::FDNuRecoVsTrue(), NuMMHelperCPT::FDNuTauEfficiency(), NuMMHelperCPT::FDTauBarEfficiency(), NuMMHelperCPT::FDWrongSignBarEfficiency(), NuMMHelperCPT::FDWrongSignNuEfficiency(), NuMMRunNuBar::ffdBarData, NuMMRun::fFDFidMass, NuMMRunNuBar::ffdNuData, NuMMRunNuBar::fHelper, NuMMRunNuBar::fndBarData, NuMMRun::fNDFidMass, NuMMRunNuBar::fndNuData, NuMMRun::fQuietMode, NuMatrix::GetPOT(), Msg::kInfo, Msg::kWarning, MAXMSG, MSG, NuMatrixSpectrum::Multiply(), NuMMHelperCPT::NDBarEfficiency(), NuMMHelperCPT::NDBarPurity(), NuMMHelperCPT::NDBarRecoVsTrue(), NuMMHelperCPT::NDNuEfficiency(), NuMMHelperCPT::NDNuPurity(), NuMMHelperCPT::NDNuRecoVsTrue(), NuMMHelperCPT::NuBeamMatrix(), NuMMHelperCPT::NuBeamMatrixNuMuCCXSec(), NuMMHelperCPT::NuBeamMatrixTauCCXSec(), NuMMHelperCPT::NuXSecGraph(), PreCalcbarAppeared, PreCalcbarNCBackground, PreCalcbarPrediction, PreCalcnuAppeared, PreCalcnuNCBackground, PreCalcnuPrediction, PreCalcpotentialNuTaus, PreCalcpotentialNuTaus1, PreCalcpotentialNuTaus2, PreCalcpotentialTauBars, PreCalcpotentialTauBars1, PreCalcpotentialTauBars2, PreCalcwrongSignNuBars, PreCalcwrongSignNuBars_NuBar, PreCalcwrongSignNuBars_NuMu, PreCalcwrongSignNuMus, PreCalcwrongSignNuMus_NuBar, PreCalcwrongSignNuMus_NuMu, NuMatrixSpectrum::RecoToTrue(), ResetCache(), NuMatrix::ResetPOT(), NuMatrixSpectrum::SetName(), NuMatrixSpectrum::TrueToReco(), NuMMHelperCPT::XSecGraphNuTaus(), and NuMMHelperCPT::XSecGraphTauBars().

Referenced by MakeFDPred(), MakeFDPredNuBar(), MakeFDPredNuMu(), and TrueComponents().

00147 {
00148   // Only cache once.
00149   if (fCached) return;
00150   
00151   ResetCache();
00152   
00153   if (!fQuietMode) {
00154     MSG("NuMMRunTransitions",Msg::kInfo) << "Pre-calculating static extrapolation variables..." << endl;
00155   }
00156  
00157 //  cout<<"cache working"<<endl; 
00158   // Caching variables, to reduce the workload of each step of the fitter
00159   // Initialise these the way they would be before
00160   PreCalcnuPrediction = new NuMatrixSpectrum(*fndNuData);
00161   PreCalcbarPrediction = new NuMatrixSpectrum(*fndBarData);
00162   
00163   // Be lazy and just alias the pointers to references so can keep
00164   // the code the same
00165   NuMatrixSpectrum &nuPrediction = *PreCalcnuPrediction;
00166   NuMatrixSpectrum &barPrediction = *PreCalcbarPrediction;
00167 
00168   
00169   // Now we have prepared the extrapolation variables, calculate them
00170   // Get the neutrinos to the FD
00171   nuPrediction.Multiply(fHelper->NDNuPurity());
00172   nuPrediction.RecoToTrue(fHelper->NDNuRecoVsTrue());
00173   nuPrediction.Divide(fHelper->NDNuEfficiency());
00174   nuPrediction.Divide(fndNuData->GetPOT());
00175   nuPrediction.Divide(fNDFidMass);
00176   
00177   //Get the antineutrinos to the FD
00178   barPrediction.Multiply(fHelper->NDBarPurity());
00179   barPrediction.RecoToTrue(fHelper->NDBarRecoVsTrue());
00180   barPrediction.Divide(fHelper->NDBarEfficiency());
00181   barPrediction.Divide(fndBarData->GetPOT());
00182   barPrediction.Divide(fNDFidMass);
00183 
00184   // Create and link to the Tau components, before we include the
00185     // cross section in the predictions (tau Xsec is different)
00186     PreCalcpotentialNuTaus = new NuMatrixSpectrum(nuPrediction);
00187     PreCalcpotentialTauBars = new NuMatrixSpectrum(barPrediction);
00188     NuMatrixSpectrum &potentialNuTaus = *PreCalcpotentialNuTaus;
00189     NuMatrixSpectrum &potentialTauBars = *PreCalcpotentialTauBars;
00190   
00191   
00192   // nutaus can come from 'numu to nutau' or 'numubar to nutau'
00193   PreCalcpotentialNuTaus1 = new NuMatrixSpectrum(nuPrediction);  // 'numu to nutau' 
00194   PreCalcpotentialNuTaus2 = new NuMatrixSpectrum(barPrediction); // 'numubar to nutau'
00195   // nutaubars can come from 'numubar to nutaubar' or 'numu to nutaubar'
00196   PreCalcpotentialTauBars1 = new NuMatrixSpectrum(barPrediction); // 'numubar to nutaubar' 
00197   PreCalcpotentialTauBars2 = new NuMatrixSpectrum(nuPrediction); // 'numu to nutaubar'
00198   // Create and link to the Tau components, before we include the
00199   // cross section in the predictions (tau Xsec is different)
00200   NuMatrixSpectrum &potentialNuTaus1 = *PreCalcpotentialNuTaus1;
00201   NuMatrixSpectrum &potentialNuTaus2 = *PreCalcpotentialNuTaus2;
00202   
00203   NuMatrixSpectrum &potentialTauBars1 = *PreCalcpotentialTauBars1;
00204   NuMatrixSpectrum &potentialTauBars2 = *PreCalcpotentialTauBars2;
00205   
00206 
00207   nuPrediction.ExtrapolateNDToFD(fHelper->NuBeamMatrixNuMuCCXSec());
00208   nuPrediction.Multiply(fFDFidMass);
00209   nuPrediction.Multiply(ffdNuData->GetPOT());
00210   // Reset the spectrum to the new PoT
00211   nuPrediction.ResetPOT(ffdNuData->GetPOT());
00212 
00213   barPrediction.ExtrapolateNDToFD(fHelper->BarBeamMatrixNuMuCCXSec());
00214   barPrediction.Multiply(fFDFidMass);
00215   barPrediction.Multiply(ffdBarData->GetPOT());
00216   // Reset the spectrum to the new PoT
00217   barPrediction.ResetPOT(ffdBarData->GetPOT());
00218   
00219 
00220   if (fHelper->DoTaus()) {
00221         // Get the Neutrinos to the FD
00222         if (fHelper->NuBeamMatrixTauCCXSec()) {
00223             potentialNuTaus.ExtrapolateNDToFD(fHelper->NuBeamMatrixTauCCXSec());
00224             potentialTauBars.ExtrapolateNDToFD(fHelper->BarBeamMatrixTauCCXSec());
00225 
00226             potentialNuTaus1.ExtrapolateNDToFD(fHelper->NuBeamMatrixTauCCXSec());
00227             potentialNuTaus2.ExtrapolateNDToFD(fHelper->NuBeamMatrixTauCCXSec());
00228 
00229             potentialTauBars1.ExtrapolateNDToFD(fHelper->BarBeamMatrixTauCCXSec());
00230             potentialTauBars2.ExtrapolateNDToFD(fHelper->BarBeamMatrixTauCCXSec());
00231         }
00232         else {
00233             MAXMSG("NuMMRunTransitions",Msg::kWarning,10) << "Doing tau extrapolation the old way." << endl;
00234             potentialNuTaus.Divide(fHelper->NuXSecGraph());
00235             potentialNuTaus.ExtrapolateNDToFD(fHelper->NuBeamMatrix());
00236             potentialNuTaus.Multiply(fHelper->XSecGraphNuTaus());
00237             potentialTauBars.Divide(fHelper->BarXSecGraph());
00238             potentialTauBars.ExtrapolateNDToFD(fHelper->BarBeamMatrix());
00239             potentialTauBars.Multiply(fHelper->XSecGraphTauBars());
00240 
00241             potentialNuTaus1.Divide(fHelper->NuXSecGraph());
00242             potentialNuTaus1.ExtrapolateNDToFD(fHelper->NuBeamMatrix());
00243             potentialNuTaus1.Multiply(fHelper->XSecGraphNuTaus());
00244 
00245             potentialNuTaus2.Divide(fHelper->NuXSecGraph());
00246             potentialNuTaus2.ExtrapolateNDToFD(fHelper->NuBeamMatrix());
00247             potentialNuTaus2.Multiply(fHelper->XSecGraphNuTaus());
00248 
00249             potentialTauBars1.Divide(fHelper->BarXSecGraph());
00250             potentialTauBars1.ExtrapolateNDToFD(fHelper->BarBeamMatrix());
00251             potentialTauBars1.Multiply(fHelper->XSecGraphTauBars());
00252 
00253             potentialTauBars2.Divide(fHelper->BarXSecGraph());
00254             potentialTauBars2.ExtrapolateNDToFD(fHelper->BarBeamMatrix());
00255             potentialTauBars2.Multiply(fHelper->XSecGraphTauBars());
00256         }
00257 
00258         potentialNuTaus.Multiply(fFDFidMass);
00259         potentialNuTaus.Multiply(ffdNuData->GetPOT());
00260         // Reset the spectrum to the new PoT
00261         potentialNuTaus.ResetPOT(ffdNuData->GetPOT());
00262         potentialNuTaus.Multiply(fHelper->FDNuTauEfficiency());
00263 
00264         potentialNuTaus1.Multiply(fFDFidMass);
00265         potentialNuTaus1.Multiply(ffdNuData->GetPOT());
00266         // Reset the spectrum to the new PoT
00267         potentialNuTaus1.ResetPOT(ffdNuData->GetPOT());
00268         potentialNuTaus1.Multiply(fHelper->FDNuTauEfficiency());
00269 
00270         potentialNuTaus2.Multiply(fFDFidMass);
00271         potentialNuTaus2.Multiply(ffdNuData->GetPOT());
00272         // Reset the spectrum to the new PoT
00273         potentialNuTaus2.ResetPOT(ffdNuData->GetPOT());
00274         potentialNuTaus2.Multiply(fHelper->FDNuTauEfficiency());
00275 
00276         // Get the Antineutrinos to the FD
00277         potentialTauBars.Multiply(fFDFidMass);
00278         potentialTauBars.Multiply(ffdBarData->GetPOT());
00279         // Reset the spectrum to the new PoT
00280         potentialTauBars.ResetPOT(ffdBarData->GetPOT());
00281         potentialTauBars.Multiply(fHelper->FDTauBarEfficiency());
00282 
00283         // Get the Antineutrinos to the FD .... mine
00284         potentialTauBars1.Multiply(fFDFidMass);
00285         potentialTauBars1.Multiply(ffdBarData->GetPOT());
00286         // Reset the spectrum to the new PoT
00287         potentialTauBars1.ResetPOT(ffdBarData->GetPOT());
00288         potentialTauBars1.Multiply(fHelper->FDTauBarEfficiency());
00289 
00290         // Get the Antineutrinos to the FD
00291         potentialTauBars2.Multiply(fFDFidMass);
00292         potentialTauBars2.Multiply(ffdBarData->GetPOT());
00293         // Reset the spectrum to the new PoT
00294         potentialTauBars2.ResetPOT(ffdBarData->GetPOT());
00295         potentialTauBars2.Multiply(fHelper->FDTauBarEfficiency());
00296           
00297     }
00298     else {
00299         PreCalcpotentialNuTaus->Multiply(0.0);
00300         PreCalcpotentialTauBars->Multiply(0.0);
00301 
00302         PreCalcpotentialNuTaus1->Multiply(0.0);
00303         PreCalcpotentialNuTaus2->Multiply(0.0);
00304 
00305         PreCalcpotentialTauBars1->Multiply(0.0);
00306         PreCalcpotentialTauBars2->Multiply(0.0);
00307 
00308     }
00309 
00310 
00311 
00312   // Create and link to the appeared components, before we include the
00313   // cross section in the predictions (Xsec is different)
00314   PreCalcnuAppeared = new NuMatrixSpectrum(barPrediction);
00315   PreCalcbarAppeared = new NuMatrixSpectrum(nuPrediction);
00316   NuMatrixSpectrum &appearedNuMus = *PreCalcnuAppeared;
00317   NuMatrixSpectrum &appearedNuBars = *PreCalcbarAppeared;
00318   appearedNuMus.SetName("Appeared_NuMus");
00319   appearedNuBars.SetName("Appeared_NuMuBars");
00320   
00321   // Apply the appeared cross-sectinos, efficiencies, true to reco
00322   appearedNuMus.Divide(fHelper->BarXSecGraph());
00323   appearedNuMus.Multiply(fHelper->NuXSecGraph());
00324   appearedNuMus.Multiply(fHelper->FDNuEfficiency());
00325  
00326   appearedNuBars.Divide(fHelper->NuXSecGraph());
00327   appearedNuBars.Multiply(fHelper->BarXSecGraph());
00328   appearedNuBars.Multiply(fHelper->FDBarEfficiency());
00329   
00330   // Now build the NC Spectrum, after the cross section calculation!
00331   PreCalcnuNCBackground = new NuMatrixSpectrum(*PreCalcnuPrediction);
00332   PreCalcbarNCBackground = new NuMatrixSpectrum(*PreCalcbarPrediction);
00333   NuMatrixSpectrum &nuNCBackground = *PreCalcnuNCBackground;
00334   NuMatrixSpectrum &barNCBackground = *PreCalcbarNCBackground;  
00335 
00336   //Get the neutrino NC background
00337   nuNCBackground.Multiply(fHelper->FDNuEfficiency());
00338   nuNCBackground.TrueToReco(fHelper->FDNuRecoVsTrue());
00339   nuNCBackground.Divide(fHelper->FDNuPurity());
00340   nuNCBackground.Multiply(fHelper->FDNuNCContamination());
00341   
00342   //Get the antineutrino NC background
00343   barNCBackground.Multiply(fHelper->FDBarEfficiency());
00344   barNCBackground.TrueToReco(fHelper->FDBarRecoVsTrue());
00345   barNCBackground.Divide(fHelper->FDBarPurity());
00346   barNCBackground.Multiply(fHelper->FDBarNCContamination());
00347   
00348   // Now build the wrong sign spectrum
00349   PreCalcwrongSignNuMus = new NuMatrixSpectrum(*PreCalcnuPrediction);
00350   PreCalcwrongSignNuBars = new NuMatrixSpectrum(*PreCalcbarPrediction);
00351   NuMatrixSpectrum &wrongSignNuMus = *PreCalcwrongSignNuMus;
00352   NuMatrixSpectrum &wrongSignNuBars = *PreCalcwrongSignNuBars;
00353 
00354   PreCalcwrongSignNuMus_NuMu = new NuMatrixSpectrum(*PreCalcnuPrediction);
00355   PreCalcwrongSignNuMus_NuBar = new NuMatrixSpectrum(*PreCalcbarPrediction);
00356   NuMatrixSpectrum &wrongSignNuMus_NuMu = *PreCalcwrongSignNuMus_NuMu;
00357   NuMatrixSpectrum &wrongSignNuMus_NuBar = *PreCalcwrongSignNuMus_NuBar;
00358 
00359   PreCalcwrongSignNuBars_NuMu = new NuMatrixSpectrum(*PreCalcnuPrediction);
00360   PreCalcwrongSignNuBars_NuBar = new NuMatrixSpectrum(*PreCalcbarPrediction);
00361   NuMatrixSpectrum &wrongSignNuBars_NuMu = *PreCalcwrongSignNuBars_NuMu; 
00362   NuMatrixSpectrum &wrongSignNuBars_NuBar = *PreCalcwrongSignNuBars_NuBar; 
00363 
00364   // Do the Neutrino background for the antineutrino spectrum
00365   wrongSignNuMus.Multiply(fHelper->FDWrongSignNuEfficiency());
00366   wrongSignNuMus_NuMu.Multiply(fHelper->FDWrongSignNuEfficiency());  // wrong sign nus coming from nus in ND
00367   wrongSignNuMus_NuBar.Multiply(fHelper->FDWrongSignNuEfficiency());  // wrong sign nus coming from nubars in ND
00368   
00369   // Do the Antineutrino background for the neutrino spectrum
00370   wrongSignNuBars.Multiply(fHelper->FDWrongSignBarEfficiency());
00371   wrongSignNuBars_NuMu.Multiply(fHelper->FDWrongSignBarEfficiency()); // wrong sign nubars coming from nus in ND
00372   wrongSignNuBars_NuBar.Multiply(fHelper->FDWrongSignBarEfficiency());  // wrong sign nubars coming from nubars in ND
00373 
00374 
00375   // Now apply oscillations, efficeicnies to the main spectrum
00376   nuPrediction.Multiply(fHelper->FDNuEfficiency());
00377   barPrediction.Multiply(fHelper->FDBarEfficiency());
00378   
00379 
00380   // If we have precalculated the wrong sign oscillations,
00381   // then we need to recalculate them next time we are given
00382   // the oscilaltion parameters. Reset the flags used to
00383   // control this
00384   //fDoneWsNeutrinos = fDoneWsAntineutrinos = false;
00385 
00386   nuPrediction.SetName("NuMus");
00387   barPrediction.SetName("NuMuBars");
00388 
00389   
00390   if (fHelper->DoTaus()) {
00391     potentialNuTaus.SetName("NuTaus");
00392     potentialTauBars.SetName("NuTauBars");
00393        
00394     potentialNuTaus1.SetName("NuTaus1");
00395     potentialNuTaus2.SetName("NuTaus2");
00396     potentialTauBars1.SetName("NuTauBars1");
00397     potentialTauBars2.SetName("NuTauBars2");
00398 
00399     }
00400 
00401   wrongSignNuMus.SetName("WrongSign_NuMus");
00402   wrongSignNuBars.SetName("WrongSign_NuMuBars");
00403 
00404   wrongSignNuMus_NuMu.SetName("WrongSignNuMus_FromNuMu");
00405   wrongSignNuMus_NuBar.SetName("WrongSignNuMus_FromNuBar");
00406   wrongSignNuBars_NuMu.SetName("WrongSignNuMuBars_FromNuMu");
00407   wrongSignNuBars_NuBar.SetName("WrongSignNuMuBars_FromNuBar");
00408 
00409   
00410   nuNCBackground.SetName("NC_NuMus");
00411   barNCBackground.SetName("NC_NuMuBars");
00412   
00413   // Now perform all the Reco To True's for the Reco prediction
00414   fCached = true;
00415     
00416   // Disable this output for now - this process should be fast anyway!
00417   if (!fQuietMode)
00418     cout << "Extrapolation variables cached..." << endl;
00419 }

Double_t NuMMRunTransSME::ComparePredWithData ( const NuMMParameters pars  )  [virtual]

Implements NuMMRun.

Definition at line 1212 of file NuMMRunTransSME.cxx.

References NuMMRunNuBar::ffdBarData, NuMMRunNuBar::ffdNuData, NuMMRunNuBar::fPredictAntiNeutrinos, NuMMRunNuBar::fPredictNeutrinos, MakeFDPred(), NuMatrixSpectrum::Spectrum(), and NuMMRun::StatsLikelihood().

01213 {
01214   // Generate a prediction
01215   const pair<NuMatrixSpectrum,NuMatrixSpectrum> predictions =
01216     this->MakeFDPred(pars);
01217   Double_t like = 0;
01218   // Now, only count a likelihood for components we want - this is just
01219   // a constant offset, but should help speed things up a little.
01220   if (fPredictNeutrinos) {
01221     like += this->StatsLikelihood(predictions.first.Spectrum(),
01222                                 ffdNuData->Spectrum());
01223     //  cout <<"Stats likelihood neutrinos:" << like; 
01224   }
01225   if (fPredictAntiNeutrinos) {
01226     like += this->StatsLikelihood(predictions.second.Spectrum(),
01227                                 ffdBarData->Spectrum());
01228     //  cout <<"Stats likelihood antineutrinos:" << like;
01229   }
01230   return like;
01231 }

virtual void NuMMRunTransSME::FitSysts ( Bool_t  fitSysts  )  [inline, virtual]

Definition at line 58 of file NuMMRunTransSME.h.

References fFitSysts.

00058 {fFitSysts = fitSysts;}

void NuMMRunTransSME::InverseOscillateNuTau ( NuMatrixSpectrum nutau,
const NuMMParameters pars 
) const

Helper. Inverse oscillates a true energy nutau spectrum.

Does nothing if the model is decay or decoherence

Definition at line 1256 of file NuMMRunTransSME.cxx.

References NuMMParameters::Alpha(), NuMatrix1D::DecayMuToTau(), NuMMParameters::Dm2(), NuMMParameters::Epsilon(), NuMMRun::fDisappearanceModel, NuMatrix1D::InverseDecohere(), NuMatrixSpectrum::InverseOscillate(), NuMatrixSpectrum::InverseOscillateNSI(), NuMMParameters::Mu2(), and NuMMParameters::Sn2().

01258 {
01259   if(fDisappearanceModel == 0) {
01260     nutau.InverseOscillate(pars.Dm2(), pars.Sn2());
01261   }
01262   if(fDisappearanceModel == 1)
01263     nutau.DecayMuToTau(pars.Dm2(), pars.Sn2(), pars.Alpha());
01264   if(fDisappearanceModel == 2)
01265     nutau.InverseDecohere(pars.Dm2(), pars.Sn2(), pars.Mu2());
01266   if(fDisappearanceModel == 3)
01267     nutau.InverseOscillateNSI(pars.Dm2(), pars.Sn2(), pars.Epsilon(), 1.0);
01268 }

void NuMMRunTransSME::InverseOscillateTauBar ( NuMatrixSpectrum nutau,
const NuMMParameters pars 
) const

Helper. Inverse oscillates a true energy taubar spectrum.

Does nothing if the model is decay or decoherence

Definition at line 1295 of file NuMMRunTransSME.cxx.

References NuMMParameters::Alpha(), NuMatrix1D::DecayMuToTau(), NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), NuMMParameters::Epsilon(), NuMMRun::fDisappearanceModel, NuMatrix1D::InverseDecohere(), NuMatrixSpectrum::InverseOscillate(), NuMatrixSpectrum::InverseOscillateNSI(), NuMMParameters::Mu2(), NuMMParameters::Sn2(), and NuMMParameters::Sn2Bar().

01297 {
01298   if(fDisappearanceModel == 0) {
01299     nutau.InverseOscillate(pars.Dm2Bar(), pars.Sn2Bar());
01300   }
01301   if(fDisappearanceModel == 1)
01302     nutau.DecayMuToTau(pars.Dm2Bar(), pars.Sn2Bar(), pars.Alpha());
01303   if(fDisappearanceModel == 2)
01304     nutau.InverseDecohere(pars.Dm2Bar(), pars.Sn2Bar(), pars.Mu2());
01305   if(fDisappearanceModel == 3)
01306     nutau.InverseOscillateNSI(pars.Dm2(), pars.Sn2(), pars.Epsilon(), -1.0);
01307 }

pair< NuMatrixSpectrum, NuMatrixSpectrum > NuMMRunTransSME::MakeFDPred ( const NuMMParameters pars  )  [virtual]

Return a pair (NQ,PQ) of NuMatrixSpectra of the far detector prediction with oscillation parameters described in pars. Pure virtual.

Implements NuMMRunNuBar.

Definition at line 745 of file NuMMRunTransSME.cxx.

References NuMatrix1D::Add(), CacheExtrapolation(), NuMMParameters::CMuMu(), NuMMParameters::CTauTau(), NuMMParameters::DecayPipe(), NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), NuMMHelperCPT::DoTaus(), NuMMHelperCPT::FDBarRecoVsTrue(), NuMMHelperCPT::FDNuRecoVsTrue(), NuMMHelperCPT::FDNuTauRecoVsTrue(), NuMMHelperCPT::FDTauBarRecoVsTrue(), NuMMHelperCPT::FDWrongSignBarRecoVsTrue(), NuMMHelperCPT::FDWrongSignNuRecoVsTrue(), fFitSysts, NuMMRunNuBar::fHelper, fhOne, NuMMRunNuBar::fPredictAntiNeutrinos, NuMMRunNuBar::fPredictNeutrinos, NuMMRun::fQuietMode, fSystMinusNQ, fSystMinusPQ, fSystPlusNQ, fSystPlusPQ, NuMMParameters::GetParameterByIndex(), NuMMParameters::GMu(), NuMMParameters::GTau(), NuMMParameters::kDecayPipe, kDecayPipe, NuMMParameters::kNormTrans, kNumSyst, kShwScaleTrans, NuMMParameters::kShwScaleTrans, kTrkScaleTrans, NuMMParameters::kTrkScaleTrans, NuMatrixSpectrum::Multiply(), NuMMParameters::NormalisationTrans(), NuMMParameters::NuBarXSecDIS(), NuMatrixSpectrum::OscillateMu2Mu(), NuMatrixSpectrum::OscillateMu2MuBar(), NuMatrixSpectrum::OscillateMu2Tau(), NuMatrixSpectrum::OscillateMu2TauBar(), PreCalcbarAppeared, PreCalcbarNCBackground, PreCalcbarPrediction, PreCalcnuAppeared, PreCalcnuNCBackground, PreCalcnuPrediction, PreCalcpotentialNuTaus, PreCalcpotentialNuTaus1, PreCalcpotentialNuTaus2, PreCalcpotentialTauBars, PreCalcpotentialTauBars1, PreCalcpotentialTauBars2, PreCalcwrongSignNuBars, PreCalcwrongSignNuBars_NuBar, PreCalcwrongSignNuBars_NuMu, PreCalcwrongSignNuMus, PreCalcwrongSignNuMus_NuBar, PreCalcwrongSignNuMus_NuMu, NuMatrixSpectrum::RemoveErrors(), NuMatrix1D::Scale(), NuMMParameters::ShwEnScaleTrans(), NuMMParameters::Sn2(), NuMMParameters::Sn2Bar(), NuMatrixSpectrum::Spectrum(), NuMMParameters::TrkEnScaleTrans(), and NuMatrixSpectrum::TrueToReco().

Referenced by ComparePredWithData().

00746 {
00747 /*  if (!fQuietMode) {
00748     cout << "sn2: " << pars.Sn2() << "; dm2: " << pars.Dm2() 
00749      << "; sn2bar: " << pars.Sn2Bar() << "; dm2bar: " << pars.Dm2Bar()
00750      << "; alpha: " << pars.TransitionProb()
00751      << endl;
00752   }
00753   */
00754   CacheExtrapolation(pars);
00755   if (!fQuietMode){
00756   cout<<",dm2 = "<<pars.Dm2()<<", sn2 = "<<pars.Sn2()<<endl;
00757   cout<<",dm2bar= "<<pars.Dm2Bar()<<", sn2bar = "<<pars.Sn2Bar()<<", bmu = "<<pars.GMu()<<endl;
00758   cout<<",normTrans = "<<pars.NormalisationTrans()<<",trkEn = "<<pars.TrkEnScaleTrans()<<", shwEn = "<<pars.ShwEnScaleTrans()<< ", decayP = "<<pars.DecayPipe()<< ", nubar = "<<pars.NuBarXSecDIS()<<endl;
00759   }
00761   // For all neutrinos
00762   
00763   // Make copies of the prepared extrapolation variables, so that
00764   // we can alter them locally with the oscillation parameters
00765 
00766   NuMatrixSpectrum nuPrediction(*PreCalcnuPrediction);
00767   NuMatrixSpectrum barPrediction(*PreCalcbarPrediction);
00768  
00769   NuMatrixSpectrum potentialNuTaus(*PreCalcpotentialNuTaus);
00770   NuMatrixSpectrum potentialTauBars(*PreCalcpotentialTauBars);
00771  
00772   NuMatrixSpectrum potentialNuTaus1(*PreCalcpotentialNuTaus1);
00773   NuMatrixSpectrum potentialNuTaus2(*PreCalcpotentialNuTaus2);
00774   
00775   NuMatrixSpectrum potentialTauBars1(*PreCalcpotentialTauBars1);
00776   NuMatrixSpectrum potentialTauBars2(*PreCalcpotentialTauBars2);
00777 
00778   NuMatrixSpectrum appearedNuMus(*PreCalcnuAppeared);
00779   NuMatrixSpectrum appearedNuBars(*PreCalcbarAppeared);
00780   
00781   NuMatrixSpectrum wrongSignNuMus(*PreCalcwrongSignNuMus);
00782   NuMatrixSpectrum wrongSignNuBars(*PreCalcwrongSignNuBars);
00783   
00784   NuMatrixSpectrum wrongSignNuMus_NuMu(*PreCalcwrongSignNuMus_NuMu);
00785   NuMatrixSpectrum wrongSignNuMus_NuBar(*PreCalcwrongSignNuMus_NuBar);
00786 
00787   NuMatrixSpectrum wrongSignNuBars_NuMu(*PreCalcwrongSignNuBars_NuMu);
00788   NuMatrixSpectrum wrongSignNuBars_NuBar(*PreCalcwrongSignNuBars_NuBar);
00789 
00790   NuMatrixSpectrum nuNCBackground(*PreCalcnuNCBackground);
00791   NuMatrixSpectrum barNCBackground(*PreCalcbarNCBackground);  
00792   
00794   // Processes for neutrinos
00795   
00796   if (fPredictNeutrinos) {  
00797 
00798     if (fHelper->DoTaus()){
00799       potentialNuTaus1.OscillateMu2Tau(pars.Dm2(),pars.Sn2(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00800       potentialNuTaus2.OscillateMu2TauBar(pars.Dm2(),pars.Sn2(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00801       potentialNuTaus1.TrueToReco(fHelper->FDNuTauRecoVsTrue());
00802       potentialNuTaus2.TrueToReco(fHelper->FDNuTauRecoVsTrue());
00803     }
00804 
00805     // wrong sign: 1) numu -> nubar
00806     //             2) nubar -> nubar
00807  
00808     wrongSignNuBars_NuMu.OscillateMu2MuBar(pars.Dm2(),pars.Sn2(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau()); 
00809     wrongSignNuBars_NuMu.TrueToReco(fHelper->FDWrongSignBarRecoVsTrue());
00810 
00811     wrongSignNuBars_NuBar.OscillateMu2Mu(pars.Dm2(),pars.Sn2(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau()); 
00812     wrongSignNuBars_NuBar.TrueToReco(fHelper->FDWrongSignBarRecoVsTrue());
00813 
00814     appearedNuMus.OscillateMu2MuBar(pars.Dm2(),pars.Sn2(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00815     nuPrediction.OscillateMu2Mu(pars.Dm2(),pars.Sn2(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00816     appearedNuMus.TrueToReco(fHelper->FDNuRecoVsTrue());
00817     nuPrediction.TrueToReco(fHelper->FDNuRecoVsTrue());
00818  
00819     if (fFitSysts){
00820     //systematics... joao's code
00821     // Apply normalization systematics
00822            
00823            nuPrediction.Scale(1+(0.016*pars.GetParameterByIndex(NuMMParameters::kNormTrans)));
00824            // Apply energy dependent systematics
00825            for(int kPar=0; kPar<kNumSyst; kPar++){
00826 
00827               double parvalue;
00828               // Only apply shifts to spectra that do change
00829               if(kPar==kShwScaleTrans) parvalue = pars.GetParameterByIndex(NuMMParameters::kShwScaleTrans);
00830               else if(kPar==kDecayPipe)  parvalue = pars.GetParameterByIndex(NuMMParameters::kDecayPipe);
00831               else if(kPar==kTrkScaleTrans)  parvalue = pars.GetParameterByIndex(NuMMParameters::kTrkScaleTrans);
00832               //else if(kPar==kNuBarXSecDIS)  parvalue = pars.GetParameterByIndex(NuMMParameters::kNuBarXSecDIS);
00833 
00834               // Get the shift value in sigmas
00835               // double parvalue = pars.GetParameterByIndex(fParMap.at(kPar));
00836 
00837               // Skip parameters that do nothing
00838               if(parvalue==0) continue;
00839 
00840               // Get the shift histogram by sign
00841               TH1D *ScaledSystNQ;
00842               //if (kPar==kShwScaleTrans||kPar==kNuBarXSecDIS||kPar==kTrkScaleTrans){  
00843               if (kPar==kShwScaleTrans||kPar==kDecayPipe||kPar==kTrkScaleTrans){  
00844                  if(parvalue<0)
00845                     {
00846                     ScaledSystNQ = (TH1D*)fSystMinusNQ[kPar]->Clone();
00847                     }
00848                     else
00849                     {
00850                     ScaledSystNQ = (TH1D*)fSystPlusNQ[kPar]->Clone();
00851                     }
00852                  }
00853 
00854               ScaledSystNQ->Scale(TMath::Abs(parvalue)); // Scale to shift value
00855               ScaledSystNQ->Add(fhOne[kPar]);                  // Add one
00856               nuPrediction.Multiply(ScaledSystNQ);         // Apply systematic to NQ spectra
00857               
00858           }//loop number of parameters
00859        }//if fitsysts
00860 
00861  
00862     //Add everything together
00863     nuPrediction.Add(appearedNuMus);
00864     nuPrediction.Add(wrongSignNuBars_NuMu);
00865     nuPrediction.Add(wrongSignNuBars_NuBar);
00866     nuPrediction.Add(nuNCBackground);
00867 
00868     if (fHelper->DoTaus()){
00869       nuPrediction.Add(potentialNuTaus1);    
00870       nuPrediction.Add(potentialNuTaus2);  
00871       }  
00872     } else {
00873       // Scale the histogram to zero, so as not to confuse people
00874       // if they look at the spectrum with neutrinos turned off
00875       nuPrediction.Spectrum()->Scale(0);
00876     }
00877   
00879   // Processes for Antineutrinos
00880   if (fPredictAntiNeutrinos)
00881     {
00882       if (fHelper->DoTaus()) {
00883       // Oscillate the taubars
00884       potentialTauBars1.OscillateMu2Tau(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00885       potentialTauBars2.OscillateMu2TauBar(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau()); 
00886       potentialTauBars1.TrueToReco(fHelper->FDTauBarRecoVsTrue());
00887       potentialTauBars2.TrueToReco(fHelper->FDTauBarRecoVsTrue()); 
00888       
00889       }
00890 
00891     //Oscillate wrong sign with the nubar parameters (not the nu parameters)
00892     //wrong sign: 1) numu -> numu
00893     //            2) nubar -> numu
00894 
00895     wrongSignNuMus_NuMu.OscillateMu2Mu(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00896     wrongSignNuMus_NuMu.TrueToReco(fHelper->FDWrongSignNuRecoVsTrue());
00897 
00898     wrongSignNuMus_NuBar.OscillateMu2MuBar(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00899     wrongSignNuMus_NuBar.TrueToReco(fHelper->FDWrongSignNuRecoVsTrue());
00900 
00901 
00902     //Oscillate appeared nubars and nubar prediction
00903     appearedNuBars.OscillateMu2MuBar(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());    
00904     barPrediction.OscillateMu2Mu(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00905     // reco to true
00906     appearedNuBars.TrueToReco(fHelper->FDBarRecoVsTrue());
00907     barPrediction.TrueToReco(fHelper->FDBarRecoVsTrue());
00908 
00909 
00910     if (fFitSysts){
00911     // Apply normalization systematics
00912            barPrediction.Scale(1+(0.04*pars.GetParameterByIndex(NuMMParameters::kNormTrans)));
00913  
00914            // Apply energy dependent systematics
00915            for(int kPar=0; kPar<kNumSyst; kPar++){
00916               double parvalue;
00917               // Only apply shifts to spectra that do change
00918 
00919               if(kPar==kShwScaleTrans) parvalue = pars.GetParameterByIndex(NuMMParameters::kShwScaleTrans);
00920               else if(kPar==kTrkScaleTrans) parvalue = pars.GetParameterByIndex(NuMMParameters::kTrkScaleTrans);
00921               else if(kPar==kDecayPipe)  parvalue = pars.GetParameterByIndex(NuMMParameters::kDecayPipe);
00922               //else if(kPar==kNuBarXSecDIS)  parvalue = pars.GetParameterByIndex(NuMMParameters::kNuBarXSecDIS);
00923 
00924  
00925              // Skip parameters that do nothing
00926              if(parvalue==0) continue;
00927  
00928              // Get the shift histogram by sign
00929              TH1D *ScaledSystPQ;
00930              //if (kPar==kShwScaleTrans||kPar==kNuBarXSecDIS||kPar==kTrkScaleTrans){
00931              if (kPar==kShwScaleTrans||kPar==kDecayPipe||kPar==kTrkScaleTrans){
00932                 if(parvalue<0)
00933                    {ScaledSystPQ = (TH1D*)fSystMinusPQ[kPar]->Clone();
00934                    }
00935                    else
00936                    {ScaledSystPQ = (TH1D*)fSystPlusPQ[kPar]->Clone();
00937                    }
00938                 }
00939  
00940               ScaledSystPQ->Scale(TMath::Abs(parvalue)); // Scale to shift value
00941               ScaledSystPQ->Add(fhOne[kPar]);                  // Add one
00942               barPrediction.Multiply(ScaledSystPQ);         // Apply systematic to PQ
00943               
00944           }//loop over systematics
00945        }//if fitsysts
00946 
00947 //.................
00948 
00949 
00950     // Add in backgrounds
00951     barPrediction.Add(appearedNuBars);
00952     barPrediction.Add(wrongSignNuMus_NuMu);
00953     barPrediction.Add(wrongSignNuMus_NuBar);
00954     barPrediction.Add(barNCBackground);
00955 
00956    if (fHelper->DoTaus()){
00957        barPrediction.Add(potentialTauBars1);
00958        barPrediction.Add(potentialTauBars2);
00959        }
00960   } //if predict antineutrinos
00961 
00962   else {
00963     // Scale the histogram to zero, so as not to confuse people
00964     barPrediction.Spectrum()->Scale(0);
00965   }
00966   nuPrediction.RemoveErrors(); 
00967   barPrediction.RemoveErrors();
00968   pair<NuMatrixSpectrum,NuMatrixSpectrum> predictions(nuPrediction,barPrediction);
00969   return predictions;
00970 /*
00971   // hacked for the time being!
00972   NuMatrixSpectrum wrongSignNuBarsTotal = wrongSignNuBars_NuMu;
00973   wrongSignNuBarsTotal.Add(wrongSignNuBars_NuBar);
00974   wrongSignNuBarsTotal.SetName("WS_NuMus");
00975   NuMatrixSpectrum wrongSignNuMusTotal = wrongSignNuMus_NuMu;
00976   wrongSignNuMusTotal.Add(wrongSignNuMus_NuBar);
00977   wrongSignNuMusTotal.SetName("WS_NuMuBars");
00978   pair<NuMatrixSpectrum,NuMatrixSpectrum>
00979     predictionsWS(wrongSignNuBarsTotal,wrongSignNuMusTotal);
00980   return predictionsWS;
00981 */
00982   //hacked for NC
00983   //Return the NC predictions, as a pair
00984   //pair<NuMatrixSpectrum,NuMatrixSpectrum>
00985   //  predictionsNC(nuNCBackground,barNCBackground);
00986   //return predictionsNC;
00987 
00988 }

NuMatrixSpectrum * NuMMRunTransSME::MakeFDPredNuBar ( const NuMMParameters pars  )  [virtual]

Get just the NuMuBar component of the prediction with oscillations parameters pars.

Reimplemented from NuMMRunNuBar.

Definition at line 1041 of file NuMMRunTransSME.cxx.

References NuMatrix1D::Add(), CacheExtrapolation(), NuMMParameters::CMuMu(), NuMMParameters::CTauTau(), NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), NuMMHelperCPT::DoTaus(), NuMMRunNuBar::fHelper, NuMMRunNuBar::fPredictAntiNeutrinos, NuMMParameters::GMu(), NuMMParameters::GTau(), NuMatrixSpectrum::OscillateMu2Mu(), NuMatrixSpectrum::OscillateMu2MuBar(), NuMatrixSpectrum::OscillateMu2Tau(), NuMatrixSpectrum::OscillateMu2TauBar(), PreCalcbarAppearedReco, PreCalcbarNCBackground, PreCalcbarPredictionReco, PreCalcpotentialTauBarsReco1, PreCalcpotentialTauBarsReco2, PreCalcwrongSignNuMusReco_NuBar, PreCalcwrongSignNuMusReco_NuMu, NuMMParameters::Sn2Bar(), and NuMatrixSpectrum::Spectrum().

01042 {
01043 /*  if (!fQuietMode) {
01044     cout << "sn2: " << pars.Sn2() << "; dm2: " << pars.Dm2() 
01045     << "; sn2bar: " << pars.Sn2Bar() << "; dm2bar: " << pars.Dm2Bar()
01046     << "; alpha: " << pars.TransitionProb()
01047     << endl;
01048   }
01049   */
01050   CacheExtrapolation(pars);
01051   cout<<"MakeFDPredNuBar.."<<endl;
01052   cout<<"bmu = "<<pars.GMu()<<", btau =  "<<pars.GTau()<<", cmumu = "<<pars.CMuMu()<<", ctautau = "<<pars.CTauTau()<<", dm2 ==="<<pars.Dm2()<<endl;
01053   //pars.PrintStatus();
01054   
01056   // For all neutrinos
01057   
01058   // Make copies of the prepared extrapolation variables, so that
01059   // we can alter them locally with the oscillation parameters
01060   NuMatrixSpectrum *barPrediction = new NuMatrixSpectrum(*PreCalcbarPredictionReco);
01061   
01062   NuMatrixSpectrum potentialTauBars1(*PreCalcpotentialTauBarsReco1);
01063   NuMatrixSpectrum potentialTauBars2(*PreCalcpotentialTauBarsReco2);
01064   NuMatrixSpectrum appearedNuBars(*PreCalcbarAppearedReco);  
01065   NuMatrixSpectrum wrongSignNuMus_NuMu(*PreCalcwrongSignNuMusReco_NuMu);  
01066   NuMatrixSpectrum wrongSignNuMus_NuBar(*PreCalcwrongSignNuMusReco_NuBar);  
01067   NuMatrixSpectrum barNCBackground(*PreCalcbarNCBackground);  
01068   
01070   // Processes for Antineutrinos
01071   
01072   if (fPredictAntiNeutrinos) {
01073 
01074   if (fHelper->DoTaus()){
01075       potentialTauBars1.OscillateMu2Tau(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
01076       potentialTauBars2.OscillateMu2TauBar(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau()); 
01077     } 
01078     appearedNuBars.OscillateMu2MuBar(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());    
01079     barPrediction->OscillateMu2Mu(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());  
01080     //Add everything together
01081     barPrediction->Add(appearedNuBars);
01082     if (fHelper->DoTaus()){
01083       barPrediction->Add(potentialTauBars1);
01084       barPrediction->Add(potentialTauBars2);
01085       }
01086 
01087     wrongSignNuMus_NuMu.OscillateMu2Mu(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
01088     wrongSignNuMus_NuBar.OscillateMu2MuBar(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
01089     barPrediction->Add(wrongSignNuMus_NuMu);
01090     barPrediction->Add(wrongSignNuMus_NuBar);
01091     barPrediction->Add(barNCBackground);
01092   } else {
01093     // Scale the histogram to zero, so as not to confuse people
01094     barPrediction->Spectrum()->Scale(0);
01095   }
01096 
01097   return barPrediction;
01098 }

NuMatrixSpectrum * NuMMRunTransSME::MakeFDPredNuMu ( const NuMMParameters pars  )  [virtual]

Get just the NuMu component of the prediction with oscillations parameters pars.

Reimplemented from NuMMRunNuBar.

Definition at line 991 of file NuMMRunTransSME.cxx.

References NuMatrix1D::Add(), CacheExtrapolation(), NuMMHelperCPT::DoTaus(), NuMMRunNuBar::fHelper, NuMMRunNuBar::fPredictNeutrinos, PreCalcnuAppearedReco, PreCalcnuNCBackground, PreCalcnuPredictionReco, PreCalcpotentialNuTausReco1, PreCalcpotentialNuTausReco2, PreCalcwrongSignNuBarsReco_NuBar, PreCalcwrongSignNuBarsReco_NuMu, and NuMatrixSpectrum::Spectrum().

00992 {
00993 /*  if (!fQuietMode) {
00994     cout << "sn2: " << pars.Sn2() << "; dm2: " << pars.Dm2() 
00995     << "; sn2bar: " << pars.Sn2Bar() << "; dm2bar: " << pars.Dm2Bar()
00996     << "; alpha: " << pars.TransitionProb()
00997     << endl;
00998   }
00999   */
01000   CacheExtrapolation(pars);
01002   // For all neutrinos
01003   
01004   // Make copies of the prepared extrapolation variables, so that
01005   // we can alter them locally with the oscillation parameters
01006   NuMatrixSpectrum *nuPrediction = new NuMatrixSpectrum(*PreCalcnuPredictionReco);
01007   
01008   NuMatrixSpectrum potentialNuTaus1(*PreCalcpotentialNuTausReco1);  
01009   NuMatrixSpectrum potentialNuTaus2(*PreCalcpotentialNuTausReco2);  
01010   NuMatrixSpectrum appearedNuMus(*PreCalcnuAppearedReco);  
01011   NuMatrixSpectrum wrongSignNuBars_NuMu(*PreCalcwrongSignNuBarsReco_NuMu);  
01012   NuMatrixSpectrum wrongSignNuBars_NuBar(*PreCalcwrongSignNuBarsReco_NuBar);  
01013   NuMatrixSpectrum nuNCBackground(*PreCalcnuNCBackground);
01014   
01015   
01017   // Processes for neutrinos
01018   
01019   if (fPredictNeutrinos) {  
01020     
01021     //Add everything together
01022     nuPrediction->Add(appearedNuMus);
01023     if (fHelper->DoTaus()) {
01024       nuPrediction->Add(potentialNuTaus1);    
01025       nuPrediction->Add(potentialNuTaus2);    
01026       }
01027     nuPrediction->Add(wrongSignNuBars_NuMu);
01028     nuPrediction->Add(wrongSignNuBars_NuBar);
01029     nuPrediction->Add(nuNCBackground);
01030   } else {
01031     // Scale the histogram to zero, so as not to confuse people
01032     // if they look at the spectrum with neutrinos turned off
01033     nuPrediction->Spectrum()->Scale(0);
01034   }
01035   
01036   return nuPrediction;
01037 }

void NuMMRunTransSME::OscillateNuBar ( NuMatrixSpectrum numu,
const NuMMParameters pars 
) const

Helper. Oscillates, decays or decoheres a true energy nubar spectrum.

Definition at line 1272 of file NuMMRunTransSME.cxx.

References NuMMParameters::Alpha(), NuMatrix1D::DecayCC(), NuMatrix1D::Decohere(), NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), NuMMParameters::Epsilon(), NuMMRun::fDisappearanceModel, NuMMParameters::Mu2(), NuMatrixSpectrum::Oscillate(), NuMatrixSpectrum::OscillateNSI(), NuMMParameters::Sn2(), and NuMMParameters::Sn2Bar().

01274 {
01275   switch(fDisappearanceModel){
01276   case 0:
01277     numu.Oscillate(pars.Dm2Bar(), pars.Sn2Bar());
01278     return;
01279   case 1:
01280     numu.DecayCC(pars.Dm2Bar(), pars.Sn2Bar(), pars.Alpha());
01281     return;
01282   case 2:
01283     numu.Decohere(pars.Dm2Bar(), pars.Sn2Bar(), pars.Mu2());
01284     return;
01285   case 3:
01286     numu.OscillateNSI(pars.Dm2(), pars.Sn2(), pars.Epsilon(), -1.0);
01287     return;
01288   default:
01289     assert(0 && "Badly configured disappearance model");
01290   }
01291 }

void NuMMRunTransSME::OscillateNuMu ( NuMatrixSpectrum numu,
const NuMMParameters pars 
) const

Definition at line 1234 of file NuMMRunTransSME.cxx.

References NuMMParameters::Alpha(), NuMatrix1D::DecayCC(), NuMatrix1D::Decohere(), NuMMParameters::Dm2(), NuMMParameters::Epsilon(), NuMMRun::fDisappearanceModel, NuMMParameters::Mu2(), NuMatrixSpectrum::Oscillate(), NuMatrixSpectrum::OscillateNSI(), and NuMMParameters::Sn2().

01235                                                                {
01236   switch(fDisappearanceModel){
01237   case 0:
01238     numu.Oscillate(pars.Dm2(), pars.Sn2());
01239     return;
01240   case 1:
01241     numu.DecayCC(pars.Dm2(), pars.Sn2(), pars.Alpha());
01242     return;
01243   case 2:
01244     numu.Decohere(pars.Dm2(), pars.Sn2(), pars.Mu2());
01245     return;
01246   case 3:
01247     numu.OscillateNSI(pars.Dm2(), pars.Sn2(), pars.Epsilon(), 1.0);
01248     return;
01249   default:
01250     assert(0 && "Badly configured disappearance model");
01251   }
01252 }

void NuMMRunTransSME::ResetCache (  )  [private]

Definition at line 100 of file NuMMRunTransSME.cxx.

References NuMMRunNuBar::fCached, PreCalcbarAppeared, PreCalcbarNCBackground, PreCalcbarPrediction, PreCalcnuAppeared, PreCalcnuNCBackground, PreCalcnuPrediction, PreCalcpotentialNuTaus, PreCalcpotentialNuTaus1, PreCalcpotentialNuTaus2, PreCalcpotentialTauBars, PreCalcpotentialTauBars1, PreCalcpotentialTauBars2, PreCalcwrongSignNuBars, PreCalcwrongSignNuBars_NuBar, PreCalcwrongSignNuBars_NuMu, PreCalcwrongSignNuMus, PreCalcwrongSignNuMus_NuBar, and PreCalcwrongSignNuMus_NuMu.

Referenced by CacheExtrapolation(), and ~NuMMRunTransSME().

00101 {
00102   // delete the precalculated MM spectra we previously created
00103   if (PreCalcnuPrediction)    delete PreCalcnuPrediction;
00104   if (PreCalcbarPrediction)   delete PreCalcbarPrediction;
00105   
00106   if (PreCalcnuAppeared)      delete PreCalcnuAppeared;
00107   if (PreCalcbarAppeared)     delete PreCalcbarAppeared; 
00108   
00109   if (PreCalcpotentialNuTaus) delete PreCalcpotentialNuTaus; 
00110   if (PreCalcpotentialTauBars)delete PreCalcpotentialTauBars; 
00111 
00112   if (PreCalcpotentialNuTaus1) delete PreCalcpotentialNuTaus1;
00113   if (PreCalcpotentialTauBars1)delete PreCalcpotentialTauBars1;
00114 
00115   if (PreCalcpotentialNuTaus2) delete PreCalcpotentialNuTaus2;
00116   if (PreCalcpotentialTauBars2)delete PreCalcpotentialTauBars2;
00117   
00118   if (PreCalcnuNCBackground)  delete PreCalcnuNCBackground; 
00119   if (PreCalcbarNCBackground) delete PreCalcbarNCBackground; 
00120   
00121   if (PreCalcwrongSignNuMus)  delete PreCalcwrongSignNuMus; 
00122   if (PreCalcwrongSignNuBars) delete PreCalcwrongSignNuBars; 
00123   
00124   if (PreCalcwrongSignNuMus_NuMu)  delete PreCalcwrongSignNuMus_NuMu; 
00125   if (PreCalcwrongSignNuBars_NuMu) delete PreCalcwrongSignNuBars_NuMu; 
00126   
00127   if (PreCalcwrongSignNuMus_NuBar)  delete PreCalcwrongSignNuMus_NuBar; 
00128   if (PreCalcwrongSignNuBars_NuBar) delete PreCalcwrongSignNuBars_NuBar; 
00129 
00130   PreCalcnuPrediction   = PreCalcbarPrediction      = 0;
00131   PreCalcnuAppeared     = PreCalcbarAppeared        = 0;
00132   PreCalcpotentialNuTaus= PreCalcpotentialTauBars   = 0;
00133   PreCalcpotentialNuTaus1= PreCalcpotentialTauBars1   = 0;
00134   PreCalcpotentialNuTaus2= PreCalcpotentialTauBars2   = 0;
00135   PreCalcnuNCBackground = PreCalcbarNCBackground    = 0;
00136   PreCalcwrongSignNuMus = PreCalcwrongSignNuBars    = 0;  
00137   PreCalcwrongSignNuMus_NuMu = PreCalcwrongSignNuBars_NuMu    = 0;  
00138   PreCalcwrongSignNuMus_NuBar = PreCalcwrongSignNuBars_NuBar    = 0;  
00139   
00140 
00141   fCached = false;
00142 }

void NuMMRunTransSME::SetSystShifts ( TFile *  systFilePQ,
TFile *  systFileNQ 
)

Definition at line 708 of file NuMMRunTransSME.cxx.

References FD_data, fhOne, NuMMRunNuBar::fPredictAntiNeutrinos, NuMMRunNuBar::fPredictNeutrinos, fSystMinusNQ, fSystMinusPQ, fSystPlusNQ, fSystPlusPQ, kDecayPipe, kNumSyst, kShwScaleTrans, and kTrkScaleTrans.

00709 {
00710    string SystName[kNumSyst];
00711    SystName[kShwScaleTrans] = "ShowerEnergyScaleFunctionBoth";
00712    SystName[kTrkScaleTrans] = "TrackEnergyOverall";
00713    SystName[kDecayPipe] = "DecayPipe";
00714    //SystName[kNuBarXSecDIS] = "NuMuBarXSecDISMultip2";
00715 
00716    FD_data = new TFile("/minos/data/users/richa/transitions/Include_FHCNu_FHCNuBar/helpers/FDRun1LEMCSummaryNorm.root");
00717    assert(FD_data);
00718    TH1D *binTemplate = (TH1D*)FD_data->Get("RecoEnergyPQ_FD");
00719 
00720    assert(systFilePQ);
00721    assert(systFileNQ);
00722 
00723    for(int kPar=0; kPar<kNumSyst; kPar++){
00724 
00725      if (fPredictNeutrinos) {
00726        //get NQ bands
00727        systFileNQ->GetObject(("NuMus_"+SystName[kPar]+"Plus").c_str(),fSystPlusNQ[kPar]);
00728        systFileNQ->GetObject(("NuMus_"+SystName[kPar]+"Minus").c_str(),fSystMinusNQ[kPar]);
00729      }
00730      if (fPredictAntiNeutrinos) {
00731        //get PQ bands
00732        systFilePQ->GetObject(("NuMuBars_"+SystName[kPar]+"Plus").c_str(),fSystPlusPQ[kPar]);
00733        systFilePQ->GetObject(("NuMuBars_"+SystName[kPar]+"Minus").c_str(),fSystMinusPQ[kPar]);
00734      }
00735 
00736        fhOne[kPar] = (TH1D*)binTemplate->Clone();
00737        for(int i=0; i<=fhOne[kPar]->GetNbinsX()+1; i++){
00738        fhOne[kPar]->SetBinContent(i,1);
00739        fhOne[kPar]->SetBinError(i,0);
00740   }
00741  }
00742 }

NuMatrixSpectrum NuMMRunTransSME::TrueComponents ( const NuMMParameters pars,
Sample_t  s 
)

Definition at line 1103 of file NuMMRunTransSME.cxx.

References CacheExtrapolation(), and TrueComponents().

01104 {
01105   CacheExtrapolation(pars);
01106   
01107   // Return the const version
01108   return static_cast<const NuMMRunTransSME *>(this)->TrueComponents(pars, s);
01109 }

NuMatrixSpectrum NuMMRunTransSME::TrueComponents ( const NuMMParameters pars,
Sample_t  s 
) const [virtual]

Return the individual true energy prediction of the far detector sub sample s with oscillation parameters described in pars. Default implementation just drops an error.

Reimplemented from NuMMRun.

Definition at line 1112 of file NuMMRunTransSME.cxx.

References NuMMRunNuBar::fCached, NuMMRunNuBar::fPredictAntiNeutrinos, NuMMRunNuBar::fPredictNeutrinos, NuMMRun::fQuietMode, Sample::kAppearedNQ, Sample::kAppearedPQ, Msg::kError, Sample::kNCNQ, Sample::kNCPQ, Sample::kSignalNQ, Sample::kSignalPQ, Sample::kTauNQ, Sample::kTauPQ, Msg::kWarning, Sample::kWrongSignNQ, Sample::kWrongSignPQ, MSG, NuMatrixSpectrum::Multiply(), PreCalcbarAppeared, PreCalcbarNCBackground, PreCalcbarPrediction, PreCalcnuAppeared, PreCalcnuNCBackground, PreCalcnuPrediction, PreCalcpotentialNuTaus, PreCalcpotentialTauBars, PreCalcwrongSignNuBars, PreCalcwrongSignNuMus, NuMatrixSpectrum::SetName(), NuMatrixSpectrum::Spectrum(), and NuMMParameters::TransitionProb().

Referenced by TrueComponents().

01113 {
01114 /*  if (!fQuietMode) {
01115     cout << "sn2: " << pars.Sn2() << "; dm2: " << pars.Dm2() 
01116     << "; sn2bar: " << pars.Sn2Bar() << "; dm2bar: " << pars.Dm2Bar()
01117     << "; alpha: " << pars.TransitionProb()
01118     << endl;
01119   }
01120   */
01121   // Request NQ spectrum, but no NuBar Predictions OR
01122   // Request PQ spectrum, but no NuMu Predictions
01123   if ( (s < 0 && !fPredictNeutrinos) || 
01124        (s > 0 && !fPredictAntiNeutrinos) ) {
01125     NuMatrixSpectrum blank(*PreCalcnuPrediction);
01126     blank.SetName("Blank");
01127     blank.Spectrum()->Scale(0);
01128     return blank;
01129   }
01130   
01131   if (!fCached) {
01132     // CacheExtrapolation(pars);
01133     MSG("NuMMRunTransitions",Msg::kError) << "Calling TrueComponents without previously caching the extrapolation" << endl;
01134     assert(fCached);
01135   }
01136   
01137 
01138   if (s == kSignalPQ) {
01139     if (!fQuietMode) cout << "Getting kSignalPQ" << endl;
01140     NuMatrixSpectrum barPrediction(*PreCalcbarPrediction);
01141     if (!fQuietMode) cout << "Returning kSignalPQ" << endl;
01142     return barPrediction;
01143   }
01144   else if (s == kSignalNQ) {
01145     if (!fQuietMode) cout << "Getting kSignalNQ" << endl;
01146     NuMatrixSpectrum nuPrediction(*PreCalcnuPrediction);  
01147     if (!fQuietMode) cout << "Returning kSignalNQ" << endl;
01148     return nuPrediction;
01149   }
01150   else if (s == kWrongSignPQ) {
01151     if (!fQuietMode) cout << "Getting kWrongSignPQ" << endl;
01152     NuMatrixSpectrum wrongSignNuMus(*PreCalcwrongSignNuMus);
01153     if (!fQuietMode) cout << "Returning kWrongSignPQ" << endl;
01154     return wrongSignNuMus;
01155   }
01156   else if (s == kWrongSignNQ) {
01157     if (!fQuietMode) cout << "Getting kWrongSignNQ" << endl;
01158     NuMatrixSpectrum wrongSignNuBars(*PreCalcwrongSignNuBars);
01159     if (!fQuietMode) cout << "Returning kWrongSignNQ" << endl;
01160     return wrongSignNuBars;
01161   }  
01162   else if (s == kNCPQ) {
01163     if (!fQuietMode) cout << "Getting kNCPQ" << endl;
01164     NuMatrixSpectrum barNCBackground(*PreCalcbarNCBackground);
01165     if (!fQuietMode) cout << "Returning kNCPQ" << endl;
01166     return barNCBackground;
01167   }
01168   else if (s == kNCNQ) {
01169     if (!fQuietMode) cout << "Getting kNCNQ" << endl;
01170     NuMatrixSpectrum nuNCBackground(*PreCalcnuNCBackground);
01171     if (!fQuietMode) cout << "Returning kNCNQ" << endl;
01172     return nuNCBackground;
01173   }
01174   else if (s == kTauPQ) {
01175     if (!fQuietMode) cout << "Getting kTauPQ" << endl;
01176     NuMatrixSpectrum potentialTauBars(*PreCalcpotentialTauBars);
01177     potentialTauBars.Multiply(1.-pars.TransitionProb());
01178     if (!fQuietMode) cout << "Returning kTauPQ" << endl;
01179     return potentialTauBars;
01180   }
01181   else if (s == kTauNQ) {
01182     if (!fQuietMode) cout << "Getting kTauNQ" << endl;
01183     NuMatrixSpectrum potentialNuTaus(*PreCalcpotentialNuTaus);
01184     potentialNuTaus.Multiply(1.-pars.TransitionProb());
01185     if (!fQuietMode) cout << "Returning kTauNQ" << endl;
01186     return potentialNuTaus;
01187   }
01188   else if (s == kAppearedPQ) {
01189     if (!fQuietMode) cout << "Getting kAppearedPQ" << endl;
01190     NuMatrixSpectrum appearedNuBars(*PreCalcbarAppeared);
01191     appearedNuBars.Multiply(pars.TransitionProb());
01192     if (!fQuietMode) cout << "Returning kAppearedPQ" << endl;
01193     return appearedNuBars;
01194   }
01195   else if (s == kAppearedNQ) {
01196     if (!fQuietMode) cout << "Getting kAppearedNQ" << endl;
01197     NuMatrixSpectrum appearedNuMus(*PreCalcnuAppeared);
01198     appearedNuMus.Multiply(pars.TransitionProb());
01199     if (!fQuietMode) cout << "Returning kAppearedNQ" << endl;
01200     return appearedNuMus;
01201   }
01202   else {
01203     MSG("NuMMRunTransitions",Msg::kWarning) << "Requested unknown prediction #" << s << ", returning blank." << endl;
01204     NuMatrixSpectrum blank(*PreCalcnuPrediction);
01205     blank.SetName("Blank");
01206     blank.Spectrum()->Scale(0);
01207     return blank;
01208   }
01209 }

std::vector< TH1D > NuMMRunTransSME::WriteFDPredHistos ( const NuMMParameters pars  )  const [virtual]

Implements NuMMRun.

Definition at line 424 of file NuMMRunTransSME.cxx.

References NuMatrix1D::Add(), NuMatrixSpectrum::Divide(), NuMatrixSpectrum::ExtrapolateNDToFD(), NuMatrixSpectrum::Multiply(), NuMatrixSpectrum::OscillateMu2Mu(), NuMatrixSpectrum::OscillateMu2MuBar(), NuMatrixSpectrum::OscillateMu2Tau(), NuMatrixSpectrum::OscillateMu2TauBar(), NuMatrixSpectrum::RecoToTrue(), NuMatrixSpectrum::Spectrum(), and NuMatrixSpectrum::TrueToReco().

00425 {
00426   //Set up a vector to push the histograms into.
00427   vector<TH1D> vHistos;
00428 
00429   //Put the ND data in the vector:
00430   TH1D hNDNuData(*(fndNuData->Spectrum()));
00431   hNDNuData.SetName("ndDataNQ");
00432   hNDNuData.SetTitle("ndDataNQ");
00433   vHistos.push_back(*(new TH1D(hNDNuData)));
00434   TH1D hNDBarData(*(fndBarData->Spectrum()));
00435   hNDBarData.SetName("ndDataPQ");
00436   hNDBarData.SetTitle("ndDataPQ");
00437   vHistos.push_back(*(new TH1D(hNDBarData)));
00438 
00439 //cout << "sn2: " << pars.Sn2() << "; dm2: " << pars.Dm2() 
00440 //     << "; sn2bar: " << pars.Sn2Bar() << "; dm2bar: " << pars.Dm2Bar()
00441 //     << "; Prob: " << pars.TransitionProb()
00442 //     << endl;
00443   //Get the neutrinos to the FD
00444   NuMatrixSpectrum nuPrediction(*fndNuData);
00445   nuPrediction.Multiply(fHelper->NDNuPurity());
00446   nuPrediction.RecoToTrue(fHelper->NDNuRecoVsTrue());
00447   nuPrediction.Divide(fHelper->NDNuEfficiency());
00448   nuPrediction.Divide(fHelper->NuXSecGraph());
00449   nuPrediction.Divide(fndNuData->GetPOT());
00450   nuPrediction.Divide(fNDFidMass);
00451   nuPrediction.ExtrapolateNDToFD(fHelper->NuBeamMatrix());
00452   nuPrediction.Multiply(fFDFidMass);
00453   nuPrediction.Multiply(ffdNuData->GetPOT());
00454 
00455   //Get the antineutrinos to the FD
00456   NuMatrixSpectrum barPrediction(*fndBarData);
00457   barPrediction.Multiply(fHelper->NDBarPurity());
00458   barPrediction.RecoToTrue(fHelper->NDBarRecoVsTrue());
00459   barPrediction.Divide(fHelper->NDBarEfficiency());
00460   barPrediction.Divide(fHelper->BarXSecGraph());
00461   barPrediction.Divide(fndBarData->GetPOT());
00462   barPrediction.Divide(fNDFidMass);
00463   barPrediction.ExtrapolateNDToFD(fHelper->BarBeamMatrix());
00464   barPrediction.Multiply(fFDFidMass);
00465   barPrediction.Multiply(ffdBarData->GetPOT());
00466 
00467   NuMatrixSpectrum potentialNuTaus1(nuPrediction);
00468   NuMatrixSpectrum potentialNuTaus2(barPrediction); 
00469   NuMatrixSpectrum potentialTauBars1(barPrediction);
00470   NuMatrixSpectrum potentialTauBars2(nuPrediction);
00471 
00472   if (fHelper->DoTaus())  {
00473     //Get the taus (ignoring any wrong-sign taubars)          // numu to nutau
00474     potentialNuTaus1.Multiply(fHelper->XSecGraphNuTaus());
00475     potentialNuTaus1.Multiply(fHelper->FDNuTauEfficiency());
00476     potentialNuTaus1.OscillateMu2Tau(pars.Dm2(),pars.Sn2(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00477     potentialNuTaus1.TrueToReco(fHelper->FDNuTauRecoVsTrue());
00478 
00479     potentialNuTaus2.Multiply(fHelper->XSecGraphNuTaus());     // numubar to nutau
00480     potentialNuTaus2.Multiply(fHelper->FDNuTauEfficiency());
00481     potentialNuTaus2.OscillateMu2TauBar(pars.Dm2(),pars.Sn2(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00482     potentialNuTaus2.TrueToReco(fHelper->FDNuTauRecoVsTrue());
00483 
00484     //Get the taubars (ignoring any wrong-sign taus)      // numubar to nutaubar
00485     potentialTauBars1.Multiply(fHelper->XSecGraphTauBars());
00486     potentialTauBars1.Multiply(fHelper->FDTauBarEfficiency());
00487     potentialTauBars1.OscillateMu2Tau(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00488     potentialTauBars1.TrueToReco(fHelper->FDTauBarRecoVsTrue());
00489  
00490     potentialTauBars2.Multiply(fHelper->XSecGraphTauBars());  // numu to nutaubar
00491     potentialTauBars2.Multiply(fHelper->FDTauBarEfficiency());
00492     potentialTauBars2.OscillateMu2TauBar(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00493     potentialTauBars2.TrueToReco(fHelper->FDTauBarRecoVsTrue());
00494 
00495   }
00496   else {
00497     potentialNuTaus1.Multiply(0.0);
00498     potentialTauBars1.Multiply(0.0);
00499     potentialNuTaus2.Multiply(0.0);
00500     potentialTauBars2.Multiply(0.0);
00501   }
00502   
00503   //Get the transitioned neutrinos (they're now antineutrinos)
00504   NuMatrixSpectrum appearedNuBars(nuPrediction);
00505   TH1D hRawAppearPQ(*(appearedNuBars.Spectrum()));
00506   hRawAppearPQ.SetName("hRawAppearPQ");
00507   vHistos.push_back(hRawAppearPQ); 
00508   
00509   appearedNuBars.Multiply(fHelper->BarXSecGraph());
00510   TH1D hTrueEffAppearPQ(*(appearedNuBars.Spectrum()));
00511   hTrueEffAppearPQ.SetName("hTrueEffAppearPQ");
00512   vHistos.push_back(hTrueEffAppearPQ); 
00513   
00514   appearedNuBars.Multiply(fHelper->FDBarEfficiency());
00515   TH1D hTrueAppearPQ(*(appearedNuBars.Spectrum()));
00516   hTrueAppearPQ.SetName("hTrueAppearPQ");
00517   vHistos.push_back(hTrueAppearPQ); 
00518   
00519 //  appearedNuBars.InverseOscillateTrans(pars.Dm2(),pars.Sn2());
00520 //  appearedNuBars.Multiply(pars.TransitionProb());
00521 
00522   appearedNuBars.OscillateMu2MuBar(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00523   TH1D hTrueTransitionPQ(*(appearedNuBars.Spectrum()));
00524   hTrueTransitionPQ.SetName("hTrueTransitionPQ");
00525   vHistos.push_back(hTrueTransitionPQ); 
00526     
00527   appearedNuBars.TrueToReco(fHelper->FDBarRecoVsTrue());
00528   TH1D hRecoAppearPQ(*(appearedNuBars.Spectrum()));
00529   hRecoAppearPQ.SetName("hRecoAppearPQ");
00530   vHistos.push_back(hRecoAppearPQ);
00531   
00532   //Get the transitioned antineutrinos (they're now neutrinos)
00533   NuMatrixSpectrum appearedNuMus(barPrediction);
00534   TH1D hRawAppearNQ(*(appearedNuMus.Spectrum()));
00535   hRawAppearNQ.SetName("hRawAppearNQ");
00536   vHistos.push_back(hRawAppearNQ); 
00537   
00538   appearedNuMus.Multiply(fHelper->NuXSecGraph());
00539   TH1D hTrueEffAppearNQ(*(appearedNuMus.Spectrum()));
00540   hTrueEffAppearNQ.SetName("hTrueEffAppearNQ");
00541   vHistos.push_back(hTrueEffAppearNQ); 
00542   
00543   appearedNuMus.Multiply(fHelper->FDNuEfficiency());
00544   TH1D hTrueAppearNQ(*(appearedNuMus.Spectrum()));
00545   hTrueAppearNQ.SetName("hTrueAppearNQ");
00546   vHistos.push_back(hTrueAppearNQ); 
00547   
00548   appearedNuMus.OscillateMu2MuBar(pars.Dm2(),pars.Sn2(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00549 //  appearedNuMus.Multiply(pars.TransitionProb());
00550   TH1D hTrueTransitionNQ(*(appearedNuMus.Spectrum()));
00551   hTrueTransitionNQ.SetName("hTrueTransitionNQ");
00552   vHistos.push_back(hTrueTransitionNQ); 
00553   
00554   appearedNuMus.TrueToReco(fHelper->FDNuRecoVsTrue());
00555   TH1D hRecoAppearNQ(*(appearedNuMus.Spectrum()));
00556   hRecoAppearNQ.SetName("hRecoAppearNQ");
00557   vHistos.push_back(hRecoAppearNQ);
00558   
00559 
00560   //Cross sections
00561   nuPrediction.Multiply(fHelper->NuXSecGraph());
00562   barPrediction.Multiply(fHelper->BarXSecGraph());
00563   
00564   // True NuMus (flux x sigma)
00565   TH1D hTrueNuMus(*(nuPrediction.Spectrum()));
00566   hTrueNuMus.SetName("hTrueNuMus");
00567   vHistos.push_back(hTrueNuMus); 
00568   
00569   // True NuBars (flux x sigma)
00570   TH1D hTrueNuBars(*(barPrediction.Spectrum()));
00571   hTrueNuBars.SetName("hTrueNuBars");
00572   vHistos.push_back(hTrueNuBars); 
00573 
00574   //Get the neutrino NC background
00575   NuMatrixSpectrum nuNCBackground(nuPrediction);
00576   nuNCBackground.Multiply(fHelper->FDNuEfficiency());
00577   nuNCBackground.TrueToReco(fHelper->FDNuRecoVsTrue());
00578   nuNCBackground.Divide(fHelper->FDNuPurity());
00579   nuNCBackground.Multiply(fHelper->FDNuNCContamination());
00580   
00581   //Get the antineutrino NC background
00582   NuMatrixSpectrum barNCBackground(barPrediction);
00583   barNCBackground.Multiply(fHelper->FDBarEfficiency());
00584   barNCBackground.TrueToReco(fHelper->FDBarRecoVsTrue());
00585   barNCBackground.Divide(fHelper->FDBarPurity());
00586   barNCBackground.Multiply(fHelper->FDBarNCContamination());
00587   
00588   //Get the wrong-sign neutrino background for the antineutrino
00589   //prediction
00590   NuMatrixSpectrum wrongSignNuMus(nuPrediction);
00591   wrongSignNuMus.Multiply(fHelper->FDWrongSignNuEfficiency());
00592   wrongSignNuMus.OscillateMu2Mu(pars.Dm2(),pars.Sn2(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00593   
00594   //Get the wrong-sign antineutrino background for the neutrino
00595   //prediction
00596   NuMatrixSpectrum wrongSignNuBars(barPrediction);
00597   wrongSignNuBars.Multiply(fHelper->FDWrongSignBarEfficiency());
00598   wrongSignNuBars.OscillateMu2Mu(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());  
00599 
00600 
00601   TH1D hFDWrongSignPQ(*(wrongSignNuBars.Spectrum()));
00602   hFDWrongSignPQ.SetName("fdWrongSignPQ");
00603   vHistos.push_back(hFDWrongSignPQ);
00604 
00605   TH1D hFDWrongSignNQ(*(wrongSignNuMus.Spectrum()));
00606   hFDWrongSignNQ.SetName("fdWrongSignNQ");
00607   vHistos.push_back(hFDWrongSignNQ);
00608 
00609   // Do the truetoreco seperately so that we can write out the histos first
00610   wrongSignNuBars.TrueToReco(fHelper->FDWrongSignBarRecoVsTrue());
00611   wrongSignNuMus.TrueToReco(fHelper->FDWrongSignNuRecoVsTrue()); 
00612 
00613   //Oscillations
00614   nuPrediction.OscillateMu2Mu(pars.Dm2(),pars.Sn2(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00615   barPrediction.OscillateMu2Mu(pars.Dm2Bar(),pars.Sn2Bar(),pars.GMu(),pars.GTau(),pars.CMuMu(),pars.CTauTau());
00616   //Efficiencies
00617   nuPrediction.Multiply(fHelper->FDNuEfficiency());
00618   barPrediction.Multiply(fHelper->FDBarEfficiency());
00619 
00620   //True to reco
00621   nuPrediction.TrueToReco(fHelper->FDNuRecoVsTrue());
00622   barPrediction.TrueToReco(fHelper->FDBarRecoVsTrue());
00623   
00624   // Dump out the prediction spectrum before adding in the backgrounds
00625   TH1D hFDNoBackNQ(*(nuPrediction.Spectrum()));
00626   hFDNoBackNQ.SetName("fdBasePredictionNQ");
00627   vHistos.push_back(hFDNoBackNQ);
00628   
00629   TH1D hFDNoBackPQ(*(barPrediction.Spectrum()));
00630   hFDNoBackPQ.SetName("fdBasePredictionPQ");
00631   vHistos.push_back(hFDNoBackPQ);
00632   
00633   //Add in backgrounds
00634   nuPrediction.Add(wrongSignNuBars);
00635   barPrediction.Add(wrongSignNuMus);
00636   nuPrediction.Add(nuNCBackground);
00637 
00638   if (fHelper->DoTaus()) { nuPrediction.Add(potentialNuTaus1);
00639                            nuPrediction.Add(potentialNuTaus2);
00640                          }
00641   barPrediction.Add(barNCBackground);
00642   if (fHelper->DoTaus()) { barPrediction.Add(potentialTauBars1);
00643                            barPrediction.Add(potentialTauBars2);
00644                          }
00645 
00646   
00647   TH1D hFDPredictionNoTransNQ(*(nuPrediction.Spectrum()));
00648   hFDPredictionNoTransNQ.SetName("fdPredictionNoTransNQ");
00649   vHistos.push_back(hFDPredictionNoTransNQ);
00650   
00651   TH1D hFDPredictionNoTransPQ(*(barPrediction.Spectrum()));
00652   hFDPredictionNoTransPQ.SetName("fdPredictionNoTransPQ");
00653   vHistos.push_back(hFDPredictionNoTransPQ);
00654   
00655   
00656   
00657   //Add in transitions
00658   nuPrediction.Add(appearedNuMus);
00659   barPrediction.Add(appearedNuBars);
00660   
00661   //Put the componants of the predictions in the vector
00662 
00663   TH1D hFDPredictionNQ(*(nuPrediction.Spectrum()));
00664   hFDPredictionNQ.SetName("fdPredictionNQ");
00665   vHistos.push_back(hFDPredictionNQ);
00666 
00667   TH1D hFDPredictionPQ(*(barPrediction.Spectrum()));
00668   hFDPredictionPQ.SetName("fdPredictionPQ");
00669   vHistos.push_back(hFDPredictionPQ);
00670 
00671   TH1D hFDNCNQ(*(nuNCBackground.Spectrum()));
00672   hFDNCNQ.SetName("fdNCNQ");
00673   vHistos.push_back(hFDNCNQ);
00674 
00675   TH1D hFDNCPQ(*(barNCBackground.Spectrum()));
00676   hFDNCPQ.SetName("fdNCPQ");
00677   vHistos.push_back(hFDNCPQ);
00678   
00679   if (fHelper->DoTaus()) {
00680     TH1D hFDTausNQ1(*(potentialNuTaus1.Spectrum()));
00681     hFDTausNQ1.SetName("fdTausNQ_fromNuMu");
00682     vHistos.push_back(hFDTausNQ1);
00683     
00684     TH1D hFDTausPQ1(*(potentialTauBars1.Spectrum()));
00685     hFDTausPQ1.SetName("fdTausPQ_fromNuMuBar");
00686     vHistos.push_back(hFDTausPQ1);
00687 
00688     TH1D hFDTausNQ2(*(potentialNuTaus2.Spectrum()));
00689     hFDTausNQ2.SetName("fdTausNQ_fromNuMuBar");
00690     vHistos.push_back(hFDTausNQ2);
00691 
00692     TH1D hFDTausPQ2(*(potentialTauBars2.Spectrum()));
00693     hFDTausPQ2.SetName("fdTausPQ_fromNuMu");
00694     vHistos.push_back(hFDTausPQ2);
00695   }
00696   
00697   TH1D hFDDataNQ(*(ffdNuData->Spectrum()));
00698   hFDDataNQ.SetName("fdDataNQ");
00699   vHistos.push_back(hFDDataNQ);
00700 
00701   TH1D hFDDataPQ(*(ffdBarData->Spectrum()));
00702   hFDDataPQ.SetName("fdDataPQ");
00703   vHistos.push_back(hFDDataPQ);
00704 
00705   return vHistos;
00706 }


Member Data Documentation

TFile * NuMMRunTransSME::FD_data [private]

Definition at line 97 of file NuMMRunTransSME.h.

Referenced by SetSystShifts().

TH1D* NuMMRunTransSME::FD_dataCC [private]

Definition at line 98 of file NuMMRunTransSME.h.

TFile* NuMMRunTransSME::FD_MC [private]

Definition at line 97 of file NuMMRunTransSME.h.

TFile * NuMMRunTransSME::FD_NuTau [private]

Definition at line 97 of file NuMMRunTransSME.h.

Definition at line 59 of file NuMMRunTransSME.h.

Referenced by FitSysts(), MakeFDPred(), and NuMMRunTransSME().

TH1D* NuMMRunTransSME::fhOne[kNumSyst]

Definition at line 55 of file NuMMRunTransSME.h.

Referenced by MakeFDPred(), and SetSystShifts().

TH1D* NuMMRunTransSME::fhOneNQ[kNumSyst]

Definition at line 54 of file NuMMRunTransSME.h.

TH1D* NuMMRunTransSME::fhOnePQ[kNumSyst]

Definition at line 53 of file NuMMRunTransSME.h.

Definition at line 52 of file NuMMRunTransSME.h.

Referenced by MakeFDPred(), and SetSystShifts().

Definition at line 51 of file NuMMRunTransSME.h.

Referenced by MakeFDPred(), and SetSystShifts().

Definition at line 50 of file NuMMRunTransSME.h.

Referenced by MakeFDPred(), and SetSystShifts().

Definition at line 49 of file NuMMRunTransSME.h.

Referenced by MakeFDPred(), and SetSystShifts().

Definition at line 112 of file NuMMRunTransSME.h.

Referenced by MakeFDPredNuBar().

Definition at line 106 of file NuMMRunTransSME.h.

Referenced by MakeFDPredNuBar().

Definition at line 111 of file NuMMRunTransSME.h.

Referenced by MakeFDPredNuMu().

Definition at line 105 of file NuMMRunTransSME.h.

Referenced by MakeFDPredNuMu().

Definition at line 120 of file NuMMRunTransSME.h.

Referenced by CacheExtrapolation(), MakeFDPred(), NuMMRunTransSME(), and ResetCache().

Definition at line 121 of file NuMMRunTransSME.h.

Referenced by CacheExtrapolation(), MakeFDPred(), NuMMRunTransSME(), and ResetCache().

Definition at line 118 of file NuMMRunTransSME.h.

Definition at line 124 of file NuMMRunTransSME.h.

Referenced by MakeFDPredNuMu().

Definition at line 125 of file NuMMRunTransSME.h.

Referenced by MakeFDPredNuMu().

Definition at line 122 of file NuMMRunTransSME.h.

Referenced by CacheExtrapolation(), MakeFDPred(), NuMMRunTransSME(), and ResetCache().

Definition at line 123 of file NuMMRunTransSME.h.

Referenced by CacheExtrapolation(), MakeFDPred(), NuMMRunTransSME(), and ResetCache().

Definition at line 119 of file NuMMRunTransSME.h.

Definition at line 126 of file NuMMRunTransSME.h.

Referenced by MakeFDPredNuBar().

Definition at line 127 of file NuMMRunTransSME.h.

Referenced by MakeFDPredNuBar().

Definition at line 147 of file NuMMRunTransSME.h.

Referenced by CacheExtrapolation(), MakeFDPred(), NuMMRunTransSME(), and ResetCache().

Definition at line 142 of file NuMMRunTransSME.h.

Referenced by CacheExtrapolation(), MakeFDPred(), NuMMRunTransSME(), and ResetCache().

Definition at line 139 of file NuMMRunTransSME.h.

Definition at line 149 of file NuMMRunTransSME.h.

Referenced by MakeFDPredNuMu().

Definition at line 144 of file NuMMRunTransSME.h.

Referenced by MakeFDPredNuMu().

Definition at line 146 of file NuMMRunTransSME.h.

Referenced by CacheExtrapolation(), MakeFDPred(), NuMMRunTransSME(), and ResetCache().

Definition at line 141 of file NuMMRunTransSME.h.

Referenced by CacheExtrapolation(), MakeFDPred(), NuMMRunTransSME(), and ResetCache().

Definition at line 138 of file NuMMRunTransSME.h.

Definition at line 148 of file NuMMRunTransSME.h.

Referenced by MakeFDPredNuBar().

Definition at line 143 of file NuMMRunTransSME.h.

Referenced by MakeFDPredNuBar().


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

Generated on 22 Nov 2017 for loon by  doxygen 1.6.1