NuMatrix1D Class Reference

#include <NuMatrix1D.h>

Inheritance diagram for NuMatrix1D:
NuMatrix NuMatrixSpectrum Rock1DMatrix

List of all members.

Public Member Functions

 NuMatrix1D ()
 NuMatrix1D (const TH1D &spectrum)
 The spectrum-only constructor.
 NuMatrix1D (Double_t POT)
 The POT-only constructor.
 NuMatrix1D (const TH1D &spectrum, Double_t POT)
 The spectrum *and* POT constructor.
 NuMatrix1D (const NuMatrix1D &original)
NuMatrix1Doperator= (const NuMatrix1D &orig)
 NuMatrix1D (const NuMatrixSpectrum &original)
 Copies a NuMatrixSpectrum.
NuMatrix1Doperator= (const NuMatrixSpectrum &orig)
 Assignment from a NuMatrixSpectrum.
 NuMatrix1D (const TString &filename, const TString &histname, Double_t POT)
 NuMatrix1D (const TString &filename, const TString &histname, const TString &pothistname="hTotalPot")
virtual ~NuMatrix1D ()
virtual NuMatrix1DCopy (void) const
 Makes a deep-object copy.
virtual Double_t StatsLikelihood (const NuMatrix *with) const
 Checks the passed object is a NuMatrix1D then returns the log-likelihood difference.
virtual Double_t StatsLikelihood (const NuMatrix1D *with) const
 Compares two 2D histograms and returns the log-likelihood difference.
NuMatrixSpectrum GetNuMatrixSpectrum (void) const
 Returns a pointer to the built in data storage.
Bool_t Complete () const
 Returns true if the spectrum is properly configured.
virtual void Divide (const NuMatrix &by)
virtual void Divide (const NuMatrix1D &by)
 Make a correction by Dividing.
virtual void Divide (const TH1 &correction, const Option_t *option="")
virtual void Divide (const TGraph &correction)
virtual void Multiply (const NuMatrix &by)
virtual void Multiply (const NuMatrix1D &by)
 Make a correction by Multiplying.
virtual void Multiply (const TH1 &correction, const Option_t *option="")
virtual void Multiply (const TGraph &correction)
virtual void Add (const TH1 *correction)
 Functions for addition and subtraction.
virtual void Add (const NuMatrix &correction, Bool_t addpot=false)
virtual void Add (const NuMatrix1D &correction, Bool_t addpot=false)
 Functions for addition and subtraction.
virtual void Subtract (const TH1 *correction)
virtual void Subtract (const NuMatrix &correction)
virtual void Subtract (const NuMatrix1D &correction)
virtual void RecoToTrue (const TH2D &correction)
virtual void TrueToReco (const TH2D &correction)
virtual void ExtrapolateNDToFD (const TH2D &beammatrix)
virtual void OscillatedNDRecoToTrue (const TH2D &unnormalizedCorrection, const NuMMParameters &oscpars)
virtual void ScaleToPOT (Double_t new_POT)
virtual void Scale (Double_t value)
 Scales the histogram.
virtual void SetValue (Double_t val)
 Sets the entire contents of the spectrum to a single value.
virtual void Oscillate (const Double_t dm2, const Double_t sn2)
 Oscillation function.
virtual void DecayCC (const Double_t dm2, const Double_t sn2, const Double_t alpha)
virtual void DecayNC (const Double_t sn2, const Double_t alpha)
virtual void DecayMuToTau (const Double_t dm2, const Double_t sn2, const Double_t alpha)
virtual void Decohere (const Double_t dm2, const Double_t sn2, const Double_t mu2)
virtual void InverseDecohere (const Double_t dm2, const Double_t sn2, const Double_t mu2)
 Decohered neutrinos aren't gone, they have to turn up somewhere, as taus.
virtual void OscillateNSI (const Double_t dm2, const Double_t sn2, const Double_t epsilon, const Double_t sign)
virtual void InverseOscillateNSI (const Double_t dm2, const Double_t sn2, const Double_t epsilon, const Double_t sign)
virtual Int_t Write (const TString &name) const
 Writes to the current file.
virtual void Read (TFile *f, const TString &name)
virtual TH1 * Spectrum (void) const
 Returns a TH1 object representing the data the spectrum holds. Or zero.
virtual void SetSpectrum (const TH1D &spectrum)
virtual void Draw (Option_t *option="")
virtual void Fill (const NuEvent &evt)
 Stores a NuEvent in the spectrum.
virtual void Interpolate ()
virtual Double_t GetEventDensity (Double_t x)
virtual Double_t GetXmin ()
virtual Double_t GetXmax ()
virtual TH1D * Rebin (Int_t nbins, Double_t *edges)
virtual Double_t GetIntegral (Double_t x)

Protected Member Functions

virtual void XToY (const TH2D &correction)
virtual void YToX (const TH2D &correction)
void FillOscCache (double dm2) const
 ClassDef (NuMatrix1D, 0)

Protected Attributes

TH1D * fSpectrum
 The storage for the 1D spectrum.
TSpline5 * fCumulativeDist

Static Protected Attributes

static double fgCachedDmsq = -123456789
 The value of dmsq that fgCachedOsc was generated at.
static NuMatrix1DfgCachedOsc = 0
 The oscillation dip calculated at dmsq = fgCachedDmsq, sinsq = 1.

Detailed Description

1D Data representation class. This class represents a 1D histogram (TH1D) and an associated POT.

Definition at line 18 of file NuMatrix1D.h.


Constructor & Destructor Documentation

NuMatrix1D::NuMatrix1D (  ) 

Definition at line 34 of file NuMatrix1D.cxx.

Referenced by Copy(), and FillOscCache().

00034                        : NuMatrix(), fSpectrum(0), fCumulativeDist(0)
00035 {
00036 
00037 }

NuMatrix1D::NuMatrix1D ( const TH1D &  spectrum  ) 

The spectrum-only constructor.

Definition at line 46 of file NuMatrix1D.cxx.

References fSpectrum.

00046                                            : NuMatrix(0), fSpectrum(0), fCumulativeDist(0)
00047 {
00048   fSpectrum = new TH1D(spectrum);
00049   fSpectrum->SetDirectory(0);
00050 }

NuMatrix1D::NuMatrix1D ( Double_t  POT  ) 

The POT-only constructor.

Definition at line 40 of file NuMatrix1D.cxx.

00040                                    : NuMatrix(POT), fSpectrum(0), fCumulativeDist(0)
00041 {
00042 
00043 }

NuMatrix1D::NuMatrix1D ( const TH1D &  spectrum,
Double_t  POT 
)

The spectrum *and* POT constructor.

Definition at line 53 of file NuMatrix1D.cxx.

References fSpectrum.

00053                                                          : NuMatrix(POT), fSpectrum(0), fCumulativeDist(0)
00054 {
00055   fSpectrum = new TH1D(spectrum);
00056   fSpectrum->SetDirectory(0);
00057 }

NuMatrix1D::NuMatrix1D ( const NuMatrix1D original  ) 

Definition at line 60 of file NuMatrix1D.cxx.

References fCumulativeDist, and fSpectrum.

00060                                                  : NuMatrix(original), fSpectrum(0), fCumulativeDist(0)
00061 {
00062   if (original.fSpectrum) {
00063     fSpectrum = new TH1D(*original.fSpectrum);
00064     fSpectrum->SetDirectory(0);
00065   }
00066   if (original.fCumulativeDist) {
00067     fCumulativeDist = new TSpline5(*original.fCumulativeDist);
00068   }
00069 }

NuMatrix1D::NuMatrix1D ( const NuMatrixSpectrum original  ) 

Copies a NuMatrixSpectrum.

Definition at line 96 of file NuMatrix1D.cxx.

References fCumulativeDist, fSpectrum, NuMatrix::GetPOT(), NuMatrix::ResetPOT(), and NuMatrixSpectrum::Spectrum().

00096                                                        : fSpectrum(0), fCumulativeDist(0)
00097 {
00098   // Copy the POT from the NuMatrixSpectrum
00099   ResetPOT(original.GetPOT());
00100 
00101   // Copy the spectrum from a NuMatrixSpectrum!
00102   fSpectrum = new TH1D(*original.Spectrum());
00103 
00104   //Copy the Cumulative Distribution from a NuMatrixSpectrum
00105   fCumulativeDist = new TSpline5(*original.fCumulativeDist);
00106 
00107 }

NuMatrix1D::NuMatrix1D ( const TString &  filename,
const TString &  histname,
Double_t  POT 
)

Load a histogram from file. This version accepts a double giving the POT associated with this histogram. If there is an error reading anything, an error is printed and the spectrum and POT are left blank.

Parameters:
filename The file to read in
histname The histogram to read from the file
POT The POT of this histogram

Definition at line 170 of file NuMatrix1D.cxx.

References exit(), fSpectrum, Msg::kError, MSG, and NuMatrix::ResetPOT().

