NuMMRunNSINu Class Reference

#include <NuMMRunNSINu.h>

Inheritance diagram for NuMMRunNSINu:
NuMMRun

List of all members.

Public Member Functions

 NuMMRunNSINu (NuMMHelperPRL *helper, NuMatrixSpectrum *ndCCData, NuMatrixSpectrum *fdCCData)
virtual Double_t ComparePredWithData (const NuMMParameters &pars)
virtual vector< NuMatrixSpectrumMakeFDPred (const NuMMParameters &pars)
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.
virtual std::vector< TH1D > WriteFDPredHistos (const NuMMParameters &pars) const
virtual const NuMatrixSpectrumGetNDDataNSI () const
 Get a const pointer to the ND data spectrum.
virtual const NuMatrixSpectrumGetFDDataNSI () const
 Get a const pointer to the FD data spectrum.
virtual void PrecalculateWSFlux (NuMMHelperPRL *helper, const NuMatrixSpectrum *ndOtherData)
virtual void SetFDData (NuMatrixSpectrum *fdData)

Private Member Functions

virtual vector< NuMatrixSpectrumCalculateFDFlux (const Double_t ncScale=1.0) const
virtual vector< NuMatrixSpectrumCalculateFDTauFlux (const Double_t ncScale=1.0) const
virtual vector< NuMatrixSpectrumCalculateFDNCBackground (const Double_t ncScale=1.0) const
void OscillateNuMu (NuMatrixSpectrum &numu, const NuMMParameters &pars) const
 Helper. Oscillates, decays or decoheres a true energy numu spectrum.
void OscillateNuBar (NuMatrixSpectrum &numu, const NuMMParameters &pars) const
void InverseOscillateNuTau (NuMatrixSpectrum &nutau, const NuMMParameters &pars) const
 Helper. Inverse oscillates a true energy nutau spectrum.
void InverseOscillateTauBar (NuMatrixSpectrum &nutau, const NuMMParameters &pars) const
virtual void SanityCheckSpecialWSFiles (NuMMHelperPRL *helper, const NuMatrixSpectrum *ndOtherData) const
virtual void WriteToJess (NuMatrixSpectrum spectrum, const std::string name, bool osc)

Private Attributes

Bool_t fSpecialWSMode
NuMatrixSpectrumfNDCCData
NuMatrixSpectrumfFDCCData
vector< NuMatrixSpectrumfFDFlux
vector< NuMatrixSpectrumfFDTauFlux
vector< NuMatrixSpectrumfFDNCBackground
vector< NuMatrixSpectrumfSpecialWSBackgroundFlux
NuMMHelperPRLfHelper

Detailed Description

Definition at line 16 of file NuMMRunNSINu.h.


Constructor & Destructor Documentation

NuMMRunNSINu::NuMMRunNSINu ( NuMMHelperPRL helper,
NuMatrixSpectrum ndCCData,
NuMatrixSpectrum fdCCData 
)

Definition at line 21 of file NuMMRunNSINu.cxx.

00024   : NuMMRun(),
00025     fSpecialWSMode(false),
00026     fNDCCData(ndCCData),
00027     fFDCCData(fdCCData),
00028     fHelper(helper)
00029 {
00030   assert(helper);
00031   assert(ndCCData);
00032   assert(fdCCData);
00033   
00034   //Independent of fit parameters for a stats-only fit
00035   fFDFlux = this->CalculateFDFlux();
00036   fFDTauFlux = this->CalculateFDTauFlux();
00037   fFDNCBackground = this->CalculateFDNCBackground();
00038 }


Member Function Documentation

vector< NuMatrixSpectrum > NuMMRunNSINu::CalculateFDFlux ( const Double_t  ncScale = 1.0  )  const [private, virtual]

Definition at line 118 of file NuMMRunNSINu.cxx.

References NuMMHelperPRL::BeamMatrixNuMuCCXSec(), NuMatrixSpectrum::Divide(), NuMatrixSpectrum::ExtrapolateNDToFD(), fFDCCData, NuMMRun::fFDFidMass, fHelper, fNDCCData, NuMMRun::fNDFidMass, NuMatrix::GetPOT(), NuMatrixSpectrum::Multiply(), NuMMHelperPRL::NDCCContamination(), NuMMHelperPRL::NDEfficiency(), NuMMHelperPRL::NDPurity(), NuMMHelperPRL::NDRecoVsTrue(), NuMatrixSpectrum::RecoToTrue(), and NuMatrix1D::Subtract().

Referenced by MakeFDPred().

00119 {
00120   //  cout << "Calculating FD Flux..." << endl;
00121   NuMatrixSpectrum prediction(*fNDCCData);
00122 //   cout << "Got ND, integral: " << prediction.Spectrum()->Integral() << endl
00123 //        << "POT" << fNDCCData->GetPOT() << endl;
00124 
00125   NuMatrixSpectrum wsBackground(prediction);
00126   
00127   NuMatrixSpectrum signal(prediction);
00128   signal.Multiply(fHelper->NDPurity());
00129   NuMatrixSpectrum ncBackground(prediction);
00130   ncBackground.Subtract(signal);
00131   ncBackground.Multiply(ncScale);
00132   prediction.Subtract(ncBackground);
00133   //  cout << "Subtract NC, integral: " << prediction.Spectrum()->Integral() << endl;
00134   prediction.Multiply(fHelper->NDPurity());
00135   
00136   wsBackground.Multiply(fHelper->NDCCContamination());
00137   prediction.Subtract(wsBackground);
00138   //  cout << "Subtract WS, integral: " << prediction.Spectrum()->Integral() << endl;
00139 
00140   prediction.RecoToTrue(fHelper->NDRecoVsTrue());
00141   prediction.Divide(fHelper->NDEfficiency());
00142   prediction.Divide(fNDCCData->GetPOT());
00143   prediction.Divide(fNDFidMass);
00144   //  cout << "Eff-POT-Mass, integral: " << prediction.Spectrum()->Integral() << endl;
00145   prediction.ExtrapolateNDToFD(fHelper->BeamMatrixNuMuCCXSec());
00146   //  cout << "Extrapolated, integral: " << prediction.Spectrum()->Integral() << endl;
00147   prediction.Multiply(fFDFidMass);
00148   //  cout << "FD Fid Mass, integral: " << prediction.Spectrum()->Integral() << endl;
00149   prediction.Multiply(fFDCCData->GetPOT());
00150 //   cout << "FD POT= " << fFDCCData->GetPOT() << endl
00151 //        << "FD POT, integral final: " << prediction.Spectrum()->Integral() << endl;
00152 
00153   vector<NuMatrixSpectrum> predictions;
00154   predictions.push_back(prediction);
00155   return predictions;
00156 }

vector< NuMatrixSpectrum > NuMMRunNSINu::CalculateFDNCBackground ( const Double_t  ncScale = 1.0  )  const [private, virtual]

Definition at line 191 of file NuMMRunNSINu.cxx.

References NuMatrix1D::Add(), NuMMHelperPRL::BeamMatrixNuMuCCXSec(), NuMatrixSpectrum::Divide(), NuMatrixSpectrum::ExtrapolateNDToFD(), NuMMHelperPRL::FDCCContamination(), NuMMHelperPRL::FDCCContaminationRecoVsTrue(), NuMMHelperPRL::FDEfficiency(), NuMMHelperPRL::FDPurity(), NuMMHelperPRL::FDRecoVsTrue(), fFDCCData, NuMMRun::fFDFidMass, fHelper, fNDCCData, NuMMRun::fNDFidMass, NuMatrix::GetPOT(), NuMatrixSpectrum::Multiply(), NuMMHelperPRL::NDCCContamination(), NuMMHelperPRL::NDEfficiency(), NuMMHelperPRL::NDPurity(), NuMMHelperPRL::NDRecoVsTrue(), NuMatrixSpectrum::RecoToTrue(), NuMatrix1D::Subtract(), and NuMatrixSpectrum::TrueToReco().

Referenced by MakeFDPred().

00192 {
00193     NuMatrixSpectrum prediction(*fNDCCData);
00194     NuMatrixSpectrum wsBackground(prediction);
00195 
00196     NuMatrixSpectrum signal(prediction);
00197     signal.Multiply(fHelper->NDPurity());
00198     NuMatrixSpectrum ndncBackground(prediction);
00199     ndncBackground.Subtract(signal);
00200     ndncBackground.Multiply(ncScale);
00201     prediction.Subtract(ndncBackground);
00202     //    prediction.Multiply(fHelper->NDPurity());
00203 
00204     wsBackground.Multiply(fHelper->NDCCContamination());
00205     prediction.Subtract(wsBackground);
00206 
00207     prediction.RecoToTrue(fHelper->NDRecoVsTrue());
00208     prediction.Divide(fHelper->NDEfficiency());
00209     prediction.Divide(fNDCCData->GetPOT());
00210     prediction.Divide(fNDFidMass);
00211     prediction.ExtrapolateNDToFD(fHelper->BeamMatrixNuMuCCXSec());
00212     prediction.Multiply(fFDFidMass);
00213     prediction.Multiply(fFDCCData->GetPOT());
00214 
00215     NuMatrixSpectrum unoscWsBackground(prediction);
00216     unoscWsBackground.Multiply(fHelper->FDCCContamination());
00217     unoscWsBackground.TrueToReco(fHelper->FDCCContaminationRecoVsTrue());
00218 
00219     prediction.Multiply(fHelper->FDEfficiency());
00220   
00221     NuMatrixSpectrum unoscTrueSpectrum(prediction);
00222     unoscTrueSpectrum.TrueToReco(fHelper->FDRecoVsTrue());
00223     unoscTrueSpectrum.Add(unoscWsBackground);
00224 
00225     NuMatrixSpectrum ncBackground(unoscTrueSpectrum);
00226     ncBackground.Divide(fHelper->FDPurity());
00227     ncBackground.Subtract(unoscTrueSpectrum);
00228     
00229     vector<NuMatrixSpectrum> predictions;
00230     predictions.push_back(ncBackground);
00231     return predictions;
00232 }

vector< NuMatrixSpectrum > NuMMRunNSINu::CalculateFDTauFlux ( const Double_t  ncScale = 1.0  )  const [private, virtual]

Definition at line 160 of file NuMMRunNSINu.cxx.

References NuMMHelperPRL::BeamMatrixTauCCXSec(), NuMatrixSpectrum::Divide(), NuMatrixSpectrum::ExtrapolateNDToFD(), fFDCCData, NuMMRun::fFDFidMass, fHelper, fNDCCData, NuMMRun::fNDFidMass, NuMatrix::GetPOT(), NuMatrixSpectrum::Multiply(), NuMMHelperPRL::NDCCContamination(), NuMMHelperPRL::NDEfficiency(), NuMMHelperPRL::NDPurity(), NuMMHelperPRL::NDRecoVsTrue(), NuMatrixSpectrum::RecoToTrue(), and NuMatrix1D::Subtract().

Referenced by MakeFDPred().