00173   : NuMatrix(), fSpectrum(0), fCumulativeDist(0)
00174 {
00175   // Attempt to open the named file
00176   TFile tf(filename);
00177   if (!tf.IsOpen() || tf.IsZombie())
00178   {
00179     MSG("NuMatrix1D",Msg::kError) << "Cannot open file " << filename << endl;
00180     exit(1);
00181   }
00182 
00183   // Now we have opened, grab the histogram
00184   fSpectrum = dynamic_cast<TH1D*>(tf.Get(histname));
00185   // Detach the spectrum from the file
00186   fSpectrum->SetDirectory(0);
00187 
00188   // Check this worked
00189   if (fSpectrum == 0)
00190   {
00191     MSG("NuMatrix1D",Msg::kError) << "Cannot read histogram " << histname
00192       << " from file " << filename << endl;
00193     exit(1);
00194   }
00195 
00196   // Now we have successfully read the histogram, set the POT
00197   this->ResetPOT(POT);
00198 
00199 }

NuMatrix1D::NuMatrix1D ( const TString &  filename,
const TString &  histname,
const TString &  pothistname = "hTotalPot" 
)

Loads a histogram from a file. This version reads the POT from a histogram in the same file by integrating it. The default name to try is the DST framework convention of "hTotalPot". If there is an error reading anything, an error is printed and the spectrum and POT are left blank.

Parameters:
filename The file to read in
histname The histogram to read from the file
pothistname The histogram to integrate to get the POT

Definition at line 204 of file NuMatrix1D.cxx.

References exit(), fSpectrum, Msg::kError, MSG, and NuMatrix::ResetPOT().

00207                       : NuMatrix(), fSpectrum(0), fCumulativeDist(0)
00208 {
00209   // Attempt to open the named file
00210   TFile tf(filename);
00211   if (!tf.IsOpen() || tf.IsZombie())
00212   {
00213     MSG("NuMatrix1D",Msg::kError) << "Cannot open file " << filename << endl;
00214     exit(1);
00215   }
00216 
00217   // Only try to read the histogram if we have been given a name!
00218   if (histname != "") {
00219     // Now we have opened, grab the histogram
00220     TH1D *spect = dynamic_cast<TH1D*>(tf.Get(histname));
00221     // Check this worked
00222     if (spect == 0)
00223     {
00224       MSG("NuMatrix1D",Msg::kError) << "Cannot read histogram " << histname
00225         << " from file " << filename << endl;
00226       return;
00227     }
00228     fSpectrum = new TH1D(*spect);
00229     // Detach the spectrum from the file
00230     fSpectrum->SetDirectory(0);
00231   }
00232 
00233   // Now we have successfully read the histogram, get the POT histogram
00234   if (!pothistname.IsNull()) {
00235     TH1 *hTotalPot = dynamic_cast<TH1*>(tf.Get(pothistname));
00236     if (hTotalPot == 0)
00237     {
00238       MSG("NuMatrix1D",Msg::kError) << "Cannot read POT histogram " << pothistname
00239         << " from file " << filename << endl;
00240       exit(1);
00241     }
00242     // Set the POT from this histogram
00243     this->ResetPOT(hTotalPot->Integral());
00244   } else {
00245     // No POT histogram is provided
00246     this->ResetPOT(0.0);
00247   }
00248 }

NuMatrix1D::~NuMatrix1D (  )  [virtual]

Definition at line 128 of file NuMatrix1D.cxx.

References fCumulativeDist, and fSpectrum.

00129 {
00130   if (fSpectrum) delete fSpectrum; fSpectrum = 0;
00131   if (fCumulativeDist) delete fCumulativeDist; fCumulativeDist = 0;
00132 }


Member Function Documentation

void NuMatrix1D::Add ( const NuMatrix1D correction,
Bool_t  addpot = false 
) [virtual]

Functions for addition and subtraction.

Implements NuMatrix.

Definition at line 365 of file NuMatrix1D.cxx.

References fSpectrum, NuMatrix::GetPOT(), NuMatrix::ResetPOT(), and Spectrum().

00366 {
00367   // If we don't have a spectrum, copy it
00368   if (!fSpectrum) {
00369     fSpectrum = new TH1D(*dynamic_cast<TH1D*>(with.Spectrum()));
00370   } else {
00371     // Just add it
00372     fSpectrum->Add(with.Spectrum());
00373   }
00374   // Add the POT if this is requested
00375   if (addpot) ResetPOT(GetPOT()+with.GetPOT());
00376 }

void NuMatrix1D::Add ( const NuMatrix correction,
Bool_t  addpot = false 
) [virtual]

Definition at line 352 of file NuMatrix1D.cxx.

References Add(), Msg::kWarning, and MSG.

00353 {
00354   // Convert the NuMatrix to a NuMatrix1D
00355   const NuMatrix1D *other = dynamic_cast<const NuMatrix1D*>(&with);
00356   if (other == 0) {
00357     MSG("NuMatrix1D",Msg::kWarning) << "Attempting to add a NuMatrix1D to something other than a NuMatrix1D" << endl;
00358   }
00359 
00360   this->Add(*other, addpot);
00361 }

void NuMatrix1D::Add ( const TH1 *  correction  )  [virtual]

Functions for addition and subtraction.

Definition at line 344 of file NuMatrix1D.cxx.

References fSpectrum.

Referenced by RockMatrixRAF0720Ext::Add(), RockMatrixRAF0720Eonly::Add(), Add(), RockMatrixRAF0720Std::Add(), NuSterileRunFC::CalcNDDataMCRatio(), NuMMRunCC2010New::CalculateFDNCBackground(), NuMMRunNDOsc::CalculateFDNCBackground(), NuMMRunNSINubar::CalculateFDNCBackground(), NuMMRunNSINu::CalculateFDNCBackground(), NuEZRunsFitter::FDDataNQ(), NuEZRunsFitter::FDDataPQ(), NuMMRunTemplates::Go(), NuMMRunCCTutorial::MakeFDPred(), NuMMRunCC2010New::MakeFDPred(), NuMMRunNoChargeCut::MakeFDPred(), NuMMRunNC2010::MakeFDPred(), NuMMRunCPT::MakeFDPred(), NuMMRunCC2010::MakeFDPred(), NuMMRunTransSME::MakeFDPred(), NuMMRunPRL::MakeFDPred(), NuMMRunNDOsc::MakeFDPred(), NuMMRunFCNSINubar::MakeFDPred(), NuMMRunTransition::MakeFDPred(), NuMMRunLED::MakeFDPred(), NuMMRunNSI::MakeFDPred(), NuMMRunFCNSINu::MakeFDPred(), NuMMRunFC::MakeFDPred(), NuMMRunFC::MakeFDPredBackgrounds(), NuMMRunLED::MakeFDPredBackgrounds(), NuMMRunTransSME::MakeFDPredNuBar(), NuMMRunTransition::MakeFDPredNuBar(), NuMMRunTransition::MakeFDPredNuMu(), NuMMRunTransSME::MakeFDPredNuMu(), NuEZRunsFitter::PredictionNQ(), NuEZRunsFitter::PredictionPQ(), RockMatrixRAF0720Eonly::Subtract(), RockMatrixRAF0720Std::Subtract(), RockMatrixRAF0720Ext::Subtract(), NuEZFitterNSI::UseFDFakeData(), NuEZFitter::UseFDFakeData(), NuMMRunNSINu::WriteFDPredHistos(), NuMMRunNoChargeCut::WriteFDPredHistos(), NuMMRunCC2010::WriteFDPredHistos(), NuMMRunTransSME::WriteFDPredHistos(), NuMMRunCC2010New::WriteFDPredHistos(), NuMMRunPRL::WriteFDPredHistos(), NuMMRunFCNSINu::WriteFDPredHistos(), NuMMRunLED::WriteFDPredHistos(), NuMMRunFC::WriteFDPredHistos(), NuMMRunFCNSINubar::WriteFDPredHistos(), NuMMRunNC2010::WriteFDPredHistos(), NuMMRunTransition::WriteFDPredHistos(), NuMMRunNSINubar::WriteFDPredHistos(), NuMMRunNDOsc::WriteFDPredHistos(), NuMMRunNSI::WriteFDPredHistos(), NuMMRunCPT::WriteFDPredHistos(), and NuMMRunCCTutorial::WriteFDPredHistos().

00345 {
00346   assert(fSpectrum);
00347   fSpectrum->Add(correction);
00348 }

NuMatrix1D::ClassDef ( NuMatrix1D  ,
 
) [protected]

Reimplemented in NuMatrixSpectrum.

Bool_t NuMatrix1D::Complete (  )  const [inline, virtual]

Returns true if the spectrum is properly configured.

Implements NuMatrix.

Definition at line 70 of file NuMatrix1D.h.

References fSpectrum.

Referenced by RockMatrixRAF0720Std::Complete(), RockMatrixRAF0720Ext::Complete(), RockMatrixRAF0720Eonly::Complete(), RockFiller::RockFiller(), and TestExperimentGenerator::testAccessSummary().

00070 { return fSpectrum; }

virtual NuMatrix1D* NuMatrix1D::Copy ( void   )  const [inline, virtual]

Makes a deep-object copy.

Implements NuMatrix.

Reimplemented in NuMatrixSpectrum, and Rock1DMatrix.

Definition at line 56 of file NuMatrix1D.h.

References NuMatrix1D().

00056 { return new NuMatrix1D(*this); };

void NuMatrix1D::DecayCC ( const Double_t  dm2,
const Double_t  sn2,
const Double_t  alpha 
) [virtual]
void NuMatrix1D::DecayMuToTau ( const Double_t  dm2,
const Double_t  sn2,
const Double_t  alpha 
) [virtual]
void NuMatrix1D::DecayNC ( const Double_t  sn2,
const Double_t  alpha 
) [virtual]

Implements NuMatrix.

Definition at line 478 of file NuMatrix1D.cxx.

References NuOscProbCalc::DecayWeightNC(), and fSpectrum.

Referenced by NuMMRunNC2010::Oscillate().

00479 {
00480   for(int i=1;i<=fSpectrum->GetNbinsX();i++) {
00481     const Double_t energy = fSpectrum->GetBinCenter(i);
00482     const Double_t decayProb = NuOscProbCalc::DecayWeightNC(energy, alpha, sn2);
00483     fSpectrum->
00484       SetBinContent(i,fSpectrum->GetBinContent(i)*decayProb);
00485     fSpectrum->
00486       SetBinError(i,fSpectrum->GetBinError(i)*decayProb);
00487   }
00488 }

void NuMatrix1D::Decohere ( const Double_t  dm2,
const Double_t  sn2,
const Double_t  mu2 
) [virtual]

Implements NuMatrix.

Definition at line 492 of file NuMatrix1D.cxx.

References NuOscProbCalc::DecoherenceWeight(), fSpectrum, and GetEventDensity().

Referenced by NuMMRunCCTutorial::MakeFDPred(), NuMMRunNoChargeCut::MakeFDPred(), NuMMRunPRL::MakeFDPred(), NuMMRunNC2010::Oscillate(), NuMMRunFCNSINubar::OscillateNuBar(), NuMMRunNSINubar::OscillateNuBar(), NuMMRunLED::OscillateNuBar(), NuMMRunFCNSINu::OscillateNuBar(), NuMMRunNSINu::OscillateNuBar(), NuMMRunTransSME::OscillateNuBar(), NuMMRunFC::OscillateNuBar(), NuMMRunNSINu::OscillateNuMu(), NuMMRunFCNSINu::OscillateNuMu(), NuMMRunLED::OscillateNuMu(), NuMMRunNSINubar::OscillateNuMu(), NuMMRunCC2010New::OscillateNuMu(), NuMMRunCC2010::OscillateNuMu(), NuMMRunFC::OscillateNuMu(), NuMMRunTransSME::OscillateNuMu(), NuMMRunFCNSINubar::OscillateNuMu(), NuMMRunNoChargeCut::WriteFDPredHistos(), NuMMRunPRL::WriteFDPredHistos(), and NuMMRunCCTutorial::WriteFDPredHistos().

00495 {
00496   for(int i=1;i<=fSpectrum->GetNbinsX();i++) {
00497     const Double_t energy = fSpectrum->GetBinCenter(i);
00498     const Double_t binWidth = fSpectrum->GetBinWidth(i);
00499 
00500     Double_t decoherenceProb = 0; //NuOscProbCalc::DecoherenceWeight(energy, mu2, sn2, dm2, -1, binWidth);
00501     Double_t densityNorm = 0;
00502 
00503     for(int j=0; j<100; j++){
00504       double ej = energy + binWidth*( (j+0.5)/100 - 0.5 );
00505       double dj = this->GetEventDensity(ej);
00506       decoherenceProb += dj*NuOscProbCalc::DecoherenceWeight(ej, mu2, sn2, dm2, -1, 0*binWidth/100);
00507       densityNorm += dj;
00508     }
00509 
00510     if(densityNorm>0) decoherenceProb /= densityNorm;
00511     else decoherenceProb = 1;
00512 
00513     if(energy<3) cout << decoherenceProb << endl;
00514 
00515     fSpectrum->
00516       SetBinContent(i,fSpectrum->GetBinContent(i)*decoherenceProb);
00517     fSpectrum->
00518       SetBinError(i,fSpectrum->GetBinError(i)*decoherenceProb);
00519   }
00520 }

void NuMatrix1D::Divide ( const TGraph &  correction  )  [virtual]

Implements NuMatrix.

Definition at line 289 of file NuMatrix1D.cxx.

References fSpectrum.

00290 {
00291   // copied form NuMatrixSpectrum
00292   for(int i=1; i<=fSpectrum->GetNbinsX(); i++) {
00293     fSpectrum->SetBinContent(i,fSpectrum->GetBinContent(i) /
00294                                correction.Eval(fSpectrum->GetBinCenter(i),0,""));
00295     fSpectrum->SetBinError(i,fSpectrum->GetBinError(i) /
00296                              correction.Eval(fSpectrum->GetBinCenter(i),0,""));
00297   }
00298 
00299 }

void NuMatrix1D::Divide ( const TH1 &  correction,
const Option_t *  option = "" 
) [virtual]

Implements NuMatrix.

Definition at line 282 of file NuMatrix1D.cxx.

References fSpectrum.

00283 {
00284   fSpectrum->Divide(fSpectrum, &correction, 1, 1, option);
00285 }

void NuMatrix1D::Divide ( const NuMatrix1D correction  )  [virtual]

Make a correction by Dividing.

Implements NuMatrix.

Definition at line 273 of file NuMatrix1D.cxx.

References fSpectrum, and Spectrum().

00274 {
00275 
00276  fSpectrum->Divide(by.Spectrum());
00277 
00278 }

void NuMatrix1D::Divide ( const NuMatrix by  )  [virtual]

Definition at line 260 of file NuMatrix1D.cxx.

References Msg::kWarning, and MSG.

00261 {
00262   // Convert the NuMatrix to a NuMatrix1D
00263   const NuMatrix1D *other = dynamic_cast<const NuMatrix1D*>(&by);
00264   if (other == 0) {
00265     MSG("NuMatrix1D",Msg::kWarning) << "Attempting to divide a NuMatrix1D to something other than a NuMatrix1D" << endl;
00266   }
00267 
00268   this->Divide(*other);
00269 }

virtual void NuMatrix1D::Draw ( Option_t *  option = ""  )  [inline, virtual]

Implements NuMatrix.

Definition at line 163 of file NuMatrix1D.h.

References fSpectrum, and option.

Referenced by RatioPlot::Draw().

00163 { fSpectrum->Draw(option); };

void NuMatrix1D::ExtrapolateNDToFD ( const TH2D &  beammatrix  )  [virtual]

Implements NuMatrix.

Definition at line 763 of file NuMatrix1D.cxx.

References XToY().

00764 {
00765   // MSG("NuMatrix1D",Msg::kWarning) << "ExtrapolateNDToFD not implemented uet for NuMatrix1D" << endl;
00766   this->XToY(correction);
00767 }

void NuMatrix1D::Fill ( const NuEvent evt  )  [virtual]

Stores a NuEvent in the spectrum.

Reimplemented from NuMatrix.

Reimplemented in Rock1DMatrix.

Definition at line 852 of file NuMatrix1D.cxx.

References NuEvent::energy, NuEvent::energyMC, fSpectrum, NuMatrix::GetFillTypeEnergy(), Msg::kError, NuMatrix::kReco, NuMatrix::kTrue, MSG, and NuEvent::rw.

Referenced by NuFCRunInfo::BuildFDSummaries().

00853 {
00854   MatrixFillType en = GetFillTypeEnergy();
00855   if (en == kTrue) {
00856     fSpectrum->Fill(evt.energyMC, evt.rw);
00857   } else if (en == kReco) {
00858     fSpectrum->Fill(evt.energy, evt.rw);
00859   } else {
00860     MSG("NuMatrix1D",Msg::kError)
00861       << "Calling NuMatrix1D::Fill but energy fill state "
00862       << en << " is unknown!" << endl;
00863   }
00864 }

void NuMatrix1D::FillOscCache ( double  dm2  )  const [protected]

Ensure fgCachedOsc is filled with a usable oscillation dip

If dmsq is wrong, or the binning is wrong, regenerate the cache, else leave it be.

Reimplemented in NuMatrixSpectrum.

Definition at line 409 of file NuMatrix1D.cxx.

References NuUtilities::CheckConsistency(), fgCachedDmsq, fgCachedOsc, fSpectrum, n, NuMatrix1D(), and NuOscProbCalc::OscillationWeight().