00161 {
00162     NuMatrixSpectrum prediction(*fNDCCData);
00163     NuMatrixSpectrum wsBackground(prediction);
00164 
00165     NuMatrixSpectrum signal(prediction);
00166     signal.Multiply(fHelper->NDPurity());
00167     NuMatrixSpectrum ncBackground(prediction);
00168     ncBackground.Subtract(signal);
00169     ncBackground.Multiply(ncScale);
00170     prediction.Subtract(ncBackground);
00171     //    prediction.Multiply(fHelper->NDPurity());
00172 
00173     wsBackground.Multiply(fHelper->NDCCContamination());
00174     prediction.Subtract(wsBackground);
00175 
00176     prediction.RecoToTrue(fHelper->NDRecoVsTrue());
00177     prediction.Divide(fHelper->NDEfficiency());
00178     prediction.Divide(fNDCCData->GetPOT());
00179     prediction.Divide(fNDFidMass);
00180     prediction.ExtrapolateNDToFD(fHelper->BeamMatrixTauCCXSec());
00181     prediction.Multiply(fFDFidMass);
00182     prediction.Multiply(fFDCCData->GetPOT());
00183 
00184     vector<NuMatrixSpectrum> predictions;
00185     predictions.push_back(prediction);
00186     return predictions;
00187 }

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

Implements NuMMRun.

Definition at line 580 of file NuMMRunNSINu.cxx.

00581 {
00582   const vector<NuMatrixSpectrum> predictions =
00583     this->MakeFDPred(pars);
00584   Double_t like = 0.0;
00585   
00586   like += this->StatsLikelihood(predictions[0].Spectrum(),
00587                                 fFDCCData->Spectrum());
00588 
00589   return like;
00590 }

virtual const NuMatrixSpectrum* NuMMRunNSINu::GetFDDataNSI (  )  const [inline, virtual]

Get a const pointer to the FD data spectrum.

Definition at line 34 of file NuMMRunNSINu.h.

References fFDCCData.

00034 {return fFDCCData;};

virtual const NuMatrixSpectrum* NuMMRunNSINu::GetNDDataNSI (  )  const [inline, virtual]

Get a const pointer to the ND data spectrum.

Definition at line 32 of file NuMMRunNSINu.h.

References fNDCCData.

Referenced by NuFCExperimentFactoryNSI::GenerateND().

00032 {return fNDCCData;};

void NuMMRunNSINu::InverseOscillateNuTau ( NuMatrixSpectrum nutau,
const NuMMParameters pars 
) const [private]

Helper. Inverse oscillates a true energy nutau spectrum.

Definition at line 642 of file NuMMRunNSINu.cxx.

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

Referenced by MakeFDPred(), and TrueComponents().

00644 {
00645   switch(fDisappearanceModel){
00646   case 0:
00647     nutau.InverseOscillate(pars.Dm2(), pars.Sn2());
00648     return;
00649   case 1:
00650     nutau.DecayMuToTau(pars.Dm2(), pars.Sn2(), pars.Alpha());
00651     return;
00652   case 2:
00653     nutau.InverseDecohere(pars.Dm2(), pars.Sn2(), pars.Mu2());
00654     return; 
00655   case 3:
00656     nutau.InverseOscillateNSI(pars.Dm2(), pars.Sn2(), pars.Epsilon(), 1.0);
00657     return;
00658   default:
00659     assert(0 && "Badly configured disappearance model");
00660   }
00661 }

void NuMMRunNSINu::InverseOscillateTauBar ( NuMatrixSpectrum nutau,
const NuMMParameters pars 
) const [private]

Definition at line 663 of file NuMMRunNSINu.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().

00665 {
00666   switch(fDisappearanceModel){
00667   case 0:
00668     nutau.InverseOscillate(pars.Dm2Bar(), pars.Sn2Bar());
00669     return;
00670   case 1:
00671     nutau.DecayMuToTau(pars.Dm2Bar(), pars.Sn2Bar(), pars.Alpha());
00672     return;
00673   case 2:
00674     nutau.InverseDecohere(pars.Dm2Bar(), pars.Sn2Bar(), pars.Mu2());
00675     return; 
00676   case 3:
00677     // Oscillate nubar with dm2 and sn2 too. The only difference is the 
00678     // sign of epsilon
00679     nutau.InverseOscillateNSI(pars.Dm2(), pars.Sn2(), pars.Epsilon(), -1.0);
00680     return;
00681   default:
00682     assert(0 && "Badly configured disappearance model");
00683   }
00684 }

vector< NuMatrixSpectrum > NuMMRunNSINu::MakeFDPred ( const NuMMParameters pars  )  [virtual]

Definition at line 236 of file NuMMRunNSINu.cxx.

References CalculateFDFlux(), CalculateFDNCBackground(), CalculateFDTauFlux(), NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), NuMMParameters::Epsilon(), NuMMHelperPRL::FDCCContamination(), NuMMHelperPRL::FDCCContaminationRecoVsTrue(), NuMMHelperPRL::FDEfficiency(), NuMMHelperPRL::FDRecoVsTrue(), NuMMHelperPRL::FDTauEfficiency(), NuMMHelperPRL::FDTauRecoVsTrue(), fFDCCData, fFDFlux, fFDNCBackground, fFDTauFlux, fHelper, NuMMRun::fQuietMode, fSpecialWSBackgroundFlux, fSpecialWSMode, InverseOscillateNuTau(), NuMatrixSpectrum::Multiply(), NuMMParameters::NCBackgroundScale(), NuMMParameters::Normalisation(), OscillateNuBar(), OscillateNuMu(), NuMMParameters::PrintStatus(), NuMatrixSpectrum::RemoveErrors(), NuMatrixSpectrum::SetName(), NuMatrixSpectrum::SetTitle(), NuMMParameters::ShwEnScale(), NuMMParameters::Sn2(), NuMMParameters::Sn2Bar(), NuMMParameters::TrkEnScale(), NuMatrixSpectrum::TrueToReco(), and WriteToJess().

00237 {
00238   if (!fQuietMode) {
00239     cout << "sn2: " << pars.Sn2() << "; dm2: " << pars.Dm2() 
00240          << "; sn2bar: " << pars.Sn2Bar() << "; dm2bar: " << pars.Dm2Bar()
00241          << "; epsilon: "<< pars.Epsilon()
00242          << endl
00243          << "; norm: " << pars.Normalisation() 
00244          << "; NCBack: " << pars.NCBackgroundScale()
00245          << "; TrkEn: " << pars.TrkEnScale()
00246          << "; ShwEn: " << pars.ShwEnScale() << endl;
00247     pars.PrintStatus();
00248   }
00249 
00250   fFDFlux = this->CalculateFDFlux(pars.NCBackgroundScale());
00251   fFDTauFlux = this->CalculateFDTauFlux(pars.NCBackgroundScale());
00252   fFDNCBackground = this->CalculateFDNCBackground(pars.NCBackgroundScale());
00253   
00254   //Pick up the pre-calculated stuff
00255   // TODO CJB - the destruction of these histograms takes ages
00256   NuMatrixSpectrum nuPrediction = fFDFlux[0];  
00257   NuMatrixSpectrum potentialTaus = fFDTauFlux[0];
00258   NuMatrixSpectrum wsBackground = fFDFlux[0];
00259   NuMatrixSpectrum ncBackground = fFDNCBackground[0];
00260 
00261   //Do extrapolation
00262   potentialTaus.Multiply(fHelper->FDTauEfficiency());
00263   
00264   wsBackground.Multiply(fHelper->FDCCContamination());
00265   if (fSpecialWSMode){
00266     wsBackground = fSpecialWSBackgroundFlux[0];
00267   }
00268 
00269   WriteToJess(wsBackground,"wsBackground", false);
00270   OscillateNuBar(wsBackground, pars);
00271   WriteToJess(wsBackground,"wsBackground", true);
00272   wsBackground.TrueToReco(fHelper->FDCCContaminationRecoVsTrue());
00273 
00274   nuPrediction.Multiply(fHelper->FDEfficiency());
00275   
00276   WriteToJess(nuPrediction,"nuPrediction", false);
00277   //  cout << "nuPrediction, before: " <<nuPrediction.Spectrum()->Integral() << endl;
00278   OscillateNuMu(nuPrediction, pars);
00279   //  cout << "nuPrediction, after: " <<nuPrediction.Spectrum()->Integral() << endl;
00280   WriteToJess(nuPrediction,"nuPrediction", true);
00281 
00282   WriteToJess(potentialTaus,"potentialTaus", false);
00283   InverseOscillateNuTau(potentialTaus, pars);
00284   WriteToJess(potentialTaus,"potentialTaus", true);
00285   
00286   potentialTaus.TrueToReco(fHelper->FDTauRecoVsTrue());
00287   nuPrediction.TrueToReco(fHelper->FDRecoVsTrue());
00288 
00289   ncBackground.Multiply(pars.NCBackgroundScale());
00290 
00291   WriteToJess(ncBackground,"ncBackground", false);
00292   // We don't oscillate the NC background. Which means it doesn't get
00293   // decayed or decohered either. Probably OK since it's a tiny component,
00294   // but not exactly ideal.
00295 
00296 
00297 
00298   // Get systematic shifts
00299   // New shift (pred/nominal), the histos are already around 1.0, don't add one
00300 
00301   //  TFile *fSyst = new TFile("/minos/app/zisvan/NSI_kNN/Systematics/fhc_syst_shifts_py.root");
00302   TFile *fSyst = new TFile("$SRT_PRIVATE_CONTEXT/NuMuBar/macros/zeynep_nsi/Systematics/fhc_syst_shifts_py.root");
00303 
00304   TH1D *trackBand = (TH1D*)fSyst->Get("fdPredTrackEnergyOverall");
00305   trackBand->Scale(pars.TrkEnScale());
00306   for(Int_t i=0; i<trackBand->GetNbinsX()+1; ++i){
00307     trackBand->SetBinContent(i, trackBand->GetBinContent(i)+1);
00308   }
00309   nuPrediction.Multiply(trackBand);
00310 
00311   TH1D *showerBand = (TH1D*)fSyst->Get("fdPredShowerEnergyScaleFunctionBoth");
00312   showerBand->Scale(pars.ShwEnScale());
00313   for(Int_t i=0; i<showerBand->GetNbinsX()+1; ++i){
00314     showerBand->SetBinContent(i, showerBand->GetBinContent(i)+1);
00315   }
00316   nuPrediction.Multiply(showerBand);
00317 
00318   TH1D *ncBand = (TH1D*)fSyst->Get("fdPredNCBackground");
00319   ncBand->Scale(pars.NCBackgroundScale());
00320   for(Int_t i=0; i<ncBand->GetNbinsX()+1; ++i){
00321     ncBand->SetBinContent(i, ncBand->GetBinContent(i)+1);
00322   }
00323 
00324   nuPrediction.Multiply(ncBand);
00325 
00326   TH1D *trackBandWS = (TH1D*)fSyst->Get("wsBackgroundTrackEnergyOverall");
00327   trackBandWS->Scale(pars.TrkEnScale());
00328   for(Int_t i=0; i<trackBandWS->GetNbinsX()+1; ++i){
00329     trackBandWS->SetBinContent(i, trackBandWS->GetBinContent(i)+1);
00330   }
00331   wsBackground.Multiply(trackBandWS);
00332 
00333   TH1D *showerBandWS = (TH1D*)fSyst->Get("wsBackgroundShowerEnergyScaleFunctionBoth");
00334   showerBandWS->Scale(pars.ShwEnScale());
00335   for(Int_t i=0; i<showerBandWS->GetNbinsX()+1; ++i){
00336     showerBandWS->SetBinContent(i, showerBandWS->GetBinContent(i)+1);
00337   }
00338   wsBackground.Multiply(showerBandWS);
00339 
00340   TH1D *ncBandWS = (TH1D*)fSyst->Get("wsBackgroundNCBackground");
00341   ncBandWS->Scale(pars.NCBackgroundScale());
00342   for(Int_t i=0; i<ncBandWS->GetNbinsX()+1; ++i){
00343     ncBandWS->SetBinContent(i, ncBandWS->GetBinContent(i)+1);
00344   }
00345   wsBackground.Multiply(ncBandWS);
00346 
00347 
00348   fSyst->Close();
00349 //   // Shifting signal
00350 //   TFile *fSyst = new TFile("/minos/app/zisvan/NSI_kNN/Systematics/fhc_syst_shifts.root");
00351 
00352 //   TH1D *trackBand = (TH1D*)fSyst->Get("hLoTrk");
00353 //   trackBand->Scale(pars.TrkEnScale());
00354 //   for(Int_t i=0; i<trackBand->GetNbinsX()+1; ++i){
00355 //     trackBand->SetBinContent(i, trackBand->GetBinContent(i)+1);
00356 //   }
00357 //   nuPrediction.Multiply(trackBand);
00358 
00359 //   TH1D *showerBand = (TH1D*)fSyst->Get("hUpShw");
00360 //   showerBand->Scale(pars.ShwEnScale());
00361 //   for(Int_t i=0; i<showerBand->GetNbinsX()+1; ++i){
00362 //     showerBand->SetBinContent(i, showerBand->GetBinContent(i)+1);
00363 //   }
00364 //   //nuPrediction.Multiply(showerBand);
00365   
00366 //   TH1D *ncBand = (TH1D*)fSyst->Get("hUpNC");
00367 //   ncBand->Scale(pars.NCBackgroundScale());
00368 //   for(Int_t i=0; i<ncBand->GetNbinsX()+1; ++i){
00369 //     ncBand->SetBinContent(i, ncBand->GetBinContent(i)+1);
00370 //   }
00371 //   //nuPrediction.Multiply(ncBand);
00372 
00373 //   fSyst->Close();
00374 
00375 //   // Shifting WS background
00376 //   TFile *fSystWS = new TFile("/minos/app/zisvan/NSI_kNN/Systematics/fhc_syst_shifts_ws.root");
00377 
00378 //   TH1D *trackBandWS = (TH1D*)fSystWS->Get("hUpTrk");
00379 //   trackBandWS->Scale(pars.TrkEnScale());
00380 //   for(Int_t i=0; i<trackBandWS->GetNbinsX()+1; ++i){
00381 //     trackBandWS->SetBinContent(i, trackBandWS->GetBinContent(i)+1);
00382 //   }
00383 //   wsBackground.Multiply(trackBandWS);
00384 
00385 //   TH1D *showerBandWS = (TH1D*)fSystWS->Get("hUpShw");
00386 //   showerBandWS->Scale(pars.ShwEnScale());
00387 //   for(Int_t i=0; i<showerBandWS->GetNbinsX()+1; ++i){
00388 //     showerBandWS->SetBinContent(i, showerBandWS->GetBinContent(i)+1);
00389 //   }
00390 //   wsBackground.Multiply(showerBandWS);
00391 
00392 //   TH1D *ncBandWS = (TH1D*)fSystWS->Get("hUpNC");
00393 //   ncBandWS->Scale(pars.NCBackgroundScale());
00394 //   for(Int_t i=0; i<ncBandWS->GetNbinsX()+1; ++i){
00395 //     ncBandWS->SetBinContent(i, ncBandWS->GetBinContent(i)+1);
00396 //   }
00397 //   wsBackground.Multiply(ncBandWS);
00398 
00399 
00400 //   fSystWS->Close();
00401 
00402 //   //Done shifting//
00403 
00404 
00405   nuPrediction.Add(ncBackground);
00406   nuPrediction.Add(potentialTaus);
00407   nuPrediction.Add(wsBackground);
00408 
00409   nuPrediction.Multiply(pars.Normalisation());
00410 
00411   WriteToJess(*fFDCCData, "data", false);
00412 
00413   //Prediction is done. Now fix the stupid error bars and give it a
00414   //better name
00415   nuPrediction.RemoveErrors();
00416   nuPrediction.SetName("fdPred");  
00417   nuPrediction.SetTitle("fdPred");  
00418   vector<NuMatrixSpectrum> ccPredictions;
00419   ccPredictions.push_back(nuPrediction);
00420 
00421   //Return these as well to be used with Feldman-Cousins sampling for NSI
00422   wsBackground.RemoveErrors();
00423   wsBackground.SetName("wsBackground");  
00424   wsBackground.SetTitle("wsBackground");  
00425 
00426   ncBackground.RemoveErrors();
00427   ncBackground.SetName("ncBackground");  
00428   ncBackground.SetTitle("ncBackground");  
00429 
00430   potentialTaus.RemoveErrors();
00431   potentialTaus.SetName("potentialTaus");  
00432   potentialTaus.SetTitle("potentialTaus");  
00433 
00434   ccPredictions.push_back(wsBackground);
00435   ccPredictions.push_back(ncBackground);
00436   ccPredictions.push_back(potentialTaus);
00437 
00438   return ccPredictions;
00439 }