Referenced by Oscillate().

00410 {
00411   // Cache at wrong dmsq can't help us
00412   if(fgCachedDmsq != dm2){
00413     delete fgCachedOsc;
00414     fgCachedOsc = 0;
00415   }
00416 
00417   // Cache with wrong binning can't help us
00418   if(fgCachedOsc &&
00419      !NuUtilities::CheckConsistency(fSpectrum, fgCachedOsc->fSpectrum)){
00420     delete fgCachedOsc;
00421     fgCachedOsc = 0;
00422   }
00423 
00424   // Cache is bad one way or another. Fill it in
00425   if(!fgCachedOsc){
00426     // Get correct binning
00427     NuMatrix1D* osc = new NuMatrix1D(*this);
00428     for(int n = 0; n < osc->fSpectrum->GetNbinsX()+2; ++n)
00429       osc->fSpectrum->SetBinContent(n, 1);
00430 
00431     for(int i = 1; i <= osc->fSpectrum->GetNbinsX(); ++i){
00432       const double energy = osc->fSpectrum->GetBinCenter(i);
00433       const double oscProb = NuOscProbCalc::OscillationWeight(energy, dm2, 1);
00434 
00435       osc->fSpectrum->SetBinContent(i,osc->fSpectrum->GetBinContent(i)*oscProb);
00436       osc->fSpectrum->SetBinError(i,osc->fSpectrum->GetBinError(i)*oscProb);
00437     } // end for i
00438 
00439     fgCachedOsc = osc;
00440     fgCachedDmsq = dm2;
00441   }
00442 
00443   // We either started with a good cache, or we now have one
00444 }

Double_t NuMatrix1D::GetEventDensity ( Double_t  x  )  [virtual]

Definition at line 1009 of file NuMatrix1D.cxx.

References fCumulativeDist, GetXmax(), and GetXmin().

Referenced by Decohere(), InverseDecohere(), and MakeSterilePred_Interpolated::MCOscillated().

01009                                               {
01010 
01011   // If x is outside the range return 0
01012   if(x < this->GetXmin() || x > this->GetXmax() ){
01013     return 0;
01014   }
01015 
01016   // Values should be positive, but protect against negative events anyway
01017   return TMath::Max(fCumulativeDist->Derivative(x),0.0);
01018 
01019 }

Double_t NuMatrix1D::GetIntegral ( Double_t  x  )  [virtual]

Definition at line 1067 of file NuMatrix1D.cxx.

References fCumulativeDist, GetXmax(), and GetXmin().

Referenced by Rebin().

01067                                           {
01068 
01069   double xmin = this->GetXmin();
01070   double xmax = this->GetXmax();
01071 
01072   // If x is outside the range return 0 or the total integral
01073   if(x < xmin){
01074     return 0;
01075   }
01076   if(x > xmax){
01077     return fCumulativeDist->Eval(xmax);
01078   }
01079 
01080   return fCumulativeDist->Eval(x);
01081 
01082 }

NuMatrixSpectrum NuMatrix1D::GetNuMatrixSpectrum ( void   )  const

Returns a pointer to the built in data storage.

Returns a NuMatrixSpectrum object as a clone of this one

Definition at line 253 of file NuMatrix1D.cxx.

References fSpectrum, and NuMatrix::GetPOT().

Referenced by NuMMRunFCNSINubar::ResetFDData(), NuMMRunFCNSINu::ResetFDData(), NuMMRunFC::ResetFDData(), NuSterileRunFC::ResetFDData(), and NuMMRunLED::ResetFDData().

00254 {
00255   return NuMatrixSpectrum(*fSpectrum, GetPOT());
00256 }

Double_t NuMatrix1D::GetXmax (  )  [virtual]

Definition at line 1097 of file NuMatrix1D.cxx.

References fCumulativeDist, and Interpolate().

Referenced by GetEventDensity(), GetIntegral(), and MakeSterilePred_Interpolated::MCOscillated().

01097                             {
01098 
01099   // Interpolate if it hasn't been done yet
01100   if(!fCumulativeDist) this->Interpolate();
01101 
01102   return fCumulativeDist->GetXmax();
01103 
01104 }

Double_t NuMatrix1D::GetXmin (  )  [virtual]

Definition at line 1086 of file NuMatrix1D.cxx.

References fCumulativeDist, and Interpolate().

Referenced by GetEventDensity(), GetIntegral(), and MakeSterilePred_Interpolated::MCOscillated().

01086                             {
01087 
01088   // Interpolate if it hasn't been done yet
01089   if(!fCumulativeDist) this->Interpolate();
01090 
01091   return fCumulativeDist->GetXmin();
01092 
01093 }

void NuMatrix1D::Interpolate (  )  [virtual]

Definition at line 868 of file NuMatrix1D.cxx.

References MuELoss::e, fCumulativeDist, and fSpectrum.

Referenced by GetXmax(), and GetXmin().

00868                             {
00869 
00870   // Determine the interpolation nodes
00871   int npx = fSpectrum->GetNbinsX()+1;
00872   vector<double> xx;
00873   vector<double> yy;
00874   bool FirstBin = true;
00875   for(int i=1;i<npx;i++){
00876     if(fSpectrum->GetBinContent(i)==0) continue;
00877     if(FirstBin){
00878       xx.push_back(fSpectrum->GetBinLowEdge(i));
00879       yy.push_back(0);
00880     }
00881     xx.push_back(fSpectrum->GetBinLowEdge(i+1));
00882     yy.push_back(fSpectrum->Integral(1,i));
00883     FirstBin = false;
00884   }
00885   npx = xx.size();
00886 
00887   // Interpolate initial nodes
00888   TSpline5 *spl = new TSpline5("",&xx[0],&yy[0],npx,"b1e1b2e2",0,0,0,0);
00889 
00890   // Assume simple nth order polynomial (x-x0)^n on the boundaries
00891   // Shift initial and final nodes to accommodate this boundary condition
00892   int boundorder = 5;
00893   double xi = xx[1]-boundorder*(yy[1]-yy[0])/spl->Derivative(xx[1]); 
00894   double xf = xx[npx-2]-boundorder*(yy[npx-2]-yy[npx-1])/spl->Derivative(xx[npx-2]);
00895   
00896   // Check that new initial and final nodes are within corresponding bins
00897   if(xi>xx[0]&&xi<xx[1])  xx[0] = xi;
00898   if(xf>xx[npx-2]&&xf<xx[npx-1]) xx[npx-1] = xf;
00899 
00900   // Reinterpolate using new nodes
00901   if(fCumulativeDist) delete fCumulativeDist;
00902   fCumulativeDist = new TSpline5("",&xx[0],&yy[0],npx,"b1e1b2e2",0,0,0,0);
00903   
00904   delete spl;
00905   spl = new TSpline5("",&xx[0],&yy[0],npx,"b1e1b2e2",0,0,0,0);
00906 
00907 //  cout << fCumulativeDist->GetXmin() << " - " << fCumulativeDist->GetXmax() << endl;
00908 
00909   // Store first and second derivatives
00910   vector<double> dyy;
00911   vector<double> ddyy;
00912 
00913   for(int i=0;i<npx;i++){
00914 
00915     double b1, b2;
00916     b1 = spl->Derivative(xx[i]);
00917     if(i>0) b2 = (b1-spl->Derivative(xx[i]-0.001))/0.001;
00918     else b2 = (spl->Derivative(xx[i]+0.001)-b1)/0.001;
00919 
00920     dyy.push_back(b1);
00921     ddyy.push_back(b2);
00922 
00923   }// end of loop over nodes
00924 
00925   delete spl;
00926 
00927   // Force interpolation to be monotone
00928   // See Randall L. Dougherty et al, Mathematics of Computation, Vol. 52, No. 186. (1989), pp. 471-494
00929   for(int j=0;j<2;j++){
00930   for(int i=0;i<npx;i++){
00931 
00932     double Sp, Sm;
00933     if(i>0) Sm = (yy[i]-yy[i-1])/(xx[i]-xx[i-1]);
00934     else Sm = 0;
00935     if(i<npx-1) Sp = (yy[i+1]-yy[i])/(xx[i+1]-xx[i]);
00936     else Sp = 0;
00937 
00938     double sigma = 0;
00939     if(Sp*Sm>0) sigma = Sp/TMath::Abs(Sp);
00940     else if(Sp*Sm<0&&dyy[i]!=0) sigma = dyy[i]/TMath::Abs(dyy[i]);
00941 
00942     if(sigma>=0) dyy[i] = TMath::Min(TMath::Max(0.0,dyy[i]),5*TMath::Min(Sp,Sm));
00943     else dyy[i] = TMath::Max(TMath::Min(0.0,dyy[i]),-5*TMath::Min(Sp,Sm));
00944 
00945     double dp = dyy[i];
00946     double dm = dyy[i];
00947     if(dyy[i]*Sp<=0) dp = 0;
00948     if(dyy[i]*Sm<=0) dm = 0;
00949 
00950     double aa = 0;
00951     if(Sm!=0&&i>0) aa = TMath::Max(0.0,dyy[i-1]/Sm);
00952     double bb = 0;
00953     if(Sp!=0&&i<npx-1) bb = TMath::Max(0.0,dyy[i+1]/Sp);
00954 
00955     double dxp;
00956     if(i<npx-1) dxp = (xx[i+1]-xx[i]);
00957     else dxp = (xx[i]-xx[i-1]);
00958     double Ai = ( -7.9*dp - 0.26*dp*bb ) / dxp;
00959     double Af = ( (20-2*bb)*Sp - 8*dp - 0.48*dp*bb ) / dxp;
00960     double dxm = dxp;
00961     if(i>0) dxm = xx[i] - xx[i-1];
00962     double Bi = ( (-20+2*aa)*Sm + 8*dm + 0.48*dm*aa ) / dxm;
00963     double Bf = ( 7.9*dm + 0.26*dm*aa ) / dxm;
00964 
00965     if(Af<Bi||Ai>Bf){
00966 
00967       dyy[i] = ( (20-2*bb)*Sp/dxp + (20-2*aa)*Sm/dxm ) / 
00968                ( (8+0.48*bb)/dxp + (8+0.48*aa)/dxm );
00969 
00970       dp = dyy[i];
00971       dm = dyy[i];
00972       if(dyy[i]*Sp<=0) dp = 0;
00973       if(dyy[i]*Sm<=0) dm = 0;
00974       Ai = ( -7.9*dp - 0.26*dp*bb ) / dxp;
00975       Af = ( (20-2*bb)*Sp - 8*dp - 0.48*dp*bb ) / dxp;
00976       Bi = ( (-20+2*aa)*Sm + 8*dm + 0.48*dm*aa ) / dxm;
00977       Bf = ( 7.9*dm + 0.26*dm*aa ) / dxm;
00978 
00979     }
00980 
00981     if(j==0) continue;
00982 
00983     Af = TMath::Min(Af,Bf);
00984     Ai = TMath::Max(Ai,Bi);
00985 
00986     ddyy[i] = TMath::Max(TMath::Min(ddyy[i],Af),Ai);
00987 
00988   }}// end of loop over nodes
00989  
00990   // Apply monotone constraints to the interpolation
00991   for(int i=0;i<npx-1;i++){
00992 
00993     TSpline5* spli = new TSpline5("",&xx[i],&yy[i],2,"b1e1b2e2",dyy[i],dyy[i+1],ddyy[i],ddyy[i+1]);
00994 
00995     Double_t x, y, b, c, d, e, f;
00996     spli->GetCoeff(spli->FindX((xx[i+1]+xx[i])/2),x,y,b,c,d,e,f);
00997     fCumulativeDist->SetPointCoeff(fCumulativeDist->FindX((xx[i+1]+xx[i])/2),b,c,d,e,f);
00998 
00999     delete spli;
01000 
01001   }// end of loop over nodes
01002 
01003   return;
01004 
01005 }