void NuMMRunNSINu::OscillateNuBar ( NuMatrixSpectrum numu,
const NuMMParameters pars 
) const [private]

Definition at line 617 of file NuMMRunNSINu.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().

Referenced by MakeFDPred(), and TrueComponents().

00619 {
00620   switch(fDisappearanceModel){
00621   case 0:
00622     numu.Oscillate(pars.Dm2Bar(), pars.Sn2Bar());
00623     return;
00624   case 1:
00625     numu.DecayCC(pars.Dm2Bar(), pars.Sn2Bar(), pars.Alpha());
00626     return;
00627   case 2:
00628     numu.Decohere(pars.Dm2Bar(), pars.Sn2Bar(), pars.Mu2());
00629     return;
00630   case 3: 
00631     // Oscillate nubar with dm2 and sn2 too. The only difference is the 
00632     // sign of epsilon
00633     numu.OscillateNSI(pars.Dm2(), pars.Sn2(), pars.Epsilon(), -1.0);
00634     return;
00635   default:
00636     assert(0 && "Badly configured disappearance model");
00637   }
00638 }

void NuMMRunNSINu::OscillateNuMu ( NuMatrixSpectrum numu,
const NuMMParameters pars 
) const [private]

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

Definition at line 594 of file NuMMRunNSINu.cxx.

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

Referenced by MakeFDPred(), and TrueComponents().

00597 {
00598   switch(fDisappearanceModel){
00599   case 0:
00600     numu.Oscillate(pars.Dm2(), pars.Sn2());
00601     return;
00602   case 1:
00603     numu.DecayCC(pars.Dm2(), pars.Sn2(), pars.Alpha());
00604     return;
00605   case 2:
00606     numu.Decohere(pars.Dm2(), pars.Sn2(), pars.Mu2());
00607     return;
00608   case 3: 
00609     numu.OscillateNSI(pars.Dm2(), pars.Sn2(), pars.Epsilon(), 1.0);
00610     return;
00611   default:
00612     assert(0 && "Badly configured disappearance model");
00613   }
00614 }

void NuMMRunNSINu::PrecalculateWSFlux ( NuMMHelperPRL helper,
const NuMatrixSpectrum ndOtherData 
) [virtual]

Definition at line 552 of file NuMMRunNSINu.cxx.

References NuMMHelperPRL::BeamMatrixNuMuCCXSec(), NuMatrixSpectrum::Divide(), NuMatrixSpectrum::ExtrapolateNDToFD(), NuMMHelperPRL::FDEfficiencyOtherChargeWS(), fFDCCData, NuMMRun::fFDFidMass, NuMMRun::fNDFidMass, fSpecialWSBackgroundFlux, fSpecialWSMode, NuMatrix::GetPOT(), NuMatrixSpectrum::Multiply(), NuMMHelperPRL::NDCCContamination(), NuMMHelperPRL::NDEfficiency(), NuMMHelperPRL::NDPurity(), NuMMHelperPRL::NDRecoVsTrue(), NuMatrixSpectrum::RecoToTrue(), SanityCheckSpecialWSFiles(), and NuMatrix1D::Subtract().

00554 {
00555   cout << "Running PrecalculateWSFlux" << endl;
00556   this->SanityCheckSpecialWSFiles(helper,ndOtherData);
00557   fSpecialWSMode = true;
00558   fSpecialWSBackgroundFlux.clear();
00559   NuMatrixSpectrum prediction(*ndOtherData);
00560   NuMatrixSpectrum wsBackground(prediction);
00561   prediction.Multiply(helper->NDPurity());
00562   
00563   wsBackground.Multiply(helper->NDCCContamination());
00564   prediction.Subtract(wsBackground);
00565   
00566   prediction.RecoToTrue(helper->NDRecoVsTrue());
00567   prediction.Divide(helper->NDEfficiency());
00568   prediction.Divide(ndOtherData->GetPOT());
00569   prediction.Divide(fNDFidMass);
00570   prediction.ExtrapolateNDToFD(helper->BeamMatrixNuMuCCXSec());
00571   prediction.Multiply(fFDFidMass);
00572   prediction.Multiply(fFDCCData->GetPOT());
00573   prediction.Multiply(helper->FDEfficiencyOtherChargeWS());
00574   
00575   fSpecialWSBackgroundFlux.push_back(prediction);
00576 }

void NuMMRunNSINu::SanityCheckSpecialWSFiles ( NuMMHelperPRL helper,
const NuMatrixSpectrum ndOtherData 
) const [private, virtual]

Definition at line 497 of file NuMMRunNSINu.cxx.

References NuMMHelperPRL::NDEfficiency(), and NuMatrixSpectrum::Spectrum().

Referenced by PrecalculateWSFlux().

00499 {
00500 
00501   string dataName = fNDCCData->Spectrum()->GetName();
00502   string otherDataName = ndOtherData->Spectrum()->GetName();
00503   string helperEffName = helper->NDEfficiency()->GetName();
00504 
00505   if ("RecoEnergy_ND" == dataName){
00506     if ("RecoEnergyPQ_ND" == otherDataName){
00507       cout << "You are extrapolating NuMu data. "
00508            << "WS histograms sanity checked as OK."
00509            << endl;
00510     }
00511     else{
00512       cout << "You are extrapolating NuMu data. "
00513            << "WS histograms failing sanity check."
00514            << endl;
00515       assert(false);
00516     }
00517     if ("EfficiencyPQ_ND" == helperEffName){
00518       cout << "WS helper sanity checked as OK."
00519            << endl;
00520     }
00521     else{
00522       cout << "WS helper failing sanity check."
00523            << endl;
00524       assert(false);
00525     }
00526   }
00527   else if ("RecoEnergyPQ_ND" == dataName){
00528     if ("RecoEnergy_ND" == otherDataName){
00529       cout << "You are extrapolating NuMu data. "
00530            << "WS histograms sanity checked as OK."
00531            << endl;
00532     }
00533     else{
00534       cout << "You are extrapolating NuMu data. "
00535            << "WS histograms failing sanity check."
00536            << endl;
00537       assert(false);
00538     }
00539     if ("Efficiency_ND" == helperEffName){
00540       cout << "WS helper sanity checked as OK."
00541            << endl;
00542     }
00543     else{
00544       cout << "WS helper failing sanity check."
00545            << endl;
00546       assert(false);
00547     }
00548   }
00549 }