void NuMatrix1D::InverseDecohere ( const Double_t  dm2,
const Double_t  sn2,
const Double_t  mu2 
) [virtual]

Decohered neutrinos aren't gone, they have to turn up somewhere, as taus.

Why not as nues? We start as c*(nu1+nu2)+s*nu3 - the phase relationship between the cos part and the sin part gets scrambled by the decoherence, so we end up with some mixture of numu and nutau. This doesn't happen between nu1 and nu2 (which would give us nue appearance) because (by assumption) the mass-splitting is too small compared to our L/E so very little decoherence occurs here.

Definition at line 523 of file NuMatrix1D.cxx.

References NuOscProbCalc::DecoherenceWeight(), fSpectrum, and GetEventDensity().

Referenced by NuMMRunNSINu::InverseOscillateNuTau(), NuMMRunFC::InverseOscillateNuTau(), NuMMRunCC2010New::InverseOscillateNuTau(), NuMMRunFCNSINu::InverseOscillateNuTau(), NuMMRunCC2010::InverseOscillateNuTau(), NuMMRunTransSME::InverseOscillateNuTau(), NuMMRunNSINubar::InverseOscillateNuTau(), NuMMRunLED::InverseOscillateNuTau(), NuMMRunFCNSINubar::InverseOscillateNuTau(), NuMMRunTransSME::InverseOscillateTauBar(), NuMMRunFCNSINubar::InverseOscillateTauBar(), NuMMRunFC::InverseOscillateTauBar(), NuMMRunNSINu::InverseOscillateTauBar(), NuMMRunNSINubar::InverseOscillateTauBar(), NuMMRunLED::InverseOscillateTauBar(), and NuMMRunFCNSINu::InverseOscillateTauBar().

00526 {
00527   for(int i=1;i<=fSpectrum->GetNbinsX();i++) {
00528     const Double_t energy = fSpectrum->GetBinCenter(i);
00529     const Double_t binWidth = fSpectrum->GetBinWidth(i);
00530 
00531 //    const Double_t decoherenceProb = 1-NuOscProbCalc::DecoherenceWeight(energy, mu2, sn2, dm2, -1, binWidth);
00532 
00533     Double_t decoherenceProb = 0; //NuOscProbCalc::DecoherenceWeight(energy, mu2, sn2, dm2, -1, binWidth);
00534     Double_t densityNorm = 0;
00535 
00536     for(int j=0; j<100; j++){
00537       double ej = energy + binWidth*( (j+0.5)/100 - 0.5 );
00538       double dj = this->GetEventDensity(ej);
00539       decoherenceProb += dj*(1-NuOscProbCalc::DecoherenceWeight(ej, mu2, sn2, dm2, -1, binWidth/100));
00540       densityNorm += dj;
00541     }
00542 
00543     if(densityNorm>0) decoherenceProb /= densityNorm;
00544     else decoherenceProb = 1;
00545 
00546     fSpectrum->
00547       SetBinContent(i,fSpectrum->GetBinContent(i)*decoherenceProb);
00548     fSpectrum->
00549       SetBinError(i,fSpectrum->GetBinError(i)*decoherenceProb);
00550   }
00551 }

void NuMatrix1D::InverseOscillateNSI ( const Double_t  dm2,
const Double_t  sn2,
const Double_t  epsilon,
const Double_t  sign 
) [virtual]

Reimplemented in NuMatrixSpectrum.

Definition at line 589 of file NuMatrix1D.cxx.

References fSpectrum, and NuOscProbCalc::NSIWeight().

00593 {
00594   for(int i=1;i<=fSpectrum->GetNbinsX();i++) {
00595     const Double_t energy = fSpectrum->GetBinCenter(i);
00596     const Double_t nsiProb = 1-NuOscProbCalc::NSIWeight(energy, dm2, sn2, epsilon, sign);
00597     fSpectrum->
00598       SetBinContent(i,fSpectrum->GetBinContent(i)*nsiProb);
00599     fSpectrum->
00600       SetBinError(i,fSpectrum->GetBinError(i)*nsiProb);
00601   }
00602 
00603 }

void NuMatrix1D::Multiply ( const TGraph &  correction  )  [virtual]

Implements NuMatrix.

Definition at line 330 of file NuMatrix1D.cxx.

References fSpectrum.

00331 {
00332   // copied form NuMatrixSpectrum
00333   for(int i=1; i<=fSpectrum->GetNbinsX(); i++) {
00334     fSpectrum->SetBinContent(i,fSpectrum->GetBinContent(i) *
00335                                correction.Eval(fSpectrum->GetBinCenter(i),0,""));
00336     fSpectrum->SetBinError(i,fSpectrum->GetBinError(i) *
00337                              correction.Eval(fSpectrum->GetBinCenter(i),0,""));
00338   }
00339 
00340 }

void NuMatrix1D::Multiply ( const TH1 &  correction,
const Option_t *  option = "" 
) [virtual]

Implements NuMatrix.

Definition at line 323 of file NuMatrix1D.cxx.

References fSpectrum.

00324 {
00325   fSpectrum->Multiply(fSpectrum, &correction, 1, 1, option);
00326 }

void NuMatrix1D::Multiply ( const NuMatrix1D correction  )  [virtual]

Make a correction by Multiplying.

Implements NuMatrix.

Definition at line 316 of file NuMatrix1D.cxx.

References fSpectrum, and Spectrum().

00317 {
00318   fSpectrum->Multiply(by.Spectrum());
00319 }

void NuMatrix1D::Multiply ( const NuMatrix by  )  [virtual]

Definition at line 303 of file NuMatrix1D.cxx.

References Msg::kWarning, and MSG.

00304 {
00305   // Convert the NuMatrix to a NuMatrix1D
00306   const NuMatrix1D *other = dynamic_cast<const NuMatrix1D*>(&by);
00307   if (other == 0) {
00308     MSG("NuMatrix1D",Msg::kWarning) << "Attempting to multiply a NuMatrix1D to something other than a NuMatrix1D" << endl;
00309   }
00310 
00311   this->Multiply(*other);
00312 }