virtual void NuMMRunNSINu::SetFDData ( NuMatrixSpectrum fdData  )  [inline, virtual]

Definition at line 40 of file NuMMRunNSINu.h.

References fFDCCData.

00040 {fFDCCData = fdData;}

NuMatrixSpectrum NuMMRunNSINu::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 441 of file NuMMRunNSINu.cxx.

References NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), NuMMParameters::Epsilon(), NuMMHelperPRL::FDCCContamination(), NuMMHelperPRL::FDEfficiency(), fFDFlux, fFDNCBackground, fFDTauFlux, fHelper, NuMMRun::fQuietMode, InverseOscillateNuTau(), Sample::kNCBackNSI, Sample::kSignalNSI, Sample::kTauNSI, Msg::kWarning, Sample::kWSBackNSI, MSG, NuMatrixSpectrum::Multiply(), NuMMParameters::NCBackgroundScale(), NuMMParameters::Normalisation(), OscillateNuBar(), OscillateNuMu(), NuMatrixSpectrum::SetName(), NuMMParameters::ShwEnScale(), NuMMParameters::Sn2(), NuMMParameters::Sn2Bar(), and NuMatrixSpectrum::Spectrum().

Referenced by NuFCExperimentFactoryNSI::GenerateNewExperiment().

00442 {
00443   // Output the oscillation settings!
00444    if (!fQuietMode) {
00445     cout << "sn2: " << pars.Sn2() << "; dm2: " << pars.Dm2() 
00446          << "; sn2bar: " << pars.Sn2Bar() << "; dm2bar: " << pars.Dm2Bar()
00447          << "; epsilon: " << pars.Epsilon() << endl
00448          << "; norm: " << pars.Normalisation() 
00449          << "; NCBack: " << pars.NCBackgroundScale()
00450          << "; ShwEn: " << pars.ShwEnScale() << endl;
00451 
00452     }
00453 
00454   if (s == kSignalNSI) {
00455     NuMatrixSpectrum prediction = fFDFlux[0];
00456     prediction.Multiply(fHelper->FDEfficiency());
00457     //    cout << s << ", before osc. " << prediction.Spectrum()->Integral() << endl;
00458     OscillateNuMu(prediction,pars);
00459     //    cout << s << ", after osc. " << prediction.Spectrum()->Integral() << endl;
00460     return prediction;
00461   }
00462 
00463   else if (s == kWSBackNSI) {
00464     NuMatrixSpectrum wrongSign = fFDFlux[0];
00465     wrongSign.Multiply(fHelper->FDCCContamination());
00466     //    cout << s << ", before osc. " << wrongSign.Spectrum()->Integral() << endl;
00467     OscillateNuBar(wrongSign,pars);
00468     //    cout << s << ", after osc. " << wrongSign.Spectrum()->Integral() << endl;
00469     return wrongSign;
00470   }
00471   else if (s == kNCBackNSI) {
00472     NuMatrixSpectrum ncBackground = fFDNCBackground[0];
00473     //    cout << s << ", not oscillating " << ncBackground.Spectrum()->Integral() << endl;
00474     return ncBackground;
00475   }
00476 
00477   else if (s == kTauNSI) {
00478     NuMatrixSpectrum potentialTaus = fFDTauFlux[0];
00479     //    cout << s << ", before osc. " << potentialTaus.Spectrum()->Integral() << endl;    
00480     InverseOscillateNuTau(potentialTaus,pars);
00481     //    cout << s << ", after osc. " << potentialTaus.Spectrum()->Integral() << endl;    
00482     return potentialTaus;
00483   }
00484 
00485   else {
00486     MSG("NuMMRunNSINu",Msg::kWarning) << "Requested unknown prediction #" << s << ", returning blank." << endl;
00487     NuMatrixSpectrum blank = fFDFlux[0];
00488     blank.SetName("Blank");
00489     blank.Spectrum()->Scale(0);
00490     return blank;
00491   }
00492 }

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

Implements NuMMRun.

Definition at line 42 of file NuMMRunNSINu.cxx.

References NuMatrix1D::Add(), NuMatrixSpectrum::Multiply(), NuMatrixSpectrum::Spectrum(), and NuMatrixSpectrum::TrueToReco().