NuMatrix1D & NuMatrix1D::operator= ( const NuMatrixSpectrum orig  ) 

Assignment from a NuMatrixSpectrum.

Definition at line 110 of file NuMatrix1D.cxx.

References fCumulativeDist, fSpectrum, NuMatrix::GetPOT(), NuMatrix::ResetPOT(), and NuMatrixSpectrum::Spectrum().

00111 {
00112   // Copy the POT count
00113   this->ResetPOT(orig.GetPOT());
00114 
00115   // Copy the spectrum from a NuMatrixSpectrum!
00116   if (fSpectrum) delete fSpectrum;
00117   fSpectrum = new TH1D(*orig.Spectrum());
00118 
00119   //Copy the Cumulative Distribution from a NuMatrixSpectrum
00120   if (fCumulativeDist) delete fCumulativeDist;
00121   fCumulativeDist = new TSpline5(*orig.fCumulativeDist);
00122 
00123   return *this;
00124 }

NuMatrix1D & NuMatrix1D::operator= ( const NuMatrix1D orig  ) 

Reimplemented in Rock1DMatrix.

Definition at line 72 of file NuMatrix1D.cxx.

References fCumulativeDist, and fSpectrum.

00073 {
00074   if (this != &orig) {
00075     // Copy the base class stuff
00076     this->NuMatrix::operator=(orig);
00077 
00078     // Copy the histogram
00079     if (fSpectrum) delete fSpectrum; fSpectrum = 0;
00080     if (orig.fSpectrum) {
00081       fSpectrum = new TH1D(*orig.fSpectrum);
00082       fSpectrum->SetDirectory(0);
00083     }
00084 
00085     //Copy the Cumulative Distribution
00086     if (fCumulativeDist) delete fCumulativeDist; fCumulativeDist = 0;
00087     if (orig.fCumulativeDist) {
00088       fCumulativeDist = new TSpline5(*orig.fCumulativeDist);
00089     }
00090 
00091   }
00092   return *this;
00093 }

void NuMatrix1D::Oscillate ( const Double_t  dm2,
const Double_t  sn2 
) [virtual]

Oscillation function.

Implements NuMatrix.

Reimplemented in NuMatrixSpectrum.

Definition at line 448 of file NuMatrix1D.cxx.

References fgCachedOsc, FillOscCache(), fSpectrum, and n.

00449 {
00450   FillOscCache(dm2);
00451 
00452   for(int n = 0; n < fSpectrum->GetNbinsX()+2; ++n){
00453     fSpectrum->SetBinContent(n,
00454                              fSpectrum->GetBinContent(n)*
00455                              (1-sn2*(1-fgCachedOsc->fSpectrum->GetBinContent(n))));
00456   }
00457 }

void NuMatrix1D::OscillatedNDRecoToTrue ( const TH2D &  unnormalizedCorrection,
const NuMMParameters oscpars 
) [virtual]

Definition at line 694 of file NuMatrix1D.cxx.

References NuMMParameters::Delta1(), NuMMParameters::Delta2(), NuMMParameters::Delta3(), NuUtilities::DistanceTargetToND(), NuMMParameters::Dm2(), NuMMParameters::Dm221(), NuMMParameters::Dm243(), NuOscProb::FourFlavourDisappearanceWeight(), RecoToTrue(), NuMMParameters::Theta12(), NuMMParameters::Theta13(), NuMMParameters::Theta14(), NuMMParameters::Theta23(), NuMMParameters::Theta24(), and NuMMParameters::Theta34().

00696 {
00697 
00698   //Create the normalized correction
00699   TH2D correction(unnormalizedCorrection);
00700   correction.Reset();
00701   //correction.Clear();
00702 
00703   for (Int_t xBin = 1; xBin < unnormalizedCorrection.GetNbinsX()+1; ++xBin){//true
00704     for (Int_t yBin = 1; yBin < unnormalizedCorrection.GetNbinsY()+1; ++yBin){
00705 
00706       Double_t energy = unnormalizedCorrection.GetXaxis()->GetBinCenter(xBin);
00707       Double_t oscProb =
00708         NuOscProbCalc::FourFlavourDisappearanceWeight(energy,
00709                                                       oscpars.Dm2(), oscpars.Theta23(),
00710                                                       oscpars.Dm221(), oscpars.Dm243(),
00711                                                       oscpars.Theta12(),
00712                                                       oscpars.Theta13(), oscpars.Theta14(),
00713                                                       oscpars.Theta24(), oscpars.Theta34(),
00714                                                       oscpars.Delta1(), oscpars.Delta2(),
00715                                                       oscpars.Delta3(),
00716                                                       NuUtilities::DistanceTargetToND());
00717       correction.SetBinContent(xBin,yBin,
00718                                 unnormalizedCorrection.GetBinContent(xBin,yBin)*oscProb);
00719       
00720     }
00721   }
00722 
00723   //Get the totals in each reco energy row, including the overflow
00724   std::vector<Double_t> vRowTotals;
00725   for (Int_t yBin = 1; yBin <= correction.GetNbinsY()+1; ++yBin){
00726     Double_t rowTotal = 0;
00727     for (Int_t xBin = 1; xBin <= correction.GetNbinsX()+1; ++xBin){
00728       rowTotal += correction.GetBinContent(xBin,yBin);
00729     }
00730     vRowTotals.push_back(rowTotal);
00731   }
00732   
00733   for(int i=1;i<=correction.GetNbinsX();i++){//true
00734     for(int j=1;j<=correction.GetNbinsY()+1;j++){//reco
00735       if(vRowTotals[j-1]>0 &&
00736          correction.GetBinContent(i,j)>0) {
00737         correction.SetBinContent(i,j,
00738                                  correction.GetBinContent(i,j)/vRowTotals[j-1]);
00739       }
00740       else {
00741         correction.SetBinContent(i,j,0);
00742       }
00743     }
00744   }
00745 
00746 
00747   
00748   this->RecoToTrue(correction);
00749   
00750 }

void NuMatrix1D::OscillateNSI ( const Double_t  dm2,
const Double_t  sn2,
const Double_t  epsilon,
const Double_t  sign 
) [virtual]

Reimplemented in NuMatrixSpectrum.

Definition at line 572 of file NuMatrix1D.cxx.

References fSpectrum, and NuOscProbCalc::NSIWeight().

00576 {
00577   for(int i=1;i<=fSpectrum->GetNbinsX();i++) {
00578     const Double_t energy = fSpectrum->GetBinCenter(i);
00579     const Double_t nsiProb = NuOscProbCalc::NSIWeight(energy, dm2, sn2, epsilon, sign);
00580     fSpectrum->
00581       SetBinContent(i,fSpectrum->GetBinContent(i)*nsiProb);
00582     fSpectrum->
00583       SetBinError(i,fSpectrum->GetBinError(i)*nsiProb);
00584   }
00585 
00586 }

void NuMatrix1D::Read ( TFile *  f,
const TString &  name 
) [virtual]

Definition at line 656 of file NuMatrix1D.cxx.

References fSpectrum, Msg::kError, and MSG.

Referenced by RockMatrixRAF0720Eonly::Read(), RockMatrixRAF0720Ext::Read(), and RockMatrixRAF0720Std::Read().

00657 {
00658   if(!f || f->IsZombie()){
00659     MSG("NuMatrix1D",Msg::kError) << "Invalid file passed to Read()\n"
00660       "Passing back empty histogram\n";
00661     fSpectrum = new TH1D;
00662     return;
00663   }
00664 
00665   fSpectrum = dynamic_cast<TH1D *>(f->Get(name));
00666 
00667   if(!fSpectrum){
00668     MSG("NuMatrix1D",Msg::kError) << "Failed to Read() histogram named "
00669       << name << " from file!\nPassing back empty histogram\n";
00670     fSpectrum = new TH1D;
00671     return;
00672   }
00673 
00674   // Shouldn't be associated with the file it came from...
00675   fSpectrum->SetDirectory(0);
00676 }

TH1D * NuMatrix1D::Rebin ( Int_t  nbins,
Double_t *  edges 
) [virtual]

Definition at line 1023 of file NuMatrix1D.cxx.

References MuELoss::e, fSpectrum, GetIntegral(), Msg::kWarning, and MSG.