00043 {
00044   vector<TH1D> vHistos;
00045 
00046   NuMatrixSpectrum nuPrediction = fFDFlux[0];
00047   
00048   TH1D hNDData(*(fNDCCData->Spectrum()));
00049   hNDData.SetName("ndData");
00050   hNDData.SetTitle("ndData");
00051   vHistos.push_back(*(new TH1D(hNDData)));
00052   
00053   NuMatrixSpectrum potentialTaus(fFDTauFlux[0]);
00054   potentialTaus.Multiply(fHelper->FDTauEfficiency());
00055 
00056   nuPrediction.Multiply(fHelper->FDEfficiency());
00057 
00058   NuMatrixSpectrum ncBackground = fFDNCBackground[0];
00059   
00060   //oscillate right signs as nus
00061   OscillateNuMu(nuPrediction, pars);
00062   InverseOscillateNuTau(potentialTaus, pars);
00063   
00064   //Output the prediction in true energy
00065   TH1D hTrueFDPrediction(*(nuPrediction.Spectrum()));
00066   hTrueFDPrediction.SetName("fdTruePrediction");
00067   vHistos.push_back(hTrueFDPrediction);
00068   
00069   potentialTaus.TrueToReco(fHelper->FDTauRecoVsTrue());  
00070   nuPrediction.TrueToReco(fHelper->FDRecoVsTrue());
00071   
00072   //Get the prediction in reco energy
00073   TH1D hBaseFDPrediction(*(nuPrediction.Spectrum()));
00074   hBaseFDPrediction.SetName("fdBasePrediction");
00075   vHistos.push_back(hBaseFDPrediction);
00076   
00077   // Add in the backgrounds: NC
00078   nuPrediction.Add(ncBackground);
00079   // tau
00080   if (0==fDisappearanceModel || 3==fDisappearanceModel){
00081     nuPrediction.Add(potentialTaus);
00082   }
00083   // WS
00084   //oscillate ws as nubars
00085   NuMatrixSpectrum wsBackground = fFDFlux[0];
00086   wsBackground.Multiply(fHelper->FDCCContamination());
00087   if(fSpecialWSMode) wsBackground = fSpecialWSBackgroundFlux[0];
00088   OscillateNuBar(wsBackground, pars);
00089   wsBackground.TrueToReco(fHelper->FDCCContaminationRecoVsTrue());
00090   nuPrediction.Add(wsBackground);
00091   
00092   //Get the prediction with backgrounds
00093   TH1D hFDPrediction(*(nuPrediction.Spectrum()));  
00094   hFDPrediction.SetName("fdPrediction");
00095   vHistos.push_back(hFDPrediction);
00096   
00097   TH1D hFDTaus(*(potentialTaus.Spectrum()));
00098   hFDTaus.SetName("fdTaus");
00099   vHistos.push_back(hFDTaus);
00100   
00101   TH1D hFDNCBackground(*(ncBackground.Spectrum()));
00102   hFDNCBackground.SetName("fdNCBackground");
00103   vHistos.push_back(hFDNCBackground);
00104 
00105   TH1D hFDWSBackground(*(wsBackground.Spectrum()));
00106   hFDWSBackground.SetName("fdWSBackground");
00107   vHistos.push_back(hFDWSBackground);
00108   
00109   TH1D hfdData(*(fFDCCData->Spectrum()));
00110   hfdData.SetName("fdData");
00111   vHistos.push_back(hfdData);
00112 
00113   return vHistos;
00114 }

void NuMMRunNSINu::WriteToJess ( NuMatrixSpectrum  spectrum,
const std::string  name,
bool  osc 
) [private, virtual]

Definition at line 687 of file NuMMRunNSINu.cxx.

References NuMMHelperPRL::FDCCContaminationRecoVsTrue(), NuMMHelperPRL::FDRecoVsTrue(), NuMMHelperPRL::FDTauRecoVsTrue(), fHelper, NuMMRun::fJessFile, Msg::kError, NuMatrixSpectrum::MakeJessTemplate(), MSG, NuMatrixSpectrum::RemoveErrors(), and NuMatrixSpectrum::Spectrum().

Referenced by MakeFDPred().

00690 {
00691   // If they don't want the plots don't try to give them
00692   if (!fJessFile) return;
00693 
00694   const string opName = osc ? (name+"_osc") : name;
00695 
00696   fJessFile->cd();
00697   if ("wsBackground" == name){
00698     TH2D jessHist =
00699       spectrum.MakeJessTemplate(fHelper->FDCCContaminationRecoVsTrue());
00700     jessHist.Write(opName.c_str());
00701   }
00702   else if ("nuPrediction" == name){
00703     TH2D jessHist =
00704       spectrum.MakeJessTemplate(fHelper->FDRecoVsTrue());
00705     jessHist.Write(opName.c_str());
00706   }
00707   else if ("potentialTaus" == name){
00708     TH2D jessHist =
00709       spectrum.MakeJessTemplate(fHelper->FDTauRecoVsTrue());
00710     jessHist.Write(opName.c_str());
00711   }
00712   else if ("ncBackground" == name){
00713     spectrum.RemoveErrors();
00714     spectrum.Spectrum()->Write(opName.c_str());
00715   }
00716   else if ("data" == name){
00717     spectrum.Spectrum()->Write(opName.c_str());
00718   }
00719   else{
00720     MSG("NuMatrixSpectrum", Msg::kError)
00721       << "You're trying to write a Jess template which has an "
00722       << "unrecognized name"
00723       << endl;
00724   }
00725 }


Member Data Documentation

Definition at line 70 of file NuMMRunNSINu.h.

Referenced by MakeFDPred(), and TrueComponents().

Definition at line 72 of file NuMMRunNSINu.h.

Referenced by MakeFDPred(), and TrueComponents().

Definition at line 71 of file NuMMRunNSINu.h.

Referenced by MakeFDPred(), and TrueComponents().

Definition at line 73 of file NuMMRunNSINu.h.

Referenced by MakeFDPred(), and PrecalculateWSFlux().

Bool_t NuMMRunNSINu::fSpecialWSMode [private]

Definition at line 66 of file NuMMRunNSINu.h.

Referenced by MakeFDPred(), and PrecalculateWSFlux().


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

Generated on 22 Nov 2017 for loon by  doxygen 1.6.1