01023                                                    {
01024 
01025   MSG("NuMatrix1D",Msg::kWarning) << "Error is estimated assuming scale * sqrt(NumEvents). "
01026                                   << "This will be wrong if histogram has varying weights "
01027                                   << "and bin edges don't match original histogram." << endl;
01028 
01029   // New Histogram
01030   TH1D* hOut = new TH1D("","",nbins,edges);
01031 
01032   // Loop over new bins
01033   for(int i=1; i<=nbins; i++){
01034 
01035     // Find up and low edges
01036     double lowedge = hOut->GetBinLowEdge(i);
01037     double upedge = lowedge + hOut->GetBinWidth(i);
01038 
01039     // Get first and last old bins in new bin
01040     int origlowbin = fSpectrum->FindBin(lowedge);
01041     int origupbin = fSpectrum->FindBin(upedge - 1e-10*hOut->GetBinWidth(i));
01042 
01043     // Sum contents and errors of old bins 
01044     double origcont = 0;
01045     double origerr = 0;
01046     for(int j=origlowbin; j<=origupbin; j++){
01047       origcont += fSpectrum->GetBinContent(j);
01048       origerr += pow(fSpectrum->GetBinError(j),2);
01049     }
01050     origerr = sqrt(origerr);
01051 
01052     // Set new bin content from interpolation    
01053     hOut->SetBinContent(i,this->GetIntegral(upedge) - this->GetIntegral(lowedge));
01054 
01055     // Set new bin error based on old bin error and content
01056     if(origcont!=0) hOut->SetBinError(i,origerr * sqrt(TMath::Abs(hOut->GetBinContent(i) / origcont)));
01057     else            hOut->SetBinError(i,0);
01058 
01059   }
01060 
01061   return hOut;
01062 
01063 }

void NuMatrix1D::RecoToTrue ( const TH2D &  correction  )  [virtual]

Implements NuMatrix.

Definition at line 687 of file NuMatrix1D.cxx.

References YToX().

Referenced by OscillatedNDRecoToTrue().

00688 {
00689   // MSG("NuMatrix1D",Msg::kWarning) << "RecoToTrue not implemented uet for NuMatrix1D" << endl;
00690   this->YToX(correction);
00691 }

void NuMatrix1D::Scale ( Double_t  value  )  [virtual]

Scales the histogram.

Implements NuMatrix.

Definition at line 680 of file NuMatrix1D.cxx.

References fSpectrum.

Referenced by NuFCRunInfo::BuildFDSummaries(), NuMMRunTemplates::Go(), NuMMRunTransSME::MakeFDPred(), RockMatrixRAF0720Eonly::Scale(), RockMatrixRAF0720Std::Scale(), and RockMatrixRAF0720Ext::Scale().

00681 {
00682   if (fSpectrum) fSpectrum->Scale(value);
00683 }

void NuMatrix1D::ScaleToPOT ( Double_t  new_POT  )  [virtual]

Changes the POT, by scaling the data. Beware when using this function that all the data you are storing scales with POT. Enumerations and plane numbers and such fall into this category. Re-implement the function in derived classes if it is needed.

Implements NuMatrix.

Definition at line 155 of file NuMatrix1D.cxx.

References fSpectrum, NuMatrix::GetPOT(), and NuMatrix::ResetPOT().

Referenced by NuMatrixSpectrum::BlessedFD(), NuMatrixSpectrum::BlessedFD7e20(), NuMatrixSpectrum::BlessedFDRatio(), NuMatrixSpectrum::BlessedND(), NuMMRunNC2010::NuMMRunNC2010(), NuEZFitter::ScaleToPot(), NuMatrixSpectrum::ScaleToPot(), NuEZFitterNSI::ScaleToPot(), NuEZFitterNSI::UseFDFakeData(), and NuEZFitter::UseFDFakeData().

00156 {
00157   // check we have POT and a spectrum
00158   if (GetPOT() != 0.0 && fSpectrum) {
00159     Double_t scale_factor = new_POT / GetPOT();
00160     fSpectrum->Scale(scale_factor);
00161   }
00162 
00163   // Reset the internal POT
00164   ResetPOT(new_POT);
00165 }

virtual void NuMatrix1D::SetSpectrum ( const TH1D &  spectrum  )  [inline, virtual]

Definition at line 157 of file NuMatrix1D.h.

References fSpectrum.

Referenced by RockMatrixRAF0720Eonly::RockMatrixRAF0720Eonly(), RockMatrixRAF0720Ext::RockMatrixRAF0720Ext(), and RockMatrixRAF0720Std::RockMatrixRAF0720Std().

00157                                                  {
00158     if (fSpectrum) delete fSpectrum;
00159     fSpectrum = static_cast<TH1D*>(spectrum.Clone());
00160   }

void NuMatrix1D::SetValue ( Double_t  val  )  [virtual]

Sets the entire contents of the spectrum to a single value.

Implements NuMatrix.

Definition at line 624 of file NuMatrix1D.cxx.

References fSpectrum, Msg::kWarning, and MSG.

Referenced by RockMatrixRAF0720Ext::SetValue(), RockMatrixRAF0720Eonly::SetValue(), and RockMatrixRAF0720Std::SetValue().

00625 {
00626   // No spectrum - don't set
00627   if (!fSpectrum) {
00628     MSG("NuMatrix1D",Msg::kWarning) << "Cannot set constant value: No histogram" << endl;
00629     return;
00630   }
00631 
00632   for (Int_t i = 1; i <= fSpectrum->GetNbinsX(); ++i)
00633   {
00634     fSpectrum->SetBinContent(i, val);
00635   }
00636 }

virtual TH1* NuMatrix1D::Spectrum ( void   )  const [inline, virtual]
Double_t NuMatrix1D::StatsLikelihood ( const NuMatrix1D with  )  const [virtual]

Compares two 2D histograms and returns the log-likelihood difference.

Implements NuMatrix.

Definition at line 606 of file NuMatrix1D.cxx.

References fSpectrum, NuUtilities::LogLikelihood(), and Spectrum().

00607 {
00608   Double_t like = 0;
00609 
00610   //Aim to minimise -2lnL. That's what I'm returning.
00611   for (Int_t i=1; i<= fSpectrum->GetNbinsX(); ++i)
00612   {
00613     Double_t mnu = fSpectrum->GetBinContent(i);
00614     Double_t dnu = with->Spectrum()->GetBinContent(i);
00615 
00616     like += NuUtilities::LogLikelihood(dnu,mnu);
00617   }
00618   return like;
00619 }

Double_t NuMatrix1D::StatsLikelihood ( const NuMatrix with  )  const [virtual]

Checks the passed object is a NuMatrix1D then returns the log-likelihood difference.

Definition at line 138 of file NuMatrix1D.cxx.

References Msg::kWarning, and MSG.

Referenced by RockMatrixRAF0720Eonly::StatsLikelihood(), RockMatrixRAF0720Ext::StatsLikelihood(), and RockMatrixRAF0720Std::StatsLikelihood().

00139 {
00140   // Convert the NuMatrix to a NuMatrix1D
00141   const NuMatrix1D *other = dynamic_cast<const NuMatrix1D*>(with);
00142   if (other == 0) {
00143     MSG("NuMatrix1D",Msg::kWarning) << "Attempting to compare a NuMatrix1D to something other than a NuMatrix1D" << endl;
00144     return -1;
00145   }
00146 
00147   // // Do the statistical comparison
00148   // MSG("NuMatrix1D",Msg::kError) << "Comparing statistical likelihood - but function not implemented" << endl;
00149   // return 0;
00150   return this->StatsLikelihood(other);
00151 }

void NuMatrix1D::Subtract ( const NuMatrix1D correction  )  [virtual]

Implements NuMatrix.

Definition at line 400 of file NuMatrix1D.cxx.

References fSpectrum, and Spectrum().

00401 {
00402   assert(fSpectrum);
00403   fSpectrum->Add(with.Spectrum(), -1.0);
00404 }

void NuMatrix1D::Subtract ( const NuMatrix correction  )  [virtual]

Definition at line 387 of file NuMatrix1D.cxx.

References Msg::kWarning, MSG, and Subtract().

00388 {
00389   // Convert the NuMatrix to a NuMatrix1D
00390   const NuMatrix1D *other = dynamic_cast<const NuMatrix1D*>(&with);
00391   if (other == 0) {
00392     MSG("NuMatrix1D",Msg::kWarning) << "Attempting to subtract a NuMatrix1D to something other than a NuMatrix1D" << endl;
00393   }
00394 
00395   this->Subtract(*other);
00396 }

void NuMatrix1D::Subtract ( const TH1 *  correction  )  [virtual]
void NuMatrix1D::TrueToReco ( const TH2D &  correction  )  [virtual]

Implements NuMatrix.

Definition at line 755 of file NuMatrix1D.cxx.

References XToY().

00756 {
00757   // MSG("NuMatrix1D",Msg::kWarning) << "TrueToReco not implemented uet for NuMatrix1D" << endl;
00758   this->XToY(correction);
00759 }

Int_t NuMatrix1D::Write ( const TString &  name  )  const [virtual]

Writes to the current file.

Implements NuMatrix.

Definition at line 640 of file NuMatrix1D.cxx.

References fSpectrum.

Referenced by RockMatrixRAF0720Std::Write(), RockMatrixRAF0720Eonly::Write(), and RockMatrixRAF0720Ext::Write().

00641 {
00642   if (fSpectrum) {
00643     if (name == "")
00644     {
00645       return fSpectrum->Write();
00646     } else {
00647       return fSpectrum->Write(name);
00648     }
00649   } else {
00650     return 0;
00651   }
00652 }

void NuMatrix1D::XToY ( const TH2D &  correction  )  [protected, virtual]

Definition at line 772 of file NuMatrix1D.cxx.

References fSpectrum.

Referenced by ExtrapolateNDToFD(), and TrueToReco().

00773 {  
00774 
00775   const int NX = correction.GetNbinsX();
00776   const int NY = correction.GetNbinsY();
00777   
00778   Double_t *val = new Double_t[NY+1];
00779   for(int i=1;i<=NY+1;i++){
00780     val[i-1] = 0;
00781   }
00782   for(int i=1;i<=NX;i++){ //loop over true
00783     const double Si = fSpectrum->GetBinContent(i);
00784     for(int j=1;j<=NY+1;j++){ //loop over reco
00785       // This is too slow
00786       const double Cij = correction.GetBinContent(i, j);
00787       // As is this
00788       //      const double Cij = correction.GetBinContent(i + (NX+2)*j);
00789       //const double Cij = correction.fArray[i+(NX+2)*j];
00790       val[j-1] += Si * Cij;
00791     }
00792   }
00793   
00794   if(NY != NX){
00795     //define binning for the new histogram which we will return:
00797     const  TArrayD* Bins = correction.GetYaxis()->GetXbins();
00798     string name = fSpectrum->GetName();
00799     string title = fSpectrum->GetTitle();
00800     delete fSpectrum;
00801     fSpectrum =new TH1D("","", NY,Bins->GetArray() );  
00802     fSpectrum->SetName(name.c_str());
00803     fSpectrum->SetTitle(name.c_str());
00804   }
00805  
00806 
00807   for(int i=1;i<=NY+1;i++) {
00808     fSpectrum->SetBinContent(i,val[i-1]);
00809   }
00810   
00811   delete[] val;
00812 }

void NuMatrix1D::YToX ( const TH2D &  correction  )  [protected, virtual]

Definition at line 815 of file NuMatrix1D.cxx.

References fSpectrum.

Referenced by RecoToTrue().

00816 {
00817   
00818   const int NX = correction.GetNbinsX();
00819   const int NY = correction.GetNbinsY();
00820   
00821   Double_t *val = new Double_t[NX+1];
00822   for(int i=1;i<=NX+1;i++) { val[i-1] = 0; }
00823   
00824   for(int i=1;i<=NY;i++){//reco
00825     for(int j=1;j<=NX;j++){//true
00826       val[j-1] += ( fSpectrum->GetBinContent(i) *
00827                     correction.GetBinContent(j,i) );
00828     }
00829   }
00830   
00831   if(NY != NX){
00832     //define binning for the new histogram which we will return:
00833     const  TArrayD* Bins = correction.GetXaxis()->GetXbins();
00834     string name = fSpectrum->GetName();
00835     string title = fSpectrum->GetTitle();
00836     delete fSpectrum;
00837     fSpectrum =new TH1D("","", NX,Bins->GetArray() );  
00838     fSpectrum->SetName(name.c_str());
00839     fSpectrum->SetTitle(name.c_str());
00840   }
00841 
00842   for(int i=1;i<=NX;i++) {
00843     fSpectrum->SetBinContent(i,val[i-1]);
00844   }
00845   
00846   delete[] val;
00847 }


Member Data Documentation

TSpline5* NuMatrix1D::fCumulativeDist [protected]
double NuMatrix1D::fgCachedDmsq = -123456789 [static, protected]

The value of dmsq that fgCachedOsc was generated at.

Reimplemented in NuMatrixSpectrum.

Definition at line 192 of file NuMatrix1D.h.

Referenced by FillOscCache().

NuMatrix1D * NuMatrix1D::fgCachedOsc = 0 [static, protected]

The oscillation dip calculated at dmsq = fgCachedDmsq, sinsq = 1.

Reimplemented in NuMatrixSpectrum.

Definition at line 194 of file NuMatrix1D.h.

Referenced by FillOscCache(), and Oscillate().

TH1D* NuMatrix1D::fSpectrum [protected]

The storage for the 1D spectrum.

Definition at line 178 of file NuMatrix1D.h.

Referenced by Add(), NuMatrixSpectrum::BinWidthNormalize(), NuMatrixSpectrum::BlessedFD(), NuMatrixSpectrum::BlessedFD7e20(), NuMatrixSpectrum::BlessedFDRatio(), NuMatrixSpectrum::BlessedForNC(), NuMatrixSpectrum::BlessedND(), NuMatrixSpectrum::CalcPoisson(), NuMatrixSpectrum::CC2010CompressX(), NuMatrixSpectrum::CC2014CompressX(), NuMatrixSpectrum::CenterTitles(), Complete(), NuMatrixSpectrum::CompressTopBins(), NuMatrixSpectrum::CutLastBin(), DecayCC(), DecayMuToTau(), DecayNC(), Decohere(), NuMatrixSpectrum::Divide(), Divide(), NuMatrixSpectrum::DividePurity(), Draw(), NuMatrixSpectrum::DrawPoisson(), Fill(), Rock1DMatrix::Fill(), NuMatrixSpectrum::FillOscCache(), FillOscCache(), NuMatrixSpectrum::FillOscCache_3flavour(), NuMatrixSpectrum::GenerateErrors(), NuMatrixSpectrum::GenerateSystErrors(), NuMatrixSpectrum::GetEntries(), NuMatrixSpectrum::GetName(), NuMatrixSpectrum::GetNbinsX(), GetNuMatrixSpectrum(), NuMatrixSpectrum::GetXaxis(), NuMatrixSpectrum::GetYaxis(), NuMatrixSpectrum::Integral(), Interpolate(), InverseDecohere(), NuMatrixSpectrum::InverseOscillate(), NuMatrixSpectrum::InverseOscillateHybrid(), NuMatrixSpectrum::InverseOscillateLinearInterp(), InverseOscillateNSI(), NuMatrixSpectrum::InverseOscillateNSI(), NuMatrixSpectrum::InverseOscillatePreAveraged(), NuMatrixSpectrum::InverseOscillateSimpleAverage(), NuMatrixSpectrum::MakeJessTemplate(), NuMatrixSpectrum::Multiply(), Multiply(), NuMatrixSpectrum::MultiplyPurity(), NuMatrixSpectrum::NCCompressX(), NuMatrixSpectrum::NuBarCompressX(), NuMatrix1D(), operator=(), Oscillate(), NuMatrixSpectrum::Oscillate(), NuMatrixSpectrum::Oscillate_3flavour(), NuMatrixSpectrum::OscillateBasic_3flavour(), NuMatrixSpectrum::OscillateFourFlavour(), NuMatrixSpectrum::OscillateFourFlavourAverage(), NuMatrixSpectrum::OscillateHybrid(), NuMatrixSpectrum::OscillateHybrid_3flavour(), NuMatrixSpectrum::OscillateLED(), NuMatrixSpectrum::OscillateLinearInterp(), NuMatrixSpectrum::OscillateMu2Mu(), NuMatrixSpectrum::OscillateMu2MuBar(), NuMatrixSpectrum::OscillateMu2Tau(), NuMatrixSpectrum::OscillateMu2TauBar(), OscillateNSI(), NuMatrixSpectrum::OscillateNSI(), NuMatrixSpectrum::OscillatePreAveraged(), NuMatrixSpectrum::OscillateSimpleAverage(), Rock1DMatrix::PrepareBins(), Read(), Rebin(), NuMatrixSpectrum::RebinToArray(), NuMatrixSpectrum::RebinToScheme(), NuMatrixSpectrum::RemoveErrors(), NuMatrixSpectrum::ResetSpectrum(), Scale(), ScaleToPOT(), NuMatrixSpectrum::SetFillColor(), NuMatrixSpectrum::SetFillStyle(), NuMatrixSpectrum::SetLineColor(), NuMatrixSpectrum::SetLineStyle(), NuMatrixSpectrum::SetLineWidth(), NuMatrixSpectrum::SetMarkerColor(), NuMatrixSpectrum::SetMarkerSize(), NuMatrixSpectrum::SetMarkerStyle(), NuMatrixSpectrum::SetName(), NuMatrixSpectrum::SetRangeUser(), SetSpectrum(), NuMatrixSpectrum::SetTitle(), SetValue(), NuMatrixSpectrum::Spectrum(), Spectrum(), StatsLikelihood(), Subtract(), NuMatrixSpectrum::UnoscillateND(), NuMatrixSpectrum::UnoscillateNDAverage(), Write(), XToY(), YToX(), and ~NuMatrix1D().


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

Generated on 2 Nov 2017 for loon by  doxygen 1.6.